微分方程数值解 第一章（1）
Euler法
对于常微分方程（ODE）问题：
$$
\begin{cases}
u’ = f(t, u), & 0 \leq t \leq T \
u(0) = u_0
\end{cases}
$$
其中 $f(t, u)$ 关于 $u$ 满足 Lipschitz 条件
欧拉法推导：
连续问题离散化
设真解为 $u(t)$。将区间 $[0, T]$ 离散化，令 $\displaystyle \frac{T}{N}=h$，其中 $h$ 为步长，且 $0 = t_0 < t_1 < \cdots < t_N = T$，$t_n=nh$
已知 $u(t_0)=u_0$，$u’(t_0)=f(t_0, u_0)$
利用 Taylor 展开推导近似值
对 $u(t_1)$ 进行 Taylor 展开：$\displaystyle u(t_1)=u(t_0 + h)=u(t_0)+u’(t_0)h+\frac{u’‘(\xi)}{2!}h^2$ ，其中 $\displaystyle \xi\in(t_0,t_1)$ ，记为 $\displaystyle u(t_0)+u’(t_0)h + R_0$ ，$\displaystyle R_0 = O(h^2)$ 是二阶小量
略去二阶小量 $R_0$ ，得到 $u(t_1)$ 的近似值 $\displaystyle u_1 = u_0+hf(t_0, u_0)\approx u(t_1)$
类似地，可得到 $u(t_2)$ 的近似值 $u_2 = u_1+hf(t_1, u_1)$
以此类推，得到递推公式：$u_{n + 1}=u_n+hf(t_n, u_n)$ ，这就是 Euler 法
显格式与隐格式
显格式：
$u_{n+1} = u_n + h f(t_n, u_n)$，变形为 $\displaystyle \frac{u_{n + 1}-u_n}{h}=f(t_n, u_n)$ ，即用差商 $\displaystyle \frac{u_{n + 1}-u_n}{h}$（向前差分）代替导数 $u’(t_n)$
隐格式：
$\displaystyle \frac{u_n - u_{n - 1}}{h}=f(t_n, u_n)$ ，也可写成 $\displaystyle u_{n + 1}=u_n+hf(t_{n + 1}, u_{n + 1})$ ，这里用向后差分 $\displaystyle \frac{u_n - u_{n - 1}}{h}$ 近似 $u’(t_n)$ ，隐格式一般需要迭代求解
迭代公式为
$$
\displaystyle u_{n + 1}^{[k + 1]}=u_n+hf(t_{n + 1}, u_{n + 1}^{[k]})
$$
数值积分法：ODE 数值解的等价积分形式
常微分方程 $u’(t)=f(t, u(t))$ ，$u(t_0)$ 已知，其等价积分形式为：
$$
u(t)=u(t_0)+\int_{t_0}^{t}f(\tau, u(\tau))d\tau
$$
通过对积分项进行不同的数值近似，可得到不同的数值解法
矩形公式
左矩形公式（显式 Euler 法）
公式：$\displaystyle u_1 = u_0 + f(t_0, u_0)h$
原理：用左端点的函数值 $f(t_0, u_0)$ 近似积分区间 $[t_0, t_0 + h]$ 上 $\displaystyle f(\tau, u(\tau))$ 的平均值，来计算积分 $\displaystyle \int_{t_0}^{t_0 + h}f(\tau, u(\tau))d\tau\approx hf(t_0, u_0)$ ，进而得到 $u(t_1)$（$t_1=t_0 + h$ ）的近似值 $u_1$
右矩形公式（隐式 Euler 法）
公式：$\displaystyle u_1 = u_0 + f(t_1, u_1)h$
原理：用右端点的函数值 $f(t_1, u_1)$ 近似积分区间 $[t_0, t_0 + h]$ 上 $f(\tau, u(\tau))$ 的平均值 ，但由于 $u_1$ 未知，该方法是隐式的，一般需迭代求解
梯形公式（改进 Euler 法）
公式：$\displaystyle u_{n + 1}=u_n+\frac{h}{2}[f(t_n, u_n)+f(t_{n + 1}, u_{n + 1})]$
原理：用梯形面积公式近似积分，即把积分区间 $[t_n, t_{n + 1}]$ 上 $f(\tau, u(\tau))$ 的平均值近似为 $\displaystyle \frac{1}{2}[f(t_n, u_n)+f(t_{n + 1}, u_{n + 1})]$ ，从而计算积分
$$
\displaystyle \int_{t_n}^{t_{n + 1}}f(\tau, u(\tau))d\tau\approx\frac{h}{2}[f(t_n, u_n)+f(t_{n + 1}, u_{n + 1})]
$$
Euler 误差分析
数值解与精确解
数值解：显式 Euler 法的迭代公式为 $\displaystyle u_{n + 1}=u_n+hf(t_n, u_n)$
精确解：对 $u(t_{n + 1})$ 在 $t_n$ 处进行 Taylor 展开，$\displaystyle u(t_{n + 1})=u(t_n)+hf(t_n, u(t_n))+\frac{u’‘(t_n)}{2!}h^2+\frac{u’‘’(\xi)}{3!}h^3$ ，$\xi\in(t_n, t_{n + 1})$
局部截断误差
定义：把精确解 $u(t)$ 代入差分格式得到的余项。差分格式 $\displaystyle L[u_n, h]=u_{n + 1}-u_n - hf(t_n, u_n)$ ，将真解代入得 $\displaystyle L[u(t_n), h]=u(t_{n + 1})-u(t_n)-hf(t_n, u(t_n))\neq0$ ，余项 $\displaystyle R_n=\frac{u’‘(t_n)}{2!}h^2+\frac{u’‘’(\xi)}{3!}h^3$ 。当 $h\to0$ 时，Euler 法的局部截断误差为 $O(h^2)$ ，意味着随着步长 $h$ 趋于 0 ，误差主要由 $h^2$ 阶项决定
改进 Euler 法局部截断误差
利用积分中值定理相关结论
$$
\displaystyle \int_{a}^{b}v(x)dx=\frac{1}{2}(b - a)(v(a)+v(b))-\frac{(b - a)^3}{12}v’‘(\xi)
$$
将真解代入改进 Euler 法格式 $\displaystyle u_{n + 1}-u_n-\frac{h}{2}[f(t_n, u_n)+f(t_{n + 1}, u_{n + 1})]=0$ ，得到余项为 $\displaystyle -\frac{h^3}{12}u’‘’(\xi)$ ，$\xi\in(t_n, t_{n + 1})$ ，所以改进 Euler 法的局部截断误差为 $O(h^3)$ ，比 Euler 法精度更高
整体误差
记 $\displaystyle e_n = u(t_n)-u_n$
Euler 法：
$$
\begin{align*}
u_{n + 1}&=u_n+hf(t_n,u_n)\
u(t_{n+1})&=u(t_n)+hf(t_n,u(t_n)) + R_n
\end{align*}
$$
作差：
$$
\displaystyle e_{n + 1}=e_n+h[f(t_n,u_n)-f(t_n,u(t_n))]+R_n
$$
记 $R=\max|R_n|$，由于 $f$ 满足 Lipschitz 条件得：
$$
\begin{align*}
|e_{n+1}|&\leq|e_n|+L\cdot h\cdot|e_n| + R\
|e_n|&\leq(1 + Lh)|e_{n - 1}|+R\
&\leq(1 + Lh)^2|e_{n - 2}|+(1 + Lh)R+R\
&\vdots\
&\leq(1 + Lh)^n|e_0|+R\sum_{j = 0}^{n - 1}(1 + Lh)^j\
&=(1 + Lh)^n|e_0|+R\cdot\frac{(1 + Lh)^n-1}{Lh}
\end{align*}
$$
由 $\ln(1 + Lh)\leq Lh$，$\Rightarrow n\ln(1 + Lh)\leq n\cdot Lh\leq L(T - t_0)$，$\Rightarrow(1 + Lh)^n\leq e^{L(T - t_0)}$
于是
$$
|e_n|\leq e^{L(T - t_0)}|e_0|+\frac{R}{Lh}(e^{L(T - t_0)}-1)，n = 1,2,\cdots,N
$$
对 Euler 法，可取 $e_0 = 0$，$R = ch^2$
则
$$
|e_n|\leq\frac{ch}{L}(e^{L(T - t_0)})\leq ch
$$
即 $e_n = o(h)$
注意比较局部截断误差与整体误差的关系
稳定性
定义：中途误差是否连续地依赖于初始值误差。初始误差包括测量误差、舍入误差
Euler 法稳定性分析：
从初值 $u_0$，$v_0$ 算出的节点值分别记为 ${u_n}$，${v_n}$。
$$
\begin{align*}
u_n&=u_{n - 1}+hf(t_{n - 1},u_{n - 1})\
v_n&=v_{n - 1}+hf(t_{n - 1},v_{n - 1})
\end{align*}
$$
两式相减，记 $e_n = u_n - v_n$，则
$$
e_n=e_{n - 1}+h[f(t_{n - 1},u_{n - 1})-f(t_{n - 1},v_{n - 1})]
$$
由 Lipschitz 条件，有
$$
|e_n|\leq|e_{n - 1}|+hL|e_{n - 1}|=(1 + Lh)|e_{n - 1}|
$$
递推可得：
$$
\begin{align*}
|e_n|&\leq(1 + Lh)^n|e_0|\
&\leq e^{L(T - t_0)}|e_0|
\end{align*}
$$
这表明 $e_n$ 连续地依赖于初始误差 $e_0$，即 Euler 法稳定
线性多步法
数值积分法
数值积分公式：
$$
u(t_{n + 1})=u(t_n)+\int_{t_n}^{t_{n + 1}}f(t,u(t))dt
$$
思想：适当地选取节点，用 Lagrange 插值多项式代替 $f(t,u(t))$
Adams 外插法（显式）
节点：$t_{n - k},\cdots,t_n$（$k + 1$ 个点），对应 $k + 1$ 步法
使用 $k$ 次 Lagrange 插值多项式
Adams 内插法（隐式）
节点：$\displaystyle t_{n - k},\cdots,t_n,t_{n + 1}$（未知 $t_{n + 1}$ 故仍为 $k + 1$ 步法）
使用 $k + 1$ 次 Lagrange 插值多项式
Adams 外插法
插值分解
将 $f(t, u(t))$ 分解为 $\displaystyle f(t, u(t)) = L_{n,k}(t) + r_{n,k}(t)$ ，其中 $\displaystyle L_{n,k}(t)$ 是插值多项式，$\displaystyle r_{n,k}(t)$ 是插值余项
根据积分的线性性质，有
$$
\displaystyle u(t_{n + 1}) = u(t_n)+\int_{t_n}^{t_{n + 1}}L_{n,k}(t)dt+\int_{t_n}^{t_{n + 1}}r_{n,k}(t)dt
$$
记为式 ①
离散格式
$$
\displaystyle u_{n + 1}=u_n+\int_{t_n}^{t_{n + 1}}L_{n,k}(t)dt
$$
是通过用插值多项式 $L_{n,k}(t)$ 近似 $f(t, u(t))$ 并积分得到的数值求解公式
局部截断误差
局部截断误差 $\displaystyle R_{n,k}=\int_{t_n}^{t_{n + 1}}r_{n,k}(t)dt$ ，它衡量了用插值多项式近似 $f(t, u(t))$ 并积分时产生的误差
Lagrange 插值公式回顾
Lagrange 插值多项式
对于给定节点 $x_0, x_1,\cdots, x_k$ ，函数 $f(x)$ 的 $k$ 次 Lagrange 插值多项式为 $\displaystyle L_k(x)=\sum_{i = 0}^{k}f(x_i)l_i(x)$ ，其中 $l_i(x)$ 是 Lagrange 基函数
Newton 型插值多项式
向前形式：$$\displaystyle f(x_0)+f[x_0,x_1](x - x_0)+\cdots$$
向后形式：$$\displaystyle f(x_k)+f[x_{k - 1},x_k](x - x_k)+\cdots + f[x_0,\cdots,x_k](x - x_k)\cdots(x - x_0)$$ ，其中 $\displaystyle f[x_i,\cdots,x_j]$ 是 $i$ 阶差商
差分
向前差分：$\displaystyle \Delta_+ f_j = f_{j + 1}-f_j$
向后差分：$\displaystyle \Delta_- f_j = f_j - f_{j - 1}$ ，且 $\displaystyle \Delta_-^k f_j=\Delta_+^k f_{-k}$
采用 Newton 向后插值公式
Newton 向后插值公式 $L_{n,k}(t)$ 是基于节点处函数值 $g(t_i)$ 及其差商，通过差商的线性组合来近似 $g(t)$ ，进而用于数值积分求解 $u(t_{n + 1})$ ，展开形式是利用差商和节点信息逐步构建，以逼近 $g(t)$ 在区间 $[t_n, t_{n + 1}]$ 上的积分
设 $t_{n - k},\cdots,t_{n - 1},t_n$ 为节点，令 $f(t, u(t)) = g(t)$ ，$t = t_n+\tau h$，$\tau\in[0,1]$
$$
\begin{align*}
L_{n,k}(t)&=L_{n,k}(t)=g(t_n)+g(t_n,t_{n-1})(t-t_n)+\cdots+g(t_n,\cdots,t_{n-k})(t-t_n)\cdots(t-t_{n-k+1})\
&=g(t_n)+\frac{g(t_n-g_{n-1})}{t_n-t_{n-1}}+\cdots\cdots
\end{align*}
$$
Newton 向后插值多项式展开
把 $L_{n,k}(t)$ 在 $t = t_n+\tau h$ 处展开，$\displaystyle g_{n - j}=f_{n - j}=f(t_{n - j},u_{n - j})$ 是节点处的函数值。利用向前差分 $\displaystyle \Delta_+^i g_{n - i}$ 来构建插值多项式，$\displaystyle \binom{s}{j}=\frac{s(s - 1)\cdots(s - j + 1)}{j!}$ ，$\displaystyle \binom{s}{0}=1$ ，得到 $L_{n,k}(t)$ 的具体形式
$$
L_{n,k}(t)=g(t_n)+\frac{\Delta_+^i g_{n - 1}}{h}\tau h+\frac{\Delta_+^i g_{n - 2}}{2!}\tau(\tau+1)+ \cdots+\frac{\tau(\tau+1)\cdots(\tau+k-1)}{k!}\Delta_+^i g_{n - k}
$$
则
$$
\begin{align*}
L_{n,k}(t)&=L_{n,k}(t_n+\tau h) \
&=\sum_{j=0}^k (-1)^j\binom{-\tau}{j}\Delta_+^i f_{n - j}\
\end{align*}
$$
代入积分并化简
将 $L_{n,k}(t)$ 代入 $\displaystyle u_{n + 1}=u_n+\int_{t_n}^{t_{n + 1}}L_{n,k}(t)dt$ ，通过换元 $t = t_n+\tau h$ ，$dt = h d\tau$ ，把积分区间从 $[t_n,t_{n + 1}]$ 转换为 $[0,1]$ 。然后对积分进行计算，提出与 $j$ 无关的项 $h$ 和 $\Delta_+^j f_{n - j}$ ，得到 $\displaystyle u_{n + 1}=u_n + h\sum_{j = 0}^{k}a_j\Delta_+^j f_{n - j}$ ，其中记 $\displaystyle a_j = (-1)^j\int_{0}^{1}\binom{-\tau}{j}d\tau$
具体步骤如下：
$$
\begin{align*}
u_{n+1}&=u_n+\int_{t_n}^{t_{n+1}} L_{n,k}(t)dt\
&=u_n+\int_{0}^{1} L_{n,k}(t_n+\tau h)hd\tau\
&=\displaystyle u_n + h\sum_{j = 0}^{k}(-1)^j\int_{0}^{1}\binom{-\tau}{j}d\tau\Delta_+^j f_{n - j}\
&=\displaystyle u_n + h\sum_{j = 0}^{k}a_j\Delta_+^j f_{n - j}
\end{align*}
$$
计算系数（生成函数法）
把 $(-1)^j(-\tau)(-\tau - 1)\cdots(-\tau - j + 1)$ 看成是 $(1 - t)^{-\tau}$ 关于 $t$ 在 $t = 0$ 处的导数
利用泰勒级数展开
$$
(1 - t)^{-\tau}=\sum_{j = 0}^{\infty}\frac{(-1)^j(-\tau)(-\tau - 1)\cdots(-\tau - j + 1)}{j!}t^j
$$
对上述展开式的左右两端在 $[0,1]$ 上积分并进行一系列变形
得
$$
\displaystyle -\frac{t}{(1 - t)\ln(1 - t)}=\sum_{j = 0}^{\infty}a_j t^{j + 1}
$$
变形得
$$
-\frac{\ln(1 - t)}{t}\sum_{j = 0}^{\infty}a_j t^{j + 1}=\frac{1}{1-t}=1+t+t^2+\cdots+t^n
$$
由
$$
-\frac{\ln(1 - t)}{t}=1+\frac{t}{2}+\frac{t^2}{3}+\cdots+\frac{t^{n-1}}{n}
$$
代入得
$$
(1+\frac{t}{2}+\frac{t^2}{3}+\cdots+\frac{t^{n-1}}{n})(a_0+a_1+\cdots+a_nt^n)=(1+t+t^2+\cdots+t^n)
$$
比较系数得到 $a_j$
$$
\begin{equation}
\begin{cases}
a_0=1 \
a_1+\frac{1}{2}a_0= 1 \
\quad \quad \vdots \
a_n + \frac{1}{2}a_{n-1}+\frac{1}{3}a_{n-2}+\cdots+\frac{1}{n+1}a_{0} = 1
\end{cases}
\end{equation}
$$
计算 $\Delta_+^j f_{n-j}$
原理：利用向前差分算子 $\Delta_+$ 的性质，$\Delta_+=(S - I)$（$S$ 为移位算子，$\displaystyle Sf_n = f_{n + 1}$ ，$I$ 为恒等算子 ），通过二项式定理 $\displaystyle (a + b)^j=\sum_{l = 0}^{j}C_{j}^{l}a^{j - l}b^{l}$ 展开 $(E - I)^j$ 来计算 $\Delta_+^j f_{n - j}$
具体实现步骤：
$$
\Delta_+ f_n=f_{n+1}-f_{n}=(S-I)f_n \
$$
$$
\begin{align*}
\Delta_+ f_n &=(S-I)^jf_{n-j} \
&=\sum_{l = 0}^{j}C_{j}^{l}S^{J - l}(-I)^{l}f_{n-j}\
&=\sum_{l = 0}^{j}(-1)^l\binom{j}{l}f_{n-l}
\end{align*}
$$
因而
$$
\begin{align*}
u_{n + 1}&=u_n + h\sum_{j = 0}^{k}a_j\Delta_+^j f_{n - j} \
&=u_n + h\sum_{j = 0}^{k}a_j\sum_{l=0}^{j}(-1)^l\binom{j}{l}f_{n - l}\
&=u_n + h\sum_{l=0}^{j}\sum_{j = 0}^{k}(-1)^la_j\binom{j}{l}f_{n - l}\
&=u_n+h\sum_{j = 0}^{k}b_{kl}f_{n - l}
\end{align*}
$$
其中记 $\displaystyle b_{kl}=\sum_{j = 0}^{k}(-1)^la_j\binom{j}{l}$
$k=0$：此时公式为 $\displaystyle u_{n + 1}=u_n+hf(t_n,u_n)$ ，这就是显式 Euler 法，它是 Adams 外插法的最简单形式，基于当前点的函数值来预估下一点的值
$k=1$：$\displaystyle u_{n + 1}=u_n+\frac{h}{2}(3f_n - f_{n - 1})$ ，利用前两个节点的函数值来计算下一个节点的数值解，相比 Euler 法考虑了更多历史信息
$k=2$：$\displaystyle u_{n + 1}=u_n+\frac{h}{12}(23f_n - 16f_{n - 1}+5f_{n - 2})$ ，使用了前三个节点的函数值，精度进一步提高
$k=3$：$\displaystyle u_{n + 1}=u_n+\frac{h}{24}(55f_n - 59f_{n - 1}+37f_{n - 2}-9f_{n - 3})$ ，借助更多历史节点信息来逼近真实解
Adams 外插法局部截断误差
插值余项公式：对于函数 $f(x)$ 的 $n$ 次插值多项式 $P_n(x)$
插值余项
$$
\displaystyle E(f,x)=f(x)-P_n(x)=\frac{(x - x_0)(x - x_1)\cdots(x - x_n)}{(n + 1)!}f^{(n + 1)}(\xi)\quad(\xi 介于 x 与节点之间)
$$
应用到 Adams 外插法：这里 $g(t)=f(t,u(t))$ ，$L_{n,k}(t)$ 是插值多项式，$\displaystyle r_{n,k}(t)=g(t)-L_{n,k}(t)$ 为插值余项
根据插值余项公式
$$
\displaystyle r_{n,k}(t)=\frac{(t - t_{n - k})(t - t_{n - k + 1})\cdots(t - t_n)}{(k + 1)!}g^{(k + 1)}(\xi)
$$
$\xi$ 与 $t$ 有关，局部截断误差就是对 $\displaystyle r_{n,k}(t)$ 在积分区间上积分产生的误差
将 $t = t_n+\tau h$ 代入插值余项表达式 $\displaystyle r_{n,k}(t)$ 。根据插值余项公式 $\displaystyle r_{n,k}(t)=\frac{(t - t_{n - k})(t - t_{n - k + 1})\cdots(t - t_n)}{(k + 1)!}u^{(k + 2)}(\xi)$ ，当 $t = t_n+\tau h$ 时，$\displaystyle t - t_{n - i}=\tau h + ih = (\tau + i)h$ （$i = 0,1,\cdots,k$ ），
代入后经过整理得到
$$
\displaystyle r_{n,k}(t_n+\tau h)=\frac{(\tau + k)(\tau + k - 1)\cdots(\tau + 1)\tau h^{k + 1}}{(k + 1)!}u^{(k + 2)}(\xi)
$$
进一步变形为
$$
\displaystyle r_{n,k}(t_n+\tau h)=(-1)^{k + 1}\binom{-\tau}{k + 1}h^{k + 1}u^{(k + 2)}(\xi)
$$
积分变换：已知 $\displaystyle R_{n,k}=\int_{t_n}^{t_{n + 1}}r_{n,k}(t)dt$ ，通过换元 $t = t_n+\tau h$ ，$dt = h d\tau$ ，积分区间从 $[t_n,t_{n + 1}]$ 变为 $[0,1]$
则
$$
R_{n,k}=\int_{0}^{1}r_{n,k}(t_n+\tau h)\cdot h d\tau
$$
代入化简：将 $\displaystyle r_{n,k}(t_n+\tau h)=(-1)^{k + 1}\binom{-\tau}{k + 1}h^{k + 1}u^{(k + 2)}(\xi)$ 代入上式
得到
$$
R_{n,k}=\int_{0}^{1}(-1)^{k + 1}\binom{-\tau}{k + 1}u^{(k + 2)}(\xi)d\tau h^{k + 2}
$$
结果分析：由于 $\displaystyle u^{(k + 2)}(\xi)$ 在积分区间内有界（在一定条件下 ），令 $\displaystyle \alpha_{k + 1}=\int_{0}^{1}(-1)^{k + 1}\binom{-\tau}{k + 1}d\tau$
则
$$
\displaystyle R_{n,k}=\alpha_{k + 1}u^{(k + 2)}(\overline{\xi})h^{k + 2}
$$
可知 $\displaystyle R_{n,k}=O(h^{k + 2})$ ，这表明 Adams 外插法的局部截断误差阶为 $h^{k + 2}$ ，即随着步长 $h$ 减小，误差以 $\displaystyle h^{k + 2}$ 的速度趋近于 0
Adams 内插法
Adams 内插法为隐格式
局部截断误差比外插法高一阶，为 $O(h^{k + 3})$
$k = 0, 1, 2, 3$的内插公式：
$k = 0$：$\displaystyle u_{n + 1} = u_n + hf_{n + 1}$ → 隐式 Euler
$k = 1$：$\displaystyle u_{n + 1} = u_n + \frac{h}{2}(f_{n + 1} + f_n)$ → 改进 Euler
$k = 2$：$\displaystyle u_{n + 1} = u_n + \frac{h}{12}(5f_{n + 1} + 8f_n - f_{n - 1})$
$k = 3$：$\displaystyle u_{n + 1} = u_n + \frac{h}{24}(9f_{n + 1} + 19f_n - 5f_{n - 1} + f_{n - 2})$
多步法一般形式：
$$
\sum_{j = 0}^{k} \alpha_j u_{n + j} = h\sum_{j = 0}^{k} \beta_j f_{n + j}
$$
由 $\displaystyle u_n, u_{n + 1}, \cdots, u_{n + k - 1}$ → $u_{n + k}$ （共 $k$ 个，即$k$步法）
Adams 方法：
$$
\alpha_k = 1，\alpha_{k - 1} = -1，\displaystyle \alpha_{k - 2} = \cdots = \alpha_0 = 0
$$
注：$\alpha_k$ 总可取成 1，显方法 $\beta_k = 0$，隐方法 $\beta_k \neq 0$
实际计算时，初始 $u_0, u_1, \cdots, u_{k - 1}$，可由单步法算
待定系数法
将真解 $u(t)$ 代入格式：$\displaystyle L[u(t_n); h] = \sum_{j = 0}^{k} [\alpha_j u(t_{n + j}) - h\beta_j f(t_{n + j}, u(t_{n + j}))]$
为简便且讨论一般，在 $t_n$ 处作 Taylor 展开：
$$
\begin{align}
&u(t_{n + j}) = u(t_n + jh) = u(t_n) + u’(t_n)jh + \frac{u’‘(t_n)}{2!}(jh)^2 + \cdots + \frac{u^{(q)}(t_n)}{q!}(jh)^q + \cdots \
&f(t_{n + j}, u(t_{n + j})) = u’(t_{n + j}) = u’(t_n + jh) = u’(t_n) + u’‘(t_n)jh + \frac{u’‘’(t_n)}{2!}(jh)^2 + \cdots + \frac{u^{(q + 1)}(t_n)}{q!}(jh)^q + \cdots
\end{align}
$$
则：
$$
\begin{align*}
&\sum_{j = 0}^{k} [\alpha_j u(t_{n + j}) - h\beta_j f(t_{n + j}, u(t_{n + j}))] \
&=\sum_{j = 0}^{k} [\alpha_j (u(t_n) + u’(t_n)jh + \frac{u’‘(t_n)}{2!}(jh)^2 + \cdots + \frac{u^{(q)}(t_n)}{q!}(jh)^q + \cdots) - h\beta_j (u’(t_n) + u’‘(t_n)jh + \frac{u’‘’(t_n)}{2!}(jh)^2 + \cdots + \frac{u^{(q + 1)}(t_n)}{q!}(jh)^q)] \
&=\sum_{j = 0}^{k} \alpha_j u(t_n) + \sum_{j = 0}^{k} (j\alpha_j - \beta_j)u’(t_n) \cdot h + \cdots + \sum_{j = 0}^{k} (\alpha_j \frac{j^q}{q!} - \beta_j \frac{j^{q - 1}}{(q - 1)!}) \cdot h^q \
&= C_0 u(t_n) + C_1 u’(t_n) + \cdots + C_q u^{(q)}(t_n) \cdot h^q + \cdots
\end{align*}
$$
其中：
$$
\begin{equation}
\begin{cases}
C_0 = \alpha_0 + \alpha_1 + \cdots + \alpha_k \
C_1 = \alpha_1 + 2\alpha_2 + \cdots + k\alpha_k - (\beta_0 + \beta_1 + \cdots + \beta_k) \
\quad \quad \vdots \
C_q = \frac{1}{q!}(\alpha_1 + 2^q\alpha_2 + \cdots + k^q\alpha_k) - \frac{1}{(q - 1)!}(\beta_1 + 2^{q - 1}\beta_2 + \cdots + k^{q - 1}\beta_k)
\end{cases}
\end{equation}
$$
若 $u(t)$ 有 $p + 2$ 次连续导数，则可选 $\alpha_j, \beta_j$
使得：
$$
C_0 = C_1 = \cdots = C_p = 0，C_{p + 1} \neq 0
$$
此时：
$$
L[u(t_n); h] = C_{p + 1} u^{(p + 1)}(t_n)h^{p + 1} + O(h^{p + 2})
$$
得到了 $p$ 阶 $k$ 步法：
$$
\sum_{j = 0}^{k} \alpha_j u_{n + j} = h\sum_{j = 0}^{k} \beta_j f_{n + j}
$$
局部截断误差为 $O(h^{p + 1})$
整体误差为 $O(h^p)$
相容性、稳定性和误差估计
k 步法：
$$
\sum_{j = 0}^{k} \alpha_j u_{n + j} = h\sum_{j = 0}^{k} \beta_j f_{n + j}
$$
相容性
相容性：差分算子逼近微分算子
以 Euler 法为例：
$$
u_{n + 1} = u_n + hf(t_n, u_n)，变形为 \frac{u_{n + 1} - u_n}{h} = f(t_n, u_n)
$$
差分算子：
$$
L[u(t_n); h] = u(t_{n + 1}) - u(t_n) - hf(t_n, u(t_n))
$$
差商算子：
$$
\frac{1}{h}L[u(t_n); h] = \frac{u(t_{n + 1}) - u(t_n)}{h} - f(t_n, u(t_n))
$$
微分算子：
$$
u’(t) - f(t, u(t))
$$
Euler 法相容性：
$$
\frac{1}{h}L[u(t_n); h] - [u’(t_n) - f(t_n, u(t_n))] \to 0，当 h \to 0
$$
k 步法相容性：
$$
\frac{1}{h} [\sum_{j = 0}^{k} \alpha_j u(t_{n + j}) - h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u(t_{n + j}))] - [u’(t_n) - f(t_n, u(t_n))] = o(1)，当 h \to 0
$$
即
$$
\sum_{j = 0}^{k} \alpha_j u(t_{n + j}) - h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u(t_{n + j})) = o(h) \ \Leftrightarrow (C_0 = 0, C_1 = 0)
$$
其中
$$
\sum_{j = 0}^{k} \alpha_j u(t_{n + j}) - h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u(t_{n + j})) = C_0u(t) + C_1u’(t)h + C_2u’'(t)h^2 + \cdots
$$
特征多项式
第一特征多项式：
$$
\rho(\lambda) = \sum_{j = 0}^{k} \alpha_j \lambda^j
$$
第二特征多项式：
$$
\sigma(\lambda) = \sum_{j = 0}^{k} \beta_j \lambda^j
$$
定理3.1
$$
k步法相容 \Leftrightarrow C_0 = 0, C_1 = 0 \Leftrightarrow \rho(1) = 0, \rho’(1) = \sigma(1)
$$
其中：
$$
\begin{align*}
&C_0 = \alpha_0 + \alpha_1 + \cdots + \alpha_k = \rho(1)\
&C_1 = \alpha_1 + 2\alpha_2 + \cdots + k\alpha_k - (\beta_0 + \beta_1 + \cdots + \beta_k) = \rho’(1) - \sigma(1)
\end{align*}
$$
稳定性
稳定性：解对初值的连续依赖性
k 步法：由初值 $u_0, u_1, \cdots, u_{k - 1}$ 得到 $u_n$ 。若 $u_0, u_1, \cdots, u_{k - 1}$ 有误差，$u_n$ 的误差是否会无限增长？
定义 3.2：若 $\exists \tau > 0, h_0 > 0$ ，使得当网格步长 $h$ 满足 $0 < h < h_0$ 时，对 $k$ 步法的任意两个解 ${u_n}$， ${v_n}$
都有
$$
\max_{n \leq \frac{T}{h}} |u_n - v_n| \leq \tau \max_{0 \leq j \leq k} |u_j - v_j|
$$
则称 $k$ 步法稳定
具体方法的误差递推关系
显式 Euler：
$$
|e_{n + 1}| \leq |e_n| + Lh|e_n| + R
$$
此为误差递推关系
隐式 Euler：
$$
\begin{align}
&u_{n + 1} = u_n + hf(t_{n + 1}, u_{n + 1}) \
&v_{n + 1} = v_n + hf(t_{n + 1}, v_{n + 1}) \
&\rightarrow|e_{n + 1}| \leq |e_n| + Lh|e_{n + 1}| \
&\rightarrow|e_{n + 1}| \leq \frac{1}{1 - Lh}|e_n|
\end{align}
$$
k 步法：
已知：由 $u_0, u_1, \cdots, u_{k - 1}$ 得到 $u_n$ ，满足
$$
\sum_{j = 0}^{k} \alpha_j u_{n + j} = h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u_{n + j})
$$
由 $v_0, v_1, \cdots, v_{k - 1}$ 得到 $v_n$ ，满足
$$
\sum_{j = 0}^{k} \alpha_j v_{n + j} = h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, v_{n + j})
$$
记 $e_n = u_n - v_n$ ，则
$$
\sum_{j = 0}^{k} \alpha_j e_{n + j} = h\sum_{j = 0}^{k} \beta_j [f(t_{n + j}, u_{n + j}) - f(t_{n + j}, v_{n + j})]
$$
记 $hb_n$ 为：
$$
hb_n=h\sum_{j = 0}^{k} \beta_j [f(t_{n + j}, u_{n + j}) - f(t_{n + j}, v_{n + j})]
$$
进一步得到
$$
e_{n + k} = \frac{1}{\alpha_k}(-\alpha_{k - 1}e_{n + k - 1} - \cdots - \alpha_0 e_n + hb_n)
$$
将 $k$ 个连续的 $e_i$ 构成向量形式：
$$
\begin{bmatrix}e_{n + k}\e_{n + k - 1}\\vdots\e_{n + 1}\end{bmatrix}=
\begin{bmatrix}
-\frac{\alpha_{k - 1}}{\alpha_k} &-\frac{\alpha_{k - 2}}{\alpha_k} &\cdots& -\frac{\alpha_0}{\alpha_k} \ 1&0&\cdots&0\ \vdots&\ddots&\ddots&\vdots\0&\cdots&1&0\end{bmatrix}
\begin{bmatrix}e_{n + k - 1}\e_{n + k - 2}\\vdots\e_{n}\end{bmatrix}+\begin{bmatrix}h\frac{b_n}{\alpha_k}\0\\vdots\0\end{bmatrix}
$$
即
$$
E_{n + 1} = CE_n + B_n
$$
类似数值代数中迭代法
且
$$
E_n = C^nE_0 + \sum_{l = 0}^{n - 1} C^l B_{n - 1 - l}
$$
关于 $C^n$ 有界性分析
特征方程
计算 $\det(C - \lambda I_k)$ ，其结果为
$$
(\alpha_k\lambda^k+\cdots+\alpha_1\lambda+\alpha_0)\cdot(-1)^k\alpha_k^{-1}
$$
记
$$
\rho(\lambda)=\alpha_k\lambda^k+\cdots+\alpha_1\lambda+\alpha_0
$$
利用 Jordan 标准型计算
设 $\lambda_j$ 是 $C$ 的特征值，$\displaystyle x_j=(d_{k - 1},\cdots,d_0)^T$ 是相应特征向量 ，则有 $Cx_j=\lambda_jx_j$ ，由此推导出一系列等式：
$$
\begin{align}
-\alpha_k^{-1}&(\alpha_{k - 1}d_{k - 1}+\cdots+\alpha_1d_1+\alpha_0d_0)=\lambda_jd_{k - 1} \
&d_{k - 2}=\lambda_jd_{k - 3} ，\cdots ，d_1=\lambda_jd_0
\end{align}
$$
进而得到
$$
d_{k - 1}=\lambda_j^{k - 1}d_0\
x_j = d_0(\lambda_j^{k - 1},\cdots,\lambda_j,1)
$$
因而 $C$ 的每个特征值的特征子空间维数都是 1，即每个特征值几何重数为 1，每个特征值的 Jordan 块只有 1 块。记 $C$ 的 Jordan 标准型为 $J$
则
$$
C = PJP^{-1}
$$
Jordan 块形式：$J$ 的 Jordan 块要么是 $\lambda$ ，要么是 $J_r(\lambda)$，$r$ 是 $\lambda$ 的重数
$$
J_r(\lambda)=\begin{bmatrix}
\lambda & 1 & \cdots & 0 \
0 & \lambda & \ddots & \vdots \
\vdots & \vdots & \ddots & 1 \
0 & 0 & \cdots & \lambda
\end{bmatrix}
$$
$\displaystyle C^n = PJ^nP^{-1}$ ，且 $J^n$ 也是分块矩阵，每一块形式为 $\lambda^n$ 或 $J_r^n$
$$
J_r^n(\lambda)=\begin{bmatrix}
\lambda & c_n^1\lambda^{n-1}-\cdots c_n^{r-1}\lambda^{n-r+1} & \cdots & 0 \
0 & \lambda & \ddots & \vdots \
\vdots & \vdots & \ddots & c_n^1\lambda^{n-1} \
0 & 0 & \cdots & \lambda
\end{bmatrix}
$$
对于 $J_r=\lambda I + S$ ，（其中 $S$ 是幂零矩阵）
$$
S=\begin{bmatrix}
0& 1 & \cdots & 0 \
0 & 0 & \ddots & \vdots \
\vdots & \vdots & \ddots & 1 \
0 & 0 & \cdots & 0
\end{bmatrix}
\quad
S^2=\begin{bmatrix}
0 & 0 &1 & \cdots & 0 \
0 & 0& \ddots & \ddots &\vdots\
0 & 0 &\ddots &\ddots &1\
\vdots & \vdots & \ddots & 0 &0\
0 & 0 & \cdots & 0&0
\end{bmatrix} \quad
S^r=0
$$
利用二项式定理计算
$$
\begin{align*}
J_r^n &= (\lambda I+S)^n \
&=\sum_{j = 0}^{n}C_{n}^{j}\lambda^{n - j}S^j \
&=\sum_{j = 0}^{r - 1}C_{n}^{j}\lambda^{n - j}S^j
\end{align*}
$$
这里 $\displaystyle C_{n}^{j}=\frac{n!}{j!(n - j)!}$
引理 1.1
$C^n$有界 $\Leftrightarrow |\lambda_j| \leq 1$ ，且若 $|\lambda_j| = 1$ ，则 $\lambda_j$ 是单根
$C^n$有极限 $\Leftrightarrow |\lambda_j| \leq 1$ ，且若 $|\lambda_j| = 1$ ，则 $\lambda_j$ 不仅是单根且 $\lambda_j = 1$
$C^n \to 0 \Leftrightarrow |\lambda_j| < 1$
定理：线性多步法稳定的充要条件
线性 $k$ 步法稳定 $\Leftrightarrow$ 第一特征多项式 $\rho(\lambda)$ 的根满足 $|\lambda| \leq 1$ ，且位于单位圆周上的根都是单根（即 $\rho(\lambda)$ 满足根条件
证明
必要性证明
将 $k$ 步法用于方程 $u’ = 0$ ，此时 $f(t, u) \equiv 0$ ，则 $\displaystyle \sum_{j = 0}^{k} \alpha_j u_{n + j} = 0$ ，${u_n = 0}$ 是上述方程的平凡解。
因为 $k$ 步法稳定，所以 $\exists \tau, h_0$ ，使得当 $0 < h < h_0$ 时，$\displaystyle \max_{n \leq \frac{T}{h}} |u_n| \leq \tau \max_{0 \leq j \leq k} |u_j|$ ，即 ${u_n}$ 关于 $n, h$ 一致有界。
又因为 $u_n = C^n u_0$ （这里 $\displaystyle u_0 = [u_{k - 1}, \cdots, u_0]^T$ ，$\displaystyle u_n = [u_{n + k - 1}, \cdots, u_n]^T$ ），所以 ${C^n}$ 一致有界，从而 $\rho(\lambda)$ 满足根条件。
充分性证明
设 ${u_n}, {v_n}$ 是 $k$ 步法的两个解，令 $e_n = u_n - v_n$ ，则
$$
\displaystyle \sum_{j = 0}^{k} \alpha_j e_{n + j} = h\sum_{j = 0}^{k} \beta_j [f(t_{n + j}, u_{n + j}) - f(t_{n + j}, v_{n + j})] = hb_n
$$
设 $B = \max{|\beta_0|, \cdots, |\beta_k|}$ ，由 $f$ 关于 $u$ 满足 Lipschitz 条件，可得
$$
|b_n| \leq B\cdot L \sum_{j = 0}^{k} |e_{n + j}|
$$
记 $\displaystyle E_n = (e_{n + k - 1}, \cdots, e_n)^T$ ，$\displaystyle B_n = (h\alpha_k^{-1}b_n, 0, \cdots, 0)^T$ ，则 $E_{n + 1} = CE_n + B_n$ ，进而有
$$
E_n = C^n E_0 + \sum_{l = 0}^{n - 1} C^l B_{n - 1 - l}
$$
由于 $\rho(\lambda)$ 满足根条件，所以 $|C^n| \leq M, \forall n$ ，且
$$
|B_n| \leq h|\alpha_k^{-1}| |b_n| \leq h|\alpha_k^{-1}| \cdot B\cdot L \sum{j = 0}^{k} |e{n + j}|
$$
经过推导可得
$$
\begin{align}
|E_n| &\leq M|E_0| + \sum_{l = 0}^{n - 1} M\cdot h|\alpha_k^{-1}| \cdot B\cdot L \sum_{j = 0}^{k} |e_{n - 1 - l + j}| \&\leq M|E_0| + M_1 \sum_{l = 0}^{n - 1} \sum_{j = 0}^{k} |e_{n - 1 - l + j}|
\end{align}
$$
对 $\displaystyle \sum_{l = 0}^{n - 1} \sum_{j = 0}^{k} |e_{n - l - 1 + j}|$ 进行放缩：
$$
\sum_{l = 0}^{n - 1} \sum_{j = 0}^{k} |e_{n - l - 1 + j}|=\sum_{l = 0}^{n - 1} |e_{n + k - 1 - l}|+\sum_{l = 0}^{n - 1} \sum_{j = 0}^{k} |e_{n - l - 1 + j}|
$$
利用向量范数性质，有
$$
\sum_{l = 0}^{n - 1} |e_{n + k - 1 - l}|+\sum_{l = 0}^{n - 1} \sum_{j = 0}^{k} |e_{n - l - 1 + j}| \leq \sum_{l = 0}^{n - 1} |E_{n - l - 1}|+\sum_{l = 0}^{n - 1} \sqrt{k} |E_{n - l - 1}|
$$
进一步得到
$$
\sum_{l = 0}^{n - 1} |E_{n - l - 1}|+\sum_{l = 0}^{n - 1} \sqrt{k} |E_{n - l - 1}|=|E_n|+(1 + \sqrt{k})\sum_{i = 0}^{n - 1} |E_i|
$$
结合前面 $|E_n|$ 的表达式
$$
|E_n| \leq M|E_0| + M B L|\alpha_k^{-1}| h |E_n|+M B L|\alpha_k^{-1}| h\cdot(\sqrt{k}+1)\sum_{l = 0}^{n - 1} |E_l|
$$
在 $\displaystyle MBL|\alpha_k^{-1}| h_0 < 1$ 条件下，令 $\displaystyle K_1=\frac{M}{1 - M B L|\alpha_k^{-1}| h_0}$ ，$\displaystyle K_2 = K_1 M B L|\alpha_k^{-1}|(\sqrt{k}+1)$ ，可得
$$
|E_n| \leq M|E_0| + M B L|\alpha_k^{-1}| h |E_n|+M B L|\alpha_k^{-1}| h\cdot(\sqrt{k}+1)\sum_{l = 0}^{n - 1} |E_l|
$$
再由 Gronwall 不等式，得到
$$
|E_n| \leq e^{K_2 T}(K_1 + K_2 h_0)|E_0|
$$
进而有
$$
\max_{n \leq \frac{T}{h}} |u_n - v_n| \leq \overline{c} \max_{0 \leq j \leq k} |u_j - v_j|
$$
说明方法稳定
应用：Euler 法稳定性、Adams 外插和内插法稳定性
收敛性与误差估计
定义与设定：设真解 $u(t_n)$ 和数值解 $u_n$ 满足 ：
$$
\begin{align}
\sum_{j = 0}^{k} \alpha_j u(t_{n + j}) &= h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u(t_{n + j})) + R_n \
\sum_{j = 0}^{k} \alpha_j u_{n + j} &= h\sum_{j = 0}^{k} \beta_j f(t_{n + j}, u_{n + j})
\end{align}
$$
令 $\displaystyle e_n = u(t_n) - u_n$
则
$$
\begin{align}
\sum_{j = 0}^{k} \alpha_j e_{n + j} &= h\sum_{j = 0}^{k} \beta_j [f(t_{n + j}, u(t_{n + j})) - f(t_{n + j}, u_{n + j})]+R_n \&= hb_n + R_n
\end{align}
$$
设 $\displaystyle R_n = C_{p + 1}h^{p + 1}$ ，记 $\displaystyle E_{n + 1} = CE_n + B^{\ast}_n$ ，其中 $\displaystyle B^{\ast}_n=(h\alpha_k^{-1}b^{\ast}_n, 0, \cdots, 0)^T$ ，且
$$
\displaystyle |B_n^*| \leq h B L|\alpha_k^{-1}| \sum_{j = 0}^{k} |e_{n + j}|+M_{p + 1}h^{p + 1}
$$
定理
相容性 + 稳定性 $\Rightarrow$ 收敛。且若 $E_0 = 0$ ，则 $|E_n| \leq C h^p$ 。即当数值方法满足相容性（差分算子能逼近微分算子）和稳定性（解对初值连续依赖）时，随着步长 $h$ 趋于 0 ，数值解会收敛到真解 ，且在初始误差为 0 时，误差阶为 $O(h^p)$
