第五讲非线性方程模型实验课件.ppt

第五讲非线性方程模型实验第1页，此课件共37页哦实验目的实验目的实验内容实验内容2、学会用、学会用Matlab代数方程的数值解代数方程的数值解.1、学会用、学会用Matlab求代数方程的解析解求代数方程的解析解.1 1、求代数方程的解析解求代数方程的解析解.4 4、实验作业、实验作业.2、求代数方程的数值解、求代数方程的数值解.第2页，此课件共37页哦n问题：如下是一则房产广告。不难算出，你向银行共借了25.2万，30年内共要还51.696万，约为当初借款的两倍，这个案例中贷款年利率是多少？建筑面积总价30%首付70%按揭月还款86.98m236万10.8万30年1436元第3页，此课件共37页哦分析n有人可能会这样算年利率=(51.696-25.2)/30/25.2=3.5%错的，因为你并不是等到30年后一次性还款。第4页，此课件共37页哦设xk第k个月的欠款数；a月还款数；r为月利率，我们得到迭代关系式xk+1=(1+r)xk-a (2.1)那么 xk=(1+r)xk-1-a=(1+r)2xk-2-(1+r)a-a=(1+r)kx0-a(1+r)k-1/r第5页，此课件共37页哦n根据a=0.1436,x0=25.2,x360=0得到n25.2(1+r)360-0.1436(1+r)360-1/r=0(2.2)n关于月利率r的高次代数方程。n年利率R=12r.第6页，此课件共37页哦 非线性方程（组）简介非线性方程（组）简介若方程是未知量x的多项式，称为高次代代数数方方程程；若方程包含x的超越函数，称为超越方程超越方程。一元非线性方程的一般形式为 f(x)=0(2.3)若对于数a有f(a)=0,则称a 为方程（2.3）的解解或根根，也称为函数f(x)的零零点点。方程的根可能是实数也可能是复数。相应地称为实实根根和复复根根。如果对于数a有f(a)=0,f(a)0,则 a称 为单单 根根，如 果 有 k1,f(a)=f(a)=f(k-1)(a)=0但f(k)(a)0,称为k重重根根，对于高次代数方程，其根的个数与其次数相同（包括重数），至于超越方程，其界可能是一个或几个甚至无穷多，也可能无解。第7页，此课件共37页哦常见的求解问题有如下两重要求：一种是要求定出在给定范围内的某个解，而解的粗略位置事先从问题的物理背景或应用（作图等）其他方法得知；另一种是定出方程的全部解，或者给定区域内的所有解，而解的个数未知。除少数特殊的方程可以利用公式直接求解（如4次以下代数方程），一般都没有解析求解方法，只能靠数值方法求得近似解。常见的数值方法有二分法等。n元非线性方程组的一般形式为fi(x1,x2,xn)=0,i=1,m(2.4)非线性方程组的解极少能用解析法求得。常用的数值方法是Newton法、拟Newton法和最优化方法等。第8页，此课件共37页哦解方程和方程组的MATLAB命令roots求多项式的根fsolve方程（组）数值解fzero求一元函数实根solve符号方程（组）求解第9页，此课件共37页哦1.多项式的根roots(p)多项式p的所有复根。例x3+2x2-5的根roots(120-5)ans=-1.6209+1.1826i-1.6209-1.1826i1.2419第10页，此课件共37页哦2.一元函数零点nfzero(F,X,tol)nF为字符串表示的函数或M函数名；nx为标量时，作为迭代初值；X为向量a,b时，返回F在a,b中的一个零点，这时要求F在a,b两点异号；tol为精度（缺损值1e-4).例:y=sin(x)-0.1x第11页，此课件共37页哦 fzero(sin(x)-0.1*x,6)ans=7.0682 fzero(sin(x)-0.1*x,2,6)ans=2.8523 注：fzero 只能求零点附近变号的根，试用fzero求解(x-1)2=0,看看发生了什么？第12页，此课件共37页哦3.非线性方程组求解fsolve用法与fzero类似,例：解方程组写M函数eg2_1fun.mfunctiony=fun(x)y(1)=4*x(1)-x(2)+exp(x(1)/10-1;y(2)=-x(1)+4*x(2)+x(1)2/8;第13页，此课件共37页哦然后用然后用 x,y,f=fsolve(eg2_2fun,0,0)x=0.2326 0.0565y=1.0e-006*0.0908 0.1798f=1注：注：X返回解向量返回解向量,y返回误差向量，返回误差向量，f0则解收敛。则解收敛。第14页，此课件共37页哦或直接用或直接用 x,y,f=fsolve(4*x(1)-x(2)+exp(x(1)/10-1,-x(1)+4*x(2)+x(1).2/8,0,0)x=0.2326 0.0565y=1.0e-006*0.0908 0.1798f=1注意：注意：fsolve采用最小二乘优化法，稳定性比采用最小二乘优化法，稳定性比fzero好，但好，但fsolve 可能陷入可能陷入局部极小。试用局部极小。试用fsolve解解x2+x+1=0,看会发生什么？不要完全相信计算机。看会发生什么？不要完全相信计算机。第15页，此课件共37页哦4.解析求解solve例解ax2+bx+c=0solve(a*x2+b*x+c,x)ans=1/2/a*(-b+(b2-4*a*c)(1/2)1/2/a*(-b-(b2-4*a*c)(1/2)第16页，此课件共37页哦 x,y=solve(4*x-y+exp(x)/10=1,-x+4*y+y2/8=0,x,y)x=.23297580773115396971569236570313y=.58138324907069742242891748561961e-1注意所得的解与注意所得的解与fsolve的不同。的不同。注意：虽然注意：虽然solve可用于求数值解，但速度很慢，且有很大的局限性，可用于求数值解，但速度很慢，且有很大的局限性，不提倡使用。不提倡使用。第17页，此课件共37页哦数值解法：图解法和迭代法图解法和迭代法n1.图解法例解方程sin(x)=0.1x （2.5)显然，解在-10,10内，函数y=sinx-0.1x的零点就是(2.5)的解，作出y=sinx-0.1x在-10,10范围内的图象（图2.1）,可看出根的大致位置。作图可使用如下MATLAB语句：close;fplot(sin(x)-0.1*x,-10,10);grid;第18页，此课件共37页哦可知8.5,7,3,0附近各有一解。（在figure窗口用matlab的zoom命令演示）第19页，此课件共37页哦2、迭代法（牛顿法，切线法）n求f(x)=0的解，从几何上说xk+1为用f(x)在xk处的切线代替f(x)求得的解，故也称为切线法。当初值x0与真解足够靠近，Newton迭代法敛。单根快，重根慢。迭代格式：第20页，此课件共37页哦例求如下方程的正根（要求精度=10-6）x2-3x+ex=2解：令f(x)=x2-3x+ex-2,f(0)=-12,f(x)0,f(x)0,即f(x)单调上升，根在0，2，先用图解法找初值。fplot(x2-3*x+exp(x)-2,0,2);gridon;第21页，此课件共37页哦唯一正根在1附近，取x0=1,迭代格式M脚本脚本eg2_2.mclear,e=1e-6;format long;x1=1x0=x1+2*2;%使使while成立成立while(abs(x0-x1)e)x0=x1,x1=x0-(x02-3*x0+exp(x0)-2)/(2*x0-3+exp(x0)end;format得得x1=1.44623868596643第22页，此课件共37页哦贷款利率问题求解考虑方程(2.2).常识上，r应比当时活期存款月利率略高。用活期存款月利率0.0198/12作为迭代初值，用fzero求解。（使用Matlab)r=fzero(25.2*(1+x)360-(1+x)360-1)/x*0.1436,0.0198/12),R=12*rr=0.0046R=0.0553第23页，此课件共37页哦练习n1、作出f(x)=xsin(1/x)在-0.1,0.1内的图，可见在x=0附近f(x)=0有无穷多个解，并设法求出它的解。n2、（月还款额）作为房产公司的代理人，你要迅速准确回答用户各方面的问题。现在有个客户看中了贵公司一套建筑面积为120m2,单价5200元/m2的房子。他计划首付30%，其余70%用20年按揭贷款（年利率5.58%)。请你提供下列信息：房屋总价格、首付款额、月付还款额。第24页，此课件共37页哦补充：混沌补充：混沌线性迭代要么收敛于它的不动点，要么趋于无穷大；而不收敛的非线性迭代可能会趋于无穷大，也可能趋于一个周期解，但也可能在一个有限区域内杂乱无章地动弹，由确定性运动导致的貌似随机的现象称为混沌现象。下面就Logistic迭代研究这一现象。第25页，此课件共37页哦1.昆虫数量的昆虫数量的Logistic模型模型xk表示第 k代昆虫数量（1表示最大值）。(2.7)式反映了下一代对上一代的既依赖又竞争的关系。当上一代很少，繁殖能力不够，从而后代很少；当上一代很多，会吃掉很多食物，后代难以存活，从而后代很少。a为资源系数，0a 4保证了xk在区间(0,1)上封闭。第26页，此课件共37页哦2.平衡与稳定平衡与稳定称a为映射g(x)的平衡解平衡解或不动点不动点，若g(x)=ax(1-x).解方程x=ax(1-x)得(2.7)式两个不动点0和1-1/a.若初始值恰好为不动点，迭代式(2.7)的只永不改变。如果对于不动点x0附近的初始值，(2.7)收敛与此不动点，我们称这一不动点是稳定的稳定的。当0 x1,在0,1内只有一个不动点0，且由|g(0)|=a1,不动点0不再稳定，而由|g(1-1/a)|=|2-a|1可知1a3时不动点1-1/a稳定，说明资源适当时，昆虫稳定于一定数量。第27页，此课件共37页哦3.周期解、分叉和混沌周期解、分叉和混沌称a为映射g(x)的周期k点，若gk(a)=a,而对任意j0,即a3,出现两个周期2解，可以证明3aa,(2.4)是的迭代序列几乎杂乱无章，即所谓混沌混沌。下列例子可形象地显示上述现象。例（分叉图）对a在0,4的不同值，画出Logistic迭代的极限形态图。如下M文件对于每一个a值，随机产生一个初值。文件显示前20步迭代的变化。最后用第180200步迭代值表示极限形态，最后结果见图2-3。第29页，此课件共37页哦%M脚本脚本2_3.mclear;close;a=0:0.01:4;M=length(a);K=200;X=zeros(K,M);x(1,:)=rand(1,M);for m=1:M,for k=1:K-1 x(k+1,m)=a(m)*x(k,m)*(1-x(k,m);end,endfor k=1:20,plot(a,x(k,:),.);title(k=,int2str(k);pause(2);end;plot(a,x(180:K,:),.);xlabel(a);ylabel(x);hold off;第30页，此课件共37页哦第31页，此课件共37页哦4.混沌的特征混沌是由确定性系统产生的貌似随机的现象。一般认为混沌有如下几个特征(i)初值的敏感性：两个任意近的点出发的两条轨迹迟早会分得很开；(ii)遍历性：任意点出发的轨迹总会进入0,1内任意小的开区间。例（初值的敏感性）如下M文件eg2_4.m验证了Logistic迭代序列的初值敏感性。对于靠得很近的两个初值（相差仅1e-4),画出了两个序列50步内的误差图（图2-4)。可见10步以后，差异增大，有时甚至接近1。第32页，此课件共37页哦%M脚本脚本eg2-4.mclear;close;a=4;e=1e-4;x=zeros(50,2);x(1,:)=0.4,0.4+e;for i=2:50 x(i,:)=a*x(i-1,:).*(1-x(i-1,:);endy=x(:,1)-x(:,2);plot(y)第33页，此课件共37页哦例子：（蛛网图）混沌的遍历性昆虫数量的Logistic模型(eg2_5.m)xk表示第k代昆虫的数量（1表示最大可能数量）。平面迭代式：蛛网图正好显示迭代计算蛛网图正好显示迭代计算x0,y0,x1,y1,的一系列变化过程。的一系列变化过程。如下如下M函数函数eg2_5.m是一个通用的是一个通用的logistic 蛛网图函数。作蛛网图函数。作出系数为出系数为a，初值为，初值为x0,从第从第m步到第步到第n步的迭代过程步的迭代过程第34页，此课件共37页哦function f=eg2_5(a,x0,m,n)x=0:0.01:1;y=a*x.*(1-x);plot(x,x,r,x,y,r);hold on;clear x,y;x(1)=x0;y(1)=a*x(1)*(1-x(1);x(2)=y(1);if mm,plot(x(i),x(i),x(i+1),y(i-1),y(i),y(i);end end hold off第35页，此课件共37页哦在命令窗口执行subplot(2,2,1);eg2_5(2.7,0.1,1,100);%收敛迭代subplot(2,2,2);eg2_5(3.4,0.1,50,500);%周期2subplot(2,2,3);eg2_5(3.5,0.1,50,500);%周期4subplot(2,2,4);eg2_5(4,0.1,50,500);%混沌可见混沌迭代对于初值为0.1，轨迹遍历了0,1区间（图2-5）第36页，此课件共37页哦作业作业(Henon吸引子）混沌和分形的著名例子，迭代模型为取初值x0=0,y0=0,进行3000次迭代，对于k1000,在(xk,yx)初亮一点（注意不要连线）可得所谓Henon引力线图。第37页，此课件共37页哦