数值分析_第九章_常微分方程初值问题数值解法

\qquad 常微分方程是描述连续变化的数学语言，微分方程的求解就是确定给定方程的可微方程$y(t)$。考虑一阶常微分方程的初值问题
$$\begin{gathered} y^{{}^{\prime }}=f(x,y),x\in [x_0,b],(1) \\ y(x_0)=y_0(2) \end{gathered}$$ \qquad 如果存在实数$L>0$，使得
$$|f(x,y_1)-f(x,y_2)|\le L|y_1-y_2|,\forall y_1,y_2\in \mathbb{R},$$则称$f$关于$y$满足利普西茨($Lipschite$)条件，$L$称为$f$的利普西茨常数(简称$Lisp.$常数)
\qquad 定理1 设$f$在区域 D = { ( x , y ) ∣ a ≤ x ≤ b , y ∈ R } D=\{(x,y)|a\leq x\leq b,y\in\mathbb{R}\} D={(x,y)∣a≤x≤b,y∈R}上连续，关于$y$满足利普西茨条件，则对任意$x_0\in [a,b]，y_0\in \mathbb{R}$，常微分方程初值问题满足(1)(2)式当$x\in [a,b]$时存在唯一的连续可微解$y(x)$。
\qquad 所谓数值解法，就是寻求解$y(x)$在一系列离散节点
$$x_1<x_2<\cdots <x_n<\cdots$$ \quad 上的近似值$y_1,y_2,\cdots ,y_n,\cdots$相邻两节点的间距$h_n=x_{n+1}-x_n$称为步长。求解过程首先对常微分方程离散化，建立求数值解的递推公式。当计算$y_{n+1}$时只用到前一点$y_n$，称为单步法。用到前面$k$个点的值$y_n,y_{n-1},\cdots$时称为$k$步法。
\qquad 欧拉法 设做出折线(由切线相连($x_0,x_1,\cdots x_n$处作切线))的顶点$P_n$，过$P_n(x_n,y_n)$依方向场的方向推进到$P_{n+1}(x_{n+1},y_{n+1})$，显然两个顶点$P_n,P_{n+1}$的坐标有关系
$$\frac{y_{n+1}-y_n}{x_{n+1}-x_n}=f(x_n,y_n)$$ \qquad 对微分方程式(1)中的导数用均差近似，即
$$\frac{y_{n+1}-y_n}{h}\approx y^{{}^{\prime }}(x_n)=f(x_n,y(x_n))$$ \qquad 依均差近似逐次计算。为分析计算公式的精度，使用泰勒展开将$y_{n+1}$在$x_n$处展开，则有
$$y(x_{n+1})=y(x_n+h)=y(x_n)+y^{{}^{\prime }}(x_n)h+\frac{h^2}{2}{y}^{{}^{\prime \prime }}({\xi }_n),{\xi }_n\in (x_n,x_{n+1})$$ \qquad 在$y_n=y_{x_n}$的前提下，$f(x_n,y_n)=f(x_n,y(x_n))=y^{{}^{\prime }}(x_n)$，于是欧拉法的误差
$$y(x_{n+1})-y_{n+1}=\frac{h^2}{2}{y}^{{}^{\prime \prime }}({\xi }_n)\approx \frac{h^2}{2}{y}^{{}^{\prime \prime }}(x_n)$$称为此方法的局部截断误差。显式单步法的局部截断误差一般性的表示方法为：
$$T_{n+1}=y(x_{n+1})-y(x_n)-h\phi (x_n,y(x_n),h)=O(h^{p+1})$$则称该方法具有$p$阶精度。
\qquad 如果对微分方程从$x_n$到$x_{n+1}$积分，得
$$y(x_{n+1})=y(x_n)+{\int }_{x_n}^{x_{n+1}}f(t,y(t))dt$$右端积分使用右矩形公式$hf(x_{n+1},y(x_{n+1}))$近似得到另一个公式
$$y_{n+1}=y_n+hf(x_{n+1},y_{n+1})$$通常使用迭代法求解(实现$x_n\to x_{n+1}$的运算)，设用欧拉公式
$$y_{n+1}^{(0)}=y_n+hf(x_n,y_n)$$给出迭代初值$y_{n+1}^{(0)}$带入上式得
$$y_{n+1}^{(1)}=y_n+hf(x_{n+1},y_{n+1}^{(0)})$$反复迭代得到
$$y_{n+1}^{(k+1)}=y_n+hf(x_{n+1},y_{n+1}^{(k)}),k=0,1,\cdots$$由利普西茨条件得($证明迭代法收敛$)
$$|y_{n+1}^{(k+1)}-y_{n+1}|=h|f(x_{n+1},y_{n+1}^{(k)})-f(x_{n+1},y_{n+1}^{})|\le hL|y_{n+1}^{(k)}-y_{n+1}|$$由此，只要$hL<1$使用迭代法就能收敛到解$y_{n+1}$。
\qquad 梯形方法(增加精度) 公式为$y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},y_{n+1})]$，迭代公式为：
$$\begin{gathered} y_{n+1}^{(0)}=y_n+hf(x_n,y_n) \\ {y}_{n+1}^{(k+1)}=y_n+\frac{h}{2}(f(x_n,y_n)+f(x_{n+1},y_{n+1}^{(k)}),k=0,1,2,\cdots \end{gathered}$$ \qquad 改进欧拉公式(增加了欧拉公式的精度，减小梯形算法迭代带来的计算次数)，建立预测-校正系统。

预测 ${\bar{y}}_{n+1}=y_n+hf(x_n,y_n)$

校正 $y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},{\bar{y}}_{n+1})]$
\qquad 龙格-库塔法，由公式
$$y(x_{n+1})=y(x_n)+{\int }_{x_n}^{x_{n+1}}f(t,y(t))dt$$要使近似精度提高就是提高右式的近似精度，近似公式为：
$${\int }_{x_n}^{x_{n+1}}f(t,y(t))dt\approx h\sum _{i=1}^r{c}_{i}f(x_n+{\lambda }_{i}h,y(x_n+{\lambda }_{i}h))$$分割的段越多，$r$越大，上式右端(相当于增量函数$\phi (x_n,y(x_n)h)$)的近似精度越高。由公式
$$y_{n+1}=y_n+h\phi (x_n,y_n,h)$$其中，
$$\begin{gathered} \phi (x_n,y_n,h)=\sum _{i=1}^r{c}_i{K}_i \\ {K}_1=f(x_n,y_n) \\ {K}_i=f(x_n+{\lambda }_{i}h,y_n+\sum _{j=1}^{i-1}{\mu }_{ij}{K}_j)i=2,3,\cdots ,r \end{gathered}$$ 通过给定局部截断误差精度求解未知系数。对于步长的选择，由给定的精度$\epsilon$和偏差$\Delta =|y_{n+1}^{(\frac{h}{2})}-y_{n+1}^{(h)}|$($折半前后两次的偏差$)决定。当精度大于偏差时，减小步长(对折)，反之则反。
\qquad 单步法的收敛性 假设单步法具有$p$阶精度，且增量函数$\phi (x,y,h)$关于$y$满足利普西茨条件
$$|\phi (x,y,h)-\phi (x,\bar{y},h)|\le L_{\phi }|y-\bar{y}|$$又设初值$y_0$是准确的，即$y_0=y(x_0)$，其整体的截断误差
$$y(x_n)-y_n=O(h^p)$$稳定性 由模型方程$y^{{}^{\prime }}=\lambda y$的欧拉公式为：
$$y_{n+1}=(1+h\lambda )y_n$$设在节点值$y_n$上有一个扰动值${\epsilon }_n$，它的传播使节点$y_{n+1}$产生${\epsilon }_{n+1}$的扰动值，扰动值满足
$${\epsilon }_{n+1}=(1+h\lambda ){\epsilon }_n$$为使差分方程的解非增(解稳定)，选择$h$充分小，使得$|1+h\lambda |\le 1$，即$-2<\lambda <0$。