常微分方程初值问题的线性多步法公式的改进与变异.pdf

常微分方程初值问题的线性多步法公式的改 进与变异 重庆大学硕士学位论文 学生姓名： 郑兴武 指导教师： 杨大地 副教授 专 业： 运筹学与控制论 学科门类： 理学 重庆大学数理学院 二 OO 八年四月 Improvement and Variation on Linear Multistep Method Formulas of Initial Value Problems in Ordinary Differential Equations A Thesis Submitted to Chongqing University in Partial Fulfillment of the Requirement for the Degree of Master of Science by Zheng Xingwu Supervisor: Associate Prof. Yang Dadi Major: Operational Research ()()() k p kjj j L y x hy xihhy xihO h (2.6) 其中p为正整数，则称相应数值方法是p阶相容的。 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 9 对于线性多步法（2.3） ，根据局部截断误差 ( ); k L y x h，由 Taylor 展开可 求得多步法的系数和相容性条件。 定义定义 2.3.3 35 如果线性多步法（2.3）是 p阶的，则 1 ( ); () p k L y x hO h ，它 等价于 011 ,0 pp CCC C ，于是 1(1)2 1 ( ); ( )() ppp kp L y x hChyxO h (2.7) 故（2.3）的局部截断误差为 1(1)2 1 ( )() ppp n kp TChyxO h (2.8) 我们称 1(1) 1 ( ) pp p Chyx 为局部截断误差主项， 1p C 为误差常数。 定义定义 2.3.4 36 设 n F是正整数变量n的一个已知函数，k为正整数，方程 1 (,)0 nnnn k F y yy (2.9) 称为差分方程。 定义定义 2.3.5 35 若记 00 ( ),( ) kk ii ii ii (2.10) 分别称( ) 与( ) 为差分方程（2.9）的第一和第二特征多项式。 定理定理 2.3.1 35 线性多步法（2.3）相容的充分必要条件是 (1)0,(1)(1) (2.11) 2.4 线性多步法的收敛性 求解常微分方程的方法可以分为： 显 （隐） 式线性多步法、 显 （隐） 式 Runge-Kutta 方法、配置方法、线性多级多值方法以及由这些方法衍生出来的一些方法。如何 选择合适的积分公式构造一个高效的常微分方程数值方法是个复杂的问题，需要 考虑公式的变步长问题、显式公式作为隐式公式的预报并与隐格式构成预报校正 系统等。 其中收敛性是对数值方法的一种起码要求，不收敛的数值方法没有任何 实际应用价值，因此，对数值方法的收敛性分析，其意义是很重要的。 数值解的收敛性问题就是由数值方法得到的数值解能否近似代替原连续方程 的解。 定义定义 2.4.1 35 对于初值问题 （2.1） 的任意一种数值方法在 n xx处的数值解为 n y， 这里 0 , n xxnha b为固定点， 设( )y x为 （2.1） 的精确解， 称() nnn ey xy 为整体误差，若 () 00 limlim ()0 n nnn hh ey xy 则称数值解 n y收敛于精确解( )y x，即数值方法是收敛的。 定理定理 2.4.1 36 线性多步法（2.3）是收敛的，则它一定是相容的。 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 10 定理 2.4.1 的逆定理是不成立的，但对单步法逆定理也成立。 定义定义 2.4.2 35 如果线性多步法（2.3）的特征多项式 0 ( ) k i i i 的根都在 单位圆内并且在单位圆上只有单根出现，则称多步法（2.3）满足根条件。 定理定理 2.4.2 35 线性多步法（2.3）是1p 阶相容的，则多步法（2.3）的差分方 程收敛的充分必要条件是多步法（2.3）满足根条件。 定理定理 2.4.3 35 线性多步法（2.3）收敛的充分必要条件是方法相容，并满足根 条件。 2.5 线性多步法的稳定性 在实际进行计算时，一方面出发值不一定是完全精确的，带有一定的误差， 同时，由于计算机的字长有限，在运算中一般总会产生舍入误差，不论是单步法 或多步法，在逐次计算下去的时候，初始数据的误差（或称之为摄动）以及在计 算过程中产生的舍入误差，都会传播下去，对以后的计算结果产生影响。所谓稳 定性问题是指初始数据的误差和计算过程中产生的误差的积累和传播是否受到控 制，或者说，如果计算结果对初始数据的误差以及计算过程中的舍入误差不敏感， 就说相应的计算方法是稳定的，否则就称之为不稳定的。 定义定义 2.5.136 令,0,1, n nN及 * ,0,1, n nN是任意两个扰动， 并令,0,1, n z nN及 * ,0,1, n znN是对应以上扰动结果的解，若存在 常数 0 h和C，使对 0 (0,)hh ，只要对 * 0, nn nN，就有 * ,0, nn zzCnN 则称方法是稳定的或称为零稳定的。 这种稳定性也称渐近稳定，首先由 Dahlquist 于 1956 年给出的，也称D稳定 性。 定理定理 2.5.1 37 线性多步法（2.3）是稳定的充分必要条件是它满足根条件。 定理定理 2.5.2 3 稳定的线性k步法，当k为奇数时，阶数不超过1k ，当k为偶 数时，阶数不超过2k 。 2.6 线性多步法的绝对稳定性与绝对稳定域 稳定性概念， 是在0h的情况下讨论的， 这样的稳定性称之为渐近稳定性 （或 古典的稳定性） 。然而，实际上我们是取有限的固定步长h进行计算的，它并不能 随意地缩小。因此，重要的是，在计算过程中所产生的摄动对以后的计算结果的 影响不会增长，这种稳定性概念就是通常所说的绝对稳定性。 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 11 研究数值方法是否数值稳定，不可能也不必对每个不同的右端函数( , )f x y进 行讨论，通常只对试验方程 0 ,Re( )0 ( ) yy y ay (2.12) 进行讨论，即研究数值方法用于解方程（2.12）得到的差分方程是否数值稳定。 定义定义 2.6.1 35 一个数值方法用于解试验方程（2.12） ，若对给定步长h得到的 线性差分方程解 n y，当n时，0 n y ，就称该方法对步长h是绝对稳定的。 当n时， n y无界，则称该方法数值不稳定。 定义定义 2.6.2 35 一个数值方法用于解试验方程（2.12） ，若在 h的复平面 中的某个区域R中方法都是绝对稳定的，而在区域R外，方法是不稳定的，则称 区域R是该数值方法的绝对稳定域。 定义定义 2.6.3 35 将线性多步法（2.3）用于解（2.12）得到k阶的线性差分方程 00 kk jnjjnj jj yhf (2.13) 若记h，则得（2.3）的特征多项式 ( ,)( )( ) (2.14) 其中 0 ( ) k i i i 称为第一特征多项式， 0 ( ) k i i i 称为第二特征多项式， 它是关于的k阶多项式。 对给定的h，( ) 的根都在单位圆内，则称该方法关于h是绝对稳 定的。若存在区间( ,) ，使线性多步方法对任意的( ,) 都是绝对稳定的， 则称( ,) 为绝对稳定区间。若存在区间R，对任意的R(复平面域)都是绝对 稳定的，则称区域R为绝对稳定区域。 定理定理 2.6.1 36 设特征多项式( , )( )( ) ，h的零点为( ) j ， 1,jk则线性多步法的绝对稳定区域 ,()1,1, j Rhjk 2.7 线性多步法的边界轨迹法 求线性多步法的绝对稳定区域较为复杂，边界轨迹法是求绝对稳定区域的有 效方法。设一点在域R的边R上，则此点对应了( , )0 的一个模为 1 的 根，记为 2 ,1 i ei ，则对于线性多步法可用边界轨迹法求出公式的绝对稳定 区域的边界线。 () (),0,2 () i i i e eh e 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 12 考虑到试验方程中的一般为负实数，方法的绝对稳定区域R可分为如下情 况： （1）R为一个有界封闭区域，且有负实轴上的 00 (,0,0hR h； （2）R为一个无界区域，且有Re( )0R； （3）R为一个无界区域，且有(,0R； （4） 其它情况。 情况（1）称方法有绝对稳定区间 0 (,0h， 0 h越大，绝对稳定性越好； 情况（2）称方法具有A-稳定性，此类公式用于刚性方程组求解时步长h不受 限制。但根据 Dahlquist 第二障碍定理，A-稳定的方法最高阶数只有 2 阶； 情况（3）称方法具有(0)A-稳定，其中有一部分方法是( )A-稳定方法，是求 解刚性方程组的主要方法之一； 情况（4）方法稳定性不好，实用价值有限。 2.8 刚性方程的线性多步法稳定性概念 在航天航空、电子电路设计、工程控制、化学动力等许多重要科技领域中， 有大量物理或化学过程可用常微分方程模型来描述，这些模型的解中既包含衰减 十分迅速的分量，又包含有变化相对缓慢的分量，当人们试图在解的慢变区间上 数值求解这类问题时，尽管此时快变分量的值已衰减到可以忽略不计，但这种快 变分量的干扰仍严重影响数值解的稳定性和精度，致使所得到的数值解失真，人 们称这类微分方程的定解问题为“刚性（Stiff）问题。 刚性方程在文献中也称作为病态方程或坏条件方程，或称为具有差别大的时 间常数问题或具有大的 Lipschitz 常数问题。这类问题在实际工作中经常碰到。例 如，控制部件一般反应灵敏，是快变的，具有小的时间常数，而受控物体一般惯 性大，是慢变的，具有大的时间常数。航空中的运载器，通常是通过控制部件来 控制质心运动。姿态运动是快变的，而质心运动是慢变的。在许多化学反应中， 有些反应瞬间完成，达到稳定状态，是快变的；而有些反应速度则较为缓慢，两 者差别很大。 刚性性质是数学问题本身的固有性质。不依赖于求解这个问题的数值方法。 正是由于这个性质，使得传统的常微分方程的数值积分方法遇到极大的困难。为 了克服这个困难，刚性常微分方程数值积分方法的研究成为数值方法中最为重要 的研究方向之一。 用绝对稳定域有限的数值方法求刚性方程数值解，由于对步长h的限制，导致 计算步数很大，为减少计算步数，应选择对步长h无限制的方法。下面给出刚性方 程数值方法的稳定性相关理论。 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 13 定义定义 2.8.1 1 一个数值方法称为A-稳定的， 如果它的绝对稳定域包含整个左 半平面Re()0h。 1963 年，Dahlquist 1 在引入A-稳定的同时，证明了以下约束性结论： （1）任何显式线性多步法（包括显式 Runge-Kutta 方法）不可能是A-稳定的。 （2）A-稳定的隐式线性多步法的阶不超过 2。 （3）具有最小误差常数的 2 阶A-稳定隐式线性多步法是梯形法。 判断一个公式是否具有A-稳定性，可使用如下定理： 定理定理 2.8.1 36 一个k 步法公式是A-稳定的，如果它满足： （1）1,1,2, j jk，其中 j 是公式的第二特征多项式( ) 的根； （2） 0 ( )( )0,1 k kjj j P xT xx ， 其中 0 0 k ll l (2.15) 0 () k jljlllj l (2.16) ( ) j T x是j次 Chebyshev 多项式。 由于A-稳定的方法太少，并且梯形公式的误差常数最小。为找其他A-稳定的 数值方法，一方面是考虑隐式的 Runge-Kutta 法，另一方面就是降低对稳定性的要 求，使方法的绝对稳定域为无限域，但又不包含整个左半平面。下面给出( )A- 稳定、(0)A-稳定、 0 A稳定和刚性稳定的定义。 定义定义 2.8.2 10 一个数值方法称( )A-稳定的，如果它的绝对稳定域包含无限 的楔形区域|arg(), (0,/2)Whh ；数值方法称(0)A-稳定的， 如果对一个充分小的(0,/2)，它是( )A-稳定的。如图 2.1。 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 14 图 2.1 ( )A-稳定域 Figure 2.1 ( )A-stability domain 定义定义 2.8.3 39 如果数值方法的绝对稳定域包含整个负实轴， 则称该方法是 0 A 稳定。 ( )A-稳定的方法，它的绝对稳定域比A-稳定的方法的绝对稳定域小，因此， 如果一个方法是A-稳定，则必然是( )A-稳定。由于( )A-稳定的数值方法对步长 h没有限制，可用于解刚性方程，所以说( )A-稳定性是A-稳定性的推广。 定理定理 2.8.2 36 设下列条件成立 (1) 公式（2.3）是零稳定的； (2) 0,/0 kkk ； (3) ( ) 的根 j 满足1,1,2, j jk； (4) ) Im ()0,0, i e ； (5) ) Im ()tanRe ()0,0, ii ee ； 则k步法公式是( )A-稳定的。其中( )( )/( ) 。 在此定理的条件（3）成立时，条件（4） ， （5）又分别等价于 ( )0, 1,1 k Ixx (2.17) 2 1( )tan( )0, 1,1 kk x IxR xx (2.18) 其中， 1 ( )( ) k kjj j IxUx (2.19) 0 ( )( ) k kjj j R xT x (2.20) 0 0 0 0 (),1,2, , (),1,2, , k ll l kj jljlllj l kj jljlllj l jk jk (2.21) 其中( )sin()/sin ,( )cos(),cos jj UxjT xjx为 Chebyshev 多项式。 定义定义 2.8.4 36 一个数值方法称为刚性稳定的，如果它是收敛的，并且对某正 数, ,D 使得区域 1 |Re()RhhD 是绝对稳定的，而区域 2 |Re(), Im()RhDhh 上具有高精度且是绝对稳定或相对稳定 重庆大学硕士学位论文 2 初值问题数值方法的基本理论 15 的。如图 2.2。 00.10.20.30.40.50.60.70.80.91 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 R R -D 1 2 h平 面 图 2.2 刚性稳定域 Figure 2.2 stiff stability domain 下面的定理给出方法为刚性稳定的判别准则，它仅利用多项式( ) 和( ) 的 零点条件。 定理定理 2.8.3 36 收敛的公式为刚性稳定的充分必要条件是(1)-(4)成立。 （1）公式是 0 A稳定的。 （2）多项式( )/(1) 的任何根的模均小于 1。 （3）模为 1 的( ) 的根是单根。 （4）令是( ) 的模为 1 的根，则 ' ( )/( ) 是实数，并且是正的。 重庆大学硕士学位论文 3 线性多步法的改进 16 3 线性多步法的改进 3.1 一般线性多步法公式的介绍 线性多步法是求解常微分方程初值问题十分重要的方法，但在一般的文献中 只记载了数十个常用的重要算法。按照线性多步法的构造方法，应该还有许多公 式，其他的公式历史上是否都推导出来过，它们的性能如何，有没有实用价值， 是一个没有得到回答的问题。 定义 3.1.140 能用方程组（2.4）唯一确定的线性多步法公式称为基本的线性 多步法公式。 现有的线性多步法公式几乎都是基本的线性多步法公式。如广泛使用的 Adams 公式和 Gear 公式等。 在一些文献中也报道过非基本的公式， 如李庆扬30 给 出的对 Gear 方法的一种改进公式，要求 0111 (,1),(0,0,) kkk ， 共2k 个自由参数，而只要求公式为k阶，只需列出 01 0,0,0 k CCC，共 1k 个方程的线性方程组，此方程组有无穷多个解，可得到无穷多个线性多步法 公式。就属于非基本的线性多步法公式。 以上得到的公式并不完全是收敛的，还需要从中筛选出收敛的公式。 判断收敛： 应用定理 2.4.3，考查得到的基本公式是否收敛，因为 01 0CC就是相容性 条件，只需要判断是否满足根条件就行了，若满足根条件，则公式收敛，保留此 公式，否则放弃此公式。这样就得到了k步法全部基本且收敛的公式列表。 计算误差系数： 若公式收敛，用公式（2.4）计算并计算此公式的局部截断误差主项系数 1p C ， 可用精确分数形式表示，同时也确定此公式的阶为p。 对以上公式，用如下策略检验其绝对稳定性： 对分法判断绝对稳定性：用定义 2.6.3 的稳定多项式和定理 2.7.1 的原则，采用 对分法逐步扫描，可判断出公式是否具有绝对稳定区间，分辨出如下三种情况： 1）(,0)是绝对稳定区间，此时公式是(0)A稳定，下面将进一步判别它的A- 稳定性或( )A-稳定性； 2）存在0,(,0)是绝对稳定区间，此时公式具有有限的对稳定区域； 3）不存在0，使(,0)是绝对稳定区间，此类公式不具有绝对稳定性。 判断A-稳定性 公式是(0)A稳定，且公式只是二阶，则用定理 2.8.1 判断公式的A-稳定性。 判断( )A-稳定性 重庆大学硕士学位论文 3 线性多步法的改进 17 公式是(0)A稳定， 且公式大于二阶， 则用定理 2.8.2 判断公式的( )A-稳定性， 并计算稳定幅角。 3.2 Adams-Moulton 公式的改进 Adams-Moulton公式是一类1k 阶线性k步法隐式公式，应用十分广泛。但它 的绝对稳定区域是有界的，一般不适合用于刚性方程求解。 通过改进k步1k 阶 Adams-Moulton 公式得到了一类更稳定的k阶线性k步法 隐式公式， 对其中的 2 步 3 阶 Adams-Moulton 公式， 改进后可以得到)(A稳定性； 对 3 步 4 阶 Adams-Moulton 公式，改进后可以得到 0 A稳定性；对 4，5 步的 Adams-Moulton 公式， 改进后能使有界的绝对稳定区域增大， 还用数值实验证明了 这类公式对刚性方程问题的有效性。 3.2.1 Adams-Moulton 公式的介绍 2步3阶的Adams-Moulton公式： )85( 12 111 kkkkk fff h yy (3.1) 具有绝对稳定区间)0 , 6(，误差系数为0417. 0 24 1 4 C。 3步4阶的Adams-Moulton公式： )5199( 24 2111 kkkkkk ffff h yy (3.2) 具有绝对稳定区间)0 , 3(，局部截断误差系数为0264. 0 720 19 5 C。 4步5阶的Adams-Moulton公式： )19106264646251( 720 32111 kkkkkkk fffff h yy (3.3) 具有绝对稳定区间)0 ,8367. 1(，局部截断误差系数为0187. 0 160 3 6 C。 5步6阶的Adams-Moulton公式系数为： ) 288 95 , 1440 1427 , 240 133 , 720 241 , 1440 173 , 160 3 , 0(),1 , 1, 0 , 0 , 0 , 0 , 0( (3.4) 具有绝对稳定区间)0 ,1842. 1(，局部截断误差系数为0143. 0 60480 863 7 C。 3.2.2 Adams-Moulton 改进公式的构造 3阶Adams-Moulton公式的改进 1) 构造改进公式 现增加一项 2 15 k f h ，仍保持为3阶，重新求解方程组（2.4），得到 ) 60 29 , 15 7 , 60 7 , 15 1 (),1 , 1, 0 , 0( (3.5) 重庆大学硕士学位论文 3 线性多步法的改进 18 得到一个3步3阶的Adams公式： )472829( 60 2111 kkkkkk ffff h yy (3.6) 它的阶仍为3，误差系数为1083. 0 120 13 4 C。 2) 稳定性讨论 定理定理 3.2.1 公式（3.6）是)(A稳定的公式。 证明: 对应公式（3.6），由（3.5）有 23 )(zzz，其零点为1 , 0 , 0，满足根 条件，且由0 10 CC，公式（3.6）相容，故公式（3.6）是收敛的。 下面验证定理2.8.2的条件： 1 由上有公式（3.6）满足根条件，条件（1）满足； 20 29 60 , 0 60 29 3 3 3 条件（2）满足； 3由（3.5） 15 1 60 7 15 7 60 29 )( 23 zzzz 其零点为2518. 0,4212. 06086. 0 32, 1 ziz，17402. 0max 31 j j z， 条件 （3） 满足； 4由（2.17-2.21）计算 15 1 , 60 11 , 15 2 , 60 1 3210 0, 15 1 , 60 11 , 6 5 3210 )52118( 30 1 6 1 15 1 30 11 15 4 )( 2323 3 xxxxxxxR )27118( 30 1 10 9 30 11 15 4 )( 22 3 xxxxxI 易判断，649. 2274. 1, 0)( 3 xxI，当然有1, 0)( 3 xxI成立，条件（4）满足： 5再由（2.18）得： 0)(,1 , 1|, )( )(1 )( 33 3 3 2 maxmin 3 xRxxA xR xIx tg Ax 计算得: ' maxmax 2778)8938. 4(,8938. 4)( arctgtg。 下面用边界轨迹法绘制了这个公式的绝对稳定区域图形，C是公式的边界轨 迹， 边界C所为的封闭区域以外无限大的区域 （ 1 R） 为公式的绝对稳定区域， 2 R为 不稳定区域。如图3.1。 重庆大学硕士学位论文 3 线性多步法的改进 19 -2024681012 -15 -10 -5 0 5 10 15 46 8 10 2 5 10 -5 -10 R2 R1 C Im Re 0 图3.5 改进公式的边界轨迹图和稳定区域 Figure 3.5 The figure of root locus and stability domain improved formula 4阶Adams-Moulton公式的改进 1) 构造改进公式 现增加一项 3 24 k f h ，仍保持为4阶，重新求解方程组（2.4），得到 ) 12 5 , 8 5 , 24 1 , 8 1 , 24 1 (),1 , 1, 0 , 0 , 0( (3.7) 得到一个4步4阶的Adams公式： )31510( 24 32111 kkkkkkk fffff h yy (3.8) 它的收敛阶仍为4，误差系数为0681. 0 720 49 5 C。 易得公式（3.8）的特征多项式为： ) 1510() 1( 24 1 ) 131510( 24 1 )(,)( 2223434 zzzzzzzzzzz (3.9) 称它们分别为公式的第一特征多项式和第二特征多项式。 对 34 )(zzz，显然有0, 1 4321 zzzz，故公式满足根条件，又由于 0 10 CC就是相容性条件，故公式（3.8）是收敛的。 2) 稳定性讨论 我们给出定理： 重庆大学硕士学位论文 3 线性多步法的改进 20 定理定理3.2.2 公式（3.8）是 0 A稳定的。 证明：公式的稳定多项式为： )()(),(zzz (3.10) 即 ) 1510() 1( 24 )(),( 2234 zzzzzz (3.11) 现需要证明的是对任意的0，都有) 1510() 1( 24 )(),( 2234 zzzzzz 的零点4 , 3 , 2 , 1,jzj满足： 1ma x 41 j j z (3.12) 作变换： 1 1 , 1 1 z z z (3.13) 并记 )( 4 1 )( 4 z z r (3.14) 将z平面上单位园内1z变换到平面上的开左半平面0)Re(。经计算得 8 1 8 3 24 49 16 32 16 )( 234 r (3.15) 易知对任意的0，都有)(r的各项系数大于0。下面生成Routh序列： 00 8 1 5 00 )9136(32 9612045 4 0 8 1 )32(12 9136 3 0 8 3 16 32 2 8 1 24 49 16 1 2 2 2 注意对任意的0，都有 , 09612045, 09136 , 032 , 049 22 可知Routh序列的第一列全大于0，由Routh定理39，知)(r的零点都在开左平 面内，由变换(3.13)知，),(z的零点都在单位园内，定理3.2.2得证。 考查1z的情况，令), 0(, i ez，注意( ) z在)( 1z时有一个二 重奇点, 对充分小的0，存在0，使得( ) z在 1 ,1z，即, 0 上是正则的，用（2.20）和（2.21）分别计算( ) z的实部和虚部： 重庆大学硕士学位论文 3 线性多步法的改进 21 ) 1() 1( 3 1 3 2 3 2 3 1 )( 334 4 xxxxxxR )762)(1( 6 7 6 1 3 2 3 1 )( 223 4 xxxxxxxI 注意到 1 ,1, 0)(, 0)( 44 xxIxR，再由： 0)(,1 ,1|, ) 1( )762(1 )( 44 3 22 maxmin 4 xRxxA x xxx tg Ax 在0 时， 0)762(1 22 xxx，而0) 1( 3 x，有 max ()0tg ， max 0 ，所 以公式不是)(A稳定的，也不是)0(A稳定的，只是 0 A稳定的。 5阶Adams-Moulton公式的改进 现增加一项 4 40 k f h ，仍保持为5阶，重新求解方程组（2.4），得到 ) 720 269 , 180 139 , 60 7 , 360 37 , 720 71 , 40 1 (),1 , 1, 0 , 0 , 0 , 0( 得到一个5步5阶的Adams公式： )18717484556269( 720 432111 kkkkkkkk ffffff h yy (3.16) 它的收敛阶仍为5，局部截止断误差系数为0437. 0 160 7 6 C。它的绝对稳定区 间为)0 ,9231. 6(。 6阶Adams-Moulton公式的改进 作类似的改进可得： ) 1 , 1, 0 , 0 , 0 , 0 , 0( ) 5472 1901 , 27360 24233 , 4560 1327 , 13680 221 , 27360 3913 , 3040 263 , 57 1 ( 得到一个6步6阶的Adams公式： 1112345 19012433132722139132631 () 5472547245601368027360304057 kkkkkkkkk yyhfffffff 它的收敛阶仍为6，局部截止断误差系数为0318. 0 1149120 36557 7 C。它的绝对稳 定区间为)0 ,5331. 3(。 3.2.3 改进前后公式的稳定域及误差比较 图3.2是3阶Adams-Moulton公式和改进公式的边界轨迹和稳定区域图。 1 C是3 阶Adams-Moutlon公式的边界轨迹， 2 C是改进公式的边界轨迹， 1 R是3阶 Adams-Moulton公式的稳定区域， 2 C所界定的连通域 2 R(包含 1 R的一部分)是改进 公式的稳定区域。可见改进公式的稳定区域有很大的扩大。注意 ，大于 改进前的 ，但小于常用的3阶BDF公式的 。 4 13 120 C 1 24 3 22 重庆大学硕士学位论文 3 线性多步法的改进 22 图3.2 3阶Adams-Moulton公式和改进公式的边界轨迹图和稳定区域 Figure 3.2 The figure of root locus and stability domain Adams-Moulton formula of order 3 and improved formula 图3.3是4阶Adams-Moulton公式和改进公式的边界轨迹和稳定区域图。 1 C是4 阶Adams-Moulton公式的边界轨迹， 2 C是改进公式的边界轨迹， 1 R是4阶 Adams-Moulton公式的稳定区域， 2 C所界定的连通域 2 R(包含 1 R的一部分)是改进 公式的稳定区域。可见改进公式是 0 A稳定的。可见改进公式的稳定区域有很大的 扩大。 注意 720 49 5 C，大于改进前的 19 720 ，但小于常用的4阶BDF公式的 125 12 。 图3.3 4阶Adams-Moulton公式和改进公式的边界轨迹图和稳定区域 Figure 3.3 The figure of root locus and stability domain Adams-Moulton formula of order 4 and improved formula 重庆大学硕士学位论文 3 线性多步法的改进 23 图3.4，3.5的两个图，分别是5，6阶Adams-Moulton公式和改进公式的边界轨迹 和稳定区域图。 1 C是5，6阶Adams-Moulton公式的边界轨迹， 2 C是改进公式的边 界轨迹， 1 R是5，6阶Adams-Moulton公式的稳定区域， 2 C所界定的连通域 2 R(包 含 1 R的一部分)是改进公式的稳定区域。 可见它们的改进公式的稳定区域都有明显 扩大。 注意 6 7 160 C ，大于改进前的 3 160 ，但小于常用的5阶BDF公式的 10 137 ， 7 0.0318C ，大于改进前的0.0143，小于常用的6阶BDF公式的0.0709。 图3.4 5阶Adams-Moulton公式和改进公式的边界轨迹图和稳定区域 Figure 3.4 The figure of root locus and stability domain Adams-Moulton formula of order 5 and improved formula 图3.5 6阶Adams-Moulton公式和改进公式的边界轨迹图和稳定区域 Figure 3.5 The figure of root locus and stability domain Adams-Moulton formula of order 6 and improved formula 重庆大学硕士学位论文 3 线性多步法的改进 24 3.2.4 数值试验 对试验方程：1)0(, ' yyy作数值试验，用4步4阶改进公式（3.8）计算，3 步4阶Adams-Moulton公式作对比，启动值 3210 ,yyyy采用精确解 x ey 计算，取 固定步长1 . 0h调整，当30时，3步4阶Adams-Moulton公式发生明显振荡，而 改进公式始终保持很高的精度，图3.6、3.7是30和35时的对比效果。可见 改进公式确可用于刚性方程的求解。 00.20.40.60.811.21.41.61.82 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 精 确 解 改 进 的 4阶 Adams-moulton公 式 3步 4阶 Adams-moulton公 式 图3.6 30时的对比试验图 Figure 3.6 The figure of comparison numerical experiment while 30 00.20.40.60.811.21.41.61.82 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 精 确 解 改 进 的 4阶 Adams-moulton公 式 3步 4阶 Adams-moulton公 式 图3.7 35时的对比试验图 Figure 3.7 The figure of comparison numerical experiment while 35 重庆大学硕士学位论文 3 线性