数值分析课件

第5章常微分方程数值解法5.1.2局部误差和方法的阶5.1.1Euler方法及其有关的方法第5章常微分方程数值解法5.1.2局部误差和方法的阶5.1Euler方法5.1.1Euler方法及其有关的方法考虑一阶常微分方程初值的问题：设f（x,y）是连续函数，对y满足Lipschitz条件，这样初值问题的解是存在唯一的，而且连续依赖于初始条件。为了求得离散点上的函数值，将微分方程的连续问题进行离散化。一般是引入点列{}，这里为步长，经常考虑定长的情形，即。记为初始问题的问题准确解在处的值，用均差近似代替的导数得5.1Euler方法5.1.1Euler方法及其有关令为的近似值，将上面两个近似写成等式，整理后得（5.1.2）（5.1.3）从处的初值开始，按（5.1.2）可逐步计算以后各点上的值。称（5.1.2）式为显式Euler。由于（5.1.3）式的右端隐含有待求函数值，不能逐步显式计算，称（5.1.3）式为隐式Euler公式或后退Euler公式。如果将（5.1.2）和（5.1.3）两式作算术平均，就得梯形公式。令为的近似值，将上面两个梯形公式也是隐式公式。以上公式都是由去计算，故称它们为单步法。例5.1取h=0.1，用Euler方法、隐式Euler方法和梯形方法解

解本题有如果用Euler方法，由（5.1.2）并代入h=0.1得同理，用隐式Euler方法有（5.1.4）梯形公式也是隐式公式。以上公式都是由去计算用梯形公式有三种方法及准确解的数值结果如表5-1所示。从表中看到，在处，Euler方法和隐式Euler方法的误差分别是和，而梯形方法的误差却是。

在例5.1中，由于f（x，y）对y是线性的，所以对隐式公式也可以方便地计算。但是，当f（x，y）是y的非线性函数时，如，其隐式Euler公式为。显然，它是的非线性方程，可以选择非线性方程求根的迭代求解。以梯形公式为例，可用显式Euler公式提供迭代初值，用公式用梯形公式有三种方法及准确解表5-1Euler方法隐式Euler方法梯形法准确解011110.11.0000001.0090911.0047621.0048370.21.0100001.0264461.0185491.0187310.31.0290001.0513151.0406331.0408180.41.0561001.0830131.0700961.0703200.51.0904901.1209211.1062781.106531表5-1Euler方法反复迭式，直到其中，步长h成为迭代参数，它需要满足一定的条件，才能收敛。若将（5.1.4）式减去该迭代公式，得假设f（x，y）关于y满足Lipschiz条件，则有反复迭式，直到其中，步长h成为迭代参数，它需要满足一定的条件这里，L是Lipschiz常数。当hL/2<1即h<2/L时，迭代序列收敛。

对于隐式公式，通常采用估计-校正技术，即先用显式公式计算，得到预估值，然后以预估值作为隐式公式的迭代初值，用隐式公式迭代一次得到校正值，称为预估-校正技术。例如，用显式Euler公式作预估，用梯形公式作校正，即称该公式为改进的Euler公式。它显然等价于显式公式为，（5.1.6）这里，L是Lipschiz常数。当hL/2<1即h<2/L时也可以表示为下列平均化的形式例5.2取h=0.1，用改进的Euler方法解解按（5.1.5），改进的Euler方法解也可以表示为下列平均化的形式例5.2取h=0.1，用改进由得计算结果如表5-2。该初值问题的准确解为。表5-20.10.20.30.40.50.60.70.81.09591.18411.26621.34341.41641.48601.55251.61531.09541.18321.26491.34161.41421.48321.54921.6165由得5.1.2局部误差和方法的阶

初值问题（5.1.1）的单步法可以写成如下统一形式（5.1.7）其中与有关。若中不含，则方法是显式的，否则是隐式的，所以一般显式单步法表示为（5.1.8）例如，Euler方法中，有对于不同的方法，计算值与准确解的误差各不相同。所以有必要讨论方法的截断误差。我们称为某一方法在点的整体截断误差。显然，不单与这步的计算有关，它与以前各步的计算也有关，所以误差被称为整体的。分析和估计整体截断误差是复杂的。为此，我们假设处的没有误差，即，考虑从到这一步的误差，这就是如下的局部误差的概念。5.1.2局部误差和方法的阶（5.1.7）其中定义5.1

