温度每隔一小段时间记录一次，要估计此刻升温有多快，可以拿两个读数相减，再除以时间间隔。间隔缩小，看起来更接近导数的定义；可若温度计只能读到小数点后两位，两次读数很快就会一样。这时算出的零，不一定说明温度停止变化。

积分遇到的是另一个问题：一段时间内累计了多少。它也只能使用有限个读数，却不必把两次接近的值相减后除以很小的数。本章从这一区别出发，把公式的来源、误差和取样位置一起检查。

用 Taylor 展开看清差分误差

以下差分公式取 $h>0$，差分误差统一记为“算出的近似值减去真实值”。前向差分为

$D_+(h)=\frac{f(x+h)-f(x)}h.$

若 $f$ 在这段邻域内足够光滑，展开得

$f(x+h)=f(x)+hf'(x)+\frac{h^2}{2}f''(x)+O(h^3),$

因此 $D_+(h)=f'(x)+hf''(x)/2+O(h^2)$。若只知道 且 ，Taylor 余项仍给出明确保证

$|D_+(h)-f'(x)|\le \frac{M_2h}{2}.$

后向差分 $D_-(h)=[f(x)-f(x-h)]/h$ 的主误差符号相反。把前后两种近似平均，便得到中心差分

$D_0(h)=\frac{f(x+h)-f(x-h)}{2h}.$

将两侧 Taylor 式相减，常数项、二次项相消：

$f(x+h)-f(x-h)=2hf'(x)+\frac{h^3}{3}f'''(x)+O(h^5).$

这里最后的展开可在 $f\in C^5$ 的邻域内使用，于是

$D_0(h)=f'(x)+\frac{h^2}{6}f'''(x)+O(h^4).$

仅为了得到二阶误差界，不需要五阶导数。若 $f\in C^3$ 且 $|f'''|\le M_3$，分别对两侧使用三阶余项，仍有

$|D_0(h)-f'(x)|\le \frac{M_3h^2}{6}.$

分母是两端相隔的 $2h$，不是 $h$。虽然推导时围绕三个位置 $x-h,x,x+h$ 展开，实际只需两次函数求值，中间值的系数为零。

例：把误差算成一个看得见的数。 对 $f(x)=x^3$，在 $x=1$ 的真实导数是 3。直接展开得到

$D_+(h)=3+3h+h^2,\qquad D_0(h)=3+h^2.$

取 $h=0.1$，结果分别为 3.31 和 3.01；改为 $h=0.05$，中心差分的误差从 0.01 降为 0.0025，正好缩小四倍。这个算例没有省略高阶项，因为三次多项式可以完全展开。

对同一三次函数，前向差分使用 1 与 1.1 两处读数，中心差分使用 0.9 与 1.1 两处读数；两者分母分别为 0.1 与 0.2。

在左端点没有 $x-h$ 的数据时，中心公式无从使用。可以寻找只用 $x,x+h,x+2h$ 的近似

$D_e(h)=\frac{af(x)+bf(x+h)+cf(x+2h)}h.$

代入 Taylor 式后，常数项要消失，一阶导数系数要为 1，二阶导数项也要消失。因此

$a+b+c=0,\qquad b+2c=1,\qquad b+4c=0.$

解得 $a=-3/2,b=2,c=-1/2$，所以

$D_e(h)=\frac{-3f(x)+4f(x+h)-f(x+2h)}{2h}=f'(x)-\frac{h^2}{3}f'''(x)+O(h^3).$

这里的展开可在 $f\in C^4$ 时使用。这就是二阶端点公式。它需要三次求值，也更容易放大数据噪声；“二阶”相同，不代表成本和敏感性相同。右端点可以把方向反过来重新推导，不能只把 $h$ 的符号改掉一半。

二阶导数则使用相加：奇次项消失，留下

$D_2(h)=\frac{f(x+h)-2f(x)+f(x-h)}{h^2}=f''(x)+\frac{h^2}{12}f^{(4)}(\xi),$

其中 $f\in C^4$，$\xi$ 位于取样邻域内。两个四阶余项的平均仍是某处的四阶导数值，这是连续性保证的。对 $f(x)=x^4$、$x=1$，可以直接核对 。

步长太小，数据误差会被放大

