变步长的龙格库塔法

1
y( x i1 )
y (h) i 1
16
由此可得
y(xi1 )
(h)
y2 i 1
1 15
(
(h)
y2 i 1
y
(h) i 1
)
h
这表明以
()
y2 i 1
作为
y(xi1) 的近似值，其误差可用步
长折半前后两次计算结果的偏差
来判断所选步长是否适当
y y ( h ) 2 i 1
(h) i 1
当要求的数值精度为ε时：
完成的步数就增加了。这样会引起计算量的增大，
并且会引起舍入误差的大量积累与传播。因此微分
方程数值解法也有选择步长的问题。
以经典的四阶龙格-库塔法(7.20)为例。从节点
xi出发，先以h为步长求出一个近似值，记为
y (h) i 1
，
由于局部截断误差为 O(h5 ) ，故有
y( xi1 )
y (h) i 1
j 1
(i 2, 3L , p)
确定原则是使近似公式在( xn , yn )处的Taylor展开式与
y( x)在xn处的Taylor展开式的前面项尽可能多地重合。
7.4.6 变步长的龙格-库塔法
在微分方程的数值解中，选择适当的步长是非常
重要的。单从每一步看，步长越小，截断误差就越
小；但随着步长的缩小，在一定的求解区间内所要
变步长的龙格—库塔方法
以经典四阶龙格—库塔公式为例。从节点xn出发，以h为
步长求一近似值yn(h1) ,
y( xn1 )
y(h) n 1
ch5 ,
将步长折半，即取
h 2
为步长从xn
跨两步到xn
，求一近似
1
值y ( h／2) n 1
,
每跨一步的截断误差是c
h 2
5
,因此有
y( xn1 )
y ( h／2) n 1
2c
h 2
5
,
由上两式
y( xn1 )
y ( h／2) n 1
1.
y( xn1 )
y(h) n 1
16
y( xn1 )
y ( h／2) n 1
1 15
[
y ( h／2) n 1
y(h) n 1
].
公式
Euler公式 改进Euler公式
yn1 K1
yn f ( xn ,
hK1 yn )
K4
f (xi h ,
yi
2 2
hK2Biblioteka Baidu
1
2 2
hK3 )
§2 Runge-Kutta Method
注：
龙格-库塔法的主要运算在于计算 Ki 的值，即计算 f 的
值。Butcher 于1965年给出了计算量与可达到的最高精 度阶数的关系：
每步须算Ki 的个数 2 3
4
5
6
可达到的最高精度 O(h2 ) O(h3 ) O(h4 ) O(h4 ) O(h5 )
ch5
当h值不大时，式中的系数c可近似地看作为常数。
然后将步长折半,即以为 h 2
步长,从节点xi出发,跨
两步到节点xi+1,再求得一个近似值
h ()
y2 i 1
,每跨一步的
截断误差是 c h 5 ,因此有
2
y( xi1 )
h
()
y(
xi
2 1
)
h
2c
h 5 2
()
这样
y( x i1 )
y2 i 1
（1）如果Δ>ε，反复将步长折半进行计算，直
至Δ<ε为止,并取其最后一次步长的计算结果作为y i 1
（2）如果Δ<ε，反复将步长加倍，直到Δ>ε
为止，并以上一次步长的计算结果yi作1 为
。
这种通过步长加倍或折半来处理步长的方法称为
变步长法。表面上看，为了选择步长，每一步都要
反复判断Δ，增加了计算工作量，但在方程的解y(x)
yn
1
yn
h( 1 2
K1
1 2
K2 )
K1
f ( xn , yn )
K2 f ( xn h, yn hK1 )
一般地，RK方法设近似公式为
p
yn1 yn h ci Ki
K1
i 1
f ( xn , yn )
i 1
Ki
f ( xn ai h, yn h bij K j )
7
O(h6 )
n8
O(hn2 )
由于龙格-库塔法的导出基于泰勒展开，故精度主要受
解函数的光滑性影响。对于光滑性不太好的解，最好 采用低阶算法而将步长h 取小。
➢ 变步长的Runge—Kutta Method
§2 Runge-Kutta Method
Rn1 ch p1 y( p1) xn
变化剧烈的情况下，总的计算工作量得到减少，结
果还是合算的。
➢ 高阶Runge—Kutta Method
§2 Runge-Kutta Method
y y h[ K K ... K ]
i 1
i
11
22
mm
K f(x , y )
1
ii
K f ( x h, y hK )
2
i
2
i
21
1
K f ( x h, y hK hK )
3
i
3
i
31
1
32
2
... ...
K f ( x h, y hK hK ... hK )
m
i
m
m1
1
m2
2
m m 1
m 1
其中i ( i = 1, …, m )，i ( i = 2, …, m ) 和 ij
( i = 2, …, m; j = 1, …, i1 ) 均为待定系数，确 定这些系数的步骤与前面相似。
K4 f ( xi h, yi hK3 )
➢ Gill公式：4阶经典龙格-库塔公式的一种改进
yi1
yi
h 6
K1
2
2 K2 2
2 K3 K4
K1 f (xi , yi )
K2
f
( xi
h 2
,
yi
h 2
K1)
K3
f
( xi
h 2
,
yi
2 1 2
hK1
1
2 2
hK2 )
§2 Runge-Kutta Method
➢ 最常用为四级4阶经典龙格-库塔法 /* Classical Runge-Kutta
Method */ ：
yi 1
yi
h 6
(K1
2K2
2K3
K4)
K1 f ( xi, yi )
K2
f
( xi
h2 ,
yi
h 2
K1)
K3
f
( xi
h2 ,
yi
h 2
K2)
Q: 由局部截断误差可以看出，步长 h 越小，局部截断 误差越小；但步长减小，在一定求解范围（区间）内 要完成的步数就增加了，步数增加会引起计算量增大， 导致舍入误差积累。因此要选取适当的步长。
选择步长时要考虑两个问题： 1.如何衡量和检验计算结果的精度？ 2.如何根据所获得的精度处理步长？
HW: p.201 #6-8