第5章-常微分方程的数值解
基本思想:若微分方程有初始值 \(x_0, y_0\) ,则把微分方程转化为递推公式,从而递推出每个离散点的方程解
5.1 欧拉方法
已知:
\[\left\{
\begin{array}{l}
\frac{dy}{dx} = f(x,y) \\
y(x_0) = y_0
\end{array}
\right.
\]
通过近似
\[\frac{dy}{dx} = \frac{y_{n+1}-y_{n}}{h}
\]
从而
\[\begin{align}
\Rightarrow \quad & \frac{y_{n+1}-y_{n}}{h} = f(x_n, y_n)\nonumber \\
\Rightarrow \quad & y_{n+1} = y_{n} + hf(x_n, y_n)\nonumber \\
\end{align}
\]
此为显式欧拉法。
也可以将 \(f(x_n, y_n)\) 改为 \(f(x_{n+1}, y_{n+1})\)
\[\Rightarrow \quad \frac{y_{n+1}-y_{n}}{h} = f(x_{n+1}, y_{n+1})
\]
此为隐式欧拉法。
梯形公式
把显示欧拉法和隐式欧拉法加在一起平均一下就得到梯形公式
\[y_{n+1} = y_n + \frac{h}{2}[f(x_n, y_n) + f(x_{n+1},y_{n+1})]
\]
改进欧拉法(预测-校正法)
- 先用显式欧拉公式做预测,算出 \(\bar{y}_{n+1} = y_n + hf(x_n, y_n)\)
- 再将 \(\bar{y}_{n+1}\) 代入隐式梯形公式的右边作校正,得到
\[y_{n+1} = y_n + \frac{h}{2}[f(x_n, y_n) + f(x_{n+1},\bar{y}_{n+1})]
\]
预测-校正法具有2阶精度,稳定性高于显示欧拉法
局部截断误差和方法的阶
初值问题的单步法可用一般形式表示为:
\[y_{n+1} = y_n + h\phi(x_n, y_n, y_{n+1}, h)
\]
\(\phi\)称为增量函数。
例如,对显式欧拉法有 \(\phi(x_n, y_n, h) = f(x_n, y_n)\)
局部阶段误差:
设 \(y(x)\) 是初值问题的准确解,称
\[T_{n+1} = y(x_{n+1}) - y(x_n) - h\phi(x_n, y_n, y_{n+1}, h)
\]
为显式单步法的局部截断误差。
精度:
若某算法的局部截断误差为 \(O(h^{p+1})\) ,则称该算法有p阶精度。
欧拉法的局部截断误差
\[\begin{align}
T_{n+1} & = y(x_{n+1}) -y_{n+1} \nonumber \\
& = [y(x_n) + hy'(x_n)+ \frac{h^2}{2}y''(x_n) + O(h^3)] - [y_n + hf(x_n, y_n)] \nonumber \\
& = \frac{h^2}{2}y''(x_n)+O(h^3) \nonumber
\end{align}
\]
具有1阶精度。
5.2 龙格-库塔法
(略)
5.3 收敛性与稳定性
收敛性:
若某算法对于任意固定的 \(x = x_n = x_0 + nh\) ,当 \(h\rightarrow 0\) (同时 \(n\rightarrow \infty\) 时有 \(y_n \rightarrow y(x_n)\) ,则称该算法是收敛的。
【例】:
就初值问题
\[\left\{
\begin{array}{l}
y' = f(x,y) \\
y(0) = y_0
\end{array}
\right.
\]
考察欧拉显式格式的收敛性
【解】:
该问题的精确解未 \(y(x) = y_0 e^{\lambda x}\)
欧拉公式为 \(y_{n+1} = y_n + h\lambda y_n = (1+\lambda h)y_n \Rightarrow y_{n} = (1+\lambda h)^{n}y_0\)
又 \(x=x_n=nh\),则 \(n = x_n/h\) ,有
\[\begin{align}
& y_{n} = (1+\lambda h)^{x_n/h} \ y_0 = [(1+\lambda h)^{1/\lambda h}]^{\lambda x_n} \nonumber \\
\Rightarrow \quad & y_0 e^{\lambda x_n} = y(x_n) \nonumber
\end{align}
\]
故是收敛的
稳定性:
若某算法在计算过程中任一步产生的误差在以后的计算中都逐步衰减,则称该算法是绝对稳定的。
稳定性的检验:
- 用待检验算法计算线性微分方程 \(y' = \lambda y, \ Re(\lambda)<0\)。
- 假设只在初值处产生误差 \(\epsilon_0 = y_0 - \bar{y}_0\) ( \(\bar{y}\) 假设为算法求得的y)
- 得到 \(\epsilon_i\) 的递推式,并求得若要满足 \(\epsilon_i\) 越来越小, \(\bar{h} = \lambda h\) 应该如何取值。 \(\bar{h}\) 的取值范围称为绝对稳定区域。
【例】:
对于显式欧拉法,
\[\begin{align}
& y_{i+1} = y_i + h\lambda y_i = (1+\bar{h})^{i+1}y_0 \nonumber \\
& \epsilon_0 = y_0 - \bar{y}_0 \nonumber \\
\Rightarrow \quad & \bar{y}_{i+1} = (1+\bar{h})^{i+1}\bar{y}_0 \nonumber \\
\Rightarrow \quad & \bar{\epsilon}_{i+1} = y_{i+1} - \bar{y}_{i+1} = (1+\bar{h})^{i+1} \epsilon_0 \nonumber
\end{align}
\]
所以若要 \(\epsilon\) 逐渐衰减, \(\bar{h}\) 需要满足: \(|1+\bar{h}|<1\)