现在把函数读数写为 $\widetilde f(x_j)=f(x_j)+\eta_j$，并假设每次读数的绝对误差都不超过 。它可以代表仪器误差，也可以是对函数求值舍入误差的估计；两者都需要根据实际尺度确定。

中心差分中的输入误差为 $(\eta_+-\eta_-)/(2h)$，绝对值至多 $\delta/h$。加上截断误差，得到

$|\widetilde D_0(h)-f'(x)|\le \frac{M_3h^2}{6}+\frac\delta h.$

这个界暂时只计算函数读数中的误差；程序随后相减、相除以及表示 $x\pm h$ 的误差还需另计。它已经足以解释主要矛盾：减小 $h$ 会压低第一项，却抬高第二项。

设 $C=M_3/6>0$，最小化模型 $E(h)=Ch^2+\delta/h$。由

$E'(h)=2Ch-\frac\delta{h^2}=0$

得到 $h_*=(\delta/(2C))^{1/3}=(3\delta/M_3)^{1/3}$。这是该模型、该尺度下的选择；若导数界不合适、数据误差随位置变化，最优值也会变化。

仍用三次函数的中心差分，截断误差恰为 $h^2$。若 $\delta=10^{-6}$，模型最小点约为 $0.00794$。在 $h=0.01$，界为 $$；在 ，截断项虽然只有 ，输入误差项却已达到 。

误差趋势示意：横向向右表示步长增大。红线表示输入误差影响减小，蓝线表示截断影响增大，深色线示意总量；它不对应指定系数或一次浮点实验的精确曲线。

前向差分的输入误差上限为 $2\delta/h$；三点端点公式为 $4\delta/h$；二阶导数公式为 $4\delta/h^2$。这些数来自系数绝对值之和，不能指望误差总会互相抵消。特别是二阶求导，缩小步长时会更快地放大噪声。

下面同时画出差分结果和误差随步长的变化。先用三次函数核对刚才的精确表达式，再切换到指数函数，把步长逐级缩小。实验把“人为加入的有界读数误差”和“浏览器实际浮点计算”分开显示；真实误差会有起伏，不会被强行画成平滑的 U 形。

Richardson 外推消掉的是哪一项

若一种近似确实具有展开

$A(h)=A_*+ch^p+O(h^{p+q}),\qquad q>0,$

并且两次计算面对的是同一个 $A_*$、同一个主误差系数 $c$，则

$A(h/2)=A_*+c\frac{h^p}{2^p}+O(h^{p+q}).$

把第二式乘 $2^p$，减去第一式，再除以 $2^p-1$，主误差就被消去：

$R(h)=\frac{2^pA(h/2)-A(h)}{2^p-1}=A_*+O(h^{p+q}).$

中心差分在足够光滑时只出现偶次误差，故可用

$R_D(h)=\frac{4D_0(h/2)-D_0(h)}3$

把二阶提高到四阶。对刚才的三次函数，$D_0(h)=3+h^2$，所以外推后恰好得到 3。对一般函数，只是消去主项，剩余误差并不会全部消失。

也别忘了输入噪声。若各读数误差至多为 $\delta$，$D_0(h/2)$ 的误差可能达到 $2\delta/h$；外推组合的输入误差上限于是为 $(4\times2+1)\delta/(3h)=3\delta/h$。提高截断阶数有价值，但已经被噪声主导的结果，不能靠不停外推恢复可信数字。

求积权重来自对近似函数的积分

令 $I=\int_a^b f(x)\,dx$。若以第 4 章的插值多项式 $p(x)=\sum y_jL_j(x)$ 代替 ，积分后便是

$Q=\int_a^b p(x)\,dx=\sum_j w_jy_j,\qquad w_j=\int_a^b L_j(x)\,dx.$

权重是基函数的积分，与未知的 $f$ 无关。我们称一个公式的代数精度为 $d$，是指它对所有次数不超过 $d$ 的多项式都精确，但对某个 $d+1$ 次多项式不精确。这与细分步长时的收敛阶不是同一个概念。

在一段 $[u,u+h]$ 上，用直线连接两端，积分得到梯形公式

$Q_T=\frac h2\bigl(f(u)+f(u+h)\bigr).$

若 $|f''|\le M_2$，一次插值余项给

