实验6-非线性方程求解

1、实验6 非线性方程求解化工系 分0班 05 张亚清【实验目的】掌握用MATLAB软件求解非线性方程和方程组的基本用法，并对结果作初步分析；练习用非线性方程和方程组建立实际问题的模型并进行求解。【实验内容】题目1分别用fzero和fsolve程序求方程sinx-x22=0的所有根，准确到10-10，取不同的 初值计算，输出初值、根的近似值和迭代次数，分析不同根的收敛域；自己构造某个迭代公式（如x=2sinx12等）用迭代法求解，并自己编写牛顿法的程序进行求解和比较。【问题分析】首先做定性分析，用MATLAB做出y1=sinx，y2=x2/2的图像，研究零点的取值区间。程序如下：x=-3:5;y1

2、=sin(x);y2=x.2/2;plot(x,y1,x,y2)由图可见，原函数有两个零点，分别在,和1,2内。【问题解答】用MATLAB中fzero函数求解，程序如下：format long gopt=optimset(fzero);opt=optimset(opt,tolx,1e-10);x,fv,ef,out=fzero(inline(sin(x)-x2/2),opt)x,fv,ef,out=fzero(inline(sin(x)-x2/2),1,2,opt)运行后得到结果：x = fv = ef = 1out = intervaliterations: 0 iterations: 7

3、funcCount: 9 algorithm: bisection, interpolation message: Zero found in the interval , x = fv = ef = 1out = intervaliterations: 0 iterations: 7 funcCount: 9 algorithm: bisection, interpolation message: Zero found in the interval 1, 2由此可得原方程共有两个根，x1=0，x2=。下面取不同的 初值计算，输出初值、根的近似值和迭代次数，分析不同根的收敛域。对于x1得到如

4、下结果：用fzero求解初值根的近似值迭代次数-51213-41212-3109-287-1660007经试验，fzero函数x1的收敛域为,。用fsolve求解初值根的近似值迭代次数-587-476-365-255-1440006经试验，fsolve函数x1的收敛域为-,。两种方法中，随着初值逼近x1，迭代次数均递减，偶尔有例外情况，但总体趋势不变。fsolve的迭代次数基本都小于fzero的迭代次数。对于x2得到如下结果：用fzero求解初值根的近似值迭代次数162538412513614713812910109118128137146经试验，fzero函数x2的收敛域为，。迭代次数总体变