设是初值问题（5.1.1）的准确解，则称为单步法（5.1.7）的局部截断误差。

定义5.2

如果给定方法的局部截断误差，其中为整数，则称该方法是p阶的，或具有p阶精度。若一个p阶单步法的局部截断误差为则称其第一个非零项为该方法的局部截断误差的主项。对于Euler方法，有Taylor展开有定义5.1设是初值问题（5.1.对于隐式Euler方法，其局部截断误差为所以Euler方法是一种一阶方法，其局部截断误差的主项为。梯形方法也是一种隐式单步法，类似可得其局部截断误差所以隐式Euler方法也是一种一阶方法，该方法的局部截断误差的主项为，仅与显式Euler方法的局部截断误差的主项反一个符号。可见，梯形方法是二阶精度的。对于隐式Euler方法，其局部截断误差为所以Euler方法是5.2Runge-Kutta方法5.2.2几类显式Runge-Kutta方法5.2.1Runge-Kutta方法的基本思想5.2Runge-Kutta方法5.2.2几类显5.2.1Runge-Kutta方法的基本思想

显式Euler方法是最简单的单步法，它是一阶的，它可以看作Talylor展开后取前两项。因此，得到高阶方法的一个直接想法是用Talylor展开，如果能计算的高阶导数，则可写出p阶方法的计算方法其中是的近似值，若将分别记成则对于二阶和三阶导数可表示为5.2.1Runge-Kutta方法的基本思想这个方法并不实用，因为一般情况下，求的导数相当麻烦。从计算高阶导数的公式知道，方法的截断误差提高一阶，需要增加的计算量很大。但是由此启发我们用区间上若干个点的导数，而不是高阶导数，将它们作线性组合得到平均斜率，将其与解的Taylor展开相比较，使前面若干项吻合，从而得到具有一定阶的方法。这就是Runge-Kutta方法的基本思想，其一般形式为(5.2.1)这个方法并不实用，因为一般情况下，求其中，与的区别在于：用微分方程准确解代替中的就得到。参数和待定，确定它们的原则和方法是：将（5.2.2）式中的在处作Taylor展开，将在处作二元Taylor展开，将展开式按H的幂次整理后，令中h的低次幂的系数为零，使首项中h的幂次尽量高，比如使，则称（5.2.1）式为L级p阶Runge-Kutta方法（简称R-K法）。其中，。它的局部截断误差是，（5.2.2）其中，与的区别在于：用微分方程准确解它与显式R-K公式的区别在于：显式公式中，对系数求和的上限是，从而构成的矩阵是一个严格下三角阵。而在隐式公式中，对系数求和的上限是L，从而构成的矩阵是方阵，需要用迭代法求出近似斜率推导隐式公式的思路和方法与显式公式法类似。类似于显式R-K公式（5.2.1），稍加改变，就得到隐式R-K公式它与显式R-K公式的区别在于：显式公式中，对系数5.2.2几类显式Runge-Kutta方法

对于L=2，则其局部截断误差是（5.2.3）5.2.2几类显式Runge-Kutta方法其局部截断

将中的各项作Taylor展开，并利用则有将它们代入（5.2.3）式，整理后得将中的各项作Taylor展开，并选取和，使方法的阶尽可能高，就是使h和的系数为零，因为的系数一般不为零。于是得到方程组显然，该方程组有无穷多组解，从而得到一族二级二阶R-K方法。若以为自由参数，取得中点公式

