matlabfindroot,求助:利⽤FindRoot求解⾃定义函数的零点该楼层疑似违规已被系统折叠 隐藏此楼查看此楼
我想写⼀个⾃定义函数(⽐如f(x)),然后⽤FindRoot来求这个函数的零点,但是我想不到如何实现我想要的效果。具体如下。
第⼀个问题。我想要让FindRoot在每⼀步迭代都print出⽬前的x和⽬前函数f(x)的值,这样我能清晰地看出寻根的进程究竟如何,不然的话我只能⼀直等待漫长的求解过程⽽且还不知道是不是求解遇到了什么困难。
我看到FindRoot有个选项EvaluationMonitor有类似的效果,但这个选项只能在⽅程求解完毕后才输出x的变化,所以应该满⾜不了我的需求。
同样的问题放在Matlab上,就很容易解决。只需要写⼀个函数⽂件,在这个函数⽂件内部加上display(x)和display(f(x))就可以了。但是我不清楚MMA如何实现这个功能。因为MMA给我的印象⼀直都是函数式编程,我不太清楚在MMA中如何写出灵活的程序实现我的需求。
第⼆个问题。我先贴上我的代码:
Clear["Global`*"];
Clear[Derivative];
$MaxPrecision=Infinity;
$MinPrecision=25;
$Assumptions=L>0&&L<1;
\[Alpha]=1`25;
J[L_]:={{0,-L-I Sqrt[1-L^2],0},{I/2,-\[Alpha],-(I/2)},{0,-L+I Sqrt[1-L^2],0}};
zStar[L_]:=I L-Sqrt[1-L^2];
zbarStar[L_]:=-I L-Sqrt[1-L^2];
T[L_]:=Module[{T0=Chop@Transpose[Eigensystem[J[L]][[2]]]},Dot[T0,{{1,0,0},{0,Exp[I x]/.Quiet@Solve[{Exp[I
x]T0[[3,2]]==Exp[-I x]Conjugate[T0[[1,2]]]},x][[1]],0},{0,0,1}}]];
f[z_,y_,zbar_,L_]:={I z y,L-(z-zbar)/(2I)-\[Alpha] y,-I zbar y}
\[Alpha]1=(-
0.0569346681735483100196233345904580474`25.+0.03380509290798852465990285780535192822`25. I);
\[Beta]1=1;
\[Alpha]2=1;matlab定义函数表达式
\[Beta]2=(-0.56217287423309301896928923557469160676`25.-
0.15431958634948985984144944464738474976`25. I);
\[Alpha]3=(1.12190799111714769562606230101308634098`25.-
0.6661355039134196995510221227089331722`25. I);
\[Beta]3=0.7671179743387728812498376749441521484`25.;
u0[A_,B_]:=((\[Alpha]1 A^2+\[Beta]1 B)/2-(\[Alpha]1 A^2-\[Beta]1 B)/2 (2k-1))k^2 (1-k);
v0[A_,B_]:=((\[Alpha]2 A^2+\[Beta]2 B)/2-(\[Alpha]2 A^2-\[Beta]2 B)/2 (2k-1))k (1-k)^2;
w0[A_,B_]:=((\[Alpha]3 A^2+\[Beta]3 B)/2-(\[Alpha]3 A^2-\[Beta]3 B)/2 (2k-1))k^2 (1-k)^2;
dkdt=k(1-k);
g[u_,v_,w_,L_]:=Dot[Inverse[T[L]],f[({zStar[L],0,zbarStar[L]}+Dot[T[L],{u,v,w}])[[1]],({zStar[L],0,zbarStar[L]}+Dot[T[L],{u,v,w}]) [[2]],({zStar[L],0,zbarStar[L]}+Dot[T[L],{u,v,w}])[[3]],L]]
targetN[L_,A_,B_]:={Chop@Expand@Integrate[Chop@Apart@(g[u0[A,B],v0[A,B],w0[A,B],L][[1]]/dkdt),
{k,0,1}],Chop@Expand@Integrate[Chop@Apart@(g[u0[A,B],v0[A,B],w0[A,B],L][[2]]/dkdt),
{k,0,1}],Chop@Expand@Integrate[Chop@Apart@(g[u0[A,B],v0[A,B],w0[A,B],L][[3]]/dkdt),{k,0,1}]};
FindRoot[{targetN[L,A,B]=={0,0,0},L>0&&L<1,A>0,B>0},{{L,0.9},{A,5},{B,5}}]
这段代码⽐较繁杂,不过关键的问题出在最后的
FindRoot[{targetN[L,A,B]=={0,0,0},L>0&&L<1,A>0,B>0},{{L,0.9},{A,5},{B,5}}]
这串代码。
⼤致就是我先定义了⼀些变量和函数,然后利⽤这些来定义我最终需要求解的函数targetN. L的具体取值会影响计算中⼀些运算的结果(⽐如EigenSystem),所以我不能写出targetN的显式表达式(⽽且显式表达式也极端复杂),只好每给定⼀个L(以及A, B),计算targetN[L,A,B]的值,以此来求targetN的零点。但是如果我运⾏
FindRoot[{targetN[L,A,B]=={0,0,0},L>0&&L<1,A>0,B>0},{{L,0.9},{A,5},{B,5}}]
的话,会报错:
Integral of (b^2 <<4>>)/((1.0-1.000000000000000000000000 k)^2 (-1.000000000000000000000000+L^2)) does not converge on {0,1}.
也就是说targetN[L,A,B]先把L, A, B当作没有赋予具体数值的符号来计算targetN的值,然后报错说积分不收敛,这恰恰是我想要避免的。我想要的是每代⼊⼀个具体的L, A, B的值都能给出相应的targetN[L,A,B]的值,利⽤⽜顿法⼀步步迭代下去,这样就能避免符号计算带来的⿇烦,⾄少Matlab中是这么⼀个过程。
所以第⼆个问题⼤概就是这个。
版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。
发表评论