1、数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法1第九章第九章 常微分方程初值问题常微分方程初值问题数值解法数值解法/*Numerical Method for Ordinary Differential Equations*/许多实际问题的数学模型是微分方程或微分许多实际问题的数学模型是微分方程或微分方程的初值问题方程的初值问题,如物体运动如物体运动,电路震荡电路震荡,化学反化学反映及生物群体的变化等。映及生物群体的变化等。能用解析方法求出精确解的微分方程为数不能用解析方法求出精确解的微分方程为数不多多,而且有的方程即使有解析解而且有的方程即使有解析解,也
2、可能由于解的也可能由于解的表达式非常复杂而不易计算表达式非常复杂而不易计算,因此有必要研究微因此有必要研究微分方程的数值解法。分方程的数值解法。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法2常微分方程常微分方程中介绍的微分方程主要有中介绍的微分方程主要有:(1)(1)变量可分离的方程变量可分离的方程 (2)(2)一阶线性微分方程一阶线性微分方程(贝努利方程贝努利方程)(3)(3)可降阶的一类高阶方程可降阶的一类高阶方程 (4)(4)二阶常系数齐次微分方程二阶常系数齐次微分方程 (5)(5)二阶常系数非齐次微分方程二阶常系数非齐次微分方程 (6)(6)全微
3、分方程全微分方程本章主要介绍本章主要介绍一阶一阶常微分方程初值问题的常微分方程初值问题的数值解法数值解法。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法3图形解 xyo简单的微分方程复杂、大型的微分方程解析解 y=f(x)数值解(xi,yi)欧拉方法改进欧拉方法 梯形法龙格-库塔法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法4 初值问题及其初值问题及其数值解数值解的概念的概念1 引言引言常用的一些常用的一些解析解法解析解法：常数常数变易变易法、法、Lapalace变换等变换等分离变量分离变量法、变量法、变量代换代换、
4、一阶一阶常微分方程初值问题：常微分方程初值问题：00(,);()dyf x y axbdxy xy ()数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法5Th(,)f x y对于初值问题对于初值问题 ，如果，如果 在下列区域内连续：在下列区域内连续：()（解的（解的存在唯一存在唯一性）性）;Daxb y 且关于且关于 满足满足Lipschitz条件，即存在常数条件，即存在常数 ，使，使y0L 1212|(,)(,)|;,f x yf x yL yyx yD 则初值问题则初值问题 存在唯一解，且解是存在唯一解，且解是连续可微连续可微的。的。()Def所谓所谓数
5、值解数值解是指：在解的是指：在解的存在区间存在区间上取一系列点上取一系列点012.nxxxx 逐个求出逐个求出 的近似值的近似值1 2 3(,.)iy i ()iy x0;ixxih 等距等距节点：节点：:h步长步长数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法6 初值问题初值问题 的解析解及其数值解的的解析解及其数值解的几何几何意义：意义：()oxy初值问题初值问题 的解表示过点的解表示过点 的一条的一条曲线曲线()00(,)xynx(,)nnxy),(00yx()yy x 0 x),(00yx 2x),(22yx 1x),(11yx 初值问题初值问题
6、的数值解表示一组的数值解表示一组离散点列离散点列()(,)iixy可用可用拟合拟合方法求该组数据方法求该组数据 的的近似曲线近似曲线(,)iixy积分积分曲线曲线数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法7建立微分方程数值解法建立微分方程数值解法,首先要将微分方程离散首先要将微分方程离散化化.一般采用以下几种方法一般采用以下几种方法:(1)(1)用差商近似导数用差商近似导数)(,)(:11nnnnxyyxyy令令进进一一步步dxdyhyynn111,nnnnnny xy xxxyxdydx 建立数值解法的常用方法建立数值解法的常用方法数值分析数值分析第
7、第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法8(2)(2)用数值积分近似积分用数值积分近似积分 1)(,()()(1nnxxnndxxyxfxyxy即即)(,)(:11nnnnxyyxyy 令令进进一一步步),()(,(11nnxxnnyxfhdxxyxfyynn 实际上是矩形法实际上是矩形法宽宽高高11(,)(0,1,)nnnnxxxxdydxf x y dxndx数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法9(3)(3)用用TaylorTaylor多项式近似并可估计误差多项式近似并可估计误差)(!2)()()(!2)()()()
8、(221nnnnnnnxyhxhyxyyhxhyxyhxyxy )(,)(:11nnnnxyyxyy 令令进进一一步步),(1nnnnyxhfyy )(max2)(211xyhyxybxann 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法102 简单的数值方法简单的数值方法 Euler方法的基本原理方法的基本原理212()()()()!nnnnh yy xy xhy x 将将 在点在点 处进行处进行Taylor展开展开1()ny x nx略去略去 项：项：2h然后用然后用 代替代替 ，即得，即得ny()ny x1()()(,()nnnny xy xhf
9、xy x 10 1 21(,),nnnnyyhf xynN 称上述公式为向前称上述公式为向前Euler 公式。公式。一、一、Euler方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法112112()()()()!nnnnh yy xy xhy x 若将若将 在点在点 处进行处进行Taylor展开展开()ny x1nx 略去略去 项：项：2h然后用然后用 代替代替 ，即得，即得ny()ny x111()()(,()nnnny xy xhf xy x 1110 1 21(,),nnnnyyhf xynN 称上述公式为向后称上述公式为向后Euler 公式。
10、公式。向后向后Euler 公式为公式为隐式隐式格式，需要利用格式，需要利用迭代法迭代法求解求解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法12 000000,0000 ,x y yy xf x yyxf x yx ydydx由出发取解曲线的切线（存在！），则斜率由于及已知，必有切线方程：00000000,()(,)dyyyxxyxxf xydxyx Euler方法的方法的几何意义几何意义数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法13）（：，可可由由切切线线算算出出，则则为为等等步步长长0001101,yxhfyy
11、yhxxh 2 1 0 ,)(n1n1n，）（：点点的的值值在在逐逐步步计计算算出出 nyxhfyyxxyynn注意：这是“折线法”而非“切线法”除第一个点是曲线切线外，其他点不是！Y=y(x)ab1x2x数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法14解：解：向前向前Euler公式：公式：1(,)nnnnyyhf xy 101()dyxydxy 例例1：分别利用向前和向后分别利用向前和向后Euler方法方法求解初值问题求解初值问题的的数值数值0 1.h 解解（取步长为（取步长为 ）10 10 90 1.nnnyxy 向后向后Euler公式：公式：111
12、(,)nnnnyyhf xy 1110 10 11 1(.).nnnyxy 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法15具体计算结果具体计算结果:数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法16数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法17)1()(),(:00bxayxyyxfdxdy 式式为为已已知知初初值值问问题题的的一一般般形形利用数值积分将微分方程离散化利用数值积分将微分方程离散化得梯形公式得梯形公式:),(),(211 nnnnyxfyxfh 1)(,()(
13、)(1nnxxnndxxyxfxyxy),(),(2111 nnnnnnyxfyxfhyy解决方法：有的可化为显格式，但有的不行解决方法：有的可化为显格式，但有的不行梯形方法为隐式算法梯形方法为隐式算法二、二、改进的改进的Euler方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法18 ,kyxfyxfhyyyxhfyynEulerknnnnnknnnnn210 ,2,210)(11)1(1)0(1，对对法法结结合合，形形成成迭迭代代算算法法与与梯形公式比梯形公式比EulerEuler法精度高一些法精度高一些,但计算量较大但计算量较大 实际计算中只迭代
14、一次，这样建立的预测实际计算中只迭代一次，这样建立的预测校正系统称作改进的校正系统称作改进的Euler公式。公式。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法19 00121211)(),(),()(2:yxyhKyhxfKyxfKKKhyynnnnnn作作等等价价变变换换1111(,)(,)(,)2nnnnnnnnnnyyhf xyhyyf xyf xy预测校正数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法20 1)0()10(2)1.0(yxyxyyh步步长长求求解解初初值值问问题题例例解解(,)2f x yyx
15、 y 00121211)(),(),()(2:yxyhKyhxfKyxfKKKhyynnnnnn作作等等价价变变换换xy21:解解的的表表达达式式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法21Euler近似解近似解精确解精确解0 1.0.1 1.10.2 1.191820.3 1.277440.4 1.358210.5 1.435130.6 1.50897y0-1y0.1-1.09545y0.2-1.18322y0.3-1.26491y0.4-1.34164y0.5-1.41421y0.6-1.483240 1.0.1 1.097740.2 1.187
16、570.3 1.271290.4 1.350130.5 1.424990.6 1.49657改进改进Euler近似解近似解结果比较结果比较数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法22三、三、常微分方程数值解法的常微分方程数值解法的稳定性稳定性Def设一个数值方法以定步长设一个数值方法以定步长 求解实验方程求解实验方程h得到线性差分方程的解得到线性差分方程的解 。当时。当时 ，若，若 ，ny0,Re()yy n 0ny 则称该方法对步长为则称该方法对步长为绝对稳定绝对稳定的；否则称为不稳定的。的；否则称为不稳定的。Def将数值方法应用于实验方程，若对一
17、切将数值方法应用于实验方程，若对一切都是绝对稳定的，则称区域都是绝对稳定的，则称区域 为该方法的为该方法的绝对稳定域绝对稳定域。hC 上述定义表明，若数值方法可使任何一步产生的误差在后面的上述定义表明，若数值方法可使任何一步产生的误差在后面的计算中都能逐步削弱，则该方法为绝对稳定。计算中都能逐步削弱，则该方法为绝对稳定。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法2311()()nnnnyyhyhy 例如例如，对于向前，对于向前Euler法：法：1(,)nnnnyyhf xy 将其应用于将其应用于实验实验方程方程1120hh 11h 当当 时，误差将逐步
18、减弱，故此时方法稳定。时，误差将逐步减弱，故此时方法稳定。向前向前Euler法法绝对稳定域绝对稳定域：2 0(,)当当 因有误差变为因有误差变为 时，则有时，则有nyny 11()nnyhy 111()()nnnnyyhyy 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法24四、四、单步方法单步方法的局部误差的局部误差和和阶阶单步法单步法的一般形式的一般形式11(,)nnnnnyyhxyyh 1(,)nnnnyyhxy h 隐式隐式单步法单步法(,)x y h 通常称通常称 为为增量增量函数函数显式显式单步法单步法Def称称 为某方法在点为某方法在点 的的
19、整体截断整体截断误差误差()nnney xy nx设设 是是准确准确的，用某种方法计算的，用某种方法计算 时产生的截时产生的截nyDef1ny 断误差，称为该方法的断误差，称为该方法的局部局部截断截断误差，即误差，即11()()(,(),)nnnnnTy xy xhxy xh (单步法：在计算单步法：在计算yn+1 时只利用时只利用yn)数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法25其中其中 为自然数，则称该方法是为自然数，则称该方法是 阶的或具有阶的或具有 阶精度。阶精度。Def如果给定方法的如果给定方法的局部截断局部截断误差为误差为如果一个如果一个
20、 阶单步方法的阶单步方法的局部局部截断截断误差为误差为11()pnTO h pppp121(,()()ppnnnTg xy xhO h 则称则称 为该方法的局部截断误差的为该方法的局部截断误差的主项主项。1(,()pnng xy xh 如向前如向前Euler方法的方法的局部局部截断截断误差误差11()()(,(),)nnnnnTy xy xhxy xh 22()()()!nnnhhy xyxhy x 2323()()!nnhhyxyx 2()O h 一阶一阶方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法26 Euler方法的误差分析方法的误差分析对
21、初值问题中的微分方程两端在区间对初值问题中的微分方程两端在区间 上积分上积分,x xh()()(,()x hxy xhy xf s y s ds 如果用如果用左矩形公式左矩形公式计算右端积分，并令计算右端积分，并令nxx()()(,()nnnnny xhy xhf xy xR 1(,()(,()nnxnnnxRf s y s dshf xy x 其中其中上述等式中如果用上述等式中如果用 代替代替 ，即得向前，即得向前Euler格式。格式。()ny xny其其局部截断局部截断误差为误差为1(,()(,()nnxnnnxRf s y s dshf xy x 数值分析数值分析第第9 9章章 常微分方
22、程初值问题数值解法常微分方程初值问题数值解法27设设 关于关于 和和 均满足均满足Lipschitz条件，即条件，即(,)f x yxy1212(,)(,)f x yf x yL yy 1212(,)(,)f xyf xyK xx 和和1(,()(,()nnxnnnxRf x y xf xy xdx 11(,()(,()(,()(,()nnnnxnxxnnnxf x y xf xy xdxf xy xf xy xdx 11()()nnnnxxnnxxKxx dxLy xy xdx 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法28121()2nnxnxKh
23、Lyxx dx 2()2hKLMR max()max(,()a x ba x bMy xf x y x 其中其中而而整体截断整体截断误差为误差为111()nnney xy ()(,()(,)nnnnnnny xhf xy xRyhf xy (,()(,)nnnnnneh f xy xf xyR 11()nnehL eR 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法2911()nnehL eR 21111()()nhLeRhL 10011()()nniihLeRhL 0111()()nnRhLehLhL 注意到注意到1hLehL 1()nhLnehL1()
24、()b a LnehL 01()()b a Lb a LReeehL 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法30Th(,)f x y对于初值问题对于初值问题 ，如果，如果 关于关于 满足满足()（向前向前Euler方法的方法的整体截断整体截断误差误差）Lipschitz条件，条件，为对应的为对应的Lipschitz常数，当常数，当,K L0h时，向前时，向前Euler方法的数值解方法的数值解 一致收敛一致收敛,x y于初值问题于初值问题 的精确解，且的精确解，且整体截断整体截断误差满足估计式误差满足估计式 ny()012()()()b a Lb a
25、 LnhKeeeMeL 如果如果 ，Euler方法的方法的整体整体截断截断误差为误差为 ()O h00()yy x 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法31一、一、Runge-Kutta方法的基本思想方法的基本思想3 龙格龙格-库塔（库塔（Runge-Kutta）方法）方法显式显式单步法的一般形式：单步法的一般形式：1(,)nnnnyyhxyf h R-K方法是利用一些点的线性组合方法是利用一些点的线性组合构造构造增量函数增量函数，使得相应方法的使得相应方法的局部局部截断误差的截断误差的阶数阶数尽可能尽可能高高。二阶二阶Runge-Kutta方法
26、方法12221(,)(,)(,(,)x y f hc f x yc f xh yhf x y 确定参数确定参数 ，使得，使得12221,c c 与与 在点在点 的的Taylor展开式有尽可能多的展开式有尽可能多的相同项相同项。()(,)y xhx y f h ()y xh x数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法3212221()(,)()(,)(,(,)y xhx y f hy xh c f x yc f xh yhf x y 122221()(,)(,)(,)(,)(,)()xyy xh c f x ycf x yhfx yhf x y fx
27、yO h 232()()()()()hy xhy xhy xyxO h22()(,)(,)(,)(,)()xyhy xh f x yf x yf x y f x yO h 比较两式的比较两式的相同项相同项得得22112c 121cc2212c 方程组有无穷多解方程组有无穷多解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法33若取其一组解若取其一组解122211112,cc21(,)nnKf xh yhK 1122()nnhyyKK 1()nnKf x y 则得到则得到改进改进的的Euler公式（公式（二阶二阶方法）方法）若取其另一组解若取其另一组解1222
28、1132443,cc 则得到则得到二阶二阶的的Heun（休恩）公式。（休恩）公式。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法34二、显式二、显式Runge-Kutta方法及其稳定性方法及其稳定性和和fDef设设 是一个正整数，代表使用函数值是一个正整数，代表使用函数值 的个数，的个数，m是一些特定是一些特定的权因子（均为的权因子（均为2 31 21,(,;,iijim ji 1 2(,)ic im 实数），则称下列方法（公式）实数），则称下列方法（公式）111()nnmmyyh c Kc K 为初值问题为初值问题 的的m级级显式显式RungeKutta
29、公式，公式，其中其中()11(,)mmnmnmiiiKf xh yhK 1(),nnKf x y 22211(,),nnKf xh yhK 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法35类似前面的处理方法，可以得到类似前面的处理方法，可以得到四级四级方法：方法：m=45()O h局部截断局部截断误差误差最常用的一种最常用的一种四阶四阶方法：方法：经典显式经典显式Runge-Kutta公式公式11234226()nnhyykkkk 1(,)nnkf xy 2122(,)nnhhkf xyk3222(,)nnhhkf xyk43(,)nnkf xh yhk
30、数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法36解：解：101()dyxydxy 01,x 例例2：用用经典的四阶经典的四阶Runge-Kutta方法求解下列初值问题方法求解下列初值问题 。0 1.h 经典的四阶经典的四阶Runge-Kutta公式：公式：11234226()nnhyykkkk 11;nnkxy21122;nnhhkxyk 431.nnkxhyhk 32122;nnhhkxyk 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法37数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问
31、题数值解法38数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法39注注：对于显式对于显式N级级R-K方法，最多只能得到方法，最多只能得到N阶方法。阶方法。上述方法的缺陷：计算非常复杂。上述方法的缺陷：计算非常复杂。可通过积分方法确定参数。可通过积分方法确定参数。例例2：确定如下三级三阶显式确定如下三级三阶显式Runge-Kutta公式中的参数：公式中的参数：1112233()nnyyh c Kc Kc K 12221133311322(,)(,)(,)nnnnnnKf xyKf xh yhKKf xh yhKhK 解：解：对微分方程对微分方程 两边积分得两边
32、积分得(,)yf x y 11()()(,()nnxnnxy xy xf x y x dx 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法40采用采用Simpson公式计算上式右端积分项公式计算上式右端积分项15114626()()(,()(,()(,()()nnnnnnnnhhhy xy xf xy xf xy xf xy xO h 可设参数可设参数132142663,ccc 则有则有112346()nnhyyKKK 12221133311322(,)(,)(,)nnnnnnKf xyKf xh yhKKf xh yhKhK 选择选择剩余参数剩余参数，使
33、得，使得411()()nny xyO h 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法41112211115331132223226()(,()(,)(,()(,)()nnnnnnnnnnhhhy xyf xy xhf xh yhKf xy xf xh yhKhKO h 221122(,()(,)nnnnhhf xy xf xh yhK 23221228(,()()()()(,()nnnnnnnhhhf xy xy xy xO hf xh yhy x 23822()(,()()nynnnhhhy xfxyy xO h 取取2211122,数值分析数值分析
34、第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法42113311322(,()(,)nnnnf xy xf xh yhKhK (,()nnf xh y xh 3313222(,()(,()nnnnnnhhf xh yhy xhf xyy x 232(,()()()()nnnnhf xh y xhy xyxO h 33132322(,()(,)(,)()(,)()nnnnnxnnnxnnf xh yhy xhf x yhhf x yy x f x yO h 2332331322(,()()()()nnnnhf xh yhy xy xOh 数值分析数值分析第第9 9章章 常微分
35、方程初值问题数值解法常微分方程初值问题数值解法43取取3313211,113311322(,()(,)nnnnf xy xf xh yhKhK 23223332122()()(,()()()()nynnnnhy xfxh yhy xhy xO hO h 23321222()()(,()()nynnnhhhy xfxyy xO h 利用利用Taylor展开式展开式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法4411343221222()()()(,)()nnnynnny xyhdhhffxyfO hdx 112211115331132223226()(,(
36、)(,)(,()(,)()nnnnnnnnnnhhhy xyf xy xhf xh yhKf xy xf xh yhKhKO h 代入代入当当 时，时，322 411()()nny xyO h 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法45例例3：求求经典四阶的经典四阶的R-K方法的方法的绝对稳定域。绝对稳定域。解：解：11234226()nnhyykkkk 1(,)nnkf xy 211222(,)()nnnhhhkf xykyk 322222(,)()nnnhhhkf xykyk 433(,)()nnnkf xh yhkyhk 22()nhy 23
37、 224()nhhy 3nyhk 22334411234()!nnhhhyhy ny 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法4623411234()()()!hhhh 其其绝对稳定域绝对稳定域为为三、隐式三、隐式Runge-Kutta方法方法m级级隐式隐式RK方法的一般形式方法的一般形式11mnniiiyyhc K 1(,)mininijjjKf xh yhK 1 2 3,m 其中系数的确定方法同其中系数的确定方法同显式显式RK方法完全类似方法完全类似数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法47(1)(1
38、)一一级级二阶二阶的隐式的隐式中点中点方法：方法：11nnyyhk 1122(,)nnhkhkf xy (2)(2)二二级级四阶四阶的隐式的隐式R-K方法：方法：1212()nnh kkyy 112131326446(),()nnhkkf xh yhk221131326464(),()nnhkkf xh yhk N级隐式级隐式R-K法法可以达到可以达到2N阶阶缺陷：需要求解缺陷：需要求解非线性方程（组）非线性方程（组）数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法48一、一、k步线性步线性多步法多步法4 线性线性多步法多步法/*Linear Mutiste
39、p Method and Predictor-Corrector Format*/所谓的线性所谓的线性多步法多步法，指的是某一步解的公式不仅，指的是某一步解的公式不仅与前一步的值有关，而且与前面若干步解的值有关的与前一步的值有关，而且与前面若干步解的值有关的方法。方法。对初值问题对初值问题 两边两边积分积分得得0(,),()dyf x yaxbdxy ay ()()(,()xxy xy xf s y s ds 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法49将将 换为节点换为节点,x x 1,nnpxx 11()()(,()nn pxnnpxy xy x
40、f s y s ds 0(,)()()(,()kkjnjnjjf x yLxlx f xy x 取节点取节点 ，构造，构造 的的k+1个点的个点的Lagrange插值多项式：插值多项式：(,)f x y1,nnn kxxx 00 1 2(),kn ijinjn iijxxlxjkxx 多步多步显式显式公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法50110()()()(,()nn pkxnnpjnjnjxjy xy xlx f xy xds 0()(,()knpjnjnjjy xhf xy x 其中其中11()nn pxjjxlx dxh 记记 m
41、ax,rp k 若函数值若函数值 已知，则得已知，则得r+1步步显式显式方法方法1,n rn rnyyy 10(,)knnpjnjnjjyyhf xy 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法51如如 时，可得二步时，可得二步显式显式阿达姆斯（阿达姆斯（Adams）格式）格式01,pk 1110()()()(,()nnxnnjnjnjxjy xy xlx f xy xds 10()(,()njnjnjjy xhf xy x 其中其中111111112()nnnnxxnxxnnxxl x dxdxhhxx 1110011132()nnnnxxnxxnn
42、xxlx dxdxhhxx 11132(,)(,)nnnnnnhyyf xyf xy 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法52 Adams显式公式的显式公式的局部局部截断误差：截断误差：(,()()(,)()!kkkd fyxf x yLxdxk 由由Lagrange插值余项知插值余项知其中其中1()()()()nnn kxxxxxxx 111()(,)()nnxnnkxy xyf x yLx dx 11111()()()!nnkxkkxdyx dxO hkdx 1(,()()!nnkxkxd fyxdxdxk （第（第二二积分积分中值中值定理）
43、定理）k阶方法阶方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法5311()()(,()nn pxnnpxy xy xf s y s ds 110(,)()()(,()kkjnjnjjf x yLxlx f xy x 取节点取节点 ，构造，构造 的的k+1个点的个点的Lagrange插值多项式：插值多项式：(,)f x y11,nnnkxxx 10110 1 2(),knijinjniijxxlxjkxx 多步多步隐式隐式公式公式1110()()(,()knnpjnjnjjy xy xhf xy x 数值分析数值分析第第9 9章章 常微分方程初值问题数值
44、解法常微分方程初值问题数值解法54其中其中11()nn pxjjxlx dxh 记记 1max,rp k 则得到则得到r+1步步q+1阶的阶的隐式隐式方法方法1110(,)knnpjnjnjjyyhf xy 如如 时，可得二步时，可得二步隐式隐式阿达姆斯（阿达姆斯（Adams）格式）格式01,pk 10111(,)(,)nnnnnnyyhf xyf xy 1100111112()nnnnxxnxxnnxxlx dxdxhhxx 1112(,)(,)nnnnnnhyyf xyf xy梯形梯形公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法55 常用的一
45、种常用的一种预测预测-校正校正公式：公式：四阶四阶Adams预测预测-校正公式：校正公式：0111555924()(,)(,)nnnnnnhyyf xyf xy011191924()(,)(,)nnnnnnhyyf xyf xy 2233379(,)(,)nnnnf xyf xy 11225(,)(,)nnnnf xyf xy （显式显式）（隐式隐式）初始迭代值由初始迭代值由4阶阶R-K方法计算方法计算数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法56201()dyxydxyy 01,x 例例4：用用Adams预测预测-校正公式校正公式求解下列初值问题求解
46、下列初值问题 。0 1.h 解：解：Adams预测预测-校正公式：校正公式：0111122555924()()()nnnnnnnnxxhyyyyyy 0111012291924()()()()nnnnnnnnxxhyyyyyy 23232322379()()nnnnnnxxyyyy 121212225()()nnnnnnxxyyyy 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法57R-K方法方法Adams预预-校法校法 精确精确解解 0 11.0000000000 0.11.0954461.0954451153 0.21.1832171.18321595
47、66 0.31.2649121.2649110640 0.41.34164135711.3416407864 0.51.41421383341.4142135623 0.61.48323982421.4832396974 0.71.54919338041.5491933384 0.81.61245153641.6124515496 0.91.67331999931.6733200530 1.01.73205071981.7320508075ix数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法58 5 一阶方程组一阶方程组与与高阶高阶方程的数值解法方程的数值解
48、法一、一、一阶微分方程组初值问题的一阶微分方程组初值问题的一般形式一般形式1112221212(,)(,)(,)mmmmmdyfx yyydxdyfx yyydxdyfx yyydx 初始条件初始条件:1122()()()mmy ay aya 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法59写成写成向量向量的的形式形式:12()()()()my xy xy xyx 12()()()()my xy xy xyx 12(,)(,)(,)(,)mf x yf x yf x yfx y 12m ()(,)()y xf x yy a 数值分析数值分析第第9 9章章
49、 常微分方程初值问题数值解法常微分方程初值问题数值解法60 n=2对应的对应的Runge-Kutta公式公式(,)(,)dyf x y zdxdzg x y zdx 0000()()y xyz xz 11234226()nnhyykkkk 11234226()nnhzzLLLL 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法611(,)nnnkf xyz 211222(,)nnnhhhkf xyk zL 322222(,)nnnhhhkf xykzL433(,)nnnkf xh yhkzhL 1(,)nnnLg xyz 211222(,)nnnhhhLg
50、xyk zL 322222(,)nnnhhhLg xykzL433(,)nnnLg xh yhk zhL 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法62例例 考虑考虑LorenzLorenz模型：模型：112322332123()()()()()()()()()()()()x tbx tx t x tx tcx tcx tx tx t x trx tx t 其中参数b=8/3，c=10，r=28解：12122323()0()()()0()()()1()x tbx tx tx tccx tx tx trx t数值分析数值分析第第9 9章章 常微分方程初值