（5.2.4）取c=2/3得Heun公式，（5.2.5）选取和，使方法的阶尽可能高，取c=1得改进的Euler公式（5.1.6）。对于L=3的情形，要计算三个斜率的近似值：类似于二阶方法的推导，可以得三阶的方法，所得系数应满足的方程组是该方程组的解也是不唯一的。常见的一种三级三阶方法是取c=1得改进的Euler公式（5.1.6）。对于L=3的情对于L=4的情形，可进行类似推导。最常用的四级四阶方法是如下经典R-K方法（5.2.6）对于L=4的情形，可进行类似推导。最常用的四级四阶方法是如下为了分析经典R-K公式的计算量和计算精度，将四阶经典R-K公式（5.2.6）与一阶显式Euler公式（5.1.2）及二阶改进的Euler公式相比较。一般说来，公式的级数越大，计算右端项f的次数越多，计算量越大。在同样步长的情况下，Euler方法每步只计算一个函数值，而经典方法要计算4个函数值。四阶R-K法的计算量差不多是改进的Euler公式的2倍，是显式Euler公式的4倍。下面的例子中Euler方法用步长，二阶改进的Euler法用步长，而四阶经典公式用步长。这样，从到三种方法都计算了4个函数制，计算量大体相当。例5.3考虑初值问题为了分析经典R-K公式的计算量和计算精度，将四阶经典R-其解析解为。分别用h=0.025的显式Euler方法，h=0.05改进Euler法和h=0.1的经典R-K方法计算到x=0.5。三种方法在x方向每前进0.1都要计算4个右端函数值，计算量相当。计算结果列于表5-3。从计算结果看，在工作量大致相同的情况下，还是经典方法比其他两种方法的结果好得多。在x=0.5处，三种方法的误差分别是。经典R-K法对多数好条件问题（，参见下节单步法的稳定性），能获得好的效果。其解析解为。分别用h=0.表5-3Euler法改进Euler法经典R-K法准确解h=0.025h=0.05h=0.10.10.0963120.0951230.095162500.095162580.20.1833480.1811930.181269100.181269250.30.2620010.2590850.259181580.259181780.40.3330790.3295630.329679710.329679950.50.3973120.3933370.393469060.39346934

在微分方程数值解法的实际计算中，有个如何选择步长的问题。因为单从每一步看，步长越小，截断误差越小。但随着步长的缩小，在一定求解范围内所要完成的步数就增加了。步数的增加不但引起计算量的增大，而且可能导致舍入误差的严重积累。表5-3Euler法改进E

在选择步长时，我们需要衡量和检验计算结果的精度，并依据所获得的精度处理步长。下面以经典R-K方法为例进行说明。从节点出发，先以h为步长求出一个近似值，由于公式的局部截断误差为，故有然后将步长折半，既h/2为步长，从跨两步到，再求得一个近似值，每跨一步的截断误差约为，因此有然后将步长折半，既h/2为步长，从跨两步到比较上述二式，有由此易得下列事后估计式这样，我们可以通过检查步长折半前后两次计算结果的偏差来判定所选的步长是否合适。比较上述二式，有由此易得下列事后估计式这样，我们可以通过检查

具体地说，对于给定的精度，将按两种情况处理。如果，我们反复将步长折半进行计算，直到为止，这时取最终得到的作为结果。如果，我们反复将步长加倍，直到为止，这时再将前一次步长折半的结果作为所要的结果。这种通过加倍或折半处理步长的方法称作变步长方法。虽然为了选择步长，每一步的计算量有所增加，但总体考虑是值得的。数值分析课件5.5

线性多步法5.5.2基于Taylor展开的方法5.5.1基于数值积分的方法5.5线性多步法5.5.2基于Taylor展5.5

线性多步法

