程序给出 0.7468241839,后面还可以继续打印许多位。我们该保留多少?如果别人问“它离真正要算的量差多远”,一串小数回答不了这个问题。
这一章把前面的方法放进一次完整计算里。我们会计算一个积分,检查权重和收敛,给出独立的误差依据,再把积分放到一个参数反求问题中。你会看到,什么时候应该细化网格,什么时候该收紧内层容差,以及什么时候继续计算已经碰到了输入数据的限制。
把要交付的量说清楚
本次第一个任务是
I=∫01e−t2dt,∣I−I∣≤10−6.区间、函数和目标误差都已经固定。这里要求的是绝对误差,不是“保留六位小数”,也不是“前后两次输出相差小于 10−6”。本例只是一个定义明确的数学积分,没有测量输入,因此没有另外的测量误差项。若将它解释为某个实际过程的模型,模型是否合适需要新的证据。
这一区别贯穿科学计算。计算核验关心代码是否实现了规定的方程、某次近似是否足够准确;模型验证则要看这些方程是否适用于所研究的实际对象。把已知解析解代回程序,能帮助核验算法;与独立观测作比较,才能讨论模型对现实的描述。两种工作互相支持,但不能互相替代。
比如两个程序都把 e−t2 误写成 e−t。它们可以各自表现出很整齐的收敛,甚至在许多位上相同,却都在逼近 1−e。所以在计较小数之前,还得核对输入表达式、区间和输出含义。

