本篇文章适合个人复习翻阅，不建议新手入门使用
本专栏：数值分析复习 的前置知识主要有：数学分析、高等代数、泛函分析

本节继续考虑数值积分问题

Richardson外推

命题：复合梯形公式的另一形式
设 $f\in C^{\infty }[a,b]$，记 $I={\int }_a^{b}f(x)\mathrm{d}x$ ，将复合梯形公式记为
$$T(h)=\frac{h}{2}\sum _{i=0}^{n-1}[f(x_i)+f(x_{i+1})]$$

则 $$T(h)=I+{\alpha }_1{h}^2+{\alpha }_2{h}^4+\cdots +{\alpha }_l{h}^{2l}+\cdots$$

其中 ${\alpha }_l(l=1,2,\dots )$ 与 $h$ 无关

证明
设 $x_{i+\frac{1}{2}}=\frac{x_i+x_{i+1}}{2},i=0,1,\dots ,n-1$

考虑 $f(x)$ 在 $x=x_{i+\frac{1}{2}}$ 处的Taylor展开公式
$$f(x)=f(x_{i+\frac{1}{2}})+f^{\prime }(x_{i+\frac{1}{2}})(x-x_{i+\frac{1}{2}})+\frac{f^{\prime \prime }(x_{i+\frac{1}{2}})}{2!}(x-x_{i+\frac{1}{2}}{)}^2+\cdots$$

若对上述 Taylor 公式代入 $x=x_i,x=x_{i+1}$，则得
$$f(x_{i+1})=f(x_{i+\frac{1}{2}})+f^{\prime }(x_{i+\frac{1}{2}})\frac{h}{2}+\frac{f^{\prime \prime }(x_{i+\frac{1}{2}})}{2!}(\frac{h}{2}{)}^2+\cdots$$$$f(x_i)=f(x_{i+\frac{1}{2}})+f^{\prime }(x_{i+\frac{1}{2}})(-\frac{h}{2})+\frac{f^{\prime \prime }(x_{i+\frac{1}{2}})}{2!}(-\frac{h}{2}{)}^2+\cdots$$

两式加和，得到
$$\frac{f(x_i)+f(x_{i+1})}{2}=f(x_{i+\frac{1}{2}})+\frac{h^2}{8}{f}^{\prime \prime }(x_{i+\frac{1}{2}})+\cdots$$

等式两端求和，乘以 $h$ 得到
$$T(h)=h\sum _{i=0}^{n-1}f(x_{i+\frac{1}{2}})+\frac{h^3}{8}\sum _{i=0}^{n-1}{f}^{\prime \prime }(x_{i+\frac{1}{2}})+\cdots \text{\hspace{1em}(1)}$$

另一方面，对Taylor公式从 $x_i$ 到 $x_{i+1}$ 进行积分，得到
$${\int }_{x_i}^{x_{i+1}}f(x)\mathrm{d}x=h\cdot f(x_{i+\frac{1}{2}})+\frac{f^{\prime }(x_{i+\frac{1}{2}})}{2}[(\frac{h}{2}{)}^2-(-\frac{h}{2}{)}^2]+\frac{f^{\prime \prime }(x_{i+\frac{1}{2}})}{6}[(\frac{h}{2}{)}^3-(-\frac{h}{2}{)}^3]+\cdots$$

等式两端求和得

$$I=\sum _{i=0}^{n-1}{\int }_{x_i}^{x_{i+1}}f(x)\mathrm{d}x=h\sum _{i=0}^{n-1}f(x_{i+\frac{1}{2}})+\frac{h^3}{24}\sum _{i=0}^{n-1}{f}^{\prime \prime }(x_{i+\frac{1}{2}})+\cdots \text{\hspace{1em}(2)}$$

结合（1）（2）式，可得
$$T(h)=I+\frac{h^3}{12}\sum _{i=0}^{n-1}{f}^{\prime \prime }(x_{i+\frac{1}{2}})+\cdots \text{\hspace{1em}(3)}$$

类似（2）式的推导，可得
$${\int }_a^b{f}^{\prime \prime }(x)\mathrm{d}x=h\sum _{i=0}^{n-1}{f}^{\prime \prime }(x_{i+\frac{1}{2}})+\frac{h^3}{24}\sum _{i=0}^{n-1}{f}^{(4)}(x_{i+\frac{1}{2}})+\cdots$$