常微分方程初值问题（5.1.1）的数值解法中，除了Runge-Kutta型公式等单步法之外，还有另一种类型的解法，即某一步的公式不仅与前一步解的值有关，而且与前若干步解的值有关，利用前面多步的信息预测下一步的值，这就是多步法的基本思想，可以期望获得较高的精度。构造多步法有多种途径，下面先讨论基于数值积分的方法。5.5线性多步法常微分方程初值问题（55.5.1基于数值积分的方法将（5.1.1）中的方程在区间上积分，可以得到（5.5.1）如推导Newton-Cotes求积公式一样，用等距节点的插值多项式来替代被积函数，再对插值多项式积分，这样就得到一系列求积公式。例如，用梯形方法计算积分项5.5.1基于数值积分的方法将（5.1.1）中的方程在代入（5.5.1）式有据此即可导出公式（5.1.4）。一般地，设由个数据点构造插值多项式，这里，。运用插值公式有将（5.5.1）离散化即得下列计算公式代入（5.5.1）式有据此即可导出公式（5.1.4）。将（5（5.5.2）其中由此可得（5.5.2）中的系数，其具体数值见表5-6。公式（5.5.2）是一个r+1步的显式公式，称为Adams显式公式。r=0时，即为Euler公式。（5.5.2）其中由此可得（5.5.2）中的系数，其具体数值表5-6j0123413-123-16555-5937-91901-27742616-1274251应用实例：考虑跳伞员的下落速度。自由落体运动可用牛顿第二定律描述：F=ma。实验表明，空气阻力模型为，其中，比例系数k依赖于物体的大小、形状，空气的密度和粘度。跳伞员下落的速度可描述为下列模型：表5-6j0负号表示下降。显然，当1<p<2时，适合于数值方法求解。设k/m=1.5，g=32，先用中点法提供开始值，再用下列两步而阶方法求其他需要计算的值。当p=1时，取h=0.2有负号表示下降。显然，当1<p<2时，适合于数值方法求可见，三秒末跳伞员的末速度约有21。若将模型修改为p=1.1，取h=0.2，则有计算结果：可见，三秒末跳伞员的末速度约有21可见三秒末跳伞员的末速度减慢了。计算结果如下图所示+表示p=1时的解，*表示p=1.1时的解可见三秒末跳伞员的末速度减慢了。计算结果如下图所示+表示

在上述Adams显式公式的推导中，选用了作为插值节点。这样的插值多项式在求积区间上逼近是一个外推结果。为了改善逼近效果，我们变外推为内推，即改用为插值节点，用数据点构造插值多项式，则有于是我们有如下的计算公式在上述Adams显式公式的推导中，选用了（5.5.3）其中其具体数值见表5-7。公式（5.5.3）是隐式公式，称为Adams隐式公式。r=0，1时分别为隐式Euler公式和梯形公式。（5.5.3）其中其具体数值见表5-7。公式（5.5.3）表5-7j0123411158-1919-51251646-264106-19对于隐式公式（5.5.3），需要用迭代求解。确定的迭代公式为表5-7j0迭代收敛条件为，其中的Lipschitz常数利用插值多项式的余项，可以求出Adams方法的局部截断误差。当然也可以从得到的显式和隐式Adama公式，有局部截断误差的定义来求出方法的局部截断误差。表5-8中列出了它们的局部截断误差的主项，有表5-8可以看出，Adams隐式方法的局部截断误差要小。r0123表5-8Adams显式公式Adams隐式公式迭代收敛条件为，其中5.5.2基于Taylor展开的方法

基于数值积分可以构造出一系列求解常微分方程的计算公式，下面介绍基于Taylor展开的待定系数法，它可灵活地构造出线性多步法。对固定的系数，可以选取待定系数使线性多步法的阶尽可能高。还可以根据需要，确定显式还是隐式。设构造如下具有p阶精度的线性多步公式(8.4.4)当时，则（8.4.4）为显式多步式。当时，（8.4.4）为隐式多步式。它们的局部截断误差为5.5.2基于Taylor展开的方法利用原微分方程，有（5.5.5）现利用Taylor展开定理，确定线性多步公式（5.5.4）中的待定参数，使她达到阶精度，即。对（5.5.5）式的右端各项在点处作Taylor展开有利用原微分方程，有（5.5.5）现利用Taylor展开定理，将它们代入（5.5.5）式整理后得将它们代入（5.5.5）式整理后得使的系数为零，得到关于和的线性方程组（5.5.6）而且得到线性多步法的局部截断误差使下面我们构造几个著名的四阶线性多步公式，考虑下列形式的公式（5.5.7）（5.5.8）由于r=3，p=4，由（5.5.6）得到5个方程，而（5.5.7）中有9个为知量，因此，（5.5.7）中有4个自由度。若取，由（5.5.6）式得到其他5个待定参数的方程组，解之得下面我们构造几个著名的四阶线性多步公式，考虑下列形式的公式（代入（5.5.7）和（5.5.8）式，得到常用的四步四阶显式Admas公式和它的余项：（5.5.9）（5.5.10）

若取，由（5.5.6）式得到其他5个待定参数的方程组，解之得由此构造成著名的四步四阶显式Milne公式和它的余项代入（5.5.7）和（5.5.8）式，得到常用的四步四阶显式（5.5.11）（5.5.12）

