chap常微分方程数值解实用PPT课件

1、 ( , )( , )( , )f x yyLipschitzLf x yf x yL yy只要函数适当光滑连续，且关于 满足条件，即存在常数 ，使得由常微分方程理论知，初值问题的解必存在且唯一。第1页/共119页 微分方程的数值解：设方程问题的解y(x)的存在区间是a,b，令a= x0 x1 xn =b,其中hk=xk+1-xk , 如是等距节点h=(b-a)/n , h称为步长。y(x)的解析表达式不容易得到或根本无法得到，我们用数值方法求得y(x)在每个节点xk上y(xk)的近似值，用yk表示，即yky(xk)，这样y0 , y1 ,.,yn称为微分方程的数值解。 第2页/共119页主要

2、问题v如何将微分方程离散化，并建立求其数值解的递推公式；v递推公式的局部截断误差，数值解与精确解的误差估计；v递推公式的稳定性与收敛性。第3页/共119页v用差商代替微商用差商代替微商v数值积分数值积分vTaylor展开展开微分方程离散化常用方法微分方程离散化常用方法第4页/共119页1,1A (, ()nnnnnnnny xy xdyf xy xydx xxx用差商代替微商11111 , , , 0, 1, 2,nnnnnnnnnnnnnnhyyfhxxyxyxyyxyyyxhfny用代替，则 第5页/共119页11B. : ( , ) (0,1,) nnnnxxxxdydxf x y dx

3、ndx用数值积分方法离散化1111,(), (),( , )(,) (0,1,(,)nnnnnnxnnxnnnnyyy xy xf x y dxhf xyyhxnyfy用代替对右端积分采用取左端点的矩形公式 则有第6页/共119页11111111,(), (),( , )(,) (0, 1, )(,)nnnnnnxnnxnnnnyyy xy xf x y dxhf xyyyhf xyn用代替对右端积分采用取右端点的矩形公式 则有第7页/共119页11111111,(), (),( , )(,)(,)2 (,)(,)2 (0, 1 ,)nnnnnnxnnnnxnnnnnnyyy xy xhf x

4、 y dxf xyf xyhyyf xyf xyn用代替对右端积分采用梯形公式 则有第8页/共119页22C ( ) ()()()()2 ()(, ()()2nnnnnnnnnxy xTaylorhy xhy xhy xyxhy xhf xy xyx在附近的展开：11 () () (,) 0, 1, 2, nnnnnnnhyy xy xyyhf xynTaylor取的线性部分，且得的近似值：展开法不仅可得到求数值解的公式，且容易估计截断误差。第9页/共119页1 解常微分方程初值问题的Euler方法Euler方法Euler方法的误差分析第10页/共119页v向前向前EulerEuler公式公式

5、（EulerEuler折线法或显格式）折线法或显格式）v向后Euler公式（后退Euler公式）v梯形公式（改进的Euler公式）vEuler预估校正格式一、Euler方法第11页/共119页0010, 0 1,1, nnnnny xyyyhfxyxxnhban,NhN1、向前Euler公式第12页/共119页 000000,0000 , , , , xy yy xdyf xyydx xf xyxy几何意义由出发取曲线的切线(存在!)，则 斜率由于及已知，必有切线方程。00000000, () (,)dyyyxxyxxf xyydx x由点斜式写出切线方程：第13页/共119页1011000

6、,hxxhyyyhf xy等步长 为 ，则，可由切线算出：（）11 , 0 1 2 nnnnnyy xxyyhf xyn逐步计算出（ ）在，点值 ：（）,， ，用分段的折线逼近逼近函数第14页/共119页2、向后（后退的）Euler 方法11 nnny xy xyxh用向后差商：11100, nnnnyyhfxyy xy (隐式算法)第15页/共119页 01111E , 0,1,2,1,nnnnkknnnulerhfhfknyyyxyyyx为避免解非线性方程，与法结合迭代法第16页/共119页3、梯形公式11 ,nnnnxy xy xfx y dtx由积分途径 11 ,nnnnyy xyy

7、x积分用梯形公式，令则得11100,2nnnnnnhyyfxyfxyy xy第17页/共119页(0)1(1)( )111012,2 01 2nnnnkknnnnnnEulernyyhf x yhyyf x yf xyk, , ,同样与法结合，形成迭代算法，对， ， ，第18页/共119页4、改进的尤拉公式梯形公式虽然提高了精度，但使算法复杂。而在实际计算中只迭代一次，这样建立的预测校正系统称作改进的尤拉公式。第19页/共119页1111 ( ,); ( ,)(,),2nnnnnnnnnnyyhf x yhyyf x yf xy预测校正11(,);(,);() / 2.pnnncnnpnpcy

8、yhfxyyyhfxyyyy第20页/共119页二、二、Euler方法的误差分析方法的误差分析11111) () nnnnnnnTy xyyyy xEulery局部截断误差在一步中产生的误差而非累积误差：其中是当（精确解！）时由法求出的值，即无误差！第21页/共119页 1121 () ()()(, () 2nnnnnnnnny xxTaylory xy xhy xhf xy xhyxx将在点展开： 1121111 ,(), 2nnnnnnnnnnnnnyyhf xyyy xhf xy xhTy xyyxx则第22页/共119页22212 max( ) , ( ) 2a x bnMyxy xh

9、TMO h 令充分光滑，则： 11 , , 1(1)nnnnnnnnnnnnnnnnnyyy xhf xy xyhf xyy xyh f xy xf xyLipschitzhL y xyhL e由条件第23页/共119页11000121211 1 N 0 111 1nnnNNNNNNNeThL eney xyeThL eThL ThLThLT对一切 成立,对取定，由， 则：第24页/共119页2） 总体方法误差总体方法误差 , ,nnnnnnf xyLipschitzf xy xf xyL y xy递推方法：从任意两相邻步的总体误差关系推出总体误差与步长的关系。由微分方程解的存在唯一性自然假定

10、 （ ， ）充分光滑，或满足条件：第25页/共119页1111111111111 1 , nnnnnnnnnnnnnnnnnnney xyney xyy xyyyTyyyyyy xh f x第步 的 总体 截 断 误 差 记 为则 对步：以 下 估 计其 中 ,ny x第26页/共119页221211 11 1NNNNNNTO heThL ThLThLT由局部截断误差，则11200 1 1NNNKkkkkOhLhLhT 211 11NhL O h O hhL第27页/共119页000 lim 1lim 1 NNxxNhhhLhLhLxxeh与步长无关常数第28页/共119页总体截断误差与局部截

11、断误差的关系是：1O h总体截断误差局部截断误差一般地，方法的总体截断误差阶越高，精度也越高。pO hp定义：一个方法的总体截断误差若为， 则称之为 阶方法。第29页/共119页误差分析表误差分析表EulerEuler方法方法局部截局部截断误差断误差总体截总体截断误差断误差迭代收敛迭代收敛条件条件向前向前EulerEuler方法方法O(h2)O(h)向后向后EulerEuler方法方法O(h2)O(h)0hL1梯形公式梯形公式O(h3)O(h2)0hL2(L为为Lip常数）常数）第30页/共119页向后Euler 方法收敛条件与截断误差 11111111111 , 01kkkknnnnnnkk

12、nnh ffhLhLyyyyxxyy收敛条件 21 01 nOO hhLhT局部截断误差，整体截断误差(当时)第31页/共119页梯形公式的收敛性 11111111111 ,2 2 012kkkknnnnnnkknnhffhLhLyyyyxxyy收敛条件 Euler梯形公式比法的局部与总体误差均高一阶，但每次迭代均多算一次函数值提高精度的计算代价。第32页/共119页2(01);(0)1.xyyxyy例：用尤拉公式和改进的尤拉公式解初值问题10.12().nnnnnhxEuleryyh yy解：取步长，公式为： 第33页/共119页112();2();1().2npnnnncnppnpcxyy

13、h yyxyyh yyyyy改进的尤拉公式为：第34页/共119页20.2(00.6); (0)1.hyyxyxy 例1：取步长，用欧拉法解初值问题122 , 0.2 0.80.2nnnnnnnnnnnEuleryyhfxyyyx yyx y解：格式为：第35页/共119页0122 1 0.2 0.8, 0.4 0.6144 0.6 0.461321yyyyyyy由计 算 得第36页/共119页0.283(12); (1)2.5hyyxy例2：取步长，用梯形解初值问题小数点后至少保留 位。11111 ,20.2 83832nnnnnnnnnnhyyfxyfxyyyyy解：梯形公式为： 第37页

14、/共119页 1012345716 131312, 1.22.307691.42.473371.62.562581.82.610622.02.63649nnyyyyyyyyyyyyyy故 由计 算 得 第38页/共119页0.23(12); (1)2.5hyxyxy例2 ：取步长，用梯形解初值问题小数点后至少保留 位。111111 ,2 0.1 33nnnnnnnnnnnnhyyfxyfxyyyx yxy解：梯形公式为：第39页/共119页 0111110011234511111110.30.1 332,2.6*,&,%,?nnnnkknnnnnnkyyx yyyx yxyyyyyyyyyy迭

15、代格式：第40页/共119页 021234522222,*,&,%,yyyyyy第41页/共119页20.2sin0 (1)1.1.21.45hyyyxyyy例3：取步长，用欧拉预校方法解初值问题计算及的近似值，小数点后至少保留 位。第42页/共119页111121212111,2 0.2sin0.1sin + sinnnnnnnnnnnnnnnnnnnnnnnnyyhfxyhyyfxyfxyyyyyxyyyyxyyx解：欧拉预校格式为：第43页/共119页 0112211,0.63171 1.20.7154880.476961.40.52611yyyyyyyy故 由计 算 得 第44页/共1

16、19页50 (0)1.0.10.210yyyhy例4：写出用反复迭代的欧拉预校法解初值问题的计算公式，并取步长，计算，要求迭代误差不超过。第45页/共119页 01111101111,2 0,1,2, 0,1,2, 0.90.950.05nnnnkknnnnnnnnkknnnyyhfxyhyyfxyfxyknfx yyyyyyy 解：欧拉预校的迭代格式为：取第46页/共119页 0012111341143751141101, 0.9,0.905,0.90475, 0.9047625,0.9047618756.251010, 0.10.904761875yyyyyyyyyyyy故 由计 算 得由

17、于 是 取第47页/共119页 0012232452254652252201, 0.814286,0.818809, 0.818583, 0.818595,0.8185941010, 0.20.818594yyyyyyyyyyyy故 由计 算 得由于 是 取第48页/共119页2.2.龙格龙格库塔方法库塔方法基本思想基本思想二阶二阶R-K方法方法三阶三阶R-K方法方法四阶四阶R-K方法方法变步长变步长R-K方法方法第49页/共119页1112() ()()()()()()2!( )( , ), ( )( , )( , ) ( , ),nnnppnnnnxypTaylory xyy xhhy x

18、hy xyxyxPy xf x yyxfx yfx y f x y若用 阶多项式近似函数有：其中。但由于公式中各阶偏导数计算复杂，不实用。一、基本思想第50页/共119页(0)(0)(0)(1)(1)(1)(2)(2)(2)( )(1); ;2,3,jjjjyffffyffxyffyffxyffyffjxy一般地有第51页/共119页111112121(,)11()22 (,)(,)nnnnnnnnnnEulerEulerEuleryyhKKf xyEuleryyhKKKf xyKf xh yhK如果将公式与改进公式写成下列形式：公式改进公式第52页/共119页11 ( , )()( , )(

19、 , )nnf x yy xyf x yf x y以上两组公式都使用函数在某些点上的值的线性组合来计算的近似值。Euler公式：每步计算一次的值，为一阶方法。改进Euler公式：需计算两次的值，二阶方法。第53页/共119页( , )(,)( )-nnnf x yxyTaylory xxTaylorR K于是可考虑用函数在若干点上的函数值的线性组合来构造近似公式，构造是要求近似公式在处的展开式与解在 处的展开式的前面几项重合，从而使近似公式达到所需要的阶数。即避免求偏导，又提高了方法的精度，此为方法的基本思想。第54页/共119页11111-(,)(,) (2,3,),pnniiinniini

20、nijjjiijiRKyyhc KKf xyKf xa h yhb Kipabc一般地，方法设近似公式为其中 ，都是参数.第55页/共119页,(,)( )iijinnnab cxyTaylory xxTaylor确定参数 ，的原则是:使近似公式在处的展开式与在 处的展开式的前面项尽可能多地重合。第56页/共119页二、二阶龙格库塔方法111221222112()(,)(,)nnnnnnpyyh c Kc KKf xyKf xa h yhb K当时，近似公式为 第57页/共119页112221123221( ,)( ,)(, ( ,)( ,) ( ,) ( ,)( ,) ( ,)( )nnnn

21、nnnnnnnnnnnxnnynnnnnx yTayloryyhc f x yc f xa h yhb f x yyh c f x yc f x ya hf x yhb f x y f x yO hy上式在处的展开式为12222321() ( ,)( ,) + ( ,) ( ,)( )nnxnnynnnncc f x y h c a f x yb f x y f x yhO h第58页/共119页123123()()()()()()2(,)(,)2 (,) (,)()nnnnnnnnnxnnynnnny xxTaylorhy xy xhy xy xO hhyf xy hfxyfxyf xyO

22、h在 处的展开式为第59页/共119页122 232 211 1/2 1/2 ( ),cccacbOh有 无 穷 多 组 解 ， 每 一 组 解 得 一近 似 公 式 ， 局 部 截 断 误 差 均 为这 些 方 法 统 称 二 阶 方 法 。43()()O hO h可以证明，无论这四个参数如何选择，都不能使局部截断误差达到，也即在计算两次函数值的情况下，局部截断误差的阶最高为。第60页/共119页122211121211,1,2()/2(,)(,)nnnnnnccabEuleryyh KKKf xyKf xh yhK取此为改进公式。近似公式为 第61页/共119页122211212110,1

23、,2( ,)(2,2)nnnnnnccabyyhKKf x yKf xhyhK取此为常用的二阶公式，称为中点公式 第62页/共119页三、三阶龙格库塔方法1123121312(4)6(,)(,)22(,2)nnnnnnnnRKhyyKKKKf xyhhKf xyKKf xh yhKhK常用的三阶公式为：第63页/共119页四、四阶龙格库塔方法112341213243 (22)6(,)(,)22(,)22(,)nnnnnnnnnnRKhyyKKKKKf xyhhKf xyKhhKf xyKKf xh yhK常用的四阶公式为： 第64页/共119页0.2,-83 ;(0)2.0.44 hR Kyy

24、yy例：设取步长写出用经典（标准的）四阶方法求解初值问题 的计算公式，计算的近似值，小数点后至少保留 位。12,8 - 3 , =0.2, 0.2,0.4fx yyhyyyy解：第65页/共119页112341122343(22);6,83;,5.62.1;22,6.322.37;22,4.2081.578.nnnnnnnnnnnnnnhyyKKKKKf xyyKhKfxyyKhKfxyyKf xh yKy由经典的四阶龙格库塔公式得第66页/共119页 10121.20160.549402,0.22.30040.42.4654nnyyyyyyyy由于第67页/共119页11,2,3,4,454

25、652pR KppR KpppR KR KTaylor两点说明：）当时， 公式的最高阶数恰好是 当时， 公式的最高阶数不是，如时仍为 ，时 公式的最高阶数为。） 方法的导出基于展开，故要求所求问题的解具有较高的光滑度。第68页/共119页 RKEulerRKRKEuler当解充分光滑时，四阶 方法确实优于改进法。对一般实际问题，四阶方法一般可达到精度要求。如果解的光滑性差，则用四阶 方法解的效果不如改进法。第69页/共119页五、变步长的龙格库塔方法( )1( )51115(2)1,(),2,2nhnhnnnnhnxhyy xychhxxhyc 以经典四阶龙格库塔公式为例。从节点 出发，以 为

26、步长求一近似值将步长折半，即取 为步长从 跨两步到，求一近似值每跨一步的截断误差是第70页/共119页5(2)11(2)11()11(2)(2)()1111()2,2()1.16()1().15hnnhnnhnnhhhnnnnhy xycy xyy xyy xyyy因此有由上两式 第71页/共119页 4、微分方程数值解的稳定性*111*11 nnnnnyyyy稳定性分析，对计算误差：其中是的近似计算值，误差积累会淹没真值？第72页/共119页 = , 0 12 nn kn knyyhykh定义：一种数值方法求解 试验方程其中是复常数。对给定的步长，若计算误差在计算， ， 时不产生增大的误差，

27、即，称对与这种方法是绝对稳定的。 hhh对 ， 的允许范围内是绝对稳定的，则称的全体为绝对稳定区域。第73页/共119页Euler法的绝对稳定区域1*1 nnnnnyyEuleryyhyyhy 的算法：计算值 1 nnnh误差方程：第74页/共119页11 11 nnhh 从而当是绝对稳定区域 Im(h)-2 -1 0 Re(h)hA绝对稳定区域越大， 可选大些，方法适应性越强。如果整个左半平面是绝对稳定区域称稳定的。第75页/共119页向后Euler 法的稳定性11 nnnyyEuleryyh y 对用向后法 ： 11 nnnh误差公式：仍受限制。要求稳定的。但收敛法是因此向后 , 10 h

28、hLAEuler第76页/共119页1 ()0 1nenRh只要则11/ 211 111nnhhh122122211112ehhRhh第77页/共119页梯形公式的稳定性11 2nnnyyhnyyyy 对用梯形公式）（1()0 1 A- nenRh当时，梯形公式是稳定的。第78页/共119页111122121222 2211()411()4nnnnnnhhheehhRhhR第79页/共119页R-K方法的绝对稳定区域2121233223443( , ) 11 ,22111 ()22411 ()24nnnnnnnyf x yyRKKh yKhyKyhhKhyKyhhhKhyKyhhhh 将代入公

29、式：第80页/共119页112342341 22611112624nnnyyKKKKyhhhh234111112624nnhhhh则234111 112624hhhh绝对稳定区域：第81页/共119页 2 1-3 -2 -1 0 -1 -2第82页/共119页11-r1R-K,nnnnnnyyyyyy单步法在计算时，只用到前一步的信息 。为提高精度，需重新计算多个点处的函数值，如方法，计算量较大。如何通过较多地利用前面的已知信息，如 ，来构造高精度的算法计算。基本思想4.线性多步法第83页/共119页11110111 (,), (,) ,( ,)(1, ,)00Taylornnn rnnn r

30、n rrrnin iin iiiiikkkyyyf xyf xyyyhfff x yknnn r 多步法中最常用的是线性多步法，它的计算公式中只出现， ，及的一次项，其一般形式为 其中均为常数，。若，显式；，隐式。构造线性多步公式常用展开和数值积分方法。第84页/共119页线性多步公式的导出1(),nnniiTaylorxTaylory xxTaylor 利用展开导出的基本方法是：将线性多步公式在 处进行展开，然后与在处的展开式相比较，要求它们前面的项重合，由此确定参数。第85页/共119页( )( )2( )1() (1,2,),( )( )()()2 ()()!kknnnnnnnnpppn

31、nnyyxky xxTayloryy xyy xxxxyxxO xxp记则在 处的展开为1011110111( ) ()nnnnnnry xyyyhfff以为例：设初值问题的解充分光滑，待定的两步公式为第86页/共119页231(4)(5)45(6)21111(),()( ,) (),()2!3! ()4!5!(,)()2! iiiiinnnnnnnnnnnnnnnnyy xy xf x yinyyyy xhyy hhhyyhhO hyff xyy xyy hh假设前 步计算结果都是准确的，即则有(4)(5)34(5) ()3!4!nnyyhhO h第87页/共119页1111(4)(5)23

32、4(5)(,)(,) () ()2!3!4!nnnnnnnnnnnnnff xyyff xyy xyyyyy hhhhO h第88页/共119页1(5)2561()()() 2!5!1 nnnnnnnpy xxTayloryyy xyy hhhO hp为使上式有 阶精度，只须使其与在 处的展开式的前项重合。211011101113(4)4111111(5)56111()()()2()()6222466()( )1202424nnnnnnnyyy hy hy hy hy hO h 将以上各公式代入并整理，得第89页/共119页01010111111111111112211116226111124

33、6624aaaaaa5,5,11,2,3,4)iippp 个参数只须 个条件。由推导知，如果选取参数，使其满足前个方程（，则近似公式为 阶公式。第90页/共119页11011111,0,02 ()2nnnnhyyff0如满足方程组前三个方程，故公式此为二阶公式。第91页/共119页0111011115(5)6110,1,343 (4)31 ( )90nnnnnnnhyyfffRh yOh又如：解上面方程组得,相应的线性二步四阶公式（Simpson公式)为其截断误差为由此可知，线性二步公式至多是四阶公式。第92页/共119页1010123()1( )( )1 (1,2,)nnnriirrkkii

34、iirxTaylory xxTaylorikikp一般地，线性多步公式中有个待定参数，如令其右端在 处的展开式与在处的展开式的前p+1项系数对应相等，可得方程组第93页/共119页1111(1)21 1( )(1)! (1)( ) ()prpniirpppiniphRippiyO h其解所对应的公式具有 阶精度，局部截断误差为显然，线性多步公式至多可达到2r+2阶精度。第94页/共119页二、常用的线性多步公式第95页/共119页1231010100123 (Adams)r=30,1()()1 (1,2,3 4)5559379=1,24242424riirrkkiiiiikik （一）阿达姆斯

35、公式取，并令由方程组，可解得，第96页/共119页1123153354(5)61115(5)6(5559379)24=0Adams1( )5( )()5!251()720nnnnnnniiniinhyyffffhRiiyO hh yO h相应的线性多步公式为因，此式称为显式公式，是四阶公式.局部截断误差为第97页/共119页12330101211125(5)610,91951 =1,24242424(9195)24Adams19()720nnnnnnnnhyyffffRh yO h如果令由方程组可解得，相应的线性多步公式为称其为四阶隐式公式，其局部截断误差为第98页/共119页利用数值积分方法

36、求线性多步公式111111()()( , ( )( ),( )( )( ),1nnnnxxnnxxnnn rnnn rry xy xf x y x dxF x dxxxxxxxF xrxF xr 基本思想是首先将初值问题化成等价的积分形式用过节点或的的 次插值多项式代替求积分即得阶的线性多步公式。第99页/共119页123330123303,( ) ( )( ) ()()()()()( )()() (0,1,2,3)nnnnin iinnnnin in injjj irx xxxF xL xl x F xxxxxxxxxl xxxxxi例如时，过节点的三次插值多项式为其中第100页/共119页

37、1111131301233231313233()()( )( ) ()()()()()6()()()()2()()()()2()()nnnnnnnnnnxxnnin ixxixnnnnxxnnnnxxnnnnxnny xy xL x dxl x dx F xxxxxxxF xdxhxxxxxxF xdxhxxxxxxF xdxhxxxF x1123123)()655 ()59 ()37 ()9 ()24nnxnnxnnnnxxxdxhhF xF xF xF x第101页/共119页1111233,(), (),(,)()(, () (,1,2,3),(5559379)24,nnnnkkkkkk

38、nnnnnnnnyyy xy xfxyF xf xy xkn nnnhyyffffAdamsxxAdams对上式用代替用代替则得这就是四阶显式公式。由于积分区间在插值区间外面，又称为四阶外插公式。第102页/共119页111(4)310(5)3031(5)35(5)10()()4!() ()4!,),( )251()( )4!720nnnnnnxxnnjxjxxnjxjnnxnnjxjFRxxdxyxxdxxxyRxxdxh y由插值余项公式可得其局部截断误差为由积分中值定理，存在（使得第103页/共119页11223111231,( ) ( )( ) ()()()()()( )()() (1

39、,0,1,2)nnnnin iinnnnin in in jjj ixx xxF xL xl x F xxxxxxxxxl xxxxxi 同样，如果过节点的三次插值多项式为其中代( )F x替求积分，第104页/共119页11125(5)12121 (9195)2419( )720,nnnnnnnnnnnnAdamshyyffffRh yxxxxAdamsAdams 即得四阶隐式公式其局部截断误差为由于积分区间在插值区间内，故隐式公式又称为内插公式第105页/共119页一阶常微分方程组与高阶方程 我们已介绍了一阶常微分方程初值问题的各种数值解法，对于一阶常微分方程组，可类似得到各种解法，而高阶

40、常微分方程可转化为一阶常微分方程组来求解。 一阶常微分方程组对于一阶常微分方程组的初值问题 0000( , , ), ()( , , ), ()yf x y z y xyzg x y z z xz（5.1） 可以把单个方程 中的f 和y看作向量来处理，这样就可把前面介绍的各种差分算法推广到求一阶方程组初值问题中来。 ( , )yf x y 第106页/共119页设 为节点上的近似解，则有改进的Euler格式为 0(1,2,3,);,iiixxih iy z1( ,)iiiiiyyhf x y z1( ,)iiiiizzhg x y z预报：1111( ,)(,)2iiiiiiiihyyf x

41、y zf xyz1111( ,)(,)2iiiiiiiihzzg x y zg xyz校正： （7.32） 又，相应的四阶龙格库塔格式（经典格式）为 1123411234(22)6(22)6iiiihyyKKKKhzzLLLL（7.33） 第107页/共119页112111221112312223122241334133(,)(,)(,)22(,)22(,)22(,)22(,)(,)iiiiiiiiiiiiiiiiiiiiiiiiKfxyzLg xyzhhKfxyKzLhhLg xyKzLhhKfxyKzLhhLg xyKzLKfxyhKzhLLg xyhKzhL式中 （7.34） 第108页

42、/共119页 把节点xi上的yi和zi值代入式(7.34), 依次算出 , 然 后 把 它 们 代 入 式(7.33), 算出节点xi+1上的yi+1 和zi+1值。 对于具有三个或三个以上方程的方程组的初值问题,也可用类似方法处理,只是更复杂一些而已。此外,多步方法也同样可以应用于求解方程组初值问题。 11223344,K L KL KL KL例 用改进的Euler法求解初值问题 (0)1(0)2yxyzyxyzzz 00.2x取步长h=0.1，保留六位小数。 解: 改进的Euler法公式为),(1iiiiizyxhfyy),(1iiiiizyxhgzz预报： ),(),(21111iiii

43、iiiizyxfzyxfhyy),(),(21111iiiiiiiizyxgzyxghzz校正： 将 及h=0.1代入上式,得 zyxzyxgzxyzyxfiiiiii),(,),(第109页/共119页110.1()0.1iiiiiiiiiiyyx yzxyzzz11111110.05 ()()0.05iiiiiiiiiiiiiiiiyyx yzxyzxyxyzzzz由初值 ,计算得 00(0)1,(0)2yyzz110.8000002.050000yz11(0.1)0.801500(0.1)2.046951yyzz220.6048202.090992yz22(0.2)0.604659(0.2)2.088216yyzz第110页/共119页高阶方程组 高阶微分方程(或方程组)的初值问题,原则上都可以归结为一阶方程组来求解。例如,有二阶微分方程的初值问题 0000( , ,)(),()yf x y yy xyy xy在引入新的变量 后,即化为一阶方程组初值问题:zy0000( , , ), (), ()zf x y zyz y xy z xy 式（7.36）为一个