某大学计算方法课件-第五章非线性方程数值解法.pdf

第五章非线性方程的数值解法科学工程中常常涉及求解非线性方程F(X)=0,(1)这里 X=(的,陶，Xn)T,F(X)=(力(X),这X)，,加X)T 6一般而言，其精确求解是非常困难的，例如对于高于4次的代数方程，理论分析已证明其不存在精确求根公式。为此，我们必须借助于数值算法来求解各种非线性方程。5.1几何方法在本节，我们介绍非线性标量方程的二类几何解法：二分法和弦截法。考虑非线性标量方程/=0.(2)为保证方程 在求解区间a,b内存在实根，我们假设/3GC(QM,f b.乙计算g 处的函数值/(g),若/(g)=0,则取7*=空；否则，转至 Step 2;Step 2.若/(亨)0,乙则记Ql=Q,仇=空；若/(中/0,则记。1 =喏，瓦=b。取近似根X1=1:仇 G 1,6 1；St ep 3.重复上述步骤，得根/*的下列隔离区间序列M%D M 瓦D%,蚓DM 切 一.这里，二分k 次后隔离区间山,加 的长度为bk _ CLk=近似根股=啥 G k，其有误差估计成-行|bk dk _ b a2因此，若要近似根3 达到预定精度：|即当2卜+1 .1*Xk 瞑 只 需 貂 ,k lnl)Q)/?2(2)时，可终止迭代，取 跳=中Zn 2作为欲求近似根。_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 算法5.1 二分法function x=half(a,b,tol)c=(a+b)/2;k=1;m=l+round(lo g(b a)lo g(2*t o l)/lo g(2);while k=?nif f(c)=0c=c;return;e ls e if f(a)*f(c)0b=(a+b)/2;elsea=(a+b)/2;endc=(a+b)/2;k=k+1;endk二分方法是方程求根问题的一种直接搜索方法，其优点是算法简单直观，数值解的精度易于判别。该算法的局限性是仅适用于标量方程，且事先要确定方程的根的一个隔离区间，当隔离区间较长时其计算速度较慢。例5.1用二分法求方程/=4 sin x 在区间自7 T 内的根，要求绝对误差小于10一8.解 记/(c)=一 4s in c,则/e 0(假,行 ),且其满足7T 7T/(9)=?-4 0.此外，曾 6区对有产 =1 4 cos /0,即/在区间隹句上单调递增。因此,原方程在区间诙加内有唯一根。编制函数文件f.m:function y=f(x)y=x 4：sin x运行Mat lab命令九Q/(7F/2,TT,10-8)得满足精度要求的数值解%=2.47 457 6 7 8 7 96 451,其需迭代28 次。5.1.2弦截法设方程(2)中 函 数 在 区 间a,b上除满足(3)外还满足如下条件：存在精确解7*的一个邻域SQ*,6=x:x-x d CQ,可使得f 在该邻域内二阶连续可微，且f 3 7 0,Vx G S(/*,5),(4)里叫/(叫6_ *eSQ*,8)W：=2 min If M|-WS(谈)(5)今在邻域SO*)=N*6 内 任 取 二 点 作 以(g,/(3)和Q i J(g)为端点的弦人：其与X轴相交于点力1-3 =为 一 一）一 进一步，作 以Qi JQi）和（物/（3）为端点的弦，2：y =/（6 2）+-力2）,X2 Xy其与X轴相交于点力2-3 一、磔=妆./（旬 /3）/3）；如此循环往复，可得一列逼近X*的点，/Q-（见图5.1）,其一般表达式为跳+i=八1 2 .（6）该公式所表征的求解方法称为弦截法或割线法。理论分析表明：当条件(3)-(5)成立时，由弦截法产生的序列 交 收敛于精确解*,且有如下误差估计2 min f,M max 1此 S（2*）例5 2用弦截法求方程I =4sinx在区间质河 内的根，要求xk-xk-i 10-8o解 记/(幻=1 4 sin i。例 5.1的求解过程已表明方程/(为=0在区间 f,7T内有唯一根。今取初始迭代点g=7T/2,%=7T,应用弦截法(6)到该方程得满足精度要求的数值解X=2.47457678736983,其需迭代7 次。例 5.1与例5.2相比较，显见弦截法的收敛速度要快于二分法的收敛速度。此外，若弦截法产生的迭代序列跳充分接近于精确解/*,则 以 与/(跳)均充分接近于零。因此，我们可近似地置xk-xk-i x/(跳)，从而弦截法可改造为Steffensen方法3+1=以一而八3)Jxk)Jxk Jk)一般而言，Steffensen方法的收敛速度要快于弦截法，例如：我们若取初始迭代点g=7 r/2,应 用Steffensen方 法(8)到例5.2中的方程/=4 s i n i,则 经6次迭代后，可得满足精度要求xk-ml 10 8 的数值解 7=2.47457678736983.5.2 Picard 迭代法本节将基于不动点原理给出非线性方程的Picard迭代法。5.2.1 标量方程的Picard迭代法将非线性标量方程情形等价地写成/=夕 0),(9)其中夕()为连续函数。该方程的求根问题在几何上可视为求曲线9=以7)与直线g=7的交点P*的横坐标7*。因此，我们可从这一几何观点出发来构造求解方程 的数值算法。图 5.2.Picard迭代重复上述迭代过程，可得点列外,（见 图5.2）,其横坐标满足公式叫+1 =3（4），A:=0,1,2,.（10）上迭代过程表明，若迭代函数夕及初始逼近值g 选择恰当，则点歹久只/逐 步 逼 近 广，即迭代序列以收敛于I*（见 图5.2中左图）；否则，迭代过程发散（见图5.2中右图）。由公式（10）确定的方法称为Picard迭代法。例5.3利用Picard迭代法迭代法求方程x3 2/+/2=0在x=1.8附近的近似根跳，并 使 其 满 足 血-跳_i|10-8解其方程可写成下列等价形式x 步2/2 /+2,由此得Picard迭代公式一6k+2,k 0,1,2,*,7。=：1.8.根据该迭代公式，经31次迭代后，可得其方程在x=1.8附近且满足精度要求的近似根3 1 =1.99999998890913 c若取其方程的另一等价形式x=x3+2 +2,则有迭代格式力 卜+1 +2,ZL=0,1,2,XQ 1.8.根据该迭代公式，经10次迭代后，计算机发生溢出，无法获得满足题设精度要求的近似根。5.2.2非线性方程组的Picard迭代法Picard迭代法(10)也可推广应用于非线性方程组情形(1),其有如下Picard送代格表Xk+i=D(Xk),k=0,1,2,(11)这里Xk=侬化暧)，需)T,方程组X =(X),X e C Rn(12)与方程组等价。为简单见，以下我们简称G-可微为可微。综合中值公式(1.27)与定理1.24可直接获得迭代格式(11)的收敛性判据。定理5.1 设中：。U IT 一肽”于闭凸集。QU。上连续可微，且满足(O o)u O o,q=s u p|d(X)|1,(13)XeR)这 里|是 中 给 定 向 量 范 数，则方程组(12)有唯一不动点X*D o,迭代格式(11)自任意迭代初值X。6 R)出发收敛于该不动点，且有先验误差估计|X*-Xk|S 4|X 1 Xo|,k=l,2 (14)i-q及后验误差估计|x*Xk|s 4|Xk _Xk_J|,k=L2 .(15)q由此可导出标量方程(9)的迭代格式(10)的收敛性准则。推 论 5.1若迭代格式(10)的迭代函数夕(/)在&b上连续可微，且满足(1)a b,/x E a,5 ,(2)q:=s u p|/O)|1,则方程有唯一不动点1*G Q,“，迭代格式(10)自任意迭代初值 g G a.b出发收敛于该不动点，且有误差估计式(14)及(15)成立。此外,对于Picar d迭代法，我们还有如下局部收敛性准则。定理5.2设中：。u肽，-在其不动点X*R D处及其某邻域内连续可微，且存在国，中的某向量范数I I-|使得|，(X*)|1,(16)则存在 X*的某邻域 S(X*,6)：=X e Rn：|X X*|W 6 ,使得迭代格式(11)自任意迭代初值X。G S(X*,b)出发收敛于不动点 X*。证 既然小于其不动点x*e。处及其某邻域内连续可微，且有|*(X*)|1,则存在 0 q 0 使得快(X)|一1,X-S(X*,b).从而，对一切x e S(x*,b),由中值公式有网X)X*|=|(X)(X*)|sup 惟,(x*+e(x x*)|x x*|001 q1|X X*|j即中(X)e S(X*,8).故应用定理5.1即得本结论。上述结论不仅给出了 Picard迭代法的收敛性准则，而且给出了计算程序的终止准则，即若要求计算精度满足|X*-X,|时，我们只需近似用条件Xk-X i|来控制迭代过程的终止，然后以当前迭代值Xk作为满足精度要求的数值解。例5.4利用Picard迭代法求解非线性方程组x=0.5 cos x 0.5 sin（*）y=0.5 sin x+0.5 cos g,并要求所得数值解X k满足精度要求：|Xk-Xk_i|2 10-8o解 记X =Xy(X)=I乙cos x s in ys in x+cos y 则其Jacobi矩阵(X)1=2s in/cos a;cosys in y 显然该矩阵的元素均在R2上的连续，因此中(X)在R2上连续可微。若取。()=X e肽2：l i xg si ,由于V X e股2有1 1 _I I (X)，2=-/(cosx s iny)2+(s inx+cos?/)2=+s in(y)2 v 2则中(。()c DQO又V X e肽2及人e R有det A/-仲 X)T(6(X)i 4A s in2 x cos2 y s in x cos x s in y cos y4 s in x cos x s in y cos y 4A cos2 x s in2 y=4A2 2A +|cos2(x y)其多项式的根为无,2=当 用 也 因 此=夜 (用)丁 仲 )=Jm ax J 1-疝 j 标从 而 有 su p 忡(X)|2 1。故据定理5.1,Picard迭代格式XsDoXk+1=中(Xk),k=0,1,2,；Xo=(0,0)T收敛于其方程组的解。利用该迭代格式计算方程组(1 7),经27次迭代后得满足精度要求的数值解X27=(0.22905926899460,0.54189671768053)T.5.2.3加速收敛技术下面，我们改进Picard迭代方法，使其迭代收敛过程得以加速运行。设Jacobi阵中(X)在X*的邻域内变化不大，其近似矩阵为P,则由Taylor展并式有X*孰X k)=(X*)-纵XQ=fo (Xk+WX*Xk)(X*Xk)g X P(X*Xk),即X*x(1 P)T 怦(Xk)P X J由此得加速迭代公式Xk+1=(1 P)T馋(Xk)-P Xk,k=0,1,.(18)例5.5 取迭代初始向量X0=(0,0户，利用加速迭代公式求解方程组(17),并要求所得数值解X k满足精度要求：|Xa-Xi|2 10一8。解 记P (X o)=；二.乙xy1cos x-s in yX =，(X)2 s in x+cos y0-11 0,k=0,1,则得加速迭代格式Xk+1=g；2 1卜(Xk)_ 1利用该迭代格式计算方程组(17),经1 1次迭代后得满足精度要求的数值解Xn=(0.22905926 7 518 49,0.5418 96 7 16 3496 7).例5.5表明，加速迭代格式大大加快了 Picard迭代方法的收敛速度，但该方法的缺陷是需要确定矩阵P,而这在实际问题的计算中是非常困难的。对于标量非线性方程（9）,我们有一种克服该困难的方法，即所谓的Aitken加速迭代法。记私+1=以冲），金k+i=3（通+1）,则由Taylor展开定理近似地有/*-xk+i x,(/*-xk/*-ifc+i 0(i*-由上二式消去未知常数p得*人力xk+l(诵+1 诵：+1)2金 k+1 2/卜+1+6k故得Aitken加速迭代公式!1 4+1=夕,以+1=夕(交k+1),Xk+i=诵+i 1八十 八十 xk+1-2xk+i+xk(19)例 5.6 利用Aitken加速迭代公式求方程X3 2/+/2=0在 x=1.8附近的近似根3,并使其满足xk-跳_i|10丑解取 _ _ _ _ _ _ _ _ _ _xo 1.8,#(1)=/272 7+2,应 用 Aitken加 速 迭 代 公 式(19)求解方程，经 3 次迭代后，可 得 其 方 程 在 x=L 8 附 近 且 满 足 精 度 要 求 的 近 似 根 g=2.00000000000000o若另取夕=X3+2+2,应 用 Aitken加速迭代公式(19)求解方程，经 6 次迭代后，可得其方程在x=1.8附近且满足精度要求的近似根方=2。与例5.3比较，例 5.6表明：Aitken迭代法大大加速了 Picard迭代公式的收敛速度，而且在某些情形下，Aitken迭代方法可将一个发散的Picard迭代过程改造成一个收敛的迭代过程。5.3 Newton 迭代法5.3.1标量方程的N ewt on迭代法为求解标量方程（2）,我们考虑如下迭代方法：在方程的解的隔离区间 a,b上选取适当迭代初值g,过曲线沙=f 的点（如产（3）引切线沙=/（g）+/（必）（一。o）与X轴相交于点进一步，过曲线g=与X轴相交于点_ /（g）X 1=XQ r；/（3）于 的点（g JQi）引切线沙=（的）+/丽（一i）_ f M）X2=Xi -r；/）如此循环往复，可得一列逼近X*的点力0,力1,-，跳,-（见图5.3）,其一般表达式为_ /（利-1）1 Q工k k 1 f/（.y,k,1,2,（20）该公式即称为Newton迭代公式，其相应求解方法称为Newton迭代法或切线法。由图5.3可见，当迭代初值3选择恰当时，点列，也,，将逼近标量方程 的精确解x*o图 5.3.Newton迭代法几何图示.例5.7用 Newton迭代法求方程x=4 sin x 在区间 f,7r 内的根,并使其满足|跳xk-i 10-8o解 记于=7 一4sin/,并取g=枭 应用Newton迭代法(20)于方程 x=4 sin x,经 7次迭代后，用得满足精度要求的数值解叼=2.47457678736983。将例5.1与例5.7相比较，显见Newton迭代法的收敛速度要远远快于二分法的收敛速度。5.3.2非线性方程组的Newton迭代法标量方程的Newton迭代法也可推广到非线性方程组情形。事实上，若记方程组(1)的精确解为X*,其第k 次逼近值为Xk,则由一阶Taylor展开式有9(X*)x F(XQ+尸(X 4X*-Xk).由F(X*)=0 及上式得X*X k 尸(Xk)5(X。.故得方程组情形的Newton迭代公式Xa+1=Xk F(Xk)Tp(x。，k=o,1 2 (21)为节省每步计算Jacobi阵FXk)的时间，我们有时也可采用如下简化的Newton迭代公式Xk+i=Xk k(X0)T尸(Xk),k=0,1,2,.(22)1X =,9(X)=-y，2例 5.8 取 迭 代 初 始 向 量X()=(0,0)T,利 用N ew ton迭代 法 求 解 方 程 组(1 7),并 要 求 所 得 数 值 解X k满足精度要解记2x cos/+sin y2y sin x cos y 应 用Newton迭代公式(21)计算方程组(1 7),经5次迭代后得满足精度要求的数值解=(0.22905926720286,0.54189671602062)7.由例5.8可见，Newton迭代法的收敛速度要快于Picard迭代方法及其加速迭代格式的收敛速度。但Newton迭代法有一个不足之处，即其要求选择一个恰当的迭代初值X。，方可保证迭代过程快速收敛。为尽量避免因初值X。选择不当而导致迭代过程缓慢收敛或发散，我们在Newton迭代公式中引入一个下山因子:0 A 1,从而产生下述Newton下山法Xk+1=Xk-可尸(Xk)iX k),k=0,1,2,.(23)在实际计算中，入可依次取为1/，亮 ，再 ，并且采用双精度控制：|Xk+i X疝2 勺 或|尸(Xk)，2 0,I T 上的范数|H I及 X。G DQ使得(Q)F X)-Ff(Y)aX-Y,VX、”D。、(b)I尸(X)T|S G,VX e Po,7：=IIF(Xo)T9(Xo)|,7 7 :=Q/?7 I，(d)S(X0,27):=X:|X-X o|0,使得|e f c+i|1时，称其为超线性收敛的；当=2时，称其为平方收敛的。在 适 当 条 件 下，我们可证明：弦 截 法(6)是超线性收敛的，St effen s en方法(8)是平方收敛的。对于一般Picar d迭代法，我们则有如下收敛阶判据。定 理5.4设中(X)在其不动点X*处及其某邻域内p阶连续可微，且叫X*)=0(=1,2,，0 1),(X*)#0,则Picar d迭代法(11)是夕阶收敛的。证明 由于中在其不动点X*e。处及其某邻域内连续可微，且但*)=0,则据定理5.2,存 在X*的某邻域S(X*,b),使得迭代格式(11)自任意迭代初值Xo G S(X*出发收敛于不动点X*。又据Taylor展开定理及已知条件，存 在X*的某邻域S(X*花)z)S(X*,b),使得所有Xk wS(X*,使 且有(Xk)中(X*)+吧 户(Xk-X*)。2=1+fo+0(xk-x*)(Xk X*yde中(X*)+fo (X*+3(Xk-X*)(Xk-X*)。曲.因此，画X。(X*)|s：SLIP I 即)(X)|Xk X*俨，k=0,l,P XS(X*,S)即|Xk+i X*|W sup|(X)|X,P，XeS(X*黝故迭代法(11)是月阶收敛的。对于Newton迭代法(21),我们有定理5.5在定理5.3的条件下，Newton迭代法(21)平方收敛。证明 据 定理5.3可 知Newton迭 代 法(21)收敛于方程组 在S(Xo,2)内的唯一解X*,且由微分性质(1.29)及定理5.3中条件(a)-(b),我们有|X*Xk+i|=x*Xk+F，(Xk)T X k)|F(X划T|0(X*)-F(Xk)-尸(Xk)(X*_ xfc)|拉例X*-X .因此，Newton迭代法(21)是平方收敛的。对于简化的Newton迭代法(2 2),我们有定 理5.6在定理5.3的条件下，简化的Newton迭代法(22)线性收敛。证 明 记W(X)=X 尸(Xo)5(X),X e S(Xo,27),则方程组X=W(X)与方程组F(X)=0 等价，且简化的Newton迭代法(22)可视为如下Picard迭代法Xk+i=W(Xk),A:=0,1,-.由微分性质(1.29)及 定 理 5.3 中 条 件(a)-(c),对 于任 意X,F 6S(Xo,2)有|W(X)W(Y)|=II尸(Xo)TF(y)-F(X)-尸(Xo)(Y-X)|F(x0)-1|F(y)-F(X)-尸(Xo)(y-X)|-I I X -Xo|+|y-Xo|X-y|2 叫|y -x|,|W(X)Xoll=II 尸(Xo)TF(X)F(Xo)尸(Xo)(X-Xo)+F(Xo)|II尸(Xo)T|M(X)尸(Xo)尸(Xo)(X-Xo)|+|尸(Xo)-F(Xo)|S|a/3|X Xo+7 (2a/57+1)7 2).因此，W是闭凸集S(Xo,2)上的压缩映射，且 W(S(Xo,27)US(X0,27)O故 据 定 理 1.23及 定 理 1.24可 知：在 X。的邻域S(X0,27)内方程组X =W(X)(即F(X)=0)存在唯一解X*,且简化的Newton迭代法(22)收敛于X*,其满足|X*-XM 1|2切|X*-X k|,即线性收敛。