若取由（8.4.6）式得到其他5个待定参数的方程组，从而得三步四阶隐式Admas公式及余项：（5.5.13）（5.5.14）

若取，求解（5.5.6）得著名的三步四阶隐式Hamming公式及其余项：（5.5.11）（5.5.12）若取（5.5.15）（5.5.16）

若取，求解（5.5.6）得到隐式Simpson公式及其余项：（5.5.15）（5.5.16）若取

例5.5分别取h=0.2，2，用四阶显式Milne公式和四阶隐式Hamming公式求解例5.4所给的初值问题。解我们用单步法提供多步法的初值。由4阶经典R-K公式为Milne公式提供初值，为Hamming公式提供。h=0.2和h=2时的计算结果及准确解之间的误差分别列于表8-9和表8-10。从表5-9看出，两种多步法的计算精度都很高，Hamming公式化比Milne公式更精确。这是因为Hamming公式的截断误差主项的系数比Milne公式小。从表5-10看到，当计算步长变大后，显式多步法Milne公式的计算结果误差增大，不稳定，而隐式多步法Hamming公式的计算结果仍然是稳定的，这说明隐式公式的稳定性比同阶的显式公式好。例5.5分别取h=0.2，2，用四阶显表5-9Milne方法误差Hammins方法误差2.20.942942680.942919552.41.122833491.122833862.61.306432141.306389302.81.492916251.492925823.01.681954501.68190299表5-9Milne方法表5-10Milne方法误差Hammins方法误差75.6457455.64574597.3823257.6371261110.9053169.635636134.14383111.6322611558.31071713.63224017-249.66267215.631690

经典R-K法和上述四阶线性多步法公式都是四阶精度，但每前进一步，前者要计算4次微分方程右端方程，而后者只要计算一次新的右端函数值，计算量减小了。表5-10Milne方法5.5.3预估-校正算法

显式多步法容易计算，但其精度和稳定性没有相应的隐式方法好。然而，隐式多步法需解方程，如果初值选得不当，则计算量较大。因此，设法选取好的迭代初值是必要的。初值的自然选取是采用同阶显式多步法计算得到的解作为隐式方法迭代的初值。这样，迭代次数不会多。若只迭代一次，则这样的算法就是预估-校正算法。对于线性多步法，常用的预估——校正方法有四阶Admas显隐式预估-校公式和Milne-Hamming方法。1.Adams预估-校正公式由（5.5.9）式作为预估公式，由（5.5.13）式作为校正公式，构成Adams预估-校正公式：5.5.3预估-校正算法1.Adams预估-校正公式

若需作进一步的修正，则记上式所得的，由（5.5.10）和（5.5.14）式有于是得到若需作进一步的修正，则记上式所得的由此可见，若记

则分别比更好。但注意到，的表达式中，是未知的，因此改为由此可见，若记则分别这样，得到下面的修正的Adams预估-校正公式：修正：校正：修正：预估：

在计算时，可调节计算步长h，使，其中是要求达到的计算精度。初值由同阶单步法提供，当计算时，可取。这样，得到下面的修正的Adams预估-校正公式：2.修正Hamming公式将Milne公式（5.5.11）和Hamming公式（5.5.15）结合，构成Milne-Hamming预估-校正公式：若需作进一步的修正，则记上式所得的，由（5.5.12）和（5.5.16）有2.修正Hamming公式若需作进一步的修正，则于是得到由此分别得Milne和Hamming公式的修正公式：从而构成如下的修正Hamming公式：于是得到由此分别得Milne和Hamming公式的修预估：修正：校正：修正：在计算时，可调节计算步长h，使。初值

由同阶单步法提供，当计算时，可取。预估：修正：校正：修正：在计算时，可调节计算步长h，使

