资源描述
单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,*,本资料仅供参考，不能作为科学依据。谢谢。本资料仅供参考，不能作为科学依据。本资料仅供参考，不能作为科学依据。谢谢。本资料仅供参考!,数学实验第三讲,方程求解,1/43,/10/6,1,1),掌握方程求解三种解法：解析,法、数值解法以及图形表示解方法,;,2),学会使用,MATLAB,软件求解析解、数值解,和图形解；,3,）经过范例学习怎样建立方程模型和,分析问题思想。,实验目,2/43,/10/6,2,一，建立方程,例子,【,问题背景,】,一段时间,美国原子能委员会是,按以下方式处理浓缩放射性废物,.,他们将废物装入,密封性能很好圆桶中,然后扔到水深,300,英尺海,里,.,这种做法是否会造成放射性污染,很自然地引发,了生态学家及社会各界关注,.,原子能委员会一再保,证,圆桶非常坚固,决不会破漏,这种做法是绝对安,全,.,然而一些工程师们却对此表示怀疑,他们认为,圆桶在海底相撞时有可能发生破裂,.,由此双方展开了,一场笔墨官司,.,终究谁意见正确呢,?,只能让事实,说话了,!,3/43,/10/6,3,数学建模普通步骤,模型准备,模型假设,模型组成,模型检验,模型分析,模型求解,模型应用,4/43,/10/6,4,二，方程求解,1,，解析方法,2,，图形放大法,3,，迭代方法,4,，区间方法,5/43,/10/6,5,1,，方程求解,之解析方法,主要针对一些比较简单方程以及方程组，比如多项式方程等。同学们以前对方程求解也是针对这么一些方程进行。该方法优点是能够利用纸笔得到简单有效而且准确解；缺点是能够求解方程数量太少。,Matlab,和,Maple,提供了求方程解析解函数，能够说对数学演算提供了不少方便。,6/43,/10/6,6,2,，方程求解,之,图形放大法,图形最大有点就是直观，试想假如我们有了函数准确图形，那么曲线和,x,轴交点就是我们要求方程解。所以我们能够利用图形工具得到方程解。当然，计算机上图形不可能等同于函数真实图形，因为计算机上图形是曲线上部分点轨迹而不是全部，所以经过图形不可能得到方程准确解，甚至它只是一个比较粗略解，当然，经过对图形放大能够得到更准确一些解。同时，这种方法也不适应大量数据处理。,7/43,/10/6,7,方程,f,(,x,)=0,1,）,建立坐标系，画曲线,f,(,x,),；,2,）,观察曲线,f,(,x,),与,x,轴相交交点,；,3,）,将其中一个交点进行局部放大,；,4,）,该交点横坐标值就是方程根,。,2,图形放大法,8/43,/10/6,8,例,:,求方程,x,5,+2,x,2,+,4,=0,一个根,.,画方程曲线图（,tuxfd.m,）,x=-3:0.01:3;,y=x.5+2*x.2+4;,y1=0,*,x;,plot,(x,y,x,y1),由此判断：方程一个根在区间,-2,，,2,内，所以将区间,-3,，,3,缩小至,-2,，,2,，再观察！,该方程有几个根？欲寻找其中一个实根，而且到达一定精度。,2,图形放大法,9/43,/10/6,9,逐次缩小区间，观察一个根在,-1.55-1.5,之间。,2,图形放大法,10/43,/10/6,10,3,方程求解,之迭代法,迭代法理论以及方法出现，对方程求解有着里程碑式意义。其基本思想以下,:,需要求解方程：,f,(,x,)=0,（,1,）,经过某种变形得：,x,=,j,(,x,),（,2,）,从而求解方程（,1,）转化成为求解（,2,）得不动点。（满足条件,x*=,j(,x*),点,x*,称为不动点）,为得到方程不动点，能够结构迭代过程以下：,x,n,+1,=,j,(,x,n,),，,n,=0,，,1,，,x,0,定义为迭代初值。,11/43,/10/6,11,解：,第一步,结构迭代函数：,x,=,j,(,x,),例：,用迭代方法求解方程,x,3,-,x,2,-,x,-1=0,。,12/43,/10/6,12,第二步,迭代,设定初值,x,0=1,，,x,n,+1,=,j,(,x,n,),，,n,=0,，,1,，,用,MATLAB,编程,（,died2.m,文件,),x=1;y=1;z=1;(,初始点,）,for,k=1:20,x=x3-x2-1;%,j,1,(,x,),y=(y2+y+1)(1/3);%,j,2,(,y,),z=1+1/z+1/z2;%,j,3,(,z,),end,X,，,y,，,z,13/43,/10/6,13,序号,j,3,(,x,),序号,j,2,(,x,),j,3,(,x,),1,j,2,(,x,),1.4422,3.0000,8,1.8175,1.8136,2,1.6537,1.4444,9,1.8385,1.8554,3,1.7532,2.1716,10,1.8389,1.8294,4,1.7995,1.6725,11,1.8391,1.8454,5,1.8209,1.9554,12,1.8392,1.8355,6,1.8308,1.7730,13,1.8392,1.8416,7,1.8354,1.8822,j,1,(,x,),迭代是失败,（,迭代不收敛,）。,准确解：,x,=1.8393,计算结果,14/43,/10/6,14,迭代函数,j,2,(,x,),和,j,3,(,x,),选取是成功。准确解为,x,=1.8393,。,而且选取函数,j,2,(,x,),、,j,3,(,x,),其收敛速度不一致，前者速度快些！,结论,1,、当碰到迭代不收敛时有什么处理方法？,2,、怎样提升收敛速度？,对于给定方程,f,(,x,)=0,有各种方式将它改写成等价形式,x,=,j,(,x,),。,但主要是,怎样改写使得序列收敛？,返 回,15/43,/10/6,15,当今最流行迭代法是牛顿法以及由此改进一些方法，比如拟牛顿法等。其基本思想就是结构迭代格式是利用函数导数，这类方法有收敛速度快，稳定性好等特点。对低维和高维情况都适合，也是当今一些软件均采取方法。当然，因为需要函数导数信息，所以自然对不可微问题受到制约。该方法迭代格式为：,x,k+1,=x,k,-f(x,k,)/f(x,k,),16/43,/10/6,16,例,:,利用牛顿法求方程,x,3,-,x,2,-,x,-1=0,根,.,方法,:,第一步,给出函数导函数,3x,2,-2x-1;,第二步,给出函数迭代格式,:,x,k+1,=x,k,-f(x,k,)/f(x,k,);,设置一定精度要求,到达即终止,.,定义函数,m,文件,:,function ff=mynewton(x),ff=(,x3,x2,-,x,-1,)/(,3*x2-2*x-1,);,定义命令,m,文件,:,x0=1;,x,1,=x,0,-mynewton(x,0,);,while abs(x1-x0)0.0001,x0=x1;,x,1,=x,0,-mynewton(x,0,);,end,17/43,/10/6,17,三，解方程函数格式及例子,Matlab,对方程求解提供了以下一些函数：,（,1,）多项式求根；,（,2,）线性方程组求解；,（,3,）普通非线性方程（组）求解：,18/43,/10/6,18,输出,:,-1.2131,-0.9017+0.5753i,-0.9017-0.5753i,-0.2694+0.9406i,-0.2694-0.9406i,0.4168+0.8419i,0.4168-0.8419i,0.8608+0.3344i,0.8608-0.3344i,例：求解多项式方程,x,9,+,x,8,+1=0,输入,:p=1,1,0,0,0,0,0,0,0,1;,roots(p),roots(),语句使用方法,19/43,/10/6,19,4,、线性方程组：,AX=b,其中,A,是,m,n,阶矩阵，,b,是,m,维向量。,x=A b,or x=inv(A)*b,特点：只能求出一个特解。,20/43,/10/6,20,solution to the following linear system of equations:,You can formulate and solve the problem as,A=3 11-2;1 1-2;1-1 1;,b=7;4;19;,x=Ab,x=,13.2188,-2.3438,3.4375,21/43,/10/6,21,例,:AX=b,解,:,输入,:,A=1 2 3;4 5 6;7 8 9;b=6;14;-3;,x1=Ab,x2=inv(A)*b,输出,:,警告,:,系统秩不足,.,解不唯一,.,1,、,题中,rank(A)=rank(A|b)=23,该方程组有没有穷解,。,2,、,输出结果是否一致？,3,、,怎样求方程组全部解？,思,考,Ab,和,inv(),语句使用方法,返 回,22/43,/10/6,22,函数,fzero,格式：,fzero(,函数名,初值或区间,),Example 1.Calculate by finding the zero of the sine function near 3.,x=fzero(sin,3),x=3.1416,Example 2.To find the zero of cosine between 1 and 2,x=fzero(cos,1 2),x=1.5708,23/43,/10/6,23,1,、方程,(,组,),，,f,1,(,x,)=0,f,n,(,x,)=0,x,=(,x,1,x,n,),2,、方程,(,组,),，,f,1,(,x,)=0,f,n,(,x,)=0,x,=(,x,1,x,n,),fun.m,function f=fun(x),f(1)=,f,1,(,x,);,f(n)=,f,n,(,x,),初值,1,）能够省略。,2,）,options=1,表示输出中间结果。,solve,(,f,1,(,x,),f,2,(,x,),f,n,(,x,),),X=fsolve(fun,X0,options),MATLAB,软件直接求解法,24/43,/10/6,24,输出,:1/2/a*(-b+(b2-4*a*c)(1/2),1/2/a*(-b-(b2-4*a*c)(1/2),单变量方程,solve(),语句使用方法,例,1,：求解方程,ax,2,+,bx,+,c,=0,输入,:,x=,solve(,a*x2+b*x+c,),或,solve(,a*x2+b*x+c=0,),1,）符号解,25/43,/10/6,25,例,2,：解方程：,x,3,-2x,2,=x-1,解：,s=solve(x3-2*x2=x-1),double(s),2,）数字解,该方程是否有实根？,vpa(s,10),solve(),语句使用方法,26/43,/10/6,26,例,3,求解方程：,tan(x)-sin(x)=0,3,）无穷解,输入,:,solve(tan(x)-sin(x)=0),输出：,0,(,不能给出全部解）,(tx1.m),solve(),语句使用方法,27/43,/10/6,27,输入：,x,y=solve(x2*y2,x-(y/2)-b),输出：,x=0,，,y=-2*b,0,，,-2*b,（符号解）,b,，,0,b,，,0,v=x,y,多变量方程组,例,4,或,b=2,solve(),语句使用方法,28/43,/10/6,28,例,6,：求解方程组,解 输入：,x,y,z=,solve,(sin(x)+y2+log(z)-7=0,3*x+2y-z3+1=0,x+y+z-5=0,x,y,z),x=0.59905375664056731520568183824539,y=2.3959314023778168490940003756591,z=2.0050148409816158357003177860955,输出,：,fsolve(),语句使用方法,29/43,/10/6,29,解：,1,）建立方程组,M-,函数文件,(nxxf.m),function,eq=,nxxf,(x),eq(1)=sin(x(1)+x(2)2+log(x(3)-7;,eq(2)=3*x(1)+2x(2)-x(3)3+1;,eq(3)=x(1)+x(2)+x(3)-5;,运行程序,(test4.m),y=,fsolve,(,nxxf,1,1,1),3,）运行结果,：,Optimization Terminated Successfully,y=0.5990 2.3959 2.0050,fsolve(),语句使用方法,30/43,/10/6,30,fsolve(),函数第三个输入是,options,它是一个结构型数据,能够经过函数,optimset(),进行设定,.,当不进行设定时采取缺省设置,.fsolve,函数还能够有后面参数设定,这一功效在有时候非常有用,.,fsolve(),语句使用方法,比如,求解函数,sin(ax)-x=0,最小正解,.,分析,:,已知知识我们知道,该方程没有解析形式解,也就是说极难得到这个解和,a,之间详细关系,.,数值方法能够对,a,先取确定值,这么变成一个一元方程,轻易进行求解,当我们改变,a,时,就得到很多这么解,经过这么方法,我们能够得到解和,a,之间一些大致关系,.,31/43,/10/6,31,详细方法以下,:,建立函数,m,文件,:,function ff=funpara(x,a),ff=sin(a*x)-x;,fsolve(),语句使用方法,建立对应命令,m,文件,:,B=zeros(100,1);,for a=1:100,x0=pi/(2*a)+0.01;,B(a)=fsolve(funpara,x0,a);,end,plot(B),32/43,/10/6,32,例：,某物体边缘呈圆形，经过测量边缘上是一个点坐标，数据为,（见程序）。,使用,fsolve,计算物体边缘曲线方程。,解：假定物体边缘曲线方程为（,x-a,）,2+(y-b)2=r2,为了得到,a,，,b,和,r,值。,将上述十一组数据代人方程，得到,（,xi-a,）,2+(yi-b)2=r2,（,i=1,11,）。次数方程个数是,11,个，未知量个数是,3,个。使用,fsolve,求解，其函数,m,文件为：,function y=funcircle(x),A=6.7630 5.1313 2.4713 -0.3435 -2.3887 -2.9927 -1.9572 0.3778 3.2455 5.7042 6.9465;,B=23.2879 25.6492 26.7268 26.1668 24.1531 21.3470 18.6699 17.0010 16.8883 18.3688 20.9564;,for i=1:11,y(i)=(A(i)-x(1)2+(B(i)-x(2)2-x(3)2;,end,33/43,/10/6,33,区间方法,对于一个闭区间上连续函数，我们有一个,0,点存在定理，利用这个定理，能够对不可微函数求得函数,0,点。,基本思想是经过判断函数在端点处函数值异号能够确定函数在开区间上最少有一个,0,点，然后经过缩小区间得到解近似。,优点是不需要函数导数信息，而且只要有,0,点就一定能够得到；缺点是相对牛顿法等速度较慢。,34/43,/10/6,34,问题关键在于圆桶到底能承受多大速度碰撞,?,圆桶和海底碰撞时速度有多大,?,工程师们进行了大量破坏性试验,发觉圆桶在直线速度为,40 ft/s,冲撞下会发生破裂,剩下问题就是计算圆桶沉入,300 ft,深海底时,其末速度终究有多大,?,问题分析,引例分析和求解,35/43,/10/6,35,1.,使用,55,加仑圆桶,;(1,加仑,=3.7854,升,),2.,装满放射性废物时圆桶重量为,W=527.436,磅,(1,磅,=0.4526,千克,),3.,在海水中圆桶受到浮力,B=470.327,磅,4.,圆桶下沉时受到海水阻力,D=C,v,C,为常数,经测算得,:C=0.08.,5.,建立坐标系,取垂直向下为坐标方向,y,海平面为坐标原点,.,y,0,问题假设,引例,36/43,/10/6,36,依据牛顿第二定律,圆桶下沉时应满足微分方程,:,建立模型,引例,37/43,/10/6,37,为了求出圆桶与海底碰撞速度,v(t),需要求出圆桶下沉到海底,300,英尺时时间,t,再计算,v(t),，,要做到这一点是十分困难,.,若将速度,v,看成是海水深度,y,函数,即,由复合函数求导法知,建立模型,引例,38/43,/10/6,38,借助数值方法求出,v,(300)=45.1ft/s,，显然大于,40ft/s,。,结论：放射性废物不能随意放入公海！,建立模型,引例,39/43,/10/6,39,1,、求以下方程根,1,）,xsin(x)+cos(x)-sin(x)-2x=0,在,-1,，,1,上近似解。,2,）,判定方程,x,7,+x,5,+1=0,有几个实根,?,3),找出函数,f(x)=x,3,-6x,2,-2x+12,全部零点。,2,、,求解线性方程组：,实 验 内 容,40/43,/10/6,40,3,、求以下方程组解,使用命令,solve(),或,fsolve(),。,x0=-5,-5,实 验 内 容,41/43,/10/6,41,4.,炮弹发射视为斜抛运动，已知初始速度为,200 m/s,，问要击中水平距离,360m,、垂直距离,160m,目标，当忽略空气阻力时，发射角应多大？,深入思索,：假如要考虑水平方向阻力，且设阻力与（水平方向）速度成正比，系数为,0.1,（,1/s,），结果又怎样？炮弹质量假设为,10,千克？,实 验 内 容,42/43,/10/6,42,5.,讨论,:,有这么电视节目,对一件商品要求观众猜价格,主持人对观众所猜数目标回答是高或者低,(,相对于实际价格,).,现在问题时,对于价格在,1000,元之内商品,寻找一个最好方法,确保在最少次数猜中商品价格,.,试给出竞猜,5,次得到正确答案概率位多大,.(,设商品最小价格单位为元,),实 验 内 容,结 束,提醒,:,利用二分法思想,能够比较快得到竞猜次数最大值,.,一个更有趣方法是二叉树,其深度就是对应解,.,对于平均次数,能够采取模拟方法,比如产生,10000,次随机数,看每次猜中次数,平均就得到平均次数,.,43/43,/10/6,43,
展开阅读全文
浙ICP备2021020529号-1 | 浙B2-20240490 
