欧拉方法
对于方程
y
′
(
x
)
=
f
(
x
,
y
)
y'(x)=f(x,y)
y′(x)=f(x,y)
它的数值解是
y
i
+
1
=
y
i
+
h
f
(
x
i
,
y
i
)
y_{i+1}=y_i+hf(x_i,y_i)
yi+1=yi+hf(xi,yi),给定一个
y
0
y_0
y0,可以求出一系列
y
i
y_i
yi
它来源于导数的近似解
y ′ ( x i ) ≈ y ( x i + 1 ) − y ( x i ) h y'(x_i)≈\frac{y(x_{i+1})-y(x_i)}{h} y′(xi)≈hy(xi+1)−y(xi)
中点法
y
i
+
1
=
y
i
+
h
f
(
x
i
+
h
2
,
y
i
+
h
2
f
(
x
i
,
y
i
)
)
y_{i+1}=y_i+hf(x_i+\frac{h}{2}, y_i+\frac{h}{2}f(x_i,y_i))
yi+1=yi+hf(xi+2h,yi+2hf(xi,yi))
它来源于导数的另一种近似解
y ′ ( x i ) ≈ y ( x i + h 2 ) − y ( x i − h 2 ) h y'(x_i)≈\frac{y(x_{i}+\frac{h}{2})-y(x_i-\frac{h}{2})}{h} y′(xi)≈hy(xi+2h)−y(xi−2h)
改进欧拉方法
y
i
+
1
=
y
i
+
h
2
[
f
(
x
i
,
y
i
)
+
f
(
x
i
+
1
,
y
i
+
1
)
]
y_{i+1}=y_i+\frac{h}{2}[f(x_i,y_i)+f(x_{i+1},y_{i+1})]
yi+1=yi+2h[f(xi,yi)+f(xi+1,yi+1)]
来源于梯形积分的梯形近似
y
(
x
i
+
1
)
−
y
(
x
i
)
=
∫
x
i
i
+
1
f
(
x
,
y
)
d
x
y(x_{i+1})-y(x_i)=\int^{i+1}_{x_i}f(x,y)dx
y(xi+1)−y(xi)=∫xii+1f(x,y)dx
改进欧拉方法等式右端也有待求的
y
i
+
1
y_{i+1}
yi+1项,他可以用欧拉方法先出来,然后代入等式再计算新的解
龙格-库塔方法
欧拉是泰勒展开保留零次和一次项,高阶泰勒方法可以保留更多次项以取得更高的精度
泰勒展开: y ( x i + 1 ) = y ( x i ) + h y ′ ( x i ) + h 2 2 y ′ ′ ( x i ) + . . . y(x_{i+1})=y(x_i)+hy'(x_i)+\frac{h^2}{2}y''(x_i)+... y(xi+1)=y(xi)+hy′(xi)+2h2y′′(xi)+...
数值微分可以用一些等间隔点线性组合近似各级导数;
f
(
n
)
(
x
)
≈
∑
a
k
f
(
x
k
)
f^{(n)}(x)≈\sum a_kf(x_k)
f(n)(x)≈∑akf(xk)
把泰勒展开和各级导数的近似融一起有这种形式的微分方程
y
i
+
1
=
y
i
+
h
∑
a
k
f
(
x
k
,
y
k
)
y_{i+1}=y_i+h \sum a_kf(x_k,y_k)
yi+1=yi+h∑akf(xk,yk) 这里的
a
k
a_k
ak和上面的
a
k
a_k
ak不等,它是各级导数的
a
k
a_k
ak的和
把
x
k
,
y
k
x_k,y_k
xk,yk按照等间隔h展开有
y
i
+
1
=
y
i
+
h
[
a
1
f
(
x
i
,
y
i
)
+
a
2
f
(
x
i
+
λ
2
h
,
y
i
+
u
2
h
)
+
.
.
.
+
a
k
f
(
x
i
+
λ
k
h
,
y
i
+
u
k
h
)
]
y_{i+1}=y_i+h[a_1f(x_i,y_i)+a_2f(x_i+\lambda_2h,y_i+u_2h)+...+a_kf(x_i+\lambda_kh,y_i+u_kh)]
yi+1=yi+h[a1f(xi,yi)+a2f(xi+λ2h,yi+u2h)+...+akf(xi+λkh,yi+ukh)]
另外下一个点的y通常是之前的y沿着导数方向移动而来,所以
u
i
=
∑
β
f
(
x
k
,
y
k
)
u_i=\sum \beta f(x_k,y_k)
ui=∑βf(xk,yk)
龙格-库塔方法可以写成如下形式:
{
y
i
+
1
=
y
i
+
h
[
a
1
K
1
+
a
2
K
2
+
.
.
.
]
K
1
=
f
(
x
i
,
y
i
)
K
2
=
f
(
x
i
+
λ
2
h
,
y
i
+
h
β
21
K
1
)
.
.
.
K
p
=
f
(
x
i
+
λ
p
h
,
y
i
+
h
∑
j
=
1
p
−
1
β
p
j
K
j
)
\begin{cases} y_{i+1}=y_i+h[a_1K_1+a_2K_2+...] \\ K_1=f(x_i,y_i)\\ K_2=f(x_i+\lambda_2h,y_i+h\beta_{21}K_1)\\ ...\\ K_p=f(x_i+\lambda_ph,y_i+h\sum^{p-1}_{j=1} \beta_{pj} K_j) \end{cases}
⎩
⎨
⎧yi+1=yi+h[a1K1+a2K2+...]K1=f(xi,yi)K2=f(xi+λ2h,yi+hβ21K1)...Kp=f(xi+λph,yi+h∑j=1p−1βpjKj)
这就是龙格-库塔方法的一般形式
如果
p
=
1
,
a
=
1
p=1,a=1
p=1,a=1,它就是欧拉方法
p
=
2
,
a
1
=
0
,
a
2
=
1
,
λ
2
=
0.5
,
u
2
=
0.5
f
(
x
i
,
y
i
)
p=2,a_1=0,a_2=1,\lambda_2=0.5,u_2=0.5f(x_i,y_i)
p=2,a1=0,a2=1,λ2=0.5,u2=0.5f(xi,yi)是中点法
p
=
2
,
a
1
=
0.5
,
a
2
=
0.5
,
λ
2
=
1
,
u
2
=
f
(
x
i
,
y
i
)
p=2,a_1=0.5,a_2=0.5,\lambda_2=1,u_2=f(x_i,y_i)
p=2,a1=0.5,a2=0.5,λ2=1,u2=f(xi,yi)是改进欧拉方法
龙格-库塔方法里面的系数
a
i
,
β
p
j
a_i, \beta_{pj}
ai,βpj可以用系数对照法建立关系式求出
对
K
p
K_p
Kp进行二元泰勒展开得到一个关系式,对
y
(
x
i
+
1
)
y(x_{i+1})
y(xi+1)进行泰勒展开得到另一个关系式,两相对比可以建立系数
a
i
,
β
p
j
a_i, \beta_{pj}
ai,βpj之间的关系
二阶龙格-库塔方法算出来是
{
a
1
+
a
2
=
1
a
2
λ
2
=
1
/
2
a
2
β
21
=
1
/
2
\begin{cases} a_1+a_2=1\\ a_2\lambda_2=1/2\\ a_2\beta_{21}=1/2 \end{cases}
⎩
⎨
⎧a1+a2=1a2λ2=1/2a2β21=1/2
他有无穷多解,中点法,改进欧拉方法都满足这个关系
二阶龙格-库塔方法算出来是
{
a
1
+
a
2
+
a
3
=
1
a
2
λ
2
+
a
3
λ
3
=
1
/
2
a
2
λ
2
2
+
a
3
λ
3
2
=
1
/
3
a
3
λ
2
β
32
=
1
/
6
a
2
=
β
21
a
3
=
β
31
+
β
32
\begin{cases} a_1+a_2+a_3=1\\ a_2\lambda_2+a_3\lambda_3=1/2\\ a_2\lambda_2^2+a_3\lambda_3^2=1/3\\ a_3\lambda_2\beta_{32}=1/6\\ a_2=\beta_{21}\\ a_3=\beta_{31}+\beta_{32} \end{cases}
⎩
⎨
⎧a1+a2+a3=1a2λ2+a3λ3=1/2a2λ22+a3λ32=1/3a3λ2β32=1/6a2=β21a3=β31+β32
线性多步法
阿达姆斯线性多步法从数值积分角度理解欧拉方法,并利用更多点近似积分取得更高精度
阿达姆斯线性多步法
y
i
+
1
=
y
i
+
h
∑
j
=
0
r
β
j
(
r
)
f
(
x
i
−
j
,
y
i
−
j
)
y_{i+1}=y_i+h\sum_{j=0}^r\beta_j^{(r)}f(x_{i-j},y_{i-j})
yi+1=yi+h∑j=0rβj(r)f(xi−j,yi−j)
这里面的
β
j
(
r
)
\beta_j^{(r)}
βj(r)可以查表
r\j | 0 | 1 | 2 | 3 | 4 |
---|---|---|---|---|---|
0 | 1 | ||||
1 | 3 2 \frac{3}{2} 23 | - 1 2 \frac{1}{2} 21 | |||
2 | 23 12 \frac{23}{12} 1223 | - 16 12 \frac{16}{12} 1216 | 5 12 \frac{5}{12} 125 | ||
3 | 55 24 \frac{55}{24} 2455 | - 59 24 \frac{59}{24} 2459 | 37 24 \frac{37}{24} 2437 | - 9 24 \frac{9}{24} 249 | |
4 | 1901 720 \frac{1901}{720} 7201901 | − 2774 720 -\frac{2774}{720} −7202774 | 2616 720 \frac{2616}{720} 7202616 | − 1274 720 -\frac{1274}{720} −7201274 | 251 720 \frac{251}{720} 720251 |