例5.6取h=0.2，用Milne-Hamming预估-校正公式和修正Hamming公式求解例5.4所给的初值问题。解用经典R-K法提供初值，计算结果列于表5-11。将表5-9与表5-11所示的计算结果进行比较，它们的计算精度排列次序是：修正Hamming公式的精度最好，其次是隐式Hamming公式，再次是Milne-Hamming预估-校正公式，最后是Milne公式。表5-11Milne-Hamming误差修正Hamming误差2.20.942916250.942924492.41.122828721.122839552.61.306382711.306395372.81.492918161.492931843.01.681894671.68190879例5.6取h=0.2，用Milne-Ha第5章常微分方程数值解法5.1.2局部误差和方法的阶5.1.1Euler方法及其有关的方法第5章常微分方程数值解法5.1.2局部误差和方法的阶5.1Euler方法5.1.1Euler方法及其有关的方法考虑一阶常微分方程初值的问题：设f（x,y）是连续函数，对y满足Lipschitz条件，这样初值问题的解是存在唯一的，而且连续依赖于初始条件。为了求得离散点上的函数值，将微分方程的连续问题进行离散化。一般是引入点列{}，这里为步长，经常考虑定长的情形，即。记为初始问题的问题准确解在处的值，用均差近似代替的导数得5.1Euler方法5.1.1Euler方法及其有关令为的近似值，将上面两个近似写成等式，整理后得（5.1.2）（5.1.3）从处的初值开始，按（5.1.2）可逐步计算以后各点上的值。称（5.1.2）式为显式Euler。由于（5.1.3）式的右端隐含有待求函数值，不能逐步显式计算，称（5.1.3）式为隐式Euler公式或后退Euler公式。如果将（5.1.2）和（5.1.3）两式作算术平均，就得梯形公式。令为的近似值，将上面两个梯形公式也是隐式公式。以上公式都是由去计算，故称它们为单步法。例5.1取h=0.1，用Euler方法、隐式Euler方法和梯形方法解

解本题有如果用Euler方法，由（5.1.2）并代入h=0.1得同理，用隐式Euler方法有（5.1.4）梯形公式也是隐式公式。以上公式都是由去计算用梯形公式有三种方法及准确解的数值结果如表5-1所示。从表中看到，在处，Euler方法和隐式Euler方法的误差分别是和，而梯形方法的误差却是。

在例5.1中，由于f（x，y）对y是线性的，所以对隐式公式也可以方便地计算。但是，当f（x，y）是y的非线性函数时，如，其隐式Euler公式为。显然，它是的非线性方程，可以选择非线性方程求根的迭代求解。以梯形公式为例，可用显式Euler公式提供迭代初值，用公式用梯形公式有三种方法及准确解表5-1Euler方法隐式Euler方法梯形法准确解011110.11.0000001.0090911.0047621.0048370.21.0100001.0264461.0185491.0187310.31.0290001.0513151.0406331.0408180.41.0561001.0830131.0700961.0703200.51.0904901.1209211.1062781.106531表5-1Euler方法反复迭式，直到其中，步长h成为迭代参数，它需要满足一定的条件，才能收敛。若将（5.1.4）式减去该迭代公式，得假设f（x，y）关于y满足Lipschiz条件，则有反复迭式，直到其中，步长h成为迭代参数，它需要满足一定的条件这里，L是Lipschiz常数。当hL/2<1即h<2/L时，迭代序列收敛。

对于隐式公式，通常采用估计-校正技术，即先用显式公式计算，得到预估值，然后以预估值作为隐式公式的迭代初值，用隐式公式迭代一次得到校正值，称为预估-校正技术。例如，用显式Euler公式作预估，用梯形公式作校正，即称该公式为改进的Euler公式。它显然等价于显式公式为，（5.1.6）这里，L是Lipschiz常数。当hL/2<1即h<2/L时也可以表示为下列平均化的形式例5.2取h=0.1，用改进的Euler方法解解按（5.1.5），改进的Euler方法解也可以表示为下列平均化的形式例5.2取h=0.1，用改进由得计算结果如表5-2。该初值问题的准确解为。表5-20.10.20.30.40.50.60.70.81.09591.18411.26621.34341.41641.48601.55251.61531.09541.18321.26491.34161.41421.48321.54921.6165由得5.1.2局部误差和方法的阶

初值问题（5.1.1）的单步法可以写成如下统一形式（5.1.7）其中与有关。若中不含，则方法是显式的，否则是隐式的，所以一般显式单步法表示为（5.1.8）例如，Euler方法中，有对于不同的方法，计算值与准确解的误差各不相同。所以有必要讨论方法的截断误差。我们称为某一方法在点的整体截断误差。显然，不单与这步的计算有关，它与以前各步的计算也有关，所以误差被称为整体的。分析和估计整体截断误差是复杂的。为此，我们假设处的没有误差，即，考虑从到这一步的误差，这就是如下的局部误差的概念。5.1.2局部误差和方法的阶（5.1.7）其中定义5.1