两类比较的示意:上排与已知数学结果核对,下排与观测数据比较。曲线形状仅用于说明对象不同,图上的接近本身并不能证明模型已经有效。
一份可以逐项检查的 Simpson 计算
被积函数在 [0,1] 光滑,我们采用第 5 章推导的复合 Simpson 公式。把区间分成正偶数 n 段,h=1/n,写成
Sn=这里每两个小区间构成一个面板,因而 n 必须是偶数。计算前可以用容易核对的函数检查实现:在 [0,1] 上,常数 1、t、t3 的积分分别为 1,1/2,1/4,Simpson 在精确算术下都应给出这些值。对于 ,第 5 章的误差公式给出
Sn−51=15n它应有误差,误差又可以精确算出。这个检查能发现一些“为了通过测试而把所有结果都算得像精确值”的错误。实际浮点测试则要给合适的舍入余量。
对当前的指数函数,我们还可以给出一个不依赖未知真实积分值的离散误差上界。求四次导数得到
f(4)(t)=(16t4−48t2+令 s=t2∈[0,1],多项式 16s2−48s+12 的导数为 ,因此其值从 12 降到 ,绝对值不超过 20。再用 ,得到 。
第 5 章已经证明,连续四阶导数和偶数分段保证
∣I−Sn∣≤180n420取 n=20,上界约为 6.94445×10−7。这个界控制的是精确函数值和精确算术下的 Simpson 近似。它没有自动把计算机对指数、节点、加法的误差一起包进去;这部分我们稍后用另一条独立路径检查。

四等分包含两个 Simpson 面板。共享节点的系数为 2,两处面板中点的系数为 4;所有系数乘公共因子 1/12 后,权重和等于 1。
网格加密后,数字怎样变化
同一函数、同一区间,分别计算 n=10,20,40:
| 分段数 | Simpson 近似值 | 理论离散误差界 |
|---|
| 10 | 0.746824948254443 | 约 1.11112×10⁻⁵ |
| 20 | 0.746824183875915 | 约 6.94445×10⁻⁷ |
| 40 | 0.746824136005348 | 约 4.34028×10⁻⁸ |
表中的小数经过显示舍入。用未提前截短的结果计算,相邻差值比为
S20−S40S10它接近 24=16。我们需要三张网格,才能形成两份相邻差并检查这个比值;只算两次,不能声称自己已经观测到了四阶。
按第四阶主误差模型,最细结果的误差估计为
S40−I≈15S20对应的外推值为
S40+15S40−S估计比前面的理论界小不少,并不矛盾:理论界允许整个区间都达到一个保守的导数上限。收敛比和外推依赖主误差展开,理论界依赖可验证的导数条件,它们的依据不同,报告时应分别标清。

保持函数与区间不变,步长连续减半。两个相邻差形成约 15.9676 的比值,用于检查四阶主误差模型;它与严格误差界的依据不同。
给这串小数一份独立的证据
我们不用把某个软件输出的更多位小数直接当作真值,而是从指数函数的 Taylor 余项得到上下界。
对 x≥0,e−x 的 Taylor 多项式为 ∑k=0m。Lagrange 余项的符号为 ,绝对值不超过 ,因为其中的指数因子介于 0 与 1 之间。因此,奇数次多项式在下方,偶数次多项式在上方。
代入 x=t2 并在 [0,1] 积分,定义
Pm=k=0∑mk!这是有限个有理数之和,可以用整数分子、分母精确计算。余项满足
∣I−Pm∣≤(m+1)!(2m+3)1尤其有
P11<I<P12,P两端约为 0.7468241327344955 和 0.7468241328180025。这些小数只是帮助阅读;进行严格比较时使用的是对应的精确分数。
现在把准备交付的结果固定为十位小数 q=0.7468241839,也就是一个确定的有理数。直接比较可得
P12<q,0<q−I<q−所以,这个具体报告值满足我们的绝对误差要求。此时指数函数库的每个内部舍入步骤不必再单独估计,因为我们已把最终报告值与独立的严格积分区间比较过。它证明的是这个数对这个数学积分足够准确,不是证明整段 Simpson 程序在所有输入上都正确,更没有证明任何现实模型。

独立的有理数上下界把真实积分夹在区间内,再与具体报告值比较。证书针对这个小数与这个积分,图中的框与箭头不表示数轴距离。
实验中可以逐次增加 Simpson 网格,观察近似、差值比、保守界和独立区间。切换到“错误权重”或“写错被积函数”时,注意哪些检查会报警、哪些检查仍可能看起来正常。把一个错误版本算得更细,通常并不会把它变成正确版本。
误差预算要有各自的对象
实际任务里,我们可能经历这样几次变化:真实对象用模型代替,真实数据用观测数据代替,连续模型变成离散问题,离散问题只迭代求解到某个容差,最后再用浮点数完成计算。
把关注的同一个标量输出分别记为 Q0,Q1,…,Q5,其中 Q 是真实目标, 是真实输入下的模型值, 是观测输入下的模型值, 是相应离散问题的精确值, 是有限迭代在精确算术下的输出, 是实际报告值。相邻项相减,精确地有
Q0−Q5=j=0∑若每一项都有相同输出单位下的上界,三角不等式给出
∣Q0−Q5∣≤这里没有假定误差互相独立,也没有保证它们会抵消。把平方和开根号当成总误差,通常需要额外的概率模型,不能随手用来替换最坏情形的和。若模型误差尚未评估,应写“未知”,不能填成零让预算表看起来完整。表达式写错、循环索引错等实现问题,也不能仅靠给它分配一个小预算就视为解决。
每项预算都应对应能做的检查。网格细化主要观察离散影响;保持网格不变而收紧停止条件,观察迭代影响;改变算术精度或求和顺序,检查有限精度敏感性;改动数据,检查输入敏感性。改变模型需要另外比较模型假设与适用情境。这些操作要分开进行,否则输出变化很难归因。
例如第 5 章的正权重求积公式在 [0,1] 上权重和为 1。若每个函数读数的绝对误差都不超过 δ,积分值的输入误差最多也是 δ。即便继续细化把离散误差压到 10−10,若读数允许 10 的同向偏差,仍不能据此报告 的总精度。

多种误差来源共同影响同一个输出的示意。色带宽度不表示实际大小,也不表示误差具有统计独立性;相加的是同一输出单位下已经得到的上界,未知项不能当成零。
积分放进求根之后
现在改成一个参数问题:
q(a)=∫01e−at2已知一个目标读数 b,希望从 q(a)=b 反求 a。每试一个 a,都需要算一次积分;外层求根和内层求积的容差开始互相影响。
先检查唯一性和敏感性。由于积分区间与参数区间都紧致,相关导数连续,可以用微积分基本定理写成
q(v)−q(u)=−∫uv∫因此
q′(a)=−∫01t2eq 严格递减,目标处于端点函数值之间时有唯一解。容易得到 ∣q′∣≤1/3,也有粗下界 ∣q′∣≥。为了后面使用精确有理数比较,我们再给一个可以直接算出的统一下界。
上一节已经说明,指数的三次 Taylor 多项式是下界。又因 a≤3/2,
e−at2≥e−3t2乘上 t2 积分便有
∣q′(a)∣≥31−由中值定理,对区间中任意两个参数都有
∣u−v∣≤m∣q(u)−q(v)∣.这一步把“积分值差一点”转成了“参数差多少”,正是残差与前向误差之间需要的桥梁。
内层不能含糊地交出一个符号
对一般 a,四阶导数为
dt4d4e在当前参数区间,它的绝对值不超过 16a4+48a3+12a2≤270。所以复合 Simpson 的统一离散误差界为
∣q(a)−Sn(a)∣≤n41.5n=64 时小于 9×10−8。若数值求积只告诉你 Sn(a),这个差比可能的求积误差还小,不能据此判定真实 为正。
假设内层确实提供了区间 [L(a),U(a)],保证包含 q(a)。这份保证必须覆盖内层的截断和算术误差。对递减的 q,二分法在中点 c 的处理是:
- L(c)>b:真实函数值仍偏大,根在 c 右边,可把左端移到 c。
- U(c)<b:真实函数值已经偏小,根在 c 左边,可把右端移到 。
一旦有
max{∣L(c)−b∣,∣U(c)−b∣}≤mεa,便可由中值定理保证参数误差不超过 εa。也可以保留已证实夹住根的区间,用区间宽度的一半作为其中点的误差界。两种停止依据要与最终报告的参数位置对应。
本例可以直接提供有理数区间
对 a≥0 定义
Pj(a)=k=0∑j仍由 Taylor 余项符号,有
P2r+1(a)≤q(a)≤P2r(a).当 a 是有理数时,两端可用有理数算术精确求出,不涉及指数函数的浮点求值;增加截断阶数可以缩窄区间。二分中点也始终是有理数,所以可以完整实现上面的符号判定。不要将这些精确分数提前四舍五入,再拿舍入值执行边界比较。
取 b=0.7,端点值分别约为 q(0.5)=0.85562439、q(1.5)=0.66335095,目标在二者之间。区间迭代会找到约 a=1.2646212;真正的停止依据应是所保留的参数区间或前面的残差界,而不是这几个显示数字。

在中点 1,第一份积分区间仍包含目标 0.7;增加级数项后,下界已经高于目标。利用积分函数递减,才能证实应保留右半段。
数据误差也会被反求放大
若读数 b 自身有误差,设真实参数 a∗ 满足 ∣q(a∗)−b∣≤δ,而数值结果为 。若求积近似 的误差上界为 ,数值残差满足 ,则
∣q(a)−q(a∗)∣≤ε从而在两个参数都属于上述区间时,
∣a−a∗∣≤mεεq 应包括内层全部已知误差,不能只取离散部分。若 δ=10−4、τ=10、,这个上界约为 。把外层残差再压低一百倍,改进仍很有限,主要预算来自读数。
这只是一个保证用的上界,不能把 δ/m 直接宣称为实际误差一定达到的“下限”。要证明数据本身无法支持某个精度,可以另找两组都与观测相容的参数,比较它们相隔多远。练习会让你做这件事。
在参数实验中,观察内层区间何时允许外层删去半边,何时需要增加级数项。调整读数误差后,图上的相容参数范围会改变;继续做更多外层迭代,只会缩小数值求解误差,不会自动消除数据不确定性。
让别人复算同一次计算
下面是一份可独立运行的 Python 代码。它把 Simpson 的近似计算与精确有理数检查分开。输入检查拒绝奇数分段;同一目标上的三张网格使用相同函数;用于交付的十位小数通过分数区间单独检验。
from math import exp, fsum, factorial
from fractions import Fraction
def simpson(f, n):
if not isinstance(n, int) or n < 2 or n % 2:
raise ValueError
这段代码没有随机取样,无需凭空加入随机种子。它也没有声称不同设备上 exp 的所有末位都完全相同:浮点运算的顺序、函数库和精度可能影响末位。复算时既要记录环境,也要明确判断标准是逐比特一致,还是在同一数学误差要求内一致。后者仍需要误差依据,不能任意设一个很大的比较容差。
本次报告可以写成下面这样:
目标为 ∫01e−t2dt,要求绝对误差不超过 10。采用复合 Simpson,核对正偶数分段与权重; 的相邻差比约为 15.9676。交付值为 0.7468241839。独立的有理数上下界 证明该报告值的误差小于 。本结论针对给定数学积分,不包含任何未经验证的物理模型或测量误差。
实际保存计算记录时,应连同输入、代码版本、运行环境、完整中间值和停止原因一起保留。只留下截图上的曲线,很难复现某个边界处理;只留下小数,又看不出那些位数为什么可信。
练习:下一步应该改哪里
1. 换一个指数参数。 计算 ∫01e−2t2dt,采用 Simpson。用 的保守界,求保证精确算术下离散误差不超过 所需的最小正偶数分段数。它是否也自动保证计算机最终输出的全部误差?
代入 a=2,得到导数界 256+384+48=688。要求 688/(180n4)≤,即 。 不够, 足够,所以最小正偶数为 46。这只控制离散误差,函数求值与算术误差仍要另计,或用独立的严格结果区间核对最终报告值。
2. 小数舍入也要检查。 使用正文的精确区间 [P11,P12],判断报告值 0.746824 是否满足绝对误差 10−6,是否满足 。需要把它打印成更多位才能作判断吗?
这个报告值小于两个端点,误差大约为 1.328×10−7。用分数比较可证它小于 10−6,也可证它大于 10−8。因此前一个要求满足,后一个不满足;增加显示的零不改变报告值,不能改变结论。
3. 内层还没有确定符号。 求递减函数 q(a)=b 的根,当前中点 c 的已证实函数区间为 [b−2×10−5,。能否删去某一半区间?若导数绝对值下界为 ,能否用 交付参数误差不超过 的结果?
区间包含 b,真实差值可能为正也可能为负,因此不能仅凭这份信息决定删去哪一半。但残差绝对值至多为 3×10−5,参数误差至多为 3×10−4,已经小于所要求的 。若没有数据误差等额外项,可以在 停止;不必为了确定符号而继续计算。
4. 给数值部分分配预算。 某参数反求问题有 ∣q′∣≥0.2,观测误差界为 2×10−5。希望参数总误差界不超过 2×。内层误差界 与外层残差容差 之和最多能取多少?若二者各取 ,是否满足?
需要 εq+τ+2×10−5≤0.2×,所以数值两项之和至多为 。各取 会使总界达到 ,超过目标;例如各取 才刚好满足这个预算。
5. 证明数据确实限制了精度。 假设正文的 q 在参数区间里分别有 q(a−)=b+δ、q(a,。利用 ,证明两组参数至少相距 。当观测只有 时,为什么任何一个单点报告都无法保证对这两种可能同时达到严格小于 的误差?
中值定理给出 2δ=∣q(a−)−q(a+)∣≤,故间距至少为 。任意单点到两端的距离之和不小于两端间距,所以至少一个距离不小于 。两端都与观测相容,额外迭代不能区分它们。这里得到的是由相容参数构造出的不可区分性结论,不是把某个保守上界误当作实际误差下限。
6. 两个方法一致,问题仍可能写错。 梯形与 Simpson 都使用程序中同一个函数 exp(-t),结果在细网格上接近;需求却是 exp(-t*t)。哪些检查能直接揭露问题,哪些证据不足以消除它?
核对函数表达式或在 t=1/2 检查函数值,会发现应为 e−1/4,程序却给出 e−1/2。将输出与正文对目标积分建立的独立分数区间比较,也会失败。两个方法的细网格一致、各自呈现预期收敛阶,只能说明它们可能准确求解了共同实现的另一个问题,不能核验需求表达式。
7. 检查换元是否仍在算同一个量。 若任务改成 J=∫02e−t2dt,用 t 变到 。一份代码算出 后就直接报告为 。网格加密和高精度计算能修复这个错误吗?写出正确关系。
dt=2du,所以 J=2∫01e−4。原代码遗漏了雅可比因子 2。细化和提高精度会更准确地得到 ,不能修复模型转换中的错误。需要改正表达式,重新运行相关核验。
8模型误差尚未评估时,可以在总误差表中先填 0,只要离散误差和求解残差已经很小。