数学分析研究极限、导数和积分的精确性质;数值分析进一步研究如何用有限的存储和有限次运算得到这些对象的近似值,以及所得结果有多可靠。一个公式在实数运算中成立,并不保证直接照着计算就能得到准确结果。
本篇先建立误差的基本语言,再讨论浮点运算、问题的条件性与算法稳定性。这些概念将贯穿后续的方程求根、插值、数值积分和微分方程数值解。
本文主要讨论实数标量问题,约定 x 为真值,x 为近似值,fl 表示指定浮点系统中的舍入或浮点运算结果。若没有特别说明,误差取绝对值,相对误差以真值的绝对值为分母。所有小量近似都会注明其条件。
1. 数值计算中的三类误差
1.1 数据误差
测量、参数估计或输入数据的近似表示,可能使实际参与计算的数据与数学问题中的数据不同。例如,以 1.41 代替 2 作为输入,即使后面的计算完全精确,输出仍可能存在偏差。
数据误差的影响取决于问题本身对输入扰动的敏感程度。
1.2 截断误差
定义(截断误差)
用有限过程替代无限过程,或用离散表达式替代连续对象时,在精确算术下引入的近似误差称为截断误差。
例如,若 f 在包含 x 与 x+h 的区间上二阶连续可微,h=0,则存在介于两点之间的 ξ,使
hf(x+h)−f(x)−f′(x)=2hf′′(ξ).
因此,当二阶导数在固定邻域内有界时,上式右端为 O(h)。这里的误差已经存在于精确实数运算中。
1.3 舍入误差
计算机只能保存有限精度的数。输入、运算结果或中间量无法精确表示时,需要用可表示的浮点数替代,由此产生舍入误差。
这三类误差描述不同来源,不意味着最终相对误差总能简单相加。例如,输入舍入产生的微小偏差,可能在后续运算中被放大。
一个精确的误差分解
设问题的目标值为 F(d),输入数据为 d;Ah 是由离散参数 h 决定的精确算术算法,实际计算结果为 y。只要各项有定义,就有
y−F(d)=[y−Ah(d)]+[Ah(d)−F(d)]+[F(d)−F(d)].
三项依次是计算过程中的舍入影响、离散近似误差和输入数据扰动造成的输出变化。取绝对值后可用三角不等式得到总误差上界,但不同来源的实际误差也可能相互抵消。
2. 绝对误差、相对误差与误差限
2.1 误差的度量
定义
有符号误差、绝对误差和相对误差分别为
e=x−x,eabs=∣x−x∣,
erel=∣x∣∣x−x∣,x=0.
绝对误差与原变量具有相同的量纲,相对误差没有量纲。真值为零时,上述相对误差没有定义,不能把分母强行改为一个很小的数而仍称其为这个定义下的相对误差。
例如,以 1.41 近似 2 时,
eabs≈0.00421356237309505,erel≈0.00297943852696799.
相对误差百分比约为 0.297944%。显示更多小数位只是更详细地报告结果,并不会自动减小误差。
2.2 误差限
定义(绝对误差限)
若已知 ε≥0,且
∣x−x∣≤ε,
则称 ε 为一个绝对误差限。此时
x∈[x−ε,x+ε].
误差限是上界,不一定等于实际误差。
性质(由绝对误差限控制相对误差)
若 ∣x∣>ε,则真值不为零,且
erel≤∣x∣−εε.
这是由 ∣x∣≥∣x∣−ε 得到的严格上界。只有当 ε 相比 ∣x∣ 很小时,才可用 ε/∣x∣ 作近似估计。
2.3 绝对容差与相对容差
不同尺度的问题通常需要不同的误差标准。若真值可用于验证,可采用
∣x−x∣≤atol+rtol∣x∣,
其中 atol≥0 是绝对容差,rtol≥0 是相对容差。
真值接近零时,绝对容差提供尺度;真值较大时,相对容差控制比例。实际算法通常不知道真值,必须结合误差估计或经过论证的停止指标,而不能直接检查这个含未知量的不等式。
先选择“大尺度”和“小尺度”预设,再比较同一个绝对误差在两种尺度下对应的相对误差。
图中的数轴会自动缩放,因此应结合刻度和数值读数比较,不能只凭两个点在屏幕上的距离判断误差。实验以输入的浮点数作为基准值;任意十进制输入本身也可能发生舍入。
3. 有效数字与误差阶
3.1 有效数字的误差判据
定义(本文采用的约定)
设 x=0,并写成
∣x∣=m×10k,1≤m<10,k∈Z.
若
∣x−x∣≤2110k−n+1,n≥1,
则称 x 具有至少 n 位有效数字。右端是从第一个非零数字数起,第 n 位所对应单位的一半。
例如,1.41 近似 2 的绝对误差小于 0.005,但大于 0.0005,所以按上述约定具有三位有效数字。对 x=0,应直接说明绝对误差限或精确到的小数位,而不用这个科学记数法定义。
“保留了几位数字”描述书写格式,“具有几位有效数字”需要误差依据。接近进位边界时,不能仅凭真值和近似值的小数字符串是否一致来判定。
3.2 小数位数与有效数字位数
将一个数舍入到小数点后 q 位,在采用就近舍入且没有其他误差时,绝对舍入误差不超过
2110−q.
固定小数位数对应固定的绝对尺度;有效数字位数则从第一个非零数字开始计算。因此,0.00123 有三位有效数字,但有五位小数。
3.3 误差阶
定义
若存在与 h 无关的常数 C>0 和 h0>0,使
∣E(h)∣≤C∣h∣p,0<∣h∣<h0,
则记 E(h)=O(hp),其中 p>0。
它是一个渐近上界,并不单独保证将步长减半就会使误差恰好缩小为原来的 2−p。若进一步有
E(h)=chp+o(hp),c=0,
则在 h→0+ 时,∣E(h/2)∣/∣E(h)∣→2−p。舍入误差占主导时,实际计算可能不再呈现这一比例。
4. 浮点数与舍入模型
4.1 规范化浮点数
在基数为 β≥2、有效位数为 t 的浮点系统中,非零规范化浮点数可表示为
±(d0.d1⋯dt−1)ββe,
其中 1≤d0≤β−1,0≤di≤β−1,各数字均为整数,且指数 e 在有限范围内变化。
有效位数限制精度,指数范围限制可表示数的大小。在同一个区间 [βe,βe+1) 内,相邻正规范化浮点数的间距为
βe−t+1.
因此浮点数在实轴上不是等间距分布的;数值越大,其所在区间内的绝对间距通常也越大。
4.2 单位舍入误差
性质(就近舍入的相对误差模型)
对处于规范化数范围内的实数 z=0,在就近舍入且不发生溢出等异常的条件下,可以写成
fl(z)=z(1+δ),∣δ∣≤u,
其中
u=21β1−t
称为单位舍入误差。这是统一上界,不意味着每次舍入都产生大小为 u 的相对误差。
在二进制系统中,大于 1 的下一个浮点数与 1 的间距为
εgap=21−t=2u.
“机器精度”或 machine epsilon 在不同资料中可能指 u,也可能指这个间距。引用数值前应先确认定义。
4.3 双精度中的典型数值
IEEE 754 binary64 的规范化有效位数为 53 个二进制位。因此,就近舍入下
u=2−53≈1.1102230246251565×10−16,
εgap=2−52≈2.2204460492503131×10−16.
这也解释了为何通常只能期待大约十五至十六位十进制有效数字,而不能把输出格式设为二十位就认为得到了二十位可靠结果。
十进制有限小数也未必能用有限二进制表示。例如 0.1 的二进制展开无限循环,通常需要先舍入才能存入 binary64。
4.4 上溢、下溢与非数值
- 上溢:结果的大小超出有限可表示范围;默认异常处理可能得到无穷大。
- 下溢附近的精度损失:结果接近零,进入非规范化数范围,可能无法维持统一的小相对误差。
- 非规范化数:用较少的有效位填补最小正规范化数与零之间的部分间隙,支持渐进下溢。
- NaN:表示无效运算等产生的非数值结果,例如浮点计算中的 0/0。
非规范化数不等于零,但不能无条件套用规范化范围内的相对误差界;此时可能需要绝对误差界。
5. 基本运算与误差传播
5.1 单步浮点运算
对已存储的浮点数 a,b,若运算有定义、精确结果位于规范化范围内,并采用正确的就近舍入,则
fl(a∘b)=(a∘b)(1+δ),∣δ∣≤u,
其中 ∘∈{+,−,×,/},除法还要求 b=0。精确结果为零且被准确表示时,可取 δ=0。
注意这里的 a,b 是已经存储的操作数;它们相对原始真值的输入误差,需要另行分析。单步舍入误差小,不代表整个表达式相对其数学真值的误差小。
5.2 数据误差的一阶传播
设 a=a+Δa、b=b+Δb。暂时只研究精确运算下的数据扰动,有
∣(a±b)−(a±b)∣≤∣Δa∣+∣Δb∣.
若 a±b=0,两边除以 ∣a±b∣ 就得到相对误差界。分母很小时,即使两个操作数的误差很小,结果的相对误差也可能很大。
若 a,b=0,写成 a=a(1+ηa)、b=b(1+ηb),则
abab−ab=ηa+ηb+ηaηb,
而当 1+ηb=0 时,
a/ba/b−a/b=1+ηbηa−ηb.
所以对充分小的输入相对误差,乘除法的输出相对误差在一阶近似下由两者之和控制。实际浮点运算还要计入运算本身的舍入。
5.3 一般函数的误差传播
若 f 在 x 处可微,则
f(x+Δx)−f(x)=f′(x)Δx+o(∣Δx∣).
它描述小扰动下的局部变化。若需要严格上界,可用中值定理:当 f 在连接两点的闭区间上连续、内部可微且 ∣f′(t)∣≤M 时,
∣f(x+Δx)−f(x)∣≤M∣Δx∣.
对在某个凸邻域内连续可微的多元函数 F,若该邻域内 ∣∂F/∂xi∣≤Mi,且连接 x 与 x+Δx 的线段在邻域内,则
∣F(x+Δx)−F(x)∣≤i=1∑mMi∣Δxi∣.
6. 相消与等价公式的数值表现
6.1 相近数相减
设 a,b>0,a=b,且两者的相对输入误差均不超过 η。在精确减法下,
∣a−b∣∣(a−b)−(a−b)∣≤∣a−b∣a+bη.
当 a 与 b 接近时,因子 (a+b)/∣a−b∣ 很大。减去相同的高位数字后,原来隐藏在低位的输入误差可能成为结果的主要部分,这就是相消造成精度损失的典型机制。
相消不意味着减法这一步一定发生很大的舍入误差:两个已存储浮点数的差甚至可能被精确算出,但结果仍会继承并放大操作数之前已有的误差。
6.2 平方根差的改写
考虑 x>0 时计算
g(x)=1+x−1=1+x+1x.
当 x 很小时,直接形式先算出一个很接近 1 的数,再减去 1;改写形式的分母接近 2,避免了这次相近数相减。
例如,精确数学输入 x=10−16 时,
g(x)≈4.999999999999999875×10−17.
在 binary64 就近舍入下,直接计算中的 1+x 会舍入成 1,于是直接形式得到 0。改写形式仍能保留这个非零小量。这里的问题主要来自求值方式,不能据此断言 g 在小正数处本身病态。
将下面实验的指数调至 −16,观察两种表达式的区别;再切换到十进制有效位模拟,观察精度降低后相消出现的位置如何变化。
实验的“参考值”使用改写公式在浏览器浮点精度下计算,并非任意精度真值。在浏览器模式下,改写结果与参考值由同一公式计算,差异为零是构造所致,不能解释为数学误差严格为零。十进制模式逐步舍入输入及中间量,用于演示有效位损失,不模拟完整的浮点指数范围和异常处理。
7. 问题的条件性与条件数
7.1 绝对条件数与相对条件数
定义(标量函数的局部条件数)
设 f 在 x 的某个开邻域内有定义,并在 x 处可微。绝对条件数为
κabs(f,x)=∣f′(x)∣.
当 x=0 且 f(x)=0 时,相对条件数为
κrel(f,x)=f(x)xf′(x).
由可微性,在 Δx→0 时有
∣f(x)∣∣f(x+Δx)−f(x)∣≤κrel(f,x)∣x∣∣Δx∣+o(∣x∣∣Δx∣).
条件数衡量指定问题、指定输入和指定误差度量下的局部敏感性,不依赖采用哪一种求值算法。多维情形还依赖范数和允许扰动的数据范围。
7.2 典型例子
| 问题 | 条件 | 相对条件数 |
|---|
| f(x)=xm | x>0,m∈R | ∣m∣ |
| f(x)=x | x>0 | 1/2 |
| f(x)=ex | x=0 | ∣x∣ |
| f(x)=1/(1−x) | x=0,1 | ∣x/(1−x)∣ |
一般函数的相对条件数不一定大于等于 1。绝对条件数与相对条件数也可能给出不同判断:平方根函数在 x→0+ 时的绝对条件数趋于无穷,但其相对条件数恒为 1/2。
对上一节的平方根差,直接求导并整理可得
κrel(g,x)=21+x1+x+1,x>0.
因此当 x→0+ 时,它趋于 1。问题在这种相对误差度量下并不病态,直接公式的精度损失可以通过改写缓解。
7.3 病态问题与局部近似的限制
当条件数很大时,称问题在该输入附近病态;“很大”的判断应结合输入精度和输出精度要求。条件数描述最坏局部放大尺度,并不保证每个扰动都达到同样的放大效果。
以 f(x)=1/(1−x) 为例,若 x=0.99、Δx=0.0001,则 κrel=99,输入相对变化约为 0.010101%,而输出相对变化为
∣f(0.99)∣∣f(0.9901)−f(0.99)∣=991≈1.010101%.
局部线性预测为 1%,与实际值接近但不完全相同。若扰动已接近奇点到输入点的距离,就不能继续把一阶预测当作精确结论。
选择倒数函数预设,逐渐让输入靠近 1;比较实际输出相对变化和“条件数乘输入相对变化”的线性预测。
实验将输入写成 x(1+η),并用代数化简后的相对变化公式计算曲线,减少直接相减引入的额外误差。图上的变化反映问题的敏感性,不能直接用来评价某个求值算法的稳定性。
8. 前向误差、后向误差与数值稳定性
8.1 前向误差
对问题 y=f(x),设某算法给出的结果为 y。前向绝对误差是
∣y−f(x)∣;
若 f(x)=0,前向相对误差为
∣f(x)∣∣y−f(x)∣.
前向分析直接问:计算结果距离原问题的精确答案有多远?
8.2 后向误差
后向分析问:能否把 y 看作某个邻近输入问题的精确答案?即寻找 Δx 使
y=f(x+Δx).
当 x=0 时,定义相对后向误差
ηback=inf{∣x∣∣Δx∣:x+Δx 在定义域内, f(x+Δx)=y}.
若集合为空,则取 +∞。输入为零时可改用绝对后向误差。定义还必须说明允许扰动哪些数据;方程系数、右端项或全部数据的扰动,会产生不同的后向误差问题。
例子
计算 2 得到 y=1.41,它是邻近输入 1.9881 的精确平方根。因此
Δx=1.9881−2=−0.0119,ηback=20.0119=0.00595.
这一后向误差不同于第 2 节的前向相对误差;二者回答的是不同问题。
8.3 后向稳定性
定义(后向稳定的基本含义)
若算法在规定的适用范围内,能够把计算结果解释为原始数据受到 O(u) 量级相对扰动后的精确解,则称其后向稳定。这里的常数应受控,可以依赖问题规模,但不能隐藏任意大的数据敏感性。
当这种后向扰动足够小、局部条件数分析有效时,前向误差的一阶尺度由
相对前向误差≲相对条件数×相对后向误差
描述。符号 ≲ 在这里表示忽略高阶项的局部估计,不是对任意大扰动成立的无条件上界。
所以,后向稳定算法在良态问题上通常能给出小的前向误差;在病态问题上,前向误差仍可能很大。
8.4 残差与误差
对方程 F(x)=0,近似解 x 的残差是 F(x),解误差是 x−x∗,其中 x∗ 是所考察的真根。
若 F 在连接两点的区间上连续、内部可微,并且该区间上 ∣F′(t)∣≥m>0,则中值定理给出
∣x−x∗∣≤m∣F(x)∣.
缺少这样的控制时,小残差不能单独保证小误差。例如,F(x)=10−12(x−1),取 x=2,残差仅为 10−12,解误差却是 1。把方程整体乘以小常数就能缩小残差,而根和解误差没有改变。
9. 计算表达式的组织
9.1 浮点加法不满足结合律
在 binary64 就近舍入下,令 A=254。对下列括号顺序分别逐步计算,有
fl(fl(A−A)+1)=1,
fl(A+fl(−A+1))=0.
原因是 −A+1 位于两个相邻可表示数的中间,就近舍入的偶数规则将其舍入为 −A。实数加法结合律成立,但在每步舍入的计算中不能任意改变括号后仍期待同样结果。
9.2 顺序求和的误差界
设 x1,…,xn 已经存储为浮点数,n≥2,按顺序逐项求和。假设各步满足第 5 节的误差模型,且 (n−1)u<1。记
s=i=1∑nxi,γk=1−kuku.
则计算结果 s 满足
∣s−s∣≤γn−1i=1∑n∣xi∣.
若 s=0,从而
∣s∣∣s−s∣≤γn−1∣∑ixi∣∑i∣xi∣.
后一个比值反映求和问题对分量相对扰动的敏感性。正负项强烈抵消时,比值可能很大;“每步只有很小的舍入误差”不足以保证最终相对误差小。
采用平衡的成对求和可减少误差累积的计算深度;补偿求和则尝试记录部分丢失的低位信息。按绝对值排序后求和有时有帮助,但不能把某一种顺序当作对任意数据都最准确的规则。
9.3 Horner 算法
对多项式
p(x)=anxn+an−1xn−1+⋯+a1x+a0,
Horner 算法将其写为嵌套形式。计算步骤为:
- 输入系数 a0,…,an 和求值点 x,令 bn=an。
- 对 k=n−1,n−2,…,0,依次计算 bk=ak+xbk+1。
- 输出 b0。
精确算术下 b0=p(x),共需 n 次乘法和 n 次加法。若各步单独舍入、满足标准模型且 2nu<1,则对已存储系数和求值点,有
∣p(x)−p(x)∣≤γ2nk=0∑n∣ak∣∣x∣k.
它也可解释为系数受到小的分量相对扰动后的精确求值,但接近多项式的根时,右端与 ∣p(x)∣ 之比可能很大。因此减少运算次数和具有良好的后向性质,都不自动保证在所有求值点获得很小的相对前向误差。
例如,x3−3x2+3x−1=(x−1)3。在 x 接近 1 时,展开式中的项容易相互抵消;因式形式可能改善给定输入的求值,但无法消除输入 x 自身的不确定性。
10. 概念速查与常见误区
10.1 概念对照
| 概念 | 所回答的问题 | 需要注意 |
|---|
| 绝对误差 | 结果偏离真值多少? | 与原量有相同量纲。 |
| 相对误差 | 偏差相对真值尺度有多大? | 真值为零时没有定义。 |
| 截断误差 | 数学近似本身引入多少偏差? | 在精确算术下也存在。 |
| 舍入误差 | 有限表示和运算造成多少偏差? | 需说明精度与舍入方式。 |
| 条件数 | 数据小扰动会放大多少? | 属于指定问题和输入。 |
| 前向误差 | 计算结果离原问题答案多远? | 度量输出偏差。 |
| 后向误差 | 改变多少输入才能使结果精确成立? | 需规定允许扰动的数据。 |
| 残差 | 近似解多大程度满足原方程? | 不能直接当作解误差。 |
10.2 显示更多位数不等于精度更高
显示位数、浮点格式精度和结果的有效数字是三件不同的事。只有减小误差或证明更强的误差界,才能提高可确认的准确程度。
10.3 绝对误差小不等于相对误差小
必须结合真值尺度解释误差。对接近零的结果,纯相对标准可能不合适,应明确使用绝对误差或混合容差。
10.4 问题病态不等于算法不稳定
病态性描述输入扰动的影响;算法稳定性描述计算方式是否引入了不必要的大误差。应分别分析,不能仅凭一次计算不准确就判定原因。
10.5 相消不等于减法必然算错
相近数的差可能被精确表示,但操作数已有的误差仍可能严重影响结果。应追踪完整计算过程,而不是只检查最后一步。
10.6 步长更小、迭代更多不一定更准确
减小步长可能降低截断误差,也可能放大舍入误差或数据噪声;增加迭代次数也可能遭遇舍入误差平台。停止准则需要结合算法的误差分析。
11. 本篇要点
- 数据误差、截断误差和舍入误差来源不同,应分别分析其对输出的影响。
- 绝对误差描述偏差大小,相对误差描述比例;误差限是有依据的上界。
- 有效数字需要误差保证,误差阶需要明确渐近变量和适用条件。
- 浮点数的间距随尺度变化;单位舍入误差与 1 附近的浮点间距要区分。
- 小的单步舍入误差可能经过相消或其他敏感运算而显著放大。
- 条件数属于问题,稳定性属于算法;后向稳定并不保证病态问题的前向误差小。
- 公式改写、求和组织和 Horner 算法可以改善计算,但都需要结合数据与条件性判断。
- 可视化实验展示的是有限精度计算现象,参考值和实际误差的来源必须说明。
下一篇将整理二分法与不动点迭代,进一步讨论求根过程的误差上界、收敛条件和停止准则。
参考资料
- Nicholas J. Higham, Accuracy and Stability of Numerical Algorithms, 2nd ed., SIAM, 2002:浮点误差分析、求和与多项式求值。
- David Goldberg, What Every Computer Scientist Should Know About Floating-Point Arithmetic, ACM Computing Surveys, 1991,Oracle 收录版本:浮点表示、舍入与相消。
- LAPACK Users’ Guide, How to Measure Errors:绝对误差和相对误差的度量。
- Oracle, Numerical Computation Guide,IEEE Arithmetic 与 Underflow 相关章节:双精度格式和渐进下溢。
- Nick Higham, What Is a Condition Number?;What Is Backward Error?, 2020。