设是初值问题（5.1.1）的准确解，则称为单步法（5.1.7）的局部截断误差。

定义5.2

如果给定方法的局部截断误差，其中为整数，则称该方法是p阶的，或具有p阶精度。若一个p阶单步法的局部截断误差为则称其第一个非零项为该方法的局部截断误差的主项。对于Euler方法，有Taylor展开有定义5.1设是初值问题（5.1.对于隐式Euler方法，其局部截断误差为所以Euler方法是一种一阶方法，其局部截断误差的主项为。梯形方法也是一种隐式单步法，类似可得其局部截断误差所以隐式Euler方法也是一种一阶方法，该方法的局部截断误差的主项为，仅与显式Euler方法的局部截断误差的主项反一个符号。可见，梯形方法是二阶精度的。对于隐式Euler方法，其局部截断误差为所以Euler方法是5.2Runge-Kutta方法5.2.2几类显式Runge-Kutta方法5.2.1Runge-Kutta方法的基本思想5.2Runge-Kutta方法5.2.2几类显5.2.1Runge-Kutta方法的基本思想

显式Euler方法是最简单的单步法，它是一阶的，它可以看作Talylor展开后取前两项。因此，得到高阶方法的一个直接想法是用Talylor展开，如果能计算的高阶导数，则可写出p阶方法的计算方法其中是的近似值，若将分别记成则对于二阶和三阶导数可表示为5.2.1Runge-Kutta方法的基本思想这个方法并不实用，因为一般情况下，求的导数相当麻烦。从计算高阶导数的公式知道，方法的截断误差提高一阶，需要增加的计算量很大。但是由此启发我们用区间上若干个点的导数，而不是高阶导数，将它们作线性组合得到平均斜率，将其与解的Taylor展开相比较，使前面若干项吻合，从而得到具有一定阶的方法。这就是Runge-Kutta方法的基本思想，其一般形式为(5.2.1)这个方法并不实用，因为一般情况下，求其中，与的区别在于：用微分方程准确解代替中的就得到。参数和待定，确定它们的原则和方法是：将（5.2.2）式中的在处作Taylor展开，将在处作二元Taylor展开，将展开式按H的幂次整理后，令中h的低次幂的系数为零，使首项中h的幂次尽量高，比如使，则称（5.2.1）式为L级p阶Runge-Kutta方法（简称R-K法）。其中，。它的局部截断误差是，（5.2.2）其中，与的区别在于：用微分方程准确解它与显式R-K公式的区别在于：显式公式中，对系数求和的上限是，从而构成的矩阵是一个严格下三角阵。而在隐式公式中，对系数求和的上限是L，从而构成的矩阵是方阵，需要用迭代法求出近似斜率推导隐式公式的思路和方法与显式公式法类似。类似于显式R-K公式（5.2.1），稍加改变，就得到隐式R-K公式它与显式R-K公式的区别在于：显式公式中，对系数5.2.2几类显式Runge-Kutta方法

对于L=2，则其局部截断误差是（5.2.3）5.2.2几类显式Runge-Kutta方法其局部截断

将中的各项作Taylor展开，并利用则有将它们代入（5.2.3）式，整理后得将中的各项作Taylor展开，并选取和，使方法的阶尽可能高，就是使h和的系数为零，因为的系数一般不为零。于是得到方程组显然，该方程组有无穷多组解，从而得到一族二级二阶R-K方法。若以为自由参数，取得中点公式