$|f(u+t)-p(u+t)|\le\frac{M_2}{2}t(h-t),\qquad 0\le t\le h.$

对它积分，$\int_0^h t(h-t)\,dt=h^3/6$，故单段积分误差不超过 $M_2h^3/12$。若 ，曲线位于弦下，梯形值偏大。更准确地，在 时有 。

把 $[a,b]$ 均分为 $n$ 段，$h=(b-a)/n$、$x_j=a+jh$。内部节点分别出现在左右两段，两个半权重合在一起，于是

$T_n=h\left[\frac{f(a)+f(b)}2+\sum_{j=1}^{n-1}f(x_j)\right],\qquad |I-T_n|\le\frac{(b-a)M_2h^2}{12}.$

每段误差是三次量级，段数却与 $1/h$ 成正比，所以全局是二阶。

也可以在每段只取中点，得到复合中点公式

$Q_{M,n}=h\sum_{j=0}^{n-1}f\left(a+(j+\tfrac12)h\right).$

以段中点展开，奇次的一阶项积分为零，二阶余项积分的绝对值不超过 $M_2h^3/24$。累加得 $|I-Q_{M,n}|\le(b-a)M_2h^2/24$。对凸函数，中点值偏小，与梯形恰好从两边夹住积分。

Simpson 的三个权重与偶数条件

在以 0 为中点的 $[-h,h]$ 上，设公式为 $A f(-h)+B f(0)+A f(h)$。对常数和二次函数要求精确：

$2A+B=2h,\qquad 2Ah^2=\frac{2h^3}{3}.$

得到 $A=h/3,B=4h/3$。奇次函数的积分为零，左右对称的加权值也为零，因此该公式还对三次函数精确：

$S=\frac h3\bigl[f(-h)+4f(0)+f(h)\bigr].$

它等于把三个点的二次插值多项式积分；对三次式也精确，是对称性额外消去了一项。对 $x^4$ 则不精确，代数精度恰为 3。

相邻 Simpson 面板共享端点，该节点的两个权重 1 相加为 2；面板中点各自保留权重 4。这里共有四个等宽子区间，每个面板占两段。

使用复合 Simpson 时，每个面板占两段，因此 $n$ 必须为正偶数。把相邻面板相加，公共端点的两个权重 1 变成 2，中点权重仍为 4：

$S_n=\frac h3\left[f(x_0)+4\sum_{j=1}^{n/2}f(x_{2j-1})+2\sum_{j=1}^{n/2-1}f(x_{2j})+f(x_n)\right].$

若 $n=5$，不能把最后一个系数强行改成 1，声称仍是这个公式；应重新选择偶数分段，或明确采用其他组合规则。

在 $f\in C^4[a,b]$ 时，误差为

$I-S_n=-\frac{b-a}{180}h^4 f^{(4)}(\xi),\qquad |I-S_n|\le\frac{(b-a)M_4h^4}{180}.$

这里 $M_4$ 是四阶导数绝对值的上界。这个常数不能只靠在 $x^4$ 上试算就证明；下面把余项为何成立补齐。

先在 $[-1,1]$ 记误差运算 $E(g)=\int_{-1}^1g(t)dt-[g(-1)+4g(0)+g(1)]/3$。带积分余项的 Taylor 式可写为

$g(t)=P_3(t)+\int_{-1}^1\frac{(t-s)_+^3}{6}g^{(4)}(s)\,ds,$

其中 $P_3$ 是在 $-1$ 展开的三次 Taylor 多项式，$(t-s)_+=\max(t-s,0)$。 对三次式为零，交换积分和有限个求值，便有

$E(g)=\int_{-1}^1K(s)g^{(4)}(s)ds,\qquad K(s)=E\left(\frac{(t-s)_+^3}{6}\right).$

例如 $0\le s\le1$ 时，只有右端点的求值项不为零，所以

$K(s)=\frac{(1-s)^4}{24}-\frac{(1-s)^3}{18}=-\frac{(1-s)^3(1+3s)}{72}.$

