《数学实验第三次讲稿知识讲解》由会员分享,可在线阅读,更多相关《数学实验第三次讲稿知识讲解(43页珍藏版)》请在金锄头文库上搜索。
1、数学实验第三讲,方程求解,2020/7/6,1,1) 掌握方程求解的三种解法:解析 法、数值解法以及图形表示解的方法; 2) 学会使用MATLAB软件求解析解、数值解 和图形解; 3) 通过范例学习怎样建立方程模型和 分析问题的思想。,实验目的,2020/7/6,2,一,建立方程,例子【问题背景】一段时间, 美国原子能委员会是 按以下方式处理浓缩放射性废物的. 他们将废物装入 密封性能很好的圆桶中, 然后扔到水深300英尺的海 里. 这种做法是否会造成放射性污染, 很自然地引起 了生态学家及社会各界的关注. 原子能委员会一再保 证, 圆桶非常坚固, 决不会破漏, 这种做法是绝对安 全的. 然而
2、一些工程师们却对此表示怀疑, 他们认为 圆桶在海底相撞时有可能发生破裂. 由此双方展开了 一场笔墨官司. 究竟谁的意见正确呢? 只能让事实 说话了!,2020/7/6,3,数学建模的一般步骤,2020/7/6,4,1,方程求解之解析方法,主要针对一些比较简单的方程以及方程组,比如多项式方程等。同学们以前对方程的求解也是针对这样一些方程进行的。该方法的优点是可以利用纸笔得到简单有效并且精确的解;缺点是可以求解的方程数量太少。Matlab和Maple提供了求方程解析解的函数,可以说对数学演算提供了不少方便。,2020/7/6,6,2,方程求解之图形放大法,图形最大的有点就是直观,试想如果我们有了函
3、数准确的图形,那么曲线和x轴的交点就是我们要求的方程的解。因此我们可以利用图形工具得到方程的解。当然,计算机上的图形不可能等同于函数的真实图形,因为计算机上的图形是曲线上部分点的轨迹而不是全部,因此通过图形不可能得到方程的精确解,甚至它只是一个比较粗略的解,当然,通过对图形的放大可以得到更精确一些的解。同时,这种方法也不适应大量的数据处理。,2020/7/6,7,方程 f(x)=0 1)建立坐标系,画曲线f(x); 2)观察曲线f(x)与x轴相交的交点; 3)将其中一个交点进行局部放大; 4)该交点的横坐标值就是方程的根。,2图形放大法,2020/7/6,8,例: 求方程 x5 +2x2 +
4、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图形放大法,2020/7/6,9,逐次缩小区间,观察一个根在-1.55-1.5之间。,2图形放大法,2020/7/6,10,3,方程求解之迭代法,迭代法的理论以及方法的出现,对方程求解有着里程碑式的意义。其基本思想如下:,需要求解的方程: f (x) = 0 (1) 经过某种变形得:x = j (
5、x) (2) 从而求解方程(1)转化成为求解(2)得不动点。(满足条件x*=j(x*)的点x*称为不动点) 为得到方程的不动点,可以构造迭代过程如下: xn+1 =j (xn),n =0,1, x0 定义为迭代初值。,2020/7/6,11,解: 第一步 构造迭代函数: x=j (x),例:用迭代方法求解方程 x3 -x2 -x-1 = 0。,2020/7/6,12,第二步 迭代 设定初值 x0=1, xn+1 =j (xn),n =0,1, 用 MATLAB 编程(died2.m文件),x=1;y=1;z=1;(初始点) for k=1:20 x=x3-x2-1; % j1 (x) y=(y
6、2+y+1)(1/3); % j2 (y) z=1+1/z+1/z2; % j3 (z) end X,y,z,2020/7/6,13,j1(x)的迭代是失败的(迭代不收敛 )。,精确解:x=1.8393,计算结果,2020/7/6,14,迭代函数j2(x)和j3(x)的选取是成功的。精确解为 x=1.8393。 并且选取函数j2(x)、j3(x)其收敛速度不一致,前者的速度快些!,结论,1、当遇到迭代不收敛时有什么解决办法? 2、如何提高收敛速度?,对于给定的方程 f(x) = 0, 有多种方式将它改写成等价的形式 x = j(x)。但重要的是如何改写使得序列收敛?,返 回,2020/7/6,
7、15,当今最流行的迭代法是牛顿法以及由此改进的一些方法,比如拟牛顿法等。其基本的思想就是构造迭代格式是利用函数的导数,这类方法有收敛速度快,稳定性好等特点。对低维和高维情况都适合,也是当今一些软件均采用的方法。当然,因为需要函数的导数信息,所以自然对不可微的问题受到制约。该方法的迭代格式为: xk+1=xk-f(xk)/f(xk),2020/7/6,16,例:利用牛顿法求方程x3 -x2 -x-1 = 0的根. 方法:第一步,给出函数的导函数3x2-2x-1; 第二步,给出函数的迭代格式: xk+1=xk-f(xk)/f(xk); 设置一定的精度要求,达到即终止. 定义函数m文件: funct
8、ion ff=mynewton(x) ff=(x3 x2 -x-1 )/(3*x2-2*x-1); 定义命令m文件: x0=1; x1=x0-mynewton (x0); while abs(x1-x0)0.0001 x0=x1; x1=x0-mynewton (x0); end,2020/7/6,17,三,解方程函数格式及例子,Matlab对方程的求解提供了以下的一些函数: (1)多项式求根; (2) 线性方程组求解; (3)一般的非线性方程(组)求解:,2020/7/6,18,输出: -1.2131 -0.9017 + 0.5753i -0.9017 - 0.5753i -0.2694 +
9、 0.9406i -0.2694 - 0.9406i 0.4168 + 0.8419i 0.4168 - 0.8419i 0.8608 + 0.3344i 0.8608 - 0.3344i,例:求解多项式方程 x9+x8+1=0,输入: p=1,1,0,0,0,0,0,0,0,1; roots(p),roots()语句的用法,2020/7/6,19,4、线性方程组:AX = b 其中A是mn阶矩阵,b是m维向量。,x=A b or x=inv(A)*b,特点:只能求出一个特解。,2020/7/6,20,solution to the following linear system of equ
10、ations: 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,2020/7/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
11、 和inv()语句的用法,返 回,2020/7/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,2020/7/6,23,1、方程(组), f1(x) = 0,fn(x) = 0, x = (x1,xn),2、方程(组
12、), f1(x) = 0,fn(x) = 0, x = (x1,xn),solve(f1(x),f2(x),fn(x) ),X = fsolve (fun, X0, options),MATLAB软件直接求解法,2020/7/6,24,输出: 1/2/a*(-b+(b2-4*a*c)(1/2) 1/2/a*(-b-(b2-4*a*c)(1/2),单变量方程,solve()语句的用法,例1: 求解方程 ax2+bx+c = 0,输入: x= solve(a*x2+b*x+c) 或 solve(a*x2+b*x+c=0),1)符号解,2020/7/6,25,例2: 解方程: x3-2x2=x-1
13、解: s=solve(x3-2*x2=x-1) double(s),2)数字解,该方程是否有实根? vpa(s,10),solve()语句的用法,2020/7/6,26,例3 求解方程: tan(x)-sin(x)=0,3)无穷解,输入: solve(tan(x)-sin(x)=0) 输出:0 (不能给出全部解) (tx1.m),solve()语句的用法,2020/7/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,solve()语句的用法,2020
14、/7/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),fsolve()语句的用法,2020/7/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 T
15、erminated Successfully y= 0.5990 2.3959 2.0050,fsolve()语句的用法,2020/7/6,30,fsolve()函数的第三个输入是options,它是一个结构型数据,可以通过函数optimset ()进行设定.当不进行设定时采用缺省设置. fsolve函数还可以有后面的参数设定,这一功能在有的时候非常有用.,fsolve()语句的用法,例如,求解函数sin(ax)-x=0的最小正解. 分析: 已知的知识我们知道,该方程没有解析形式的解,也就是说很难得到这个解和a之间的具体关系.数值的方法可以对a先取确定的值,这样变成一个一元方程,容易进行求解,
16、当我们变化a时,就得到很多这样的解,通过这样的方法,我们可以得到解和a之间的一些大致关系.,2020/7/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),2020/7/6,32,2020/7/6,33,例:某物体的边缘呈圆形,通过测量边沿上是一个点的坐标,数据为(见程序)。使用fsolve计算物体边缘的曲线方程。 解:假定物体边缘的曲线方程为(x-a)2+(y-b)2=