（5.2.4）取c=2/3得Heun公式，（5.2.5）选取和，使方法的阶尽可能高，取c=1得改进的Euler公式（5.1.6）。对于L=3的情形，要计算三个斜率的近似值：类似于二阶方法的推导，可以得三阶的方法，所得系数应满足的方程组是该方程组的解也是不唯一的。常见的一种三级三阶方法是取c=1得改进的Euler公式（5.1.6）。对于L=3的情对于L=4的情形，可进行类似推导。最常用的四级四阶方法是如下经典R-K方法（5.2.6）对于L=4的情形，可进行类似推导。最常用的四级四阶方法是如下为了分析经典R-K公式的计算量和计算精度，将四阶经典R-K公式（5.2.6）与一阶显式Euler公式（5.1.2）及二阶改进的Euler公式相比较。一般说来，公式的级数越大，计算右端项f的次数越多，计算量越大。在同样步长的情况下，Euler方法每步只计算一个函数值，而经典方法要计算4个函数值。四阶R-K法的计算量差不多是改进的Euler公式的2倍，是显式Euler公式的4倍。下面的例子中Euler方法用步长，二阶改进的Euler法用步长，而四阶经典公式用步长。这样，从到三种方法都计算了4个函数制，计算量大体相当。例5.3考虑初值问题为了分析经典R-K公式的计算量和计算精度，将四阶经典R-其解析解为。分别用h=0.025的显式Euler方法，h=0.05改进Euler法和h=0.1的经典R-K方法计算到x=0.5。三种方法在x方向每前进0.1都要计算4个右端函数值，计算量相当。计算结果列于表5-3。从计算结果看，在工作量大致相同的情况下，还是经典方法比其他两种方法的结果好得多。在x=0.5处，三种方法的误差分别是。经典R-K法对多数好条件问题（，参见下节单步法的稳定性），能获得好的效果。其解析解为。分别用h=0.表5-3Euler法改进Euler法经典R-K法准确解h=0.025h=0.05h=0.10.10.0963120.0951230.095162500.095162580.20.1833480.1811930.181269100.181269250.30.2620010.2590850.259181580.259181780.40.3330790.3295630.329679710.329679950.50.3973120.3933370.393469060.39346934

在微分方程数值解法的实际计算中，有个如何选择步长的问题。因为单从每一步看，步长越小，截断误差越小。但随着步长的缩小，在一定求解范围内所要完成的步数就增加了。步数的增加不但引起计算量的增大，而且可能导致舍入误差的严重积累。表5-3Euler法改进E

在选择步长时，我们需要衡量和检验计算结果的精度，并依据所获得的精度处理步长。下面以经典R-K方法为例进行说明。从节点出发，先以h为步长求出一个近似值，由于公式的局部截断误差为，故有然后将步长折半，既h/2为步长，从跨两步到，再求得一个近似值，每跨一步的截断误差约为，因此有然后将步长折半，既h/2为步长，从跨两步到比较上述二式，有由此易得下列事后估计式这样，我们可以通过检查步长折半前后两次计算结果的偏差来判定所选的步长是否合适。比较上述二式，有由此易得下列事后估计式这样，我们可以通过检查

具体地说，对于给定的精度，将按两种情况处理。如果，我们反复将步长折半进行计算，直到为止，这时取最终得到的作为结果。如果，我们反复将步长加倍，直到为止，这时再将前一次步长折半的结果作为所要的结果。这种通过加倍或折半处理步长的方法称作变步长方法。虽然为了选择步长，每一步的计算量有所增加，但总体考虑是值得的。数值分析课件5.5

线性多步法5.5.2基于Taylor展开的方法5.5.1基于数值积分的方法5.5线性多步法5.5.2基于Taylor展5.5

线性多步法

常微分方程初值问题（5.1.1）的数值解法中，除了Runge-Kutta型公式等单步法之外，还有另一种类型的解法，即某一步的公式不仅与前一步解的值有关，而且与前若干步解的值有关，利用前面多步的信息预测下一步的值，这就是多步法的基本思想，可以期望获得较高的精度。构造多步法有多种途径，下面先讨论基于数值积分的方法。5.5线性多步法常微分方程初值问题（55.5.1基于数值积分的方法将（5.1.1）中的方程在区间上积分，可以得到（5.5.1）如推导Newton-Cotes求积公式一样，用等距节点的插值多项式来替代被积函数，再对插值多项式积分，这样就得到一系列求积公式。例如，用梯形方法计算积分项5.5.1基于数值积分的方法将（5.1.1）中的方程在代入（5.5.1）式有据此即可导出公式（5.1.4）。一般地，设由个数据点构造插值多项式，这里，。运用插值公式有将（5.5.1）离散化即得下列计算公式代入（5.5.1）式有据此即可导出公式（5.1.4）。将（5（5.