第8章 常微分方程数值解法
本章主要内容：
1．欧拉法、改进欧拉法. 2．龙格-库塔法。
3．单步法的收敛性与稳定性。
重点、难点
一、微分方程的数值解法
在工程技术或自然科学中，我们会遇到的许多微分方程的问题，而我们只能对其中具有较简单形式的微分方程才能够求出它们的精确解。对于大量的微分方程问题我们需要考虑求它们的满足一定精度要求的近似解的方法，称为微分方程的数值解法。本章我们主要
?dy?dx?f(x,y)讨论常微分方程初值问题?的数值解法。
?y(x)?y0?0 数值解法的基本思想是：在常微分方程初值问题解的存在区间[a,b]内，取n+1个节点a=x0＜x1＜…＜xN=b （其中差hn= xn –xn-1称为步长，一般取h为常数，即等步长），在这些节点上把常微分方程的初值问题离散化为差分方程的相应问题，再求出这些点的上的差分方程值作为相应的微分方程的近似值（满足精度要求）。
二、欧拉法与改进欧拉法
欧拉法与改进欧拉法是用数值积分方法对微分方程进行离散化的一种方法。
n?1 将常微分方程y??f(x,y)变为y(xn?1)?y(xn)??f(t,y(t))dtxnx???
1．欧拉法（欧拉折线法）
欧拉法是求解常微分方程初值问题的一种最简单的数值解法。
欧拉法的基本思想：用左矩阵公式计算（＊）式右端积分，则得欧拉法的计算公式为：yn?1?yn?hf(xn,yn)(n?0,1,...,N?1)欧拉法局部截断误差
h?b?a NRn?1h2?y??(?n?1)2xn??n?1?xn?1或简记为O（h2）。
我们在计算时应注意欧拉法是一阶方法，计算误差较大。
欧拉法的几何意义：过点A0（x0，y0），A1（x1，y1），…，A n（x n，y n ），斜率分别为f（x0，y0），f（x1，y1），…，f（x n，y n）所连接的一条折线，所以欧拉法亦称为欧拉折线法。
例1用欧拉法解初值问题
?dy?dx??2xy??y(0)?1?(0?x?1)
在x＝0 (0.2) 1处的近似解。（计算过程保留4位小数）。
【思路】 用欧拉法求解常微分方程的初值问题时，首先熟练掌握欧拉公式的一般形式， 根据具体题目写出找出欧拉公式的迭代式，并根据初始条件和所给步长进行迭代求解。 解 ∵ f(x，y)＝－2xy ，h＝0.2，
欧拉公式为：
yn?1?yn?hf(xn,yn)?yn?0.2(?2xnyn)?(1?0.4xn)yn(n?0,1,2,3,4,5)
列表计算如下： n 0 1 2 3 4 5
xn 0 0.2 0.4 0.6 0.8 1 yn 1 1 0.92 0.7728 0.5873 0.3994 y(xn) 1 0.9608 0.8521 0.6977 0.5273 0.3679
y(xn)-yn 0 -0.0392 -0.0679 -0.0751 -0.06 -0.0315 2．改进欧拉法
改进欧拉法比欧拉法的计算准确，是对欧拉法的改进。改进欧拉法的基本思想：用梯形公式计算（＊）式右端积分，则得改进欧拉法的计算公式为：
yn?1?yn?h?f(xn,yn)?f(xn?1,yn?1)?(n?0,1,...,N?1)2h?b?aN
利用改进欧拉法计算常微分方程初值问题时，我们应注意此公式为隐式表达式，需要对它进行迭代求解。计算时可以采用一次迭代和多次迭代，因此，就有改进欧拉法预估-校正法公式和反复迭代的改进欧拉法预估-校正法公式。 改进欧拉法预估-校正法公式：
?y?n0??1?yn?hf(xn,yn)? ?h?0??yn?1?yn?f(xn,yn)?f(xn?1,yn?1)(n?0,1,...,N?1)2???
反复迭代的改进欧拉法预估-校正法公式：
?y?n0??1?yn?hf(xn,yn)? ?h?1??m??y?nm?y?f(x,y)?f(x,ynnnn?1n?1)(n?0,?,N?1,m?0,1,...,)2??1??
改进欧拉法的局部截断误差Rn?1h3?y???(?n?1)xn??n?1?xn?1或简记为O（h3）。12从局部截断误差的形式看，改进欧拉法是二阶方法，因此，它比欧拉法更精确。 例2用预估-校正法求初值问题
2??y???y?xy???y(0)?1(0?x?1)
在x=0（0.2）1的解。
【思路】掌握预估-校正法的计算公式，根据已知条件迭代求解。 解 步长h=0.2，将f(x,y)??y?xy2代入预估-校正公式，整理得
(0)2?y?0.8y?0.2xynnn?n?1 ?2(0)(0)2??yn?1?0.9yn?0.1xnyn?0.1(yn?1?xn?1(yn?1)列表计算如下：
n 0 1 2 3 4 5 x n 0 0.2 0.4 0.6 0.8 1 y n 1 0.80720 0.63690 0.49048 0.37780 0.29103 y(x n) 1 0.80463 0.63145 0.48918 0.37720 0.29100
(m)(m?1)?10?3 例3用改进欧拉法求解例1的初值问题,要求yn?yn。
【思路】掌握改进欧拉法的计算公式，根据已知条件迭代求解，并检验迭代解是否满
足精度要求，若满足则确定此解为常微分方程在某点的近似解。 解 将f(x,y)??2xy代入改进欧拉法的计算公式得：
(0)?yn?1?yn?hf(xn,yn)?yn?0.2(?2xnyn)?(1?0.4xn)yn? ?h??m?1??m?1??y?nm?y?f(x,y)?f(x,y)?y?0.2xy?xynnnn?1n?1nnnn?1n?1)2??1????
n 列表计算如下：
xn 0 0.2 0.4 0.6 0.8 1 ?0? yn?1? yn?2? yn?3?yn?yny(xn) y(xn)-yn 0 1 2 3 4 5
1 1 0.8846 0.7181 0.5337 0.3622 1 0.96 0.8523 0.7003 0.5325 0.3750 1 0.9616 0.8549 0.7025 0.5327 0.3725 1 0.9615 0.8549 0.7022 0.5327 0.3730 1 0.9608 0.8521 0.6977 0.5273 0.3679 0 -0.0007 -0.0028 -0.0045 -0.0054 -0.0051 三、龙格-库塔法 1．龙格-库塔法
龙格-库塔法具有精度高、收敛、稳定，不需要计算高阶导数等优点，是求解微分方程初值问题的一组著名的显示单步方法，广泛应用于求解常微分方程的初值问题。本章我们介绍了二、三、四阶龙格-库塔法。 龙格-库塔法的基本思想： 在计算初值问题?y(xn?1)?y(xn)?y??f(x,y)的数值解时，考虑均差，则由微分中值
h?y(x0)?y0定理可得y(xn?1)?y(xn)?hy?(xn??h)(0???1)， 由初值问题可得公式为：
y(xk?1)?y(xk)?hf(xk??h,y(xk??h))上式中f(xk??h,y(xk??h))称为区间
上的平均斜率。如果给平均斜率一种计算方法，就可得到计算y（x n+1）的近似值yn+1的公式。
如果仅取xn处的斜率值f(xn,yn)作为平均斜率的近似值，则得到的yn?1的公式为欧拉公式；
如果取xn,xn?1处的斜率值f(xn,yn)，f(xn?1,yn?1)的平均值作为平均斜率的近似值，则得到的yn?1的公式改进欧拉公式。 ㈠ 二阶龙格-库塔法的公式和局部截断误差：
?yn?1?yn?h(?1k1??2k2)???k1?f(xn,yn)???k2?f(xn?1,yn?hk1)
在上式中选择不同的参数，会得到不同的二阶龙格-库塔法公式，所以二阶龙格-库塔法公式不唯一。二阶龙格-库塔法公式的局部截断误差为Ο（h3）。
常见的二阶龙格-库塔法公式有以下两种 改进欧拉法迭代公式