第八章 刚性方程组及其数值计算
武汉大学数学与统计学院
考虑如下线性常微分方程组,
3,,
(0 ) ( 2,1,2 ),T
y M y y R
y
? ??
?
其中
0.1 49.9 0
0 50 0
0 70 30000
M
????
??
??
?? ?
??
这里矩阵 M的特征值为
1 2 30, 1,5 0,3 0 0 0 0,? ? ?? ? ? ? ? ?
上述初值问题的精确解是,
0,1 5 0
1
50
2
5 0 3 0 0 0 0
3
( ),
( ),
( ),
tt
t
tt
y t e e
y t e
y t e e
??
?
??
??
?
??
显然当 t ? ?? 时解的各个分量 ( ),1,2,3
iy t i ?
是指数衰减的,并趋于稳态解
1 2 3(,,) ( 0,0,0 ),y y y ?
1 2 3( ),( ),( )y t y t y t
趋于稳态解 (0,0,0) 的速度是
由因子
0.1te? 决定的,
假如试图利用四级 Runge-Kutta方法求解上述初
值问题,要求计算直至得到符合精度要求的稳态
解为止,我们讨论计算过程可能遇到的问题,
一,稳定性要求,
4
3
2, 7 8,1,2,3,
2, 7 8
3 0 0 0 0 1 0,
30000
ihi
h
?
? ?
??
? ? ? ?
二,为使解充分接近稳态解只需要,
0,1 0.te ? ?
0, 1 4
40.
tee
t
?? ?
??
而实际上 1t ? 后
50 30 00 0,ttee??
已经不起作用了 !!!
5
4
40
40 40
4 10
10
t
N
h ?
?
? ? ? ? ?
往后的计算我们当然希望使用大步长 !但由
于稳定性要求,仍要用小步长,从而耗费了巨
大的计算量,并且误差积累的影响也随着计
算步数的增加越来越严重,
N是计算步数
410N ?
上述例子中,系数矩阵的特征值虽然都
是负数,但绝对值相差非常悬殊,
考虑 n维非线性常微分方程组
( ) (,( ) ),0,( 1 )
( ), [0,],,[0,],n n n
y t f t y t t T
y t T R f T R R
? ? ? ?
? ? ?
设 ()yt 是 (1)定义在 [0,T]上的解,并满足
0( 0 ),yy ?
现用 ( ) ( ) ( )y t y t z t?? 表示 (1)在 ()yt 附近的解,
则 ()zt 应满足
( ) (,( ) ) ( ),0,( 2 )yz t f t y t z t t T?? ? ? ?
(1)在
()yt
处的线性化方程
记矩阵
( ) (,( ) )yJ t f t y t??
的特征值为
( ),1,2,.,,,.i t i n? ?
若条件
R e ( ) 0,1,2,.,,,,0, ( 3 )i t i n t T? ? ? ? ?
称解 ()yt 为局部稳定,否则是不稳定的,
定义,设 ( ),[ 0,]y t t T? 是方程组 (1)的一个解,
假定相应 Jacobi矩阵
()Jt
的特征值满足 (3),并且
R e ( ) R e ( ),[ 0,]m a x m i nii
ii
t t t T?? ?? ? ?
则称 (1)在 ()yt 附近为 刚性方程组,
R e ( )
R e ( )
m a x
m i n
i
i
i
i
t
t
?
?
刚性比
刚性方程组的解
快变部分 慢变部分
一般地说,隐型方法比显型方法具有更大
的绝对稳定区域,因此使用隐型方法求解
刚性方程组更为合适,
隐型 Runge-Kutta方法
Adams内插方法
武汉大学数学与统计学院
考虑如下线性常微分方程组,
3,,
(0 ) ( 2,1,2 ),T
y M y y R
y
? ??
?
其中
0.1 49.9 0
0 50 0
0 70 30000
M
????
??
??
?? ?
??
这里矩阵 M的特征值为
1 2 30, 1,5 0,3 0 0 0 0,? ? ?? ? ? ? ? ?
上述初值问题的精确解是,
0,1 5 0
1
50
2
5 0 3 0 0 0 0
3
( ),
( ),
( ),
tt
t
tt
y t e e
y t e
y t e e
??
?
??
??
?
??
显然当 t ? ?? 时解的各个分量 ( ),1,2,3
iy t i ?
是指数衰减的,并趋于稳态解
1 2 3(,,) ( 0,0,0 ),y y y ?
1 2 3( ),( ),( )y t y t y t
趋于稳态解 (0,0,0) 的速度是
由因子
0.1te? 决定的,
假如试图利用四级 Runge-Kutta方法求解上述初
值问题,要求计算直至得到符合精度要求的稳态
解为止,我们讨论计算过程可能遇到的问题,
一,稳定性要求,
4
3
2, 7 8,1,2,3,
2, 7 8
3 0 0 0 0 1 0,
30000
ihi
h
?
? ?
??
? ? ? ?
二,为使解充分接近稳态解只需要,
0,1 0.te ? ?
0, 1 4
40.
tee
t
?? ?
??
而实际上 1t ? 后
50 30 00 0,ttee??
已经不起作用了 !!!
5
4
40
40 40
4 10
10
t
N
h ?
?
? ? ? ? ?
往后的计算我们当然希望使用大步长 !但由
于稳定性要求,仍要用小步长,从而耗费了巨
大的计算量,并且误差积累的影响也随着计
算步数的增加越来越严重,
N是计算步数
410N ?
上述例子中,系数矩阵的特征值虽然都
是负数,但绝对值相差非常悬殊,
考虑 n维非线性常微分方程组
( ) (,( ) ),0,( 1 )
( ), [0,],,[0,],n n n
y t f t y t t T
y t T R f T R R
? ? ? ?
? ? ?
设 ()yt 是 (1)定义在 [0,T]上的解,并满足
0( 0 ),yy ?
现用 ( ) ( ) ( )y t y t z t?? 表示 (1)在 ()yt 附近的解,
则 ()zt 应满足
( ) (,( ) ) ( ),0,( 2 )yz t f t y t z t t T?? ? ? ?
(1)在
()yt
处的线性化方程
记矩阵
( ) (,( ) )yJ t f t y t??
的特征值为
( ),1,2,.,,,.i t i n? ?
若条件
R e ( ) 0,1,2,.,,,,0, ( 3 )i t i n t T? ? ? ? ?
称解 ()yt 为局部稳定,否则是不稳定的,
定义,设 ( ),[ 0,]y t t T? 是方程组 (1)的一个解,
假定相应 Jacobi矩阵
()Jt
的特征值满足 (3),并且
R e ( ) R e ( ),[ 0,]m a x m i nii
ii
t t t T?? ?? ? ?
则称 (1)在 ()yt 附近为 刚性方程组,
R e ( )
R e ( )
m a x
m i n
i
i
i
i
t
t
?
?
刚性比
刚性方程组的解
快变部分 慢变部分
一般地说,隐型方法比显型方法具有更大
的绝对稳定区域,因此使用隐型方法求解
刚性方程组更为合适,
隐型 Runge-Kutta方法
Adams内插方法
