1、 8 常微分方程的数值解法常微分方程的数值解法（Numerical Solution of Ordinary Differential Equations）n 本章主要内容本章主要内容n 8.1 基本概念与基本求解途径基本概念与基本求解途径 n 8.2 欧拉方法与局部截断误差欧拉方法与局部截断误差 n 8.3 龙格龙格库塔法库塔法n 8.4 单步法的收敛性与稳定性单步法的收敛性与稳定性n 8.5 线性多步法线性多步法 n 8.6 一阶常微分方程组与高阶常微分方程的数值解法一阶常微分方程组与高阶常微分方程的数值解法n 8.7 边值问题的差分方法简介边值问题的差分方法简介 n 重点：重点：欧拉法与
2、龙格欧拉法与龙格库塔法、预测库塔法、预测校正技术校正技术n 难点：难点：各种方法的收敛性与稳定性各种方法的收敛性与稳定性8.1 基本概念与基本求解途径基本概念与基本求解途径 常常微分方程的应用背景微分方程的应用背景 00)()(qtqRCtqUdtdq 桥梁扰度问题设一根长为L的两端固定的矩形截面梁，q是均匀负荷强度，E是弹性模量，S是断点作用力，I(x)是惯性矩。则描述梁的桡度y(x)的数学模型是 电容充电问题 设q是电容器上的带电量，C为电容，R为电阻，U为电源的电动势。描述电容器充电过程的数学模型是 牛顿力学定律设s是物体的位移，f是物体的受力，M是物体的质量，则是物体的速度，a是物体的
3、加速度，牛顿定律可描述为000022)(,)(vtsstsMfadtsd一阶常微分方程一阶常微分方程的初值问题的初值问题二阶常微分方程二阶常微分方程的初值问题的初值问题二阶常微分方程的二阶常微分方程的两点边值问题两点边值问题0)()0()()(2)()(22LyyLxxEIqxxyxEISdxyd8.1 基本概念与基本求解途径基本概念与基本求解途径 常常微分方程的数学模型微分方程的数学模型 一阶常微分方程初值问题的数学模型的一般形式00)(),()(yxyyxfxy二阶常微分方程初值问题的数学模型的一般形式两点边值问题（第一类边界条件）的数学模型的一般形式 0000)(,)()()()()()
4、(yxyyxyxfxqxyxpxy)(,)(),()(byayyyxfxy 数学方法及存在问题数学方法及存在问题 只有少数简单的微分方程能够用初等方法求得它们的解，多数情形只能利用近似方法求解。在常微分方程教材中的级数解法级数解法、逐步逼近法逐步逼近法等都是近似方法，这些方法可以给出解的近似表达式，统称为近似解析近似解析法法。但实际中遇到的常微分方程问题往往不能给出解析解，或者解析解不变于使用。如下面两个问题及解析解分别为0)0(21yxydxdy1)0(yydxdyxxtxeydteey022解析解要用插值法或要用插值法或数值积分法求值数值积分法求值要查指数表要查指数表8.1 基本概念与基本
5、求解途径基本概念与基本求解途径8.1 基本概念与基本求解途径基本概念与基本求解途径 微分方程的数值求解途径的数学模型微分方程的数值求解途径的数学模型 数值方法的基本思想是离散化，即把求解区间a,b分成n等分（区间离散化），同时对微分算子离散化，求函数y(x)在n+1个离散的等距节点上的近似值，而不必求y(x)的解析表达式。初初值问题的数值解法，就是寻求解函数值问题的数值解法，就是寻求解函数y(x)在在a,b的一系列等距的的一系列等距的离散节点离散节点 a=x0 x1n)的绝对误差均不超过的绝对误差均不超过 ，则称该单步，则称该单步法是法是稳定的稳定的。例例8-7 考察显欧拉方法和隐欧拉方法的绝
6、对稳定性。考察显欧拉方法和隐欧拉方法的绝对稳定性。例例8-8 考查二阶龙格考查二阶龙格-库塔公式的绝对稳定性库塔公式的绝对稳定性。例例8-9 用四阶经典用四阶经典R-K公式求解初值问题公式求解初值问题(h分别取分别取0.1和和0.2)1)0(1,020yxydxdy 表表8-3 常用单步法的绝对稳定区间常用单步法的绝对稳定区间在计算在计算yn+1时，选择前面已算出的时，选择前面已算出的yn-r,yn-r+1,yn-1,yn多个值的线性组多个值的线性组合作为初值，相应近似斜率合作为初值，相应近似斜率yn=f(xn,yn)的线性组合构造增量表达式，的线性组合构造增量表达式，以获得较高的精度。以获得
7、较高的精度。8.5 线性多步法线性多步法 线性多步法的基本思想线性多步法的基本思想 线性线性r+1步法计算公式的一般形式步法计算公式的一般形式8.5.1 线性多步法公式的构造线性多步法公式的构造101(8.14)rrnknkknkkkya yhb yb-1=0时，显式公式时，显式公式b-10时，隐式公式时，隐式公式 基于泰勒级数展开的构造方法基于泰勒级数展开的构造方法01111(8.19)()()1(1,2,)rkkrrjjkkkkakajkbjp)()()()1()(1)!1(2)1(11111pnprkrkkpkppnhOxybkpakphR p+1个方程个方程,2r+1个未知数的线个未知
8、数的线性方程组，其解解即可构造出性方程组，其解解即可构造出步步p阶精度的线性多步公式阶精度的线性多步公式 局部截断误差局部截断误差8.5 线性多步法线性多步法 8.5.1 线性多步法公式的构造线性多步法公式的构造 基于基于泰勒级数展开泰勒级数展开的构造方法的构造方法例例8-10 构造形如构造形如 的线性两步隐式公式的线性两步隐式公式)(110111101nnnnnnfbfbfbhyayay例例8-11 构造形如构造形如 的线性两步隐式公式的线性两步隐式公式)()(11011nnnnnfbfbhyyay8.5 线性多步法线性多步法 8.5.1 线性多步法公式的构造线性多步法公式的构造 基于基于数
9、值积分数值积分的构造方法的构造方法将初值问题（将初值问题（8-1）式的微分方程两端同时从）式的微分方程两端同时从xn-r到到xn+1积分，得积分，得1)(,()()(1nrnxxrnndxxyxfxyxy 对式中的被积函数对式中的被积函数f(x,y(x)用用Lagrang插值函数近似，积分后可得插值函数近似，积分后可得到不同的线性到不同的线性r+1步公式。在区间步公式。在区间xn-r,xn+1选取的插值节点越多，选取的插值节点越多，公式的精度越高，插值节点若包括公式的精度越高，插值节点若包括xn+1则为则为隐公式隐公式，否则为，否则为显公式显公式。课后作业：课后作业：使用泰勒级数法使用泰勒级数
10、法构造构造Milne公式公式8.5 线性多步法线性多步法 8.5.2 亚当姆斯亚当姆斯(Adams)公式公式 Adams显式公式显式公式（外插公式）（外插公式）r步步r阶阶Adams显式公式的一般形式为显式公式的一般形式为二阶二阶Adams显公式显公式三阶三阶Adams显公式显公式四阶四阶Adams显公式显公式101rkknknnybhyy Adams显式公式既可以通过数值积分来构造，也可以用显式公式既可以通过数值积分来构造，也可以用Taylor展开展开式来构造。现通过数值积分来构造可得式来构造。现通过数值积分来构造可得10101),(rjjnjnnrjjnjnjnnnfAyyxfAyy1)(
11、nnxxjnjndxxlA其中)3(211nnnnffhyy)51623(12211nnnnnfffhyy)9375955(243211nnnnnnffffhyy表表8-5显式显式Adams公式系数表公式系数表8.5 线性多步法线性多步法 8.5.2 亚当姆斯亚当姆斯(Adams)公式公式 Adams隐式公式隐式公式（内插公式）（内插公式）r步步r+1阶阶Adams隐式公式的一般形式为隐式公式的一般形式为211rkknknnybhyy Adams隐式公式通常隐式公式通常Taylor展开式来构造。展开式来构造。3步步4阶阶Adams隐式公式隐式公式表表8-6隐式隐式Adams公式系数表公式系数表
12、)2()()()()()()(2101111hxybhxybxybhxybhxyhxyyxyRnnnnnnnnn)()()3462241()()22261()()221()()1(54)4(211321122112101hOhxybbbhxybbbhxybbbhxybbbbnnnn 为达到四阶精度，令为达到四阶精度，令h的的14次幂的系数为零，解之得次幂的系数为零，解之得241,245,2419,832101bbbb)5199(242111nnnnnnffffhyy从而有从而有3步步4阶阶Adams隐式公式隐式公式8.5 线性多步法线性多步法 8.5.3 线性多步法预测线性多步法预测校正公式校
13、正公式 不论是单步法还是多步法，隐公式相比显公式而言其稳定性好，但在不论是单步法还是多步法，隐公式相比显公式而言其稳定性好，但在计算时会遇到两个问题：计算时会遇到两个问题：n 隐公式如何能方便的计算；隐公式如何能方便的计算；n 实际计算步长取多大。实际计算步长取多大。一般地，不难验证，如果预测公式（显公式）是一般地，不难验证，如果预测公式（显公式）是p阶或阶或p+1阶精度，校阶精度，校正公式（隐公式）是正公式（隐公式）是p+1阶精度，则用预测公式提供初值，校正公式迭代阶精度，则用预测公式提供初值，校正公式迭代一次的效果也能达到一次的效果也能达到p+1阶精度，再迭代下去的话，效果就不明显了。阶精
14、度，再迭代下去的话，效果就不明显了。预预测测校正技术校正技术既保证了计算精度，又使隐式计算显式化，克服了隐式公式既保证了计算精度，又使隐式计算显式化，克服了隐式公式需要反复迭代计算的困难。通过对预测需要反复迭代计算的困难。通过对预测校正公式使用校正公式使用外推原理，外推原理，得到误得到误差估计式，用来调整计算步长，使其达到精度要求。差估计式，用来调整计算步长，使其达到精度要求。4阶阶Adams显公式显公式预测预测，4阶阶Adams隐公式隐公式校正校正),(9),(37),(59),(55(243322111nnnnnnnnnnyxfyxfyxfyxfhyy),(),(5),(19),(9(24
15、2211111nnnnnnnnnnyxfyxfyxfyxfhyy8.5 线性多步法线性多步法 8.5.3 线性多步法预测线性多步法预测校正公式校正公式 Milne(4阶阶)公式公式预测预测，Hamming(4阶阶)公式公式校正校正),(2),(),(2(34221131nnnnnnnnyxfyxfyxfhyy),(),(2),(83)9(81111121nnnnnnnnnyxfyxfyxfhyyy 中点公式公式中点公式公式预测预测，梯形公式，梯形公式校正校正112(,)nnnnyyhf xy111(,)(,)2nnnnnnhyyf xyf xy由局部截断误差可知由局部截断误差可知 进而得到进而
16、得到3阶公式阶公式)(51),(),(2),(2111111111nnnnnnnnnnnnnnyyyyyxfyxfhyyyxhfyy改进：校正：预测：注意：注意：隐式公式的计算需要隐式公式的计算需要用单步法（如龙格用单步法（如龙格-库塔法）库塔法）若干个初值才能启动计算若干个初值才能启动计算。)()(12)()()(3)(4311143111hOxyhyxyThOxyhyxyTnnnnnnnn )(414)(4111hOTTxynnn 改进公式是一个改进公式是一个事后误事后误差估计式差估计式。通常用下式。通常用下式来选择合适的步长。来选择合适的步长。)(5111nnyy8.6 一阶常微分方程组
17、和高阶常微分方程简介一阶常微分方程组和高阶常微分方程简介8.6.1 一阶常微分方程组一阶常微分方程组 一阶常微分方程组初值问题的一般形式一阶常微分方程组初值问题的一般形式 一阶常微分方程组初值问题的解法一阶常微分方程组初值问题的解法bxacaycaycayyyyxfxyyyyxfxyyyyxfxynnnnnnn)(,)(,)(),()(),()(),()(22112121222111TnTnTnccccyxfyxfyxfyxfxyxyxyxy),(),(,),(),(),()(,),(),()(212121cayyxfxy)(),()(定义定义向量向量得到得到 类似于一阶常微分方程初值问题的数
18、值方法，但需注意求解公式也是用向量形式表示。例例8-12 二阶亚当姆斯预测二阶亚当姆斯预测校正公式求解初值问题方程组校正公式求解初值问题方程组8.6.2 高阶常微分方程高阶常微分方程 高阶常微分方程初值问题的一般形式高阶常微分方程初值问题的一般形式bxacaycaycayyyyxfxymmmm1)1(10)1()()(,)(,)(),()()1(21,mmyyyyyymjcayyyxfyyyyyyyjjmmmm,2,1)(),(1113221引入变量：转化为一阶常微分方程组的初值问题 高阶常微分方程初值问题的解法高阶常微分方程初值问题的解法求解一阶常微分方程组的初值问题例例 用欧拉法求解2阶常
19、微分方程的初值问题。)10(1)0(,0)0(2 xyyyyxy8.7 边值问题的数值解简介边值问题的数值解简介8.7.1 二阶常微分方程边值问题的二阶常微分方程边值问题的打靶法打靶法 二二 阶常微分方程边值问题的一般形式阶常微分方程边值问题的一般形式 打靶法的基本思想打靶法的基本思想)(,)(),(),(byaybaxyyxfy先猜测一个初始斜率 y(a)=s，通过解初值问题 sayaayyyxfy)()(),(找出s*使得(s*)=，即把问题转化为求方程(s)=0 的根。yx0aby x()斜率斜率=s0()s0斜率斜率=s1()s1每计算一个每计算一个(s s)都必须解一个初值问题都必须
20、解一个初值问题y(b)=(s)8.7.2 二阶常微分方程边值问题的二阶常微分方程边值问题的有限差分法有限差分法 有限差分法的基本思想有限差分法的基本思想将求解区间a,b 等分为N 份，取节点 xi=a+ih (i=0,N)，在每一个节点处将 y 和 y 离散化。有限差分法公式推导有限差分法公式推导泰勒展开泰勒展开)()()(2)(22hOhhxyxyhxy )(2)()()(2hOhhxyhxyxy NiiiiiiiyyNihyyyxfhyyy,1,.,1)2,(2011211)(12)()()()()()4(2 yhhhhxyxyhxyhxyxy 例例8-13 经验、总结与提高经验、总结与提
21、高 本章讨论求解常微分方程的一些常用数值解法。欧拉方法和龙格库塔方法都是将微分方程离散化为差分方程求解，是一种步进式的方法。欧拉方法由于精度低，在实际中很少使用，但因其公式简单、直观，对其他方法具有启示作用。四阶龙格库塔公式由于精度高、易于编程，因此得到了广泛的应用，但其缺点是对函数f(x,y)的光滑性要求较高，且每一步计算需要计算四次函数值，计算量较大。数值方法的收敛性和稳定性是判别方法好坏的两个重要指标。一般只要增量函数关于y满足李普希兹条件，单步发都是收敛的。隐式单步法是绝对稳定的，显式单步法有自己的稳定区域，在选取步长时，应格外注意。采用多步法求解常微分方程初值问题时，预测校正技术可以把显式和隐式公式结合起来，使隐式公式可以通过显示计算，同时又能保证隐公式的稳定性。边值问题比初值问题复杂得多，它有时有若干个解，有时甚至无解。打靶法是将边值问题转化为初值问题的方法之一，还有其他方法。有限差分法和有限元法是求解偏微分方程的两种最主要的方法。基于变分原理的有限元法，以较少的单元达到令人满意的结果。很容易将两点边值问题的变分原理和有限元计算过程推广到2维和3维偏微分方程（组）定解问题。
