常微分方程初值问题

第二部分 数值分析基础理论与算法
五、常微分方程初值问题
电子工程学院
1
常微分方程初值问题数值解法
5.1 基本理论
一阶常微分方程初值问题：
dy

dx
=
f
(x, y)
y(x0 ) = y0
电子工程学院
2
例：方程 xy' - 2 y = 4x ⇒ y' = 2 y + 4 x
令： f ( x, = y) 2 y + 4 且给出初值 y(1)= -3
17
欧拉方法的稳定性
分析微小误差在迭代过程中积累的情况。
问题：已知 y '( x) = f ( x, y), x ∈[a, b], y(0) = y0 + e, 其中e为初始值得误差，分析误差随迭代过程的变化情况。
讨论一：考虑简单情况，假设 f(x,y) ≡ αy，则欧拉迭代公式为：
yn+=1 yn + hα yn yn+1 = (1+ hα )n ( y0 + e) 积累误差： en+1= (1+ hα )n (e)
x
就得到一阶常微分方程的初值问题：
d=y f ( x, = y) 2 y + 4
dx
x
y(1) = −3
电子工程学院
3
解的存在及唯一性：
只要函数 f(x , y) 适当光滑连续，且关于y满足 李普希兹 ( L i p s c h i t z )条件，即存在常数 L,使得
f (x, y1) - f (x, y2 ) ≤ L y1 - y2 由常微分方程理论知，初值问题的解必存在且唯一。
y(xn+1) − h
y(xn )
两边都含有未知的
则隐式算法：
yn+1，∴隐式方程
yn+1 = yn + hf (xn+1, yn+1)

y(
x0
)
=
y0
n = 0,1,...
与前面欧拉公式的差异，
比较（前向）Euler公式： 由右矩形数值求积得到

y(x0 ) = y0 yn+1 = yn + hf (xn , yn
\欧拉折线法实际上用一系列直线组成的折线近似代
替曲线 y=y(x)，用折线交点处的yk 近似代替精确解 y(xk)
电子工程学院
10
注意：这是“折线法”而非“切线 法”除第一个点是曲线切线 外，其他点不是！
电子工程学院
11
误差分析：
近似值
欧拉公式的局部截断误差：y(xk+1) − yk+1
将y(x)在xk 处泰勒展开：
yn )
则有差分方程：
yn+1 − yn ≈ hf (xn , yn ) (n = 0,1,...)
电子工程学院
化为解差分方
程初值问题
7
C. Taylor级数法
在 xn附近将 y(x)按Taylor级数展开：
y( xn+1) = y( xn + h)
=
y(xn ) + hy′(xn ) +
h2 2
如果︱1+hα︱>1,则误差越来越大；反之，如果︱1+hα︱< 1，则误差越来越小。
讨论二：一般情况，如果对任意 ︱1+h f (x,y)︱< 1 ，则误差收敛。
电子工程学院
18
前向欧拉方法步长选择对解的影响
例：y '(x) = −15y, x ≥ 0, y(0) = 1 步长h=1/7
步长h=1/10 2.5
x n+1 f
xn
x,
y
dt
≈
h 2
(
f
(xn ,
yn ) +
f
( xn+1 ,
yn +1 ))
且令：yn+1 =y(xn+1), yn = y(xn ),
两边都含有未知的
[ ] 则得：yn+1
=yn
+
h 2
yn+1，∴隐式方程
f (xn , yn ) + f (xn+1, yn+1) ,
梯形法具有二阶精度
注意：y为x的函数。
则n阶泰勒法公式为：
y ( xn+=1 )
y ( xn
)
+
hf
( xn
,
yn
)
+
h2 2
f
'( xn
,
yn
)
+
...
+
hn n!
f (n-1) (xn , yn )
前向欧拉方法为一阶泰勒方法。 n阶泰勒方法具有n阶精度，但各阶
导数计算复杂，不实用。
y(x)= y(xk ) + y′(xk )(x
−
xk
)
+
y′′(xk ) 2!
(x
−
xk
)2
+ ...
令x=xk+1
y ( xk +1 )
= y(xk ) + hf (xk , y(xk ))+
--------------------------------------------
h2 2!
y′′( xk
(
yp
+
yc
).
0 -2
-40
1
2
3
4
5
电子工程学院
22
例 ：用欧拉公式和改进型欧拉公式解初值问题

y′

=
y-
2x y
y(0) = 1
(0 < x < 1)
解 ：取步长 h=0.1;
欧拉公式为： yn+=1
yn
+
h(
yn
-
2 xn yn
).

y=p

yn
+
h(
yn
-
2 xn yn
);
)
+
...
截断误差：
yk +1
y(xk +1) − yk +1 = o(h2 )
\精度为一阶
或：y( xk+1 ) − y= k+1
h2 y′′(ξn ) ≈
h2 2
y′′( xk )
电子工程学院
12
二、后向（后退的）Euler方法
利用向后差商：
= y′( xn +1 )
f (xn+1, yn+1) ≈
电子工程学院
4
dy

dx
=
f
(x, y)
y(x0 ) = y0
本章主要内容： （1）对常微分方程离散化； （2）建立求数值解的递推公式。
电子工程学院
5
微分方程离散化常用方法
A.用差商代替微商（导数）的方法
= dy dx ( xn , yn )
f
(xn ,
y(xn ))
≈
y(xn+1) − y(xn ) xn+1 − xn
)
n = 0,1,...
电子工程学院
13
为避免解非线性方程，与前向Euler法结合：

yn(0+=)1 y(k +1)
n+1
yn + hf (xn , yn ) = yn + hf (xn+1, yn(k+1) )
k = 0,1, 2,...
将前向欧拉折线公式的数值 解作为迭代初值，带入第2式
微分方程的数值解：
设方程解y(x)的存在区间是 [a,b]，令a= x0< x1<…< xn =b, 其中 hk=xk+1-xk , 若是等距节点 h=(b-a)/n , h称为步长。
当y(x)的解析表达式不容易得到或根本无法得到时，可
用数值方法求得y(x)在每个节点 xk 上函数值y(xk)的近似值， 用 yk 表示，即yk≈y(xk)，则称y0 , y1 ,...,yn为微分方程的数值解。
局部 误差
局部误差：初始点y(n)为精 确值条件下，估计y(n+1)时 产生的误差。
y( x(n+1) ) - y(n+1)
= x y(x(n+1) ) - y(n) + hf (x(n) , y(n) )
步长
y
积累
误差
由于单步误差在迭代过程中
产生积累，从而产生的误差为
积累误差。
x
步长
电子工程学院
几何意义：是欧拉折线法与后向欧拉法的算术平均，精度高
于前向/后向欧拉算法。
电子工程学院
20
四、改进的欧拉公式
梯形公式虽然提高了精度，但仍是隐式方法，算
法复杂，运算量大。而在实际计算中先预测右边的 yn+1，再利用梯形公式求解，这样建立的预测—校正 系统称作改进的欧拉公式。
预测 校正
yn+=1 yn + hf (xn , yn ); 代入
前向欧拉方法 的迭代公式
n = 0,1,...
几何意义：
由 (x0，y0)出发做曲线 y =y(x) 的切线（存在！），
该点切线斜率为：
dy dx (x0 , y0 ) = f (x0 , y0 )
电子工程学院
9
由于f (x0，y0) 及 (x0，y0) 已知，则切线方程为：