搜索
您的当前位置:首页正文

fsolve传递未知参量解方程

来源:哗拓教育
fsolve函数解方程

[X,FVAL,EXITFLAG,OUTPUT,JACOB]=FSOLVE(FUN,X0,...) returns the Jacobian of FUN at X.

Examples

FUN can be specified using @:

x = fsolve(@myfun,[2 3 4],optimset('Display','iter'))

where myfun is a MATLAB function such as:

function F = myfun(x)

F = sin(x);

FUN can also be an anonymous function:

x = fsolve(@(x) sin(3*x),[1 4],optimset('Display','off'))

If FUN is parameterized, you can use anonymous functions to capture the problem-dependent parameters. Suppose you want to solve the system of nonlinear equations given in the function myfun, which is parameterized by its second argument c. Here myfun is an M-file function such as

function F = myfun(x,c)

F = [ 2*x(1) - x(2) - exp(c*x(1))

-x(1) + 2*x(2) - exp(c*x(2))];

To solve the system of equations for a specific value of c, first assign the value to c. Then create a one-argument anonymous function that captures

that value of c and calls myfun with two arguments. Finally, pass this anonymous function to FSOLVE:

c = -1; % define parameter first x = fsolve(@(x) myfun(x,c),[-5;-5])

以matlab R2008a版本为例,各版本出错提示可能有所不同。有不对之处,欢迎指正。 1.solve 和 fsolve的基本含义

matlab给出的关于solve 和 fsolve的基本描述为: solve——Symbolic solution of algebraic equations

fsolve——Solve system of nonlinear equations

可见solve用于解决代数方程(组)的符号(解析)解,而fsolve用来解决非线性方程(组)的数值解。