5、化趋势是先递增，后递减。这和sin函数的基本性质有关。用fsovle求解初值根的近似值迭代次数162435465667778898108118128138148经试验，fsolve函数x2的收敛域为,+。随着初值的增大，迭代次数逐渐增加，增加速度递减。用自己构造的迭代公式：x=2sinx12 求解。简单分析可知，要使等式右边有意义，则sin x0且x0。编写程序如下：x0=1;x(1)=x0;tol=1e-10;u=1;n=100;i=1;while(abs(u)tol)&(sin(x(i)eps) x(i+1)=(2*sin(x(i); u=x(i+1)-x(i); i=i+1; if(in

6、) error(n is full); endendxi-1输出结果如下：ans = 1 ans = 12迭代次数为12。对该方法迭代次数进行统计：用x=2sinx12求解初值根的近似值迭代次数11221231471381091413141411当初值取其他值时，x=2sinx12无意义。该迭代方法收敛次数多于fzero法和fsolve法。下面用牛顿法的程序进行求解和比较。程序如下：x0=1;x(1)=x0;tol=1e-6;u=1;n=10;i=1;while(abs(u)tol) x(i+1)=x(i)-(sin(x(i)-x(i)2/2)/(cos(x(i)-x(i); u=x(i+1)

7、-x(i); i=i+1; if(in) error(n is full); endendx输出结果如下：ans = 1 迭代次数为6。对牛顿法迭代次数进行统计：用牛顿法求解初值根的近似值迭代次数162637475767788898108118128138149牛顿法迭代次数随初值渐渐远离根的真实值而逐渐增大，增大速度递减。【结果分析与讨论】以上四种迭代方法（fzero、fsolve、x=2sinx12）、牛顿），均可得到根x2的近似值。其中，牛顿法和自己构造的迭代公式无法得到根x1的值。四种迭代方法中fsolve和牛顿法迭代次数较少，收敛较快。题目6给定4种物质对应的参数ai,bi,ci和交

8、互作用矩阵Q如下：a1=a2=a3=a4=b1=b2=b3=b4=c1=c2=c3=c4=Q = 在压强p=760mmHg下，为了形成均相共沸混合物，温度和组分分别是多少请尽量找出所有可能的解。【问题分析】【问题解答】用MATLAB编程，程序如下：function f=azeofun(XT,n,P,a,b,c,Q)x(n)=1;for i=1:n-1 x(i)=XT(i); x(n)=x(n)-x(i);endT=XT(n);p=log(P);for i=1:n d(i)=x*Q(i,1:n); dd(i)=x(i)/d(i);endfor i=1:n f(i)=x(i)*(b(i)/(T+c

9、(i)+log(x*Q(i,1:n)+dd*Q(1:n,i)-a(i)-1+p);endendn=4;P=760;a=,;b=,;c=,;Q= ; ; ; ;XT0=,50;XT,Y=fsolve(azeofun,XT0,n,P,a,b,c,Q)得到结果：XT = Y = * 即四种物质组成均相共沸混合物时的比例分别为%，%，%，%，温度为。【结果分析与讨论】在上面的计算中，我们对初值XT0的取法是：四种物质各占约1/4，温度为50。如果取其他初值，还可以得到其他的均相共沸混合物，结果归纳如下表：初值解XT0 x1x2x3x4T,50,0,80 ,0,80 题目7用迭代公式 计算序列 ，分析其

10、收敛性，其中a分别取5,11,15；b（0）任意，初值 =1。观察有无混沌现象出现，并找出前几个分岔点，观察分岔点的极限趋势是否符合Feigenbaum常数揭示的规律。【问题分析】本题主要是观察及判断分岔、混沌现象。通过取三个不同的a值（b保持不变），计算出迭代50次的结果，并作出各自的图形，观察各自的分岔现象，相互之间比较可推测a变化时的整体走势；此外，还可以利用matlab作出该迭代函数的迭代序列随着参数发生变化的收敛、分岔、混沌现象图，便于观察整体走势。【问题解答】MATLAB程序如下：x0=1;x(1)=x0;n=50;i=1;a=5;b=5;for i=1:n; x(i+1)=a*x

11、(i)*exp(-b*x(i);endt=1:n+1;plot(t,x)x分别令a=5,11,15，b=5得到如下结果：a=5a=11a=15迭代次数迭代值迭代次数迭代值迭代次数迭代值1112223334445556667778889991010101111111212121313131414141515151616161717171818181919192020202121212222222323232424242525252626262727272828282929293030303131313232323333333434343535353636363737373838383939394

12、04040414141424242434343444444454545464646474747484848494949505050可作出分岔和混沌图，如下图所示：由图可以看出有混沌现象出现。下面寻找前几个分叉点：求第一个分岔点：平衡点：x=x*=ln(a)运行程序：y=diff(a*x*exp(-x) %求一阶导数表达式得到结果：y =a/exp(x) - (a*x)/exp(x)然后运行程序：syms ay1=solve(a/exp(log(a) - (a*log(a)/exp(log(a)-1)y2=solve(a/exp(log(a) - (a*log(a)/exp(log(a)+1)

13、%计算第一个分岔点时的a值得结果：y=1 y2=exp(2)= 计算第二个分岔点：function y=fenchadian2(x,a)y(1)=x(1)*x(2)*exp(-x(2)-x(3);y(2)=x(1)*x(3)*exp(-x(3)-x(2);y(3)= (x(1)/exp(x(2)-(x(1)*x(2)/exp(x(2)*(x(1)/exp(x(3)-(x(1)*x(3)/exp(x(3)-a;x0=12,1,2;x,fv,ef,out,jac=fsolve(fenchadian2,x0,-1) %分析易知a=-1而非a=1运行结果：x = 计算第三个分岔点：function y

14、=fenchadian3(x,a)y(1)=x(1)*x(2)*exp(-x(2)-x(3);y(2)=x(1)*x(3)*exp(-x(3)-x(4);y(3)=x(1)*x(4)*exp(-x(4)-x(5);y(4)=x(1)*x(5)*exp(-x(5)-x(2);y(5)=(x(1)/exp(x(2)-(x(1)*x(2)/exp(x(2)*(x(1)/exp(x(3)-(x(1)*x(3)/exp(x(3)*(x(1)/exp(x(4)-(x(1)*x(4)/exp(x(4)*(x(1)/exp(x(5)-(x(1)*x(5)/exp(x(5)-a;x0=14,1,2,3,4;x,fv,ef,out,jac=fsolve(fenchadian3,x0,-1)运行结果：x = 计算第四个分岔点（程序类似，在此省略）：得结果：x = 【结果分析与讨论】a15时，便明显产生混沌现象，x的排列变得无明显规律。下面讨论其分岔点是否符合Feigenbaum常数规律：第一个分岔点a=第二个分岔点a=第三个分岔点a=第四个分岔点a=再次计算分岔点时发现不再有好的收敛性。又有： / 由两组数据来看，前者还与Feigenbaum常数相差较远，第二组数据便与其极限值大大靠近，可以推测，继续运算下去，会得到更相近的数值因此符