数值分析常微分方程的数值解法
《计算机数学基础》数值部分第五单元辅导
14 常微分方程的数值解法
一、重点内容 1. 欧拉公式：
),...,,,()
,()(1-210=???+=+=≈01+1+n k kh x x y x hf y y x y k
k k k k k
局部截断误差是O (h 2)。
2. 改进欧拉公式：
预报－校正公式：
??
???++=+=++++)],(),([2),(1111k k k k k k k k k k y x f y x f h
y y y x hf y y 校正值
预报值
即 ))],(,(),([2
11k k k k k k k k y x hf y x f y x f h
y y +++
=++ 或表成平均的形式：
???
?
?
????
+21=+=+=1+1+)(),(),(c p k p k k c k k k p y y y y x hf y y y x hf y y
改进欧拉法的局部截断误差是O (h 3)
3. 龙格－库塔法
二阶龙格－库塔法的局部截断误差是O (h 3) 三阶龙格－库塔法的局部截断误差是O (h 4) 四阶龙格?库塔法公式： )22(6
43211κκκκ++++
=+h
y y k k
其中 κ1=f (x k ,y k )；κ2=f (x n +12h ，y k +21h κ1)；κ3=f (x k +12h ，y n +2
1
h κ2)；κ4=f (x k +h ，y k +h κ3)
四阶龙格－库塔法的局部截断误差是O (h 5)。
二、实例
例1 用欧拉法解初值问题???1
=060≤≤0--='2)()
.(y x xy y y ，取步长h =0.2。计算过程保留
4位小数。 解h =0.2, f (x )=－y －xy 2。首先建立欧拉迭代格式
)
,,)((.),(210=-420=--=+=2
1+k y x y y hx hy y y x hf y y k k k k
k k k k k k k
当k =0，x 1=0.2时，已知x 0=0,y 0=1，有
y (0.2)≈y 1=0.2×1(4－0×1)＝0.8000
当k ＝1，x 2=0.4时，已知x 1=0.2, y 1=0.8，有 y (0.4)≈y 2=0.2×0.8×(4－0.2×0.8)＝0.614 4 当k =2，x 3=0.6时，已知x 2=0.4,y 2=0.6144，有 y (0.6)≈y 3=0.2×0.6144×(4－0.4×0.4613)＝0.8000
例2 用欧拉预报－校正公式求解初值问题???1
=10
=++'2)(sin y x y y y ，取步长h =0.2,计算
y (0.2),y (0.4)的近似值，计算过程保留5位小数。 解 步长h =0.2, 此时f (x ,y )=－y －y 2sin x 欧拉预报－校正公式为：
??
?
?
?++=+=++++)]
,(),([2
),(111
1k k k k k k k k k k y x f y x f h
y y y x hf y y 校正值预报值
有迭代格式：
??
????
?
??+--=--+--+=-=--+=++++++++)sin (1.0)sin 1.09.0()]sin ()sin [(2)sin 2.08.0()
sin (12
1112112
121k k k k k k k k k k k k k k k k k k k k k k x y y x y y x y y x y y h y y x y y x y y h y y 校正值预报值 当k =0，x 0=1, y 0=1时，x 1=1.2，有 631710=11?02-80?1=20-80=0001.)sin .()sin ..(x y y y
715490=21631710+63171010-1?1?10-90?1=≈212
1.).sin ..(.)sin ..().(y y 当k =1，x 1=1.2, y 1=0.71549时，x 2=1.4，有
47697
0=21715490?02-80?715490=20-80=1112.)
.sin ..(.)sin ..(x y y y
)
.sin ..(.).sin ...(.).(41476970+47697010-21?715490?10-90?715490=≈412
2
y y
=0.52608
例3 写出用四阶龙格－库塔法求解初值问题?
??2=03-8=')(y y
y 的计算公式，取步长h =0.2计
算y (0.4)的近似值。计算过程保留4位小数。
解 此处f (x ,y )=8－3y , 四阶龙格－库塔法公式为
)22(643211κκκκ++++=+h
y y k k
其中 κ1=f (x k ,y k )；κ2=f (x n +12h ，y k +21h κ1)；κ3=f (x k +12h ，y n +2
1
h κ2)；κ4=f (x k +h ，y k +h κ3)
本例计算公式为： )(.43211++2+2+6
2
0+
=κκκκk k y y
其中 κ1=8－3 y k ；κ2=5.6－2.1 y k ；κ3=6.32－2.37y k ; κ4=4.208＋1.578y k
)
1,...,2,1,0(5494.02016.1))
578.1208.4()37.232.6(2)1.26.5(238(6
2
.01-=+=-+-+-+-+
=+n k y y y y y y y k k k k k k k 当x 0=0,y 0==2,
4654
2=30042?54940+20161=54940+20161=≈403004
2=2?54940+20161=54940+20161=≈201201......).(.....).(y y y y y y
例4 设初值问题1=00=+')(,y y y ，证明用梯形公式求解该问题的近似解为
n
n h h y ??
? ??+2-2=
证明 解初值问题的梯形公式为
)],(),([2
111+++++=k k k k k k y x f y x f h
y y (k =0,1,2,…,n －1)
y y x f -=),(
][2
11++--+=∴k k k k y y h
y y
整理成显式
k k y h h y ??
?
??+-=+221( k =0,1,2,…,n －1)
用k =n ,n －1,n －2,…,1,0反复代入上式，得到
01
23121
22...222222y h h y h h y h h y h h y n n n n n +--+??
? ??+-==??? ??+-=??? ??+-=??? ??+-=
n
n h h y y ??
?
??+2-2=∴1
=0
例5 选择填空题：
1. 取步长h =0.1, 用欧拉法求解初值问题?????
1
=1+='2)(y y
x
y y 的计算公式是 答案：1,1,...,2,1,0],)1.01(1
.01.1[02
1=-=++=+y n k k y y k k
解答：欧拉法的公式
),...,,,()
,()(1-210=??
?+=+=≈01+1+n k kh
x x y x hf y y x y k k k k k k
此处y x
y
y x f +=2),(，迭代公式为 1=210=10+11
0+11=+10+110+=02
21
+y k k y y k y y y k k k k k ,...,,,),)
.(..()).((.
2. 改进欧拉法的平均形式公式是( )
(A)?????????+21=+=+=1+)(),(),(c p k p k k c k k k p y y y y x hf y y y x hf y y (B)???
?
?
????
+21=+=+=1+1+1+)
(),(),(c p k p k k c k k k p y y y y x hf y y y x hf y y
(C)?????????+2=+=+=1+1+)(),(),(c p k p k k c k k k p y y h y y x hf y y y x hf y y (D)???
?
?
????+21=+=+=1+1+)
(),(),(c p k p k k c k k k p y y y y x hf y y y x hf y y
答案：(D)
解答：见改进欧拉法平均形式公式。 三、练习题
1.求解初值问题???=='00
y x y y x f y )()
,(欧拉法的局部截断误差是( ); 改进欧拉法的局部截断
误差是( ); 四阶龙格－库塔法的局部截断误差是( )
(A)O (h 2) (B)O (h 3) (C)O (h 4) (D)O (h 5)
2. 改进欧拉预报－校正公式是
????
?2+=+=1+1+]
[h
y y y y k k k k 校正值
预报值
改进欧拉法平均形式公式为y p = , y c = ，y k +1= 试说明它们是同一个公式。
3. 设四阶龙格－库塔法公式为
)22(643211κκκκ++++=+h
y y k k
其中 κ1=f (x k ,y k )；κ2=f (x n +12h ，y k +21h κ1)；κ3=f (x k +12h ，y n +2
1
h κ2)；κ4=f (x k +h ，y k +h κ3)
取步长h =0.3，用四阶龙格－库塔法求解初值问题?
??0=0-1=')(y y
y 的计算公式
是 。
4.取步长h =0.1, 用欧拉法求解初值问题?????1
=01≤≤02
1
=')()
(y x xy y
5. 试写出用欧拉预报－校正公式求解初值问题?
??1=00
=+')(y y y 的计算公式，并取步长
h =0.1，求y (0.2)的近似值。要求迭代误差不超过10－
5。
6. 对于初值问题???1
=0='2
)(y xy y 试用(1)欧拉法；(2)欧拉预报－校正公式；(3)四阶龙格－库
塔法分别计算y (0.2),y (0.4)的近似值。
7. 用平均形式改进欧拉法公式求解初值问题???0
=0=+')(y x
y y 在x =0.2,0.4,0.6处的近似值。
8. 证明求解初值问题的梯形公式是
y k +1=y k +)],(),([2
11+++k k k k y x f y x f h , h =x k +1－x k (k =0,1,2,…,n －1)， 四、练习题答案 1. (A), (B), (D)
2. ),(k k y x hf ; ),(),(1+1++k k k k y x f y x f y k ＋),(k k y x hf ; )(2
1
),(1c p p k k y y y x f y +++ 只需将y c , y p 的表达式代入到y k +1中，就得到预报－校正公式。 3. )1,...,2,1,0(2591625.12591625.0-=+n k y k
提示：其中 κ1=1－y k ；κ2=0.85(1－y k )；κ3=0.8725(1－y k )；κ4=0.73825(1－y k )
)
1,...,2,1,0(2591625.12591625.0))
1(73825.0)1(745.1)1(7.11(63
.01-=+=-+-+-+-+=+n k y y y y y y y k k k k k k k
4.y 1=1，y 2=1.005 000，y 3=1.010 025，y 4=1.025 175，y 5=1.045 679，
y 6=1.078 21，y 7=1.103 976，y 8=1.142 615，y 9=1.188 320，y 10=1.241 794
5. 计算公式为
)1,...,2,1,0(05.095.09.01
11-=????
?-==+++n k y y y y y k k k k
k 22211≈20???809
8180=286
8140=??
?9050=90=y y y y y y ).(,..,.. 6.欧拉法：y (0.2)≈1.000 00; y (0.4)≈1.080 00
欧拉预报－校正公式：y (0.2)≈1.020 84; y (0.4)≈1.042 40 四阶龙格－库塔法：y (0.2)≈1.002 673 ; y (0.4)≈1.021 798 7. y p =0, y c =0.04, y 1=0.02；
y p =0.056, y c =0.0888, y 2=0.0724； y p =0.13792, y c =0.164816, y 3=0.151368 8. 提示：见教材关于梯形公式的推导。