常微分方程的数值解法

常微分方程的分类

初值问题

$$\{\begin{array}{r}\frac{dy}{dx}=f(x,y)x\in [a,b]\\ y(x_0)=y_0\end{array}$$

定理只要存在正整数$L$，$f(x)$满足
$$|f(x,y_1)-f(x,y_2)|\le L|y_1-y_2|$$
李普希思条件，那么初值问题就存在唯一的$f(x)$

边值问题(非重点)

$$\{\begin{array}{r}{y}^{\prime \prime }=f(x,y,y^{\prime })a\le x\le b\\ y(a)=x_0,y(b)=y_n\end{array}$$

常微分方程数值解的基本思想

$$原表达式\{\begin{array}{r}{y}^{\prime }=f(x,y)\\ y_0=\eta \end{array}$$

化导数为差商的方法

解序列表达式
$$\{\begin{array}{r}{y}_0=\eta \\ y_{k+1}=y_k+hf(t_k,y_k)k=0,1,2,\dots n\end{array}$$
其中$t_k$表示已给定步长，第k个点

数值积分的求解思路

解序列的表达形式
$$\{\begin{array}{r}{y}_{k+1}=y_k+{\int }_{t_k}^{t_{k+1}}f(t,y)dt\\ 表达式2\end{array}$$
而积分可以用不同的数值积分的方法求得

泰勒展开法

将$y(x_{n+1})$在$x_n$点处进行泰勒展开，得到解序列
$$\{\begin{array}{r}{y}_0=\eta \\ y_{k+1}=y_k+hf(t_k,y_k)\end{array}$$

欧拉法

欧拉法的基本公式

解序列
$$\{\begin{array}{r}{y}_0=\eta \\ y_{k+1}=y_k+hf(t_k,y_k)\end{array}$$
根据解序列可以得到对应$x_1,x_2,x_3\cdots ,x_n$的函数值$y_1,y_2,\cdots ,y_n$的函数值

根据已得到的函数值，利用插值法构造原函数y(x)

向后欧拉法

$$\{\begin{array}{r}{y}_0=\eta \\ y_{n+1}=y_n+hf(x_{n+1},y_{n+1})\end{array}$$

向后欧拉法是隐式的公式，右边表达式中含有未知函数值，没有显示欧拉法计算方便

欧拉法只要满足李普希思条件就是收敛的
$$|y_{n+1}^{(k+1)}-y_{n+1}|<hL|y_{n+1}^k-y_{n+1}|$$
只要$hL<1$那么迭代式就是收敛的

梯形欧拉法

序列公式
$$y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},y_{n+1})]$$
梯形欧拉公式就是显示欧拉公式和隐式欧拉公示的算术平均

常微分方程的误差问题

由于已知常微分方程
$$y_{n+1}=y_n+hf(y_n,x_n)$$
这种形式，可知，$y_{n+1}$的值和$y_n$之间是有关联的而$y_n$,本身是有误差的故，误差会进行累积

局部误差

定义：在假设第i步精确的情况下第i+1步的误差称为局部截断误差
$$e_{i+1}(h)=y(x_{i+1})-y_i$$

累积截断误差(总体截断误差)

累积截断误差和总体截断误差是一个概念
$$E_k(h)=\sum _{i=1}^k{e}_i(h)$$
其中$e_i(h)$称为求解公式第$i$步的误差

$p$阶精度

若求解公式的局部截断误差等于
$$e_i(h)=O(h^{p+1})$$
则求解公式具有$p$阶精度

常用求解公式局部截断误差

欧拉公式的局部截断误差

$$e_{i+1}(h)=\frac{h^2}{2}{y}^{\prime \prime }({\xi }_i)$$

欧拉公式具有1阶精度

向后欧拉公式的局部截断误差

$$e_{i+1}(h)=-\frac{h^2}{2}{y}^{\prime \prime }({\xi }_i)$$

向后欧拉公式也具有1阶精度

梯形公式的截断误差

$$e_{i+1}(h)=-\frac{y^{\prime \prime \prime }({\eta }_{i+1})}{12}{h}^3$$

梯形公式具有2阶精度

梯形公式是迭代公式其收敛性与欧拉公式相似

中点欧拉公式（两步法）

$$y_{i+1}=y_{i-1}+2hf(x_i,y_i)$$

中点公式的截断误差
$$e_i=O(h^3)$$
中点公式具有2阶精度

几种求解公式直接的比较

改进的欧拉公式

先用显示的欧拉公式对$y_{n+1}$进行预测，再回代到梯形公式进行预测
$$y_{i+1}=y_i+\frac{h}{2}[f(x_i,y_i)+f(x_{i+1},y_i+hf(x_i,y_i))]$$

改进的欧拉法又称为预测-校证法，具有2阶精度

龙格库塔法

建立高精度的单步递推公式

2阶龙格-库塔法

2阶龙格库塔法满足的基本公式
$$\{\begin{array}{r}{y}_{i+1}=y_i+({\lambda }_1+{\lambda }_2)h{y}^{\prime }(x_i)+{\lambda }_2\rho h^2{y}^{\prime \prime }(x_i)\\ {\lambda }_1+{\lambda }_2=1{\lambda }_2\rho =\frac{1}{2}\end{array}$$

3阶龙格-库塔法

基本公式
$$
\left{
\begin{aligned}
y_{n+1}=y_n+\frac{h}{6}(K_1+4K_2+K_3)\
K_1=f(x_n,y_n)\
K_2=f(x_n+\frac{h}{2},y_n+\frac{h}{2}K_1)\
K_3=f(x_n+h,y_n+h(2K_2-K_1))
\end{aligned}

\right.
$$

单步法的收敛性和稳定性

收敛性

定义若某算法对任意固定的$x=x_i=x_0+ih$当$h\to 0同时i\to \infty 时候有y_i\to y(x_i)$则称该算法是收敛的

稳定性

扰动的基本定义:

$${\delta }_n=y(x_n)-y_n$$

$y_{(}{x}_n)$为数值算法得到的理想数值解

$y_n$为实际计算得到的数值解

稳定的定义

若在某一处的$x_n$点处得到的扰动为${\delta }_n$ 那么对于任意一点$m,m\ge n$都有
$$|{\delta }_m|\le |{\delta }_n|$$
则该数值方法就是稳定的

常见数值算法的稳定性

只考虑$y^{\prime }=\lambda y$的稳定性

欧拉法的稳定条件

$$0\le h\le -\frac{2}{\lambda }$$

​

隐式欧拉法稳定性

隐式欧拉法是恒稳定

线性多步法

线性多步法是利用已知$y(x)$在各点的近似值以及导数的近似值线性组合来逼近的

龙格库塔法是利用$y(x)$$在[x_i,x_{i+1}]$内若干导数预测值的线性组合来逼近的

这两种方法都是高精度的方法

Adamas方法

$$y_{n+k}=y_{n+k-1}+h\sum _{i=-1}^k{\beta }_i{f}_{n+i}$$

常用的Admans显示公式

常用的Admans隐式公式

米尔尼方法，汉明方法以及辛普森方法

米尔尼方法及其余项

$$\begin{gathered} y_{n+1}=y_{n-3}+\frac{4h}{3}(2f_n-f_{n-1}+2f_{n-2}) \\ 余项T_n=\frac{14}{45}{h}^5{y}^{(5)}({\xi }_n) \end{gathered}$$

汉明方法及其余项

$$\begin{gathered} y_{n+1}=\frac{1}{8}(9y_n-y_{n-2})+\frac{3h}{8}(f_{n+1}+2f_n-f_{n-1}) \\ 余项T_{n+1}=-\frac{1}{40}{h}^5{y}^{(5)}({\xi }_n) \end{gathered}$$

辛普森方法及其余项

$$\begin{gathered} y_{n+1}=y_{n-1}+\frac{h}{3}(f_{n+1}+4f_n+f_{n-1}) \\ 余项T_{n+1}=-\frac{1}{90}{h}^5{y}^{(5)}({\xi }_n) \end{gathered}$$

预测校证法

米尔宁汉明方法和4阶Adams 显隐式预测校证公式

$$\begin{gathered} y_{n+1}^p=y_{n-3}+\frac{4h}{3}(2f_n-f_{n-1}+2f_{n-2}) \\ {y}_{n+1}=\frac{1}{8}(9y_n-y_{n-2})+\frac{8h}{3}(f(x_{n+1},y_{n+1}^p)+2f_n-f_{n-1}) \end{gathered}$$

米尔公式和汉明公式的截断误差分别是

$$\begin{gathered} T_{n+1}^M=\frac{14}{45}{h}^5{y}^{(5)}(x_n)+O(h^6) \\ {T}_{n+1}^H=-\frac{1}{40}{h}^5{y}^{(5)}(x_n)+O(h^6) \end{gathered}$$

一阶方程组和高阶方程

常微分方程的边值问题的数值解法

边值条件的分类

• 第一边值条件 $y(a)=\alpha ,y(b)=\beta$
• 第二边值条件 $y^{\prime }(a)=\alpha ,y^{\prime }(b)=\beta$
• 第三边值条件 $y^{\prime }(a)=a_{0}y(a)+{\alpha }_1,y^{\prime }(b)=\beta y(b)+{\beta }_1$

边值的近似解，的三种基本求法

• 编制问题华为初值问题方法求解
• 差分法，用差商来替代微分方程以边值中的导数，化为代数方程求解
• 有限元法

差分法

5y^{{(5)}(x_n)+O(h}6)\
T_{n+1}^{H=-\frac{1}{40}h}5y^{{(5)}(x_n)+O(h}6)
$$

一阶方程组和高阶方程

常微分方程的边值问题的数值解法

边值条件的分类

• 第一边值条件 $y(a)=\alpha ,y(b)=\beta$
• 第二边值条件 $y^{\prime }(a)=\alpha ,y^{\prime }(b)=\beta$
• 第三边值条件 $y^{\prime }(a)=a_{0}y(a)+{\alpha }_1,y^{\prime }(b)=\beta y(b)+{\beta }_1$

边值的近似解，的三种基本求法

• 编制问题华为初值问题方法求解
• 差分法，用差商来替代微分方程以边值中的导数，化为代数方程求解
• 有限元法

差分法