下载文档原格式
⾮线性⽅程求解基于MATLAB的⾮线性⽅程的五种解法探讨摘要：本⽂利⽤matlab软件对⾮线性⽅程解法中的⼆分法、简单迭代法、⽜顿法、割线法以及Steffensen法的数值分析⽅法的算法原理及实现⽅法进⾏了探讨。
对f x x x=+-()2ln2的零点问题，分别运⽤以上五种不同的⽅法进⾏数值实验，⽐较⼏种解法的优缺点并进⾏初步分析评价。
关键词：⼆分法、简单迭代法、⽜顿法、割线法、Steffensen法1、引⾔在很多实际问题中，经常需要求⾮线性⽅程f(x) =0的根。
⽅程f(x) =0的根叫做函数f(x)的零点。
由连续函数的特性知：若f(x)在闭区间[a,b ]上连续,且()()0f a f b<.则f(x) =0在开区间(a,b)内⾄少有⼀个实根。
这时称[a,b]为⽅程f(x) =0的根的存在区间。
本⽂主要对⾮线性⽅程的数值解法进⾏分析，并介绍了⾮线性⽅程数值解法的五种⽅法。
并设=+-.f x x x()2ln2f x在[1,2]上的图形，如图1：. 显然，函数在[1,2]之间有⼀个零点。
⾸先画出()2、计算机配置操作系统Windows 7 旗舰版内存2GB处理器AMD 4核 A6-3400M APU 1.4GHz图.13、⼆分法⼆分法的基本思想是将⽅程根的区间平分为两个⼩区间,把有根的⼩区间再平分为两个更⼩的区间,进⼀步考察根在哪个更⼩的区间内。
如此继续下去,直到求出满⾜精度要求的近似值。
设函数()f x 在区间[a,b ]上连续,且f(a)·f(b) <0，则[a,b ]是⽅程f(x) =0的根的存在区间，设其内有⼀实根，记为x*。
取区间[a,b ]的中点()2k a b x +=并计算1()f x ,则必有下列三种情况之⼀成⽴： (1) 1()f x =0，x1就是⽅程的根x*;(2)()f a .1()f x <0,⽅程的根x*位于区间[a, 1x ]之中,此时令111,a a b x ==; (3)1()f x .()f b <0,⽅程的根x*位于区间[1x ,b ]之中,此时令11a x =,1b b =。
数学实验报告非线性方程求解一、实验目的1.掌握用 MATLAB 软件求解非线性方程和方程组的基本用法，并对结果作初步分析；2.练习用非线性方程和方程组建立实际问题的模型并进行求解。
二、实验内容题目1【问题描述】(Q1)小张夫妇以按揭方式贷款买了1套价值20万元的房子，首付了5万元，每月还款1000元，15年还清。
问贷款利率是多少？(Q2)某人欲贷款50 万元购房，他咨询了两家银行，第一家银行开出的条件是每月还4500元，15 年还清；第二家银行开出的条件是每年还45000 元，20 年还清。
从利率方面看，哪家银行较优惠（简单假设：年利率=月利率×12）？【分析与解】假设初始贷款金额为x0，贷款利率为p，每月还款金额为x，第i个月还完当月贷款后所欠银行的金额为x i，(i=1,2,3,......,n)。
由题意可知：x1=x0(1+p)−xx2=x0(1+p)2−x(1+p)−xx3=x0(1+p)3−x(1+p)2−x(1+p)−x……x n=x0(1+p)n−x(1+p)n−1−⋯−x(1+p)−x=x0(1+p)n−x (1+p)n−1p=0因而有：x0(1+p)n=x (1+p)n−1p (1)则可以根据上述方程描述的函数关系求解相应的变量。
(Q1)根据公式(1)，可以得到以下方程：150p(1+p)180−(1+p)180+1=0设 f(p)=150p(1+p)180−(1+p)180+1，通过计算机程序绘制f(p)的图像以判断解p的大致区间，在Matlab中编程如下：for i = 1:25t = 0.0001*i;p(i) = t;f(i) = 150*t*(1+t).^180-(1+t).^180+1;end;plot(p,f),hold on,grid on;运行以上代码得到如下图像：f(p)~p关系曲线图通过观察上图可知p∈[0.002,0.0022]。
Solution1：对于p∈[0.002,0.0022]，采用二分法求解，在Matlab 中编程如下：clear;clc;x0=150000;n=180;x=1000;p0=0.002;p1=0.0022;while (abs(p1-p0)>1e-8)f0=x0*(1+p0).^n+x*(1-(1+p0).^n)/p0;f1=x0*(1+p1).^n+x*(1-(1+p1).^n)/p1;p2=(p0+p1)/2;f2=x0*(1+p2).^n+x*(1-(1+p2).^n)/p2;if (f0*f2>0 && f1*f2<0)p0=p2;elsep1=p2;end;end;p0结果得到p0=0.00208116455078125=0.2081%.所以贷款利率是0.2081%。
实验7非线性方程求解实验目的：1． 1． 掌握用MA TLAB 软件求解非线性方程和方程组的基本用法，并对结果作出初步的分析2． 2． 练习用非线性方程组建立实际问题的模型并进行求解 实验内容：6-3．（1）小张夫妇以按揭方式贷款买了一套价值20万的房子，首付了5万元，没有还款1000元，15年还清。
问贷款利率是多少？（2）某人与贷款50万元购房，他咨询了两家银行，第一家银行开出的条件是每月还4500元，15年还清；第二家银行开出的条件是每年还45000元，20年还清。
从利率方面看，那家银行较优惠（简单的假设年利率=月利率*12）？解：（1）问题的建模：假设需付的款为0a ，每月付款b 元，经n 年后还清，付款利率为r ，于是对与按揭付款的方式不难列出以下的方程组：'00a a ='10211*(1)*(1)*(1)n n a a r b a a r b a a r b -=+-=+-=+-M于是有：223123'*122*1210*122*121*120*121*120(1)*(1)(1)*(1)*(1)(1)(1(1)(1)(1))(1)(1(1)(1)(1)(1))(1)1(1)*n n n n n n n n n n a a r b b r a r b b r b r a r b r r r a r b r r r r r a r b r----+=+--+=+--+-+=+-+++++++=+-++++++++++-=+-ML L如果对于年利率则不用考虑第一次要付的钱，所以易得：0(1)1(1)*nnn r a a r b r +-=+-按照上述的关系得到了本题的m 文件：function y=rate(x,a0,n,b,opt)if opt>=0%如果是〉0则按照月利率算y=a0*(1+x)^(n*12)-b*((1+x)^(n*12)-1)/x;%计算月利率， elsey=a0*(1+x)^n-b*((1+x)^n-1)/x;%计算年利率 end<rate.m> 主文件：clear all ;a0=150;b=1;n=15; a1=500;b1=4.5;n1=15; a2=500;b2=45;n2=20; x0=1;[x1,fv1,ef1,out1]=fzero(@rate,x0,[],a0,n,b,1);%对第一题小张夫妇买方月利率问题进行求解[x2,fv2,ef2,out2]=fzero(@rate,x0,[],a1,n1,b1,1);%对二问第一家银行的利率求解 [x3,fv3,ef3,out3]=fzero(@rate,x0,[],a2,n2,b2,-1);%对二问第二家银行的利率求解 format long ;[x1,x2,x3]%输出其中x1匙小张夫妇买房的月利率，x2是第一家银行的月利率，x3市第二家银行的年利率x2=x2*12,x3%比较第一家与第二家银行的利率大小<rate_run.m>运行MATLAB 结果： ans =0.46 0.84 0.39x2 =0.14 x3 =0.39结论：小张夫妇的贷款的银行的月利率是0.21%，第一家银行的年利率是7.02%，第二家银行的年利率是6.39%，所以第二家银行的利率较低。
6-6．给定4种物质对应的参数i i i c b a ,,和相互作用矩阵Q 如下：1a =18.607， 2a =15.841， 3a =20.443， 4a =19.293; 1b =2643.31，2b =2755.64， 3b =4628.96，4b =4117.07; 1c =239.73， 2c =219.16， 3c =252.64， 4c =227.44;Q=[1.0 0.192 2.169 1.611 0.316 1.0 0.477 0.524 0.377 0.360 1.0 0.2960.524 0.282 2.065 1.0]在压强p =760mmHg 下，为了形成均相共沸混合物，温度和组分分别是多少？请尽量找出所有可能的解。
解：均相共沸混合物的组分问题：模型：所谓的共沸混合物，使之有两种或以上物质组成的液体混合物，当在某种压力下背蒸馏后或局部气化时，在气体状态下和在液体状态下保持相同的组分。
设该混合物由n 个可能的祖坟组成，组分i 所占的比例为i x ，则：11,0.nii i xx ==≥∑（1）均相共沸混合物一概满足问题条件，即公沸混合物的每个组分在气体状态下和在液体状态下具有相同的化学势，在压强P 不大的条件下，这个条件可以表示为：,1,,.i i P P i n γ==L （2）（2）式中i P 是组分的饱和汽相压强，与温度T 有关，可以根据如下表达式确定：ln ii i i b P a T c =-+（3）其中,,i i i a b c 为常数。
（2）式中i γ是组分i 的液相活度系数，可以根据如下的表达式确定：11ln 1ln(),1,,.n n j iji j ij ni j k jk k x q x q i n x q γ==⎛⎫ ⎪⎪=--=⎪ ⎪⎝⎭∑∑∑L （4）其中ijq 表示组分i 和组分j 的交互作用参数，ijq 构成交互作用的矩阵Q ，Q 不一定是对称阵。
对（2）式两边取对数，并将（3），（4）式代入：11ln()1ln 0,1,,.n n j ij ij ij i ni j i k jk k x q b x q a P i n T c x q ==⎛⎫ ⎪⎪++--+==+⎪ ⎪⎝⎭∑∑∑L （5）只有当组分i 参与到该共沸混合物中时才需要满足（5），所以将（5）式改写为11ln()1ln 0,1,,.n n j ij i i j ij i n i j i k jk k x q b x x q a P i n T c x q ==⎡⎤⎛⎫⎢⎥ ⎪⎢⎥ ⎪++--+==+⎢⎥⎪⎪⎢⎥⎝⎭⎣⎦∑∑∑L于是从上面的分析可以对本题进行解答，对于本题参数：a=[18.607,15.841,20.443,19.293]';b=[2643.31,2755.64,4628.96,4117.07]'; c=[239.73,219.16,252.64,227.44]'; 交互矩阵Q=[1.0 0.192 2.169 1.611 0.316 1.0 0.477 0.524 0.377 0.36 1.0 0.296 0.524 0.282 2.065 1.0];于是：便可以得到m 文件：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:nd(i) = x * Q(i,1:n)';dd(i)=x(i)/d(i); endfor i=1:nf(i)=x(i)*(b(i)/(T+c(i)) + log(x*Q(i,1:n)') + dd*Q(1:n,i) - a(i) - 1 + p); end<azerofun. m> 主文件：n=4; P=760;a=[18.607,15.841,20.443,19.293]'; b=[2643.31,2755.64,4628.96,4117.07]'; c=[239.73,219.16,252.64,227.44]'; Q=[1.0 0.192 2.169 1.611 0.316 1.0 0.477 0.524 0.377 0.36 1.0 0.2960.524 0.282 2.065 1.0];%本题的参数 XT0=[0.3,0.3,0.3,58];%给定的初值[XT,Y]=fsolve(@azeofun,XT0,[],n,P,a,b,c,Q)%用fsolve 求解结果： XT =-0.0000 0.5858 0.4142 71.9657 Y =1.0e-010 *0.0002 -0.0505 0.5047 -0.3945对于上述结果：我们对初值XT0的取法：4种物质第一种0.3，第二种0.3，第三种0.3，最6-8，假设商品在t 时刻的市场价格为p （t ），需求函数为D(p(t))=c-dp(t)（c ，d≥0）。
而生产方的期望价格为q （t ），供应函数为S （q （t ））。
当供销平衡时S(q(t))=D （q(t)）。
若期望价格与市场价格不符，商品市场不均衡，生产方t+1时期的期望价格将会调整，方式为[](1)()()()(01)q t q t r p t q t r +-=-<<，以[][](())(())()c D p t c S q t p t dd --==带入：得到关q(t)的变化规律，是否有混沌现象产生？找出几个分叉点，观察分叉点是否极限趋势是否符合Feigenbaum 常数揭示的规律。
解：对于上述过程若S(x)=arctan(ux),u=4.8,r=0.3,d=0.25不难推导出q （t ）的变化规律：(1)(1)()arctan(())r rcq t r q t q t d d μ+=--+于是：可以得到本题的m 文件function y=iter_pq(x,c) u=4.8;d=0.25;r=0.3;y=(1-r)*x-r*atan(u*x)/d+r*c/d;%函数<iter_pq.m > 迭代函数：function chaos(iter_fun,x0,r,n) % 该函数没有返回值；iter_fun 是迭代函数（句柄）；x0是迭代初值；kr=0;for rr=r(1):r(3):r(2) % 输入中[r(1),r(2)]是参数变化的范围，r(3) 是步长 kr=kr+1;y(kr,1)=feval(iter_fun,x0,rr);for i=2:n(2) %输入中n(2)是迭代序列的长度，但画图时前n(1)个迭代值被舍弃y(kr,i)=feval(iter_fun,y(kr,i-1),rr); end endplot([r(1):r(3):r(2)],y(:,n(1)+1:n(2)),'k.');<chaos.m> 主函数：chaos(@iter_pq,0.5,[0.3,1.1,0.001],[100,200])%c 从0.3到1.1变化过程中，寻找分叉点输出结果：分析：从图上找到的几个分叉点的x 坐标分别是：1.079，0.949，0.907，0.897，0.8948，于是可以利用极限表达式：11limn n n n n b b b b -→∞+--可以分别的得到： 2132b b b b --=1.0790.9490.9490.907--=3.093243b b b b --=0.9490.9070.9070.897--=4.2 4354b b b b --=0.9070.8970.8970.8948--=4.54从上面的求解可以看出，当n 越大时比值越接近于4.6692，也就是Feigenbaum 常数。