【在matlab里面solve命令主要是用来求解代数方程(即多项式)的解,但是也不是说其它方程一个也不能解,不过求解非代数方程的能力相当有限,通常只能给出很特殊的实数解。从计算机的编程实现角度讲,如今的任何算法都无法准确的给出任意非代数方程的所有解,但是我们有很多成熟的算法来实现求解在某点附近的解。matlab也不例外,它也只能给出任意非代数方程在某点附近的解,函数有两个:fzero和fsolve,具体用法请用help或doc命令查询吧。如果还是不行,你还可以将问题转化为非线性最优化问题,求解非线性最优化问题的最优解,可以用的命令有:fminbnd, fminsearch, fmincon等等。】(引自:http://blog.sina.com.cn/s/blog_4c4af5c101008w9f.html,作者:ggbondg)

下面举几个例子:

     

例1:>> solve('a*x-1') ans = 1/a

例2:>> solve('exp(x)+sin(x)-2') ans =

.44867191635127271149118657202662

注:对于solve结果的显示,有时看起来比较长,可用vpa进行精度控制,如: >> vpa(solve('exp(x)+sin(x)-2'),3) ans = .449

 

例3:>> fsolve(@(x)exp(x)+sin(x)-2,0) Optimization

terminated:

first-order

optimality

is

less

than

options.TolFun.

ans =

 0.4487

2.关于solve 和 fsolve求解方程组时的书写规则

对于solve,方程可以直接书写,不需要运算符”.”;

对于fsolve,当未知量与未知量有乘除操作或未知量有开方、幂等操作时运算符”.”可写

也可不写(记得好像必须写,试了试,发现不写也行)。下面举几个例子:

 例4:>> solve('x+y.^2-1','x.^2-y-3')  ??? Error using ==> solve at 77

 ' x+y.^2-1 ' is not a valid expression or equation.  例5:>> solve('x+y^2-1','x^2-y-3')  ans =

 x: [4x1 sym]  y: [4x1 sym]

 例6:function shiyan  clc  clear  x0=[0,0];  fsolve(@mf,x0) 

 function F=mf(x)  F=[x(1)+x(2)^2-1;  x(1)^2-x(2)-3];  %%%%%Result%%%%%%

 Optimizer appears to be converging to a minimum that is not a root:  Sum of squares of the function values is > sqrt(options.TolFun).  Try again with a new starting point.  ans =

 1.6268 -0.1537

例7:把例6中的mf函数,换成如下再试试。

function F=mf(x)

F=[x(1)+x(2).^2-1; x(1).^2-x(2)-3];

例8:把例6的初值x0设为x0=[-2,2];运行结果为:

Optimization terminated: first-order optimality is less than options.TolFun. ans = -2.1875 1.7854

可见,用fsolve解非线性方程组,比较依赖处置的选择,因此建议用fsolve解方程时,能大致了解问题的求解区间,以便选择合适的初值。

3.关于solve 和 fsolve求解时,参数为多数值的求解 问题来源:

http://www.chinavib.com/forum/viewthread.php?tid=82658&page=1#pid426339 类似问题描述:k=(5.0e+4):(1e+3):(6e+4);h=1.6e-6;n1=2.2899;n0=1.5040;n2=1.000;解方程组:

p1=sqrt(k.^2.*n1.^2-b.^2;p2=sqrt(b.^2-k.^2.*n2.^2);p0=sqrt(b.^2-k.^2.*n0.^2); p1*h-pi-atan(p0./p1)-atan(p2./p1)=0。(见:http://www.chinavib.com/forum/thread-83102-1-1.html) 解决方法:这是非线性方程组,不过我们可以先用solve试试:

 例9:问题如前所述。考虑用solve是否能解。  clear  clc

 k=(5.0e+4):(1e+3):(6e+4);  h=1.6e-6;  n1=2.2899;  n0=1.5040;  n2=1.000;

 for i=1:length(k)

 y=solve(['p1=sqrt(',num2str(k(i)),'^2*',num2str(n1),'^2-b^2)'],...  ['p2=sqrt(b^2-',num2str(k(i)),'^2*',num2str(n2)','^2)'],...  ['p0=sqrt(b^2-',num2str(k(i)),'^2.*',num2str(n0),'^2)'],...  ['p1*',num2str(h),'-',num2str(pi),'-atan(p0/p1)-atan(p2/p1)=0']);  end

%%%%%Result%%%%% > In solve at 140

In shiyan at 9

y=solve(['p1=sqrt(',num2str(k(i)),'^2*',num2str(n1),'^2-b^2)'],...

['p2=sqrt(b^2-',num2str(k(i)),'^2*',num2str(n2)','^2)'],...

['p0=sqrt(b^2-',num2str(k(i)),'^2.*',num2str(n0),'^2)'],...

['p1*',num2str(h),'-',num2str(pi),'-atan(p0/p1)-atan(p2/p1)=0']); Warning: Warning, solutions may have been lost

Warning: Explicit solution could not be found. >> whos y

Name Size Bytes Class Attributes y 0x0 64 sym

由此说明利用solve并不能解决这个复杂的非线性方程组,考虑数值解法。

 例10:问题如例9,考虑fsolve求解问题。  % 主程序  clear  clc

 global k h n1 n2 n0

 k1=(5.0e+4):(1e+3):(6e+4);  h=1.6e-6;  n1=2.2899;  n0=1.5040;  n2=1.000;  x0=[0 1 0 0];

 b=zeros(1,length(k1));  for i=1:length(k1)  k=k1(i);

 y=fsolve(@myfun,x0);  b(i)=y(4);  end

 plot(k1,b);  % 子程序

 function F=myfun(x)  global k h n1 n2 n0

 % p0 p1 p2 b----x(1) x(2) x(3) x(4)  F=[sqrt(k^2*n1^2-x(4)^2)-x(2);  sqrt(x(4)^2-k^2*n2^2)-x(3);  sqrt(x(4)^2-k^2*n0^2)-x(1);

 x(2)*h-pi-atan(x(1)./x(2))-atan(x(3)./x(2))];  %%%%% Result%%%%%%%%%  y=fsolve(@myfun,x0);

 Optimizer appears to be converging to a minimum that is not a root:  Sum of squares of the function values is > sqrt(options.TolFun).  Try again with a new starting point.

这说明所选初值并不合理。

 例11 例10还可以这样改写:  % 主程序  clear  clc

 % global k h n1 n2 n0  k1=(5.0e+4):(1e+3):(6e+4);  h=1.6e-6;  n1=2.2899;  n0=1.5040;  n2=1.000;  x0=[0 1 0 0];

 b=zeros(1,length(k1));  for i=1:length(k1)  k=k1(i);

 y=fsolve(@(x)myfun(x,k, h, n1, n2, n0),x0);  b(i)=y(4);  end

 plot(k1,b);  % 子程序

 function F=myfun(x,k, h, n1, n2, n0)  % global k h n1 n2 n0

 % p0 p1 p2 b----x(1) x(2) x(3) x(4)  F=[sqrt(k^2*n1^2-x(4)^2)-x(2);  sqrt(x(4)^2-k^2*n2^2)-x(3);  sqrt(x(4)^2-k^2*n0^2)-x(1);

 x(2)*h-pi-atan(x(1)./x(2))-atan(x(3)./x(2))];

例12参见

http://www.chinavib.com/forum/viewthread.php?tid=82658&page=1#pid426339此贴,看看大侠ChaChing,dingd(1stOpt解法)等的解法。

注:利用fsolve解数值解,初值的选择十分重要。而1stOpt则对初值的选择要求比较低,不妨一试。关于这一点,请参考如下文章:

作者:dingd非线性方程组-1stOpt与fsolve的比较: http://www.chinavib.com/forum/thread-48384-1-5.html

solve,fsolve的用法

2007-06-03 23:40:42| 分类: tech | 标签:tech:技术类 |字号大中小 订阅

solve是方程,方程组的符号解法; fsolve是数值的优化方法; 两种方法各有所长吧。

第一种,幸运的话,可以得到解析解,就是那种符号解; 但是复杂的方程,往往是得不到的。

第二种的话,不出差错的话,总是可以得到一些可用的数值解; 可不要忽略了第二种哦; 举例如下,syms y;

y=solve('(5-y^2)^0.5*(10-y^2)^0.5*besselj(1,(10-y^2)^0.5)*besselj(0,(5-y^2)^0.5)+y^2*besselj(0,(10-y^2)^0.5)*besselj(1,(5-y^2)^0.5)') 是得不到隐式解的; 第二种 定义函数如下 function f=mytest(x)

f=(5-x^2)^0.5*(10-x^2)^0.5*besselj(1,(10-x^2)^0.5)*besselj(0,(5-x^2)^0.5)+x^2*besselj(0,(10-x^2)^0.5)*besselj(1,(5-x^2)^0.5) 使用

x=fsolve(mytest,1+1i) 或者

x=fsolve(@mytest,1+1i)

都是可以的;

结果如下 f =

-0.0136 + 0.1593i f =

-0.0136 + 0.1593i f =

0.0235 + 0.0353i f =

0.0235 + 0.0353i f =

0.0034 - 0.0046i f =

0.0034 - 0.0046i f =

3.7005e-005 +8.2011e-005i f =

3.6997e-005 +8.2015e-005i f =

1.2868e-008 -1.7310e-008i f =

4.4882e-009 -1.3191e-008i

Optimization terminated: first-order optimality is less than options.TolFun. x =

0.7521 + 1.1982i

======================================================================================== 初值的选择是很重要的 如果用

x=fsolve(@mytest,1)

这里用的实数初解,结果很差的 x =

9.6085e-006

============================================= 如果函数改为范数 abs的话,结果也是比较差的 function f=mytest(x)

f=abs((5-x^2)^0.5*(10-x^2)^0.5*besselj(1,(10-x^2)^0.5)*besselj(0,(5-x^2)^0.5)+x^2*besselj(0,(10-x^2)^0.5)*besselj(1,(5-x^2)^0.5)) 你可以试试,我是试过了的。

========================================= 更好的初值选择

x=fsolve(@mytest,2+2i) 结果f =

-8.3966e-009 +4.1120e-009i

Optimization terminated: first-order optimality is less than options.TolFun. x =

0.7521 + 1.1982i

因篇幅问题不能全部显示,请点此查看更多更全内容

Top