数值分析9-3 龙格-库塔方法

再看Taylor方法 一、再看 方法
h y(x + h) = y(x) + hy' (x) + y"(x) 2 h3 h4 (4) + y′′′(x) + y (x) 6 24 +⋯ 一般可取公式为如下形式
2
h h ( p) yn+1 = yn + hy'n + y"n +⋯+ yn 2! p!
2
d f ( x, y) dx 首先希望能确定系数 λ1、λ2、p，使得到的算法格式有 阶 ，使得到的算法格式有2阶 精度， 的前提假设下， ( x, y) 精度，即在 yi = y( x i ) 的前提假设下，使得 + f y ( x, y) dy = fx dx Ri = y ( x i +1 ) − y i +1 = O ( h 3 ) = f x ( x, y) + f y ( x, y) f ( x, y) y′′( x) =
第九章 常微分方程数值解 §2 龙格 - 库塔法 /* Runge-Kutta Method */
一、再看 再看Taylor方法 方法 二、Runge-Kutta 方法
dy = f ( x, y) x ∈[a, b] dx y(a) = y0
复习
h2 y(x + h) = y(x) + hy' (x) + y"(x) 2 h3 h4 (4) + y′′′(x) + y (x) +⋯ 6 24
{
}
Step 3: 将 yi+1 与 y( xi+1 ) 在 xi 点的泰勒展开作比较 点的泰勒 泰勒展开作比较
y i + 1 = y i + ( λ 1 + λ 2 ) h y ′( x i ) + λ 2 ph 2 y ′′( x i ) + O ( h 3 )
h2 y ( x i + 1 ) = y ( x i ) + h y ′( x i ) + y ′′( x i ) + O ( h 3 ) 2
由于龙格-库塔法的导出基于泰勒展开， 由于龙格 库塔法的导出基于泰勒展开， 库塔法的导出基于泰勒展开 故精度主要受解函数的光滑性影响。 故精度主要受解函数的光滑性影响。对 于光滑性不太好的解，最好采用低阶算 于光滑性不太好的解，最好采用低阶算 而将步长h 取小。 开式构造公式时，希望 近似值的泰勒展开式与精确值的泰勒展式有 尽可能多的项相同， 两者关于h 尽可能多的项相同 ， 即 ： 两者关于 h 的展式 中的系数尽可能多的项相同 尽可能多的项相同， 中的系数尽可能多的项相同，这样局部截断 误差与h 的更高次同阶， 误差与 h 的更高次同阶 ， 从而获得高精度的 近似公式。 近似公式。
= y ′ ( x i ) + ph y ′′ ( x i ) + O ( h 2 )
Step 2: 将 K2 代入第 式，得到 代入第1式
y i + 1 = y i + h λ 1 y ′( x i ) + λ 2[ y ′( x i ) + ph y ′′( x i ) + O ( h 2 )] = y i + ( λ 1 + λ 2 ) h y ′( x i ) + λ 2 ph 2 y ′′( x i ) + O ( h 3 )
yi +1 = yi + h ( K1 + 2K2 + 2K3 + K4 ) 6 K1 = f ( xi , yi ) K2 = f ( xi + h, yi + h K1 ) 2 2 K3 = f ( xi + h, yi + h K2 ) 2 2 K4 = f ( xi + h, yi + hK3 )
注： 龙格-库塔法 库塔法的运算在于计算 的值， 龙格 库塔法的运算在于计算 Ki 的值，即 的值。 计算 f 的值。Butcher 于1965年给出了计 年给出了计 算量与可达到的最高精度阶数的关系： 算量与可达到的最高精度阶数的关系：
每步须算K 每步须算 i 的 2 3 4 5 6 7 n≥8 个数 可达到的最高 4 3 4 5 6 n−2 O(h2 ) O(h ) O(h ) O(h ) O(h ) O(h ) O(h ) 精度
2
问题: 为获得更高的精度，应该如何进一步推广？ 问题 为获得更高的精度，应该如何进一步推广？
yi +1 = yi + h [ λ1 K 1 + λ2 K 2 + ... + λm K m ] K 1 = f ( xi , yi ) K 2 = f ( x i + α 2 h, yi + β 21 hK 1 ) K 3 = f ( x i + α 3 h, yi + β 31 hK 1 + β 32 hK 2 ) ... ... K m = f ( x i + α m h, y + β m1 hK 1 + β m 2 hK 2 + ... + βm m −1 hK m −1 )
3 则必须有： 要求 Ri = y( xi +1 ) − yi +1 = O( h ) ，则必须有：
1 λ1 + λ 2 = 1 , λ 2 p = 2
这里有 3 个未知 个方程。 数， 2 个方程。
存在无穷多个解。所有满足上式的格式统称为 阶龙格 存在无穷多个解。所有满足上式的格式统称为2阶龙格 - 库 无穷多个解 塔格式。注意到， 塔格式。注意到，p = 1, λ 1 = λ 2 = 1 就是改进的欧拉法。 就是改进的欧拉法。
Step 1: 将 K2 在 ( xi , yi ) 点作 Taylor 展开 用二元函数泰勒展开 展开(用二元函数泰勒展开 用二元函数泰勒展开)
K 2 = f ( xi + ph, yi + phK 1 ) = f ( xi , yi ) + phf x ( xi , yi ) + phK 1 f y ( xi , yi ) + O( h2 )
y i +1 K1 K2
= = =
1 1 yi + h K 1 + K 2 2 2 f ( xi , yi ) f ( x i + h , y i + hK 1 )
步长一定是一个h 步长一定是一个 吗？
将改进欧拉法推广为： 将改进欧拉法推广为：
yi +1 K1 K2 = = = yi + h [λ1 K 1 + λ2 K 2 ] f ( x i , yi ) f ( x i + ph, yi + phK 1 )
二元函数的Taylor展开 展开 二元函数的
一元函数的Taylor展开 展开 一元函数的
f (x, y) = f (x0 , y0 ) ′ f x (x0 , y0 )( x − x0 ) + ′(x0 , y0 )( y − y0 ) + f y f ′′(x , y )( x − x )2 0 x 0 0 1 ′ ′ + + 2 f xy (x0 , y0 )( x − x0 )( y − y0 ) +⋯ 2 2 ′ + f y′(x0 , y0 )( y − y0 )
二、Runge-Kutta 方法
建立高精度的单步递推格式。 建立高精度的单步递推格式。
单步递推法的基本思想是从 点出发， 单步递推法的基本思想是从 ( xi , yi ) 点出发，以某一斜 基本思想 率沿直线达到 ( xi+1 , yi+1 ) 点。欧拉法及其各种变形所 能达到的最高精度为2阶 能达到的最高精度为 阶。 考察改进的欧拉法，可以将其改写为： 考察改进的欧拉法，可以将其改写为： 斜率 一定取K1 K2 一定取 平均值吗 的平均值吗？
其中λi ( i = 1, …, m )，αi ( i = 2, …, m ) 和 ， βij ( i = 2, …, m; j = 1, …, i−1 ) 均为待定系 − 确定这些系数的步骤与前面相似。 数，确定这些系数的步骤与前面相似。
最常用为四级4阶经典龙格 库塔法 最常用为四级 阶经典龙格-库塔法 /* Classical Runge-Kutta Method */ ：