3

我正在尝试编写一小段代码来执行错误传播。到目前为止,我可以让 Mathematica 在函数 f(x1,x2,...,xi,...,xn) 中生成错误 delta_f 的公式,错误为 dx1,dx2,...,dxi,.. .dxn:

fError[f_, xi__, dxi__] := 
  Sum[(D[f[xi], xi[[i]]]*dxi[[i]])^2, {i, 1, Length[xi]}]^(1/2)

其中 fError 要求输入函数 f 的所有变量都被 {...} 包围。例如,

d[{mv_, Mv_, Av_}] := 10^(1/5 (mv - Mv + 5 - Av))
FullSimplify[fError[d, {mv, Mv, Av}, {dmv, dMv, dAv}]]

返回

2 Sqrt[10^(-(2/5) (Av - mv + Mv)) (dAv^2 + dmv^2 + dMv^2)] Log[10]

我的问题是,我该如何评估这个?理想情况下,我想将 fError 修改为:

fError[f_, xi__, nxi__, dxi__]

其中 nxi 是 xi 的实际值列表(分开,因为将 xi 设置为它们的数值会破坏上面的微分步骤。)这个函数应该找到误差 delta_f 的一般公式,然后在可能的情况下对其进行数值评估。我认为解决方案应该像 Hold[] 或 With[] 或类似的东西一样简单,但我似乎无法理解。

4

3 回答 3

2

我没有关注你所做的一切,因为这是两年前发布的,你很可能不再在做这件事了。我会给你我的错误传播解决方案,希望它能以某种方式帮助你或其他人。

我试图在下面链接的视频和文件中包含最好的文档。如果您打开 .cdf 文件并清除它,您应该能够看到我的代码...

文件: https ://drive.google.com/file/d/0BzKVw6gFYxk_YUk4a25ZRFpKaU0/view?pli=1

视频教程: https ://www.youtube.com/watch?v=q1aM_oSIN7w

-布赖恩

编辑: 我发布了链接,因为我无法附加文件,也不想为数学新手发布没有文档的代码。这是直接的代码。我鼓励任何发现此解决方案有助于快速查看文档的人,因为它展示了一些提高生产力的技巧。

Manipulate[
 varlist = ToExpression[variables];
 funct = ToExpression[function];
 errorFunction[variables, function]
 , {variables, "{M,m}"}, {function, "g*(M-m)/(M+m)"}, 
 DisplayAllSteps -> True, LabelStyle -> {FontSize -> 17}, 
 AutoAction -> False,
 Initialization :> (
   errorFunction[v_, f_] := (
     varlist = ToExpression[v];
     funct = ToExpression[f];
     varlength = Length[Variables[varlist]];
     theoretical = 
      Sqrt[(Total[
         Table[(D[funct, Part[varlist, n]]*
             Subscript[U, Part[varlist, n]])^2, {n, 1, 
           varlength}]])];
     Part[theoretical, 1];
     varlist;
     uncert = Table[Subscript[U, Part[varlist, n]], {n, 1, varlength}];
     uncert = DeleteCases[uncert, Alternatives @@ {0}];
     theoretical = Simplify[theoretical];
     Column[{Row[{Grid[{
           {"Variables", varlist},
           {"Uncertainties", uncert},
           {"Function", function},
           {"Uncertainty Function", theoretical}}, Alignment -> Left, 
          Spacings -> {2, 1}, Frame -> All, 
          ItemStyle -> {"Text", FontSize -> 20}, 
          Background -> {{LightGray, None}}]}],
       Row[{
         Grid[{{"Brian Gennow  March/24/2015"}}, Alignment -> Left, 
          Spacings -> {2, 1}, ItemStyle -> "Text", 
          Background -> {{None}}]
         }]}]))]
于 2015-09-14T21:27:42.827 回答
1

这个问题是 5 年前发布的,但我最近遇到了同样的问题,并认为我会分享我的解决方案(针对不相关的错误)。

我定义了一个带有errorProp两个参数的函数,func并且vars. errorProp, ,的第一个参数func是表达式的符号形式,您希望为其计算由于其参数错误而导致的其值的误差。for 的第二个参数errorProp应该是一个表单列表

{{x1,x1 value, dx1, dx1 value},{x2,x2 value, dx2, dx2 value}, ... , 
{xn,xn value, dxn, dxn value}}

其中xi's 和dxi's 是变量及其误差的符号表示,而xi valueanddxi value是变量的数值及其不确定性(参见下面的示例)。

该函数errorProp返回误差的符号形式、输入函数的值以及根据 中的输入计算得出func的误差值。这是代码:funcvars

ClearAll[errorProp];
errorProp[func_, vars_] := Module[{derivs=Table[0,{Length[vars]}], 
funcErrorForm,funcEval,funcErrorEval,rplcVals,rplcErrors},

For[ii = 1, ii <= Length[vars], ii++,
    derivs[[ii]] = D[func, vars[[ii, 1]]];
];

funcErrorForm = Sqrt[Sum[(derivs[[ii]]*vars[[ii, 3]])^2,{ii,Length[vars]}]];

SetAttributes[rplcVals, Listable];
rplcVals = Table[Evaluate[vars[[ii, 1]]] :> Evaluate[vars[[ii, 2]]], {ii, 
Length[vars]}];

SetAttributes[rplcErrors, Listable];
rplcErrors = Table[Evaluate[vars[[ii, 3]]] :> Evaluate[vars[[ii, 4]]], {ii, 
 Length[vars]}];

funcEval = func /. rplcVals;
funcErrorEval = funcErrorForm /. rplcVals /. rplcErrors;
Return[{funcErrorForm, funcEval, funcErrorEval}];
];

errorProp在这里,我展示了一个具有两个变量的相当复杂的函数的实际示例:

ClearAll[test];
test = Exp[Sqrt[1/y] - x/y];
errorProp[test, {{x, 0.3, dx, 0.005}, {y, 0.9, dy, 0.1}}]

返回

 {Sqrt[dy^2 E^(2 Sqrt[1/y] - (2 x)/y) (-(1/2) (1/y)^(3/2) + x/y^2)^2 + (
 dx^2 E^(2 Sqrt[1/y] - (2 x)/y))/y^2], 2.05599, 0.0457029}

使用错误传播公式计算返回相同的结果:

{Sqrt[(D[test, x]*dx)^2 + (D[test, y]*dy)^2], 
test /. {x :> 0.3, dx :> 0.005, y :> 0.9, dy :> 0.1}, 
Sqrt[(D[test, x]*dx)^2 + (D[test, y]*dy)^2] /. {x :> 0.3, 
dx :> 0.005, y :> 0.9, dy :> 0.1}}

返回

{Sqrt[dy^2 E^(
2 Sqrt[1/y] - (2 x)/y) (-(1/2) (1/y)^(3/2) + x/y^2)^2 + (
dx^2 E^(2 Sqrt[1/y] - (2 x)/y))/y^2], 2.05599, 0.0457029}
于 2018-04-17T01:47:07.080 回答
1

Mathematica 12 引入了使用微分方法处理错误传播的Around函数。

因此,虽然不完全符合问题所需的格式,但可能会出现这样的情况:

expression = a^2*b;
expression /. {a -> Around[aval, da], b -> Around[bval, db]}

输出:

aval^2 bval ± Sqrt[aval^4 db^2+4 bval^2 Abs[aval da]^2]

除了aval、bval、da、db,您也可以使用数值。

于 2019-04-18T02:30:10.253 回答