另一半由对称性得到相同的 $|s|$ 表达式。于是 $K\le0$，且积分为 $-1/90$；连续的 $g^{(4)}$ 对这个同号权重取加权平均，得到 $E(g)=-g^{(4}$。平移缩放到宽度为 的一个面板，多出 ；累加 个面板，恰好给出上述复合公式。四阶连续性和偶数分段各自在这里起了作用。

例：相同函数，不同方法。 对 $\int_0^2x^4dx=32/5=6.4$，取 $n=2$、。梯形为 ，中点为 ，Simpson 为 。前两者确实一高一低。Simpson 仍有误差，因为被积函数是四次式。

改为 $n=4$，$S_4=77/12$，误差为 $1/60$，是原误差的 $1/16$。本例的四阶导数恒为 24，误差公式可以精确核对，不能把这种整齐比例强求在所有函数、所有步长上。

对于这几种正权重公式，所有权重之和为 $b-a$。若每个读数误差至多 $\delta$，输入误差造成的积分变化至多为 $(b-a)\delta$，不会出现差分中的 $1/h$ 放大。但真实积分若因正负抵消而接近零，相对误差仍可能很大；“积分有平均作用”不等于任何积分都具有很小的相对条件数。

用同一批读数逐列建立 Romberg 表

梯形公式在函数足够光滑时有偶次误差展开。例如 $f\in C^6[a,b]$ 时，

$T(h)=I+\frac{h^2}{12}\bigl[f'(b)-f'(a)\bigr]-\frac{h^4}{720}\bigl[f'''(b)-f'''(a)\bigr]+O(h^6).$

可以用局部展开核对这条式子。记一段的中点为 $m$，把函数展开到足够阶数，分别计算两端平均和积分，得到

$T_{\text{段}}-I_{\text{段}}=\frac{h^3}{12}f''(m)+\frac{h^5}{480}f^{(4)}(m)+O(h^7).$

同一段左右端点记为 $u,v$，对导数也作中心展开：

$f'(v)-f'(u)=hf''(m)+\frac{h^3}{24}f^{(4)}(m)+O(h^5),\qquad f'''(v)-f'''(u)=hf^{(4)}(m)+O(h^3).$

将后一组分别乘 $h^2/12$ 与 $-h^4/720$ 后相加，四阶导数前的系数为 $h^5(1/288-1/720)=h^5/480$，恰好重现单段误差。各段相加，内部端点的导数差消去，只留下 两端； 个 累加为 ，便得到前面的展开。我们需要的是带固定系数的偶次结构；不能仅从 这个上界，就推出所有后续偶次项都存在。

第一次 Richardson 消项得到 $(4T(h/2)-T(h))/3$，正好等于步长为 $h/2$ 的 Simpson 值。再把两个这样的结果按系数 16 组合，消去四次项。以 $h_k=(b-a)/2^k$、 记行，定义

$R_{k,0}=T(h_k),\qquad R_{k,j}=R_{k,j-1}+\frac{R_{k,j-1}-R_{k-1,j-1}}{4^j-1},\quad 1\le j\le k.$

只要相应的光滑性与误差展开成立，第 $j$ 列消去前 $j$ 个偶次项，误差阶为 $2j+2$。这里从第 0 列计数，所以分母是 $4^j-1$。

沿用 $x^4$ 在 $[0,2]$ 上的例子，表格可以完整算出：

同一积分的三行 Romberg 结果：首次消去二次误差项时用分母 3，下一次用分母 15。表中每个分数都来自相邻的两项，结果 32/5 是本例的精确积分。

细分时也不必重复计算所有函数值。若已有 $n$ 段的 $T_n$，新网格只多出各段中点：

$T_{2n}=\frac12T_n+\frac h2\sum_{j=0}^{n-1}f\left(a+(j+\tfrac12)h\right),\qquad h=\frac{b-a}{n}.$

旧数据的权重随步长减半，新中点各获得权重 $h/2$，这就解释了公式。高阶列主要增加算术运算，函数求值仍能复用。不过，函数不够光滑、噪声已经主导或舍入发生严重抵消时，表格右下角不一定越来越可靠。

下面逐次增加节点，查看梯形、Simpson 和 Romberg 的实际结果。选择四次多项式时，预测哪一列会精确；换成平方根函数后，观察原来的阶数判断为何不再适用。选中一个面板，还可以检查组成积分值的具体节点和权重。

自适应求积要分配误差预算

若函数只在局部变化很快，把所有小段都细分会浪费求值。自适应方法为每段比较粗、细两次结果，再决定是否继续分裂。

在当前区间上，记单个 Simpson 面板为 $S_c$，左右半区间各用一个 Simpson 面板后的总和为 $S_f$。若已经进入四阶误差规律，细结果的误差约为粗结果的 $1/16$。设 $S_c=I+e$、，相减得

$I-S_f\approx\frac{S_f-S_c}{15}.$

因此 $\widehat e=|S_f-S_c|/15$ 是细结果绝对误差的估计，修正值为

$Q=S_f+\frac{S_f-S_c}{15}.$

符号要留意：当粗细两次都高估，细值较小，修正项应为负。对 $x^4$ 在 $[0,2]$，$S_c=20/3$、$S_f=77/12$，差为 ；修正 后得到 ，本例恰好精确。

若当前区间分得的绝对预算为 $\tau$，当 $\widehat e\le\tau$ 时接受；否则把区间一分为二，两个孩子各得 $\tau/2$。继续二分后，所有叶子预算的总和仍为原预算。如果各个误差估计可靠，逐段误差相加才有受控的依据；不能给每个孩子同样的全局预算，然后期待总误差自动不超标。

非均匀取样的示意：部分区域使用更密的节点。真正的自适应程序依据误差估计分裂区间；图中的密度不是某个指定容差下运行所得的结果。

这个算法还需要实际的停止边界：限制最大深度与求值次数，发现中点已经与端点成为同一个浮点数时停止细分，遇到非有限函数值就报告失败。达到这些限制时，应返回“预算尚未满足”，不能冒充成功。

更根本的困难是：粗细结果相近，也可能只是都没看见真正的变化。取

$f(x)=\sin^2(4\pi x),\qquad 0\le x\le1.$

初始五个取样位置 $0,1/4,1/2,3/4,1$ 的函数值在精确算术中全为零，所以 $S_c=S_f=0$，估计误差为零；但真实积分是 。函数光滑得很好，错在采样没有发现振荡。浮点求值通常留下极小的非零数，仍然会被误认为足够准确。

下面可逐步展开求积树，查看每个区间为何被接受或分裂。先对光滑的四次式检查 $1/15$ 的估计，再换成这个漏采样反例，比较“程序报告的误差”和真实误差。增加初始分段或提供已知变化位置可以帮助发现问题，但有限次取样不能对任意未知函数给出无条件保证。

Gauss 求积把节点也作为未知量

前面的等距公式固定了节点，只设计权重。如果可以自由选择取样位置，能否用更少的求值对高次多项式精确？

在 $[-1,1]$，考虑对称的两点公式 $w[f(-r)+f(r)]$。对常数与二次函数要求精确，得到

$2w=2,\qquad 2wr^2=\frac23.$

于是 $w=1,r=1/\sqrt3$。一次、三次函数因对称性也精确，得到两点 Gauss–Legendre 公式

$G_2(f)=f(-1/\sqrt3)+f(1/\sqrt3).$

它用两个内部节点达到代数精度 3。注意这些节点不是第 4 章的 Chebyshev 根节点：两类节点解决的最优化或正交问题不同，不能互换。

两点 Gauss 的位置与权重共同满足低次多项式的积分条件；这里的精确性针对三次及以下多项式，不延伸为任意函数都精确。

一般的 $m$ 点公式为何能达到 $2m-1$？把次数为 $m$、首项系数为 1 的正交多项式记为 $\pi_m$，它满足

$\int_{-1}^1\pi_m(t)q(t)\,dt=0\quad\text{对所有次数小于 }m\text{ 的 }q.$

这样的多项式可以构造：从 $t^m$ 减去一个次数较低的多项式，让它与 $1,t,\ldots,t^{m-1}$ 的积分乘积都为零。系数满足一个线性系统；其矩阵的二次型是某个非零多项式平方的积分，严格为正，所以系统可逆。

它在 $(-1,1)$ 内有 $m$ 个互异根。否则，把区间内所有真正发生变号的根各取一次，构成次数小于 $m$ 的乘积 $q$；$\pi_mq$ 就不再变号，且不恒为零，积分不可能为零，与正交性矛盾。这迫使全部 $$ 个根都在内部且为单根。

取这些根为节点，权重仍按 $w_i=\int_{-1}^1L_i(t)dt$ 定义。任意次数不超过 的多项式 ，除以 后可写为

$P=q\pi_m+r,\qquad \deg q\le m-1,\quad\deg r\le m-1.$

正交性让第一项积分为零；在所有节点上，第一项也为零。余下的 $r$ 由 $m$ 点插值精确重现，所以它的积分等于加权求和值。两边把同一项消去，就证明了 $G_m(P)=\int_{-1}^1P$。

次数不能再普遍提高：$\pi_m^2$ 是 $2m$ 次，在全部节点上为零，但积分严格为正。因此任何仅用这 $m$ 个节点函数值的公式，都不可能对所有 $2m$ 次多项式精确。

Gauss 权重也确实为正。对次数 $2m-2$ 的 $L_i^2$ 使用刚证明的精确性，得到 $w_i=\int_{-1}^1L_i^2>0$；对常数求积又得到 。

实际区间为 $[a,b]$ 时，令 $c=(a+b)/2,d=(b-a)/2$，换元 $x=c+dt$：

$\int_a^b f(x)dx\approx d\sum_{i=1}^m w_i f(c+dt_i).$

节点要移动，权重也要乘 $d$。例如 $\int_1^3(x^3-x)dx$，两点为 $2\pm1/\sqrt3$，，两次求值之和为 ，与精确积分相等。对于非多项式，“代数精度高”仍需结合光滑性、区间宽度和误差检查；它不是任意函数的准确性保证。

端点函数有限，导数也可能无界

$f(x)=\sqrt x$ 在 0 处的值是 0，闭型梯形和 Simpson 都能计算。但它的导数在 0 附近无界，前面用有界二阶或四阶导数给出的通常误差保证不适用。不能因此把 $f(0)$ 也说成无穷大。

而 $f(x)=x^{-1/2}$ 本身在 0 处无界，包含端点求值的闭型公式就无法直接使用；积分却收敛，值为 2。使用内部节点虽然避开了无穷大的读数，也没有自动修复函数在端部的非光滑性。

这两个例子都可令 $x=t^2$、$dx=2t\,dt$。于是

$\int_0^1\sqrt x\,dx=\int_0^1 2t^2dt,\qquad \int_0^1x^{-1/2}dx=\int_0^1 2dt.$

第二个换元在 $t>0$ 时计算，再把变换后常数函数连续延伸到 0；不要在程序里先算出 $0\times\infty$。现在两者分别成为二次式与常数，低阶公式已经足够。判断方法是否适合，通常从这些条件检查开始，而不是从“哪个公式阶数最高”开始。

留下几个可以独立复现的计算

1. 只有 $f(x),f(x+h),f(x+3h)$ 三个读数。推导一个二阶的单边一阶导数公式，再求它的主误差项。

2. 中心差分满足 $|f'''|\le12$，每个读数误差至多 $10^{-5}$。用本章模型选步长。若改算二阶导数，每个读数误差经过公式后最多被放大成多少？

3. 对 $f(x)=x^5$ 在 $x=1$ 的中心差分，推导 $D_0(h)$ 以及一次 Richardson 外推后的结果。外推是否精确？

4. 估计 $\int_0^1x^2dx$，取两段，分别算梯形、中点、Simpson。另对 $\int_0^1e^xdx$，只使用 ，要求 Simpson 理论误差界不超过 ，最少取哪个正偶数 ？

5. 某函数在固定区间上的梯形结果为 $T_1=2$、$T_2=1.6$、$T_4=1.5$。填出两次外推列。仅凭右下角不再变化，能否认定已经得到精确积分？

6. 一段的 Simpson 粗细值分别为 0.8 和 0.797。该段预算为 $10^{-4}$。计算误差估计与修正值，决定接受还是分裂；若分裂，两个孩子各分多少预算？

7. 三点 Gauss–Legendre 的节点为 $-\sqrt{3/5},0,\sqrt{3/5}$，权重为 。分别算 与 在 的求积结果，与真值比较；解释为什么结论不矛盾。再把三点公式换到 ，计算 ，写明新节点与新权重。

8. 对 $\int_0^1x^{-1/3}dx$，说明它是否收敛、闭型公式能否直接求值，以及换元 $x=t^3$ 后应该计算哪个积分。