结合 ${\int }_a^b{f}^{\prime \prime }(x)\mathrm{d}x=f^{\prime }(b)-f^{\prime }(a)$，可将（3）式化为
$$T(h)=I+{\alpha }_1{h}^2+h^5{c}_4\sum _{i=0}^{n-1}{f}^{(4)}(x_{i+\frac{1}{2}})+\cdots$$

重复上述操作，考虑 ${\int }_a^b{f}^{(4)}(x)\mathrm{d}x$，消去 $h^5$ 的项，得到 $h^4$ 的项，继续重复操作，可得
$$T(h)=I+{\alpha }_1{h}^2+{\alpha }_2{h}^4+\cdots +{\alpha }_l{h}^{2l}+\cdots$$

定义：Richardson外推
从低阶精度格式的截断误差的渐近展开式出发，做简单线性计算从而得到高阶精度格式的方法称为Richardson外推

例：
考虑复合梯形公式 $T(h)$ 满足的式子
$$T(h)=I+{\alpha }_1{h}^2+{\alpha }_2{h}^4+\cdots +{\alpha }_l{h}^{2l}+\cdots$$

此时截断误差量级为 $O(h^2)$

取步长为 $\frac{h}{2}$，则有
$$T(\frac{h}{2})=I+{\alpha }_1\frac{h^2}{4}+{\alpha }_2\frac{h^4}{16}+\cdots +{\alpha }_l\frac{h^{2l}}{2^{2l}}+\cdots$$

结合这两个式子，消去$h^2$项，得
$$\frac{4T(\frac{h}{2})-T(h)}{3}=I-\frac{1}{4}{\alpha }_2{h}^4+\cdots +\frac{{\alpha }_l}{3}(\frac{1}{2^{2l}}-1)h^{2l}+\cdots$$
记 $T_1(h)=\frac{4T(\frac{h}{2})-T(h)}{3}$，且
$$T_1(h)=I+{\beta }_2{h}^4+{\beta }_3{h}^6+\cdots +{\beta }^l{h}^{2l}+\cdots$$
若用 $T_1(h)$ 估计 $I$ ，则截断误差量级提高到 $O(h^4)$
类似地，可继续做……

注：只要截断误差可表示为 $h$ 的幂级数，均可使用 Richardson外推提高精度

Romberg（龙贝格）算法

在上述对复合梯形公式的截断误差进行Richardson外推的过程中，记复合梯形公式 $$T_0(h)=T(h)=\frac{h}{2}\sum _{i=0}^{n-1}[f(x_i)+f(x_{i+1})]$$

加速一次（即进行一次Richardson外推）后的估计式记为
$$T_1(h)=\frac{4T(\frac{h}{2})-T(h)}{3}$$

记加速 $n$ 次的估计式为 $T_n(h)$，则有递推式
$$T_n(h)=\frac{4^n}{4^n-1}{T}_{n-1}(\frac{h}{2})-\frac{1}{4^n-1}{T}_{n-1}(h)$$

若记 $T_m^{(k)}=T_m(\frac{h}{2^k}),k=0,1,2,\dots$，则有递推式
$$T_n^{(k)}=\frac{4^n}{4^n-1}{T}_{n-1}^{(k+1)}-\frac{1}{4^n-1}{T}_{n-1}^{(k)}$$

定理：
设被积函数 $f(x)$ 充分光滑

1. $\underset{k\to \infty }{\lim }{T}_m^{(k)}=I$
2. $\underset{m\to \infty }{\lim }{T}_m^{(k)}=I$

注：证明略去，第一个结论说明当节点数目无穷多时，$T_m^{(k)}$ 收敛于准确的积分值；第二个结论说明随着Richardson外推的进行，$T_m^{(k)}$ 也收敛于准确的积分值

上述递推式和收敛定理给出了如下的Romberg算法

定义：Romberg算法
对预先给定的精度 $\epsilon$，求 $I={\int }_a^{b}f(x)\mathrm{d}x$ 的近似值，算法如下：

初始取 $k=0,m=0,h=b-a$

1. 代入梯形公式，求 $T_0^{(k)}(k=0,1,2,\dots )$
2. 加速一次，由递推公式求 $T_1^{(k)}$
3. 直至 $|T_k^{(0)}-T_{k-1}^{(0)}|<\epsilon$，则取 $T_k^{(0)}\approx I$

注：具体求解顺序如下表

参考书籍：《数值分析》李庆扬 王能超 易大义 编