1.2.5 给出了浮点数集、机器精度与「每步运算 = = = × ( 1 + δ ) \times(1+\delta) × ( 1 + δ ) 推导 出来（它不是假设，而是 IEEE 舍入规则的推论），再给出四组常用工具的完整分析：Sterbenz 引理、求和误差界、Kahan 补偿求和、以及两个经典现象（0.1 + 0.2 ≠ 0.3 0.1+0.2\neq0.3 0.1 + 0.2 = 0.3
（一）数的舍入。 设 β = 2 \beta=2 β = 2 t t t t = 53 t=53 t = 53 [ 2 e , 2 e + 1 ) [2^e,\ 2^{e+1}) [ 2 e , 2 e + 1 )
ulp ( x ) = 2 e + 1 − t ( x ∈ [ 2 e , 2 e + 1 ) ) . \operatorname{ulp}(x)=2^{,e+1-t}\qquad(x\in[2^e,2^{e+1})) . ulp ( x ) = 2 e + 1 − t ( x ∈ [ 2 e , 2 e + 1 )) .
设实数 x x x 舍入到最近 （round to nearest, ties to even）：fl ( x ) \operatorname{fl}(x) fl ( x ) x x x
∣ fl ( x ) − x ∣ ≤ 1 2 ulp ( x ) = 2 e − t , ∣ fl ( x ) − x x ∣ ≤ 2 e − t x ≤ 2 e − t 2 e = 2 − t = u . \lvert\operatorname{fl}(x)-x\rvert\ \le\ \tfrac12\operatorname{ulp}(x)=2^{,e-t},\qquad \left\lvert\frac{\operatorname{fl}(x)-x}{x}\right\rvert\ \le\ \frac{2^{,e-t}}{x}\ \le\ \frac{2^{,e-t}}{2^{,e}}=2^{-t}=u . ∣ fl ( x ) − x ∣ ≤ 2 1 ulp ( x ) = 2 e − t , x fl ( x ) − x ≤ x 2 e − t ≤ 2 e 2 e − t = 2 − t = u .
（次正规区间内间距只会更细，同一结论保持。）
（二）运算的舍入。 IEEE 754 对 + , − , × , ÷ , +,-,\times,\div,\sqrt{\ } + , − , × , ÷ , 先把运算按实数精确算出，再把结果舍入一次 （correct rounding——「先算准、再存近」）。设 z = x o p y z=x\ \mathrm{op}\ y z = x op y
fl ( x o p y ) = z ( 1 + δ ) , ∣ δ ∣ ≤ u . \operatorname{fl}(x\ \mathrm{op}\ y)=z(1+\delta),\qquad\lvert\delta\rvert\le u . fl ( x op y ) = z ( 1 + δ ) , ∣ δ ∣ ≤ u .
这就是 1.2.5 的基本公理——每一步运算的相对误差不超过 u u u 。「只有一次」是关键：它排除了老式硬件「先舍入操作数、再算、又舍入」的双重误差（历史上的 x87 80 位扩展寄存器造成的跨平台不可复现，正是因为不同平台「何时舍入到 64 位」不一致）。
（三）哪些运算享受这待遇。 + , − , × , ÷ , +,-,\times,\div,\sqrt{\ } + , − , × , ÷ , sin , exp , ln \sin,\exp,\ln sin , exp , ln 超越函数不强制 正确舍入（正确舍入超越函数至今是研究课题），所以不同数学库的 sin(x) 可以在最后一位上合法地不同。比较浮点相等时容差 ∼ u \sim u ∼ u 0 0 0
定义 1.2.10 同时给了 e p s = β 1 − t \mathrm{eps}=\beta^{1-t} eps = β 1 − t u = 1 2 β 1 − t u=\tfrac12\beta^{1-t} u = 2 1 β 1 − t
量 值（β = 2 , t = 53 \beta=2,t=53 β = 2 , t = 53 含义
e p s \mathrm{eps} eps 2 − 52 ≈ 2.22 × 1 0 − 16 2^{-52}\approx2.22\times10^{-16} 2 − 52 ≈ 2.22 × 1 0 − 16 1.0 1.0 1.0 间距 （np.spacing(1.0)、np.finfo(float).eps）u u u 2 − 53 ≈ 1.11 × 1 0 − 16 2^{-53}\approx1.11\times10^{-16} 2 − 53 ≈ 1.11 × 1 0 − 16 舍入的最大相对误差 （单位舍入误差）
一句话：间距是 e p s = 2 u \mathrm{eps}=2u eps = 2 u u u u ——「网格上相邻两点的距离」是「落点误差」的两倍。正文「1.0 + 2 − 53 1.0+2^{-53} 1.0 + 2 − 53 1.0 1.0 1.0 1.0 + 2 − 52 1.0+2^{-52} 1.0 + 2 − 52 2 − 53 2^{-53} 2 − 53 1 1 1 平局位置 ，平局判偶（mantissa 为偶者胜）落回 1.0 1.0 1.0 e p s = 2 u \mathrm{eps}=2u eps = 2 u
顺带一个常用公式（ulp 的显式版本）：对 x ≠ 0 x\neq0 x = 0
ulp ( x ) = 2 ⌊ log 2 ∣ x ∣ ⌋ + 1 − t ( 次正规区间内取最小间距 2 e min − t + 1 ) . \operatorname{ulp}(x)=2^{,\lfloor\log_2\lvert x\rvert\rfloor+1-t}\qquad(\text{次正规区间内取最小间距}\ 2^{,e_{\min}-t+1}). ulp ( x ) = 2 ⌊ l o g 2 ∣ x ∣⌋ + 1 − t ( 次正规区间内取最小间距 2 e m i n − t + 1 ) .
引理 B.3 · Sterbenz：浮点减法的精确性
x , y x,y x , y y / 2 ≤ x ≤ 2 y y/2\le x\le 2y y /2 ≤ x ≤ 2 y fl ( x − y ) = x − y \operatorname{fl}(x-y)=x-y fl ( x − y ) = x − y 精确成立 。
证明
同号不妨设 x ≥ y > 0 x\ge y\gt0 x ≥ y > 0 y / 2 y/2 y /2 x x x y y y y ≤ x ≤ 2 y y\le x\le2y y ≤ x ≤ 2 y y = m 2 E y=m\,2^{E} y = m 2 E x = n 2 F x=n\,2^{F} x = n 2 F m , n ∈ [ 2 t − 1 , 2 t ) m,n\in[2^{t-1},\ 2^{t}) m , n ∈ [ 2 t − 1 , 2 t ) y ≤ x ≤ 2 y y\le x\le2y y ≤ x ≤ 2 y
2 E + t − 1 ≤ m 2 E ≤ n 2 F ≤ 2 m 2 E < 2 E + t + 1 ⟹ x ∈ [ 2 E + t − 1 , 2 E + t + 1 ) ⟹ F ∈ { E , E + 1 } . 2^{E+t-1}\le m\,2^E\le n\,2^F\le 2m\,2^E\lt2^{E+t+1}\ \Longrightarrow\ x\in[2^{E+t-1},\ 2^{E+t+1})\ \Longrightarrow\ F\in\{E,\ E+1\}. 2 E + t − 1 ≤ m 2 E ≤ n 2 F ≤ 2 m 2 E < 2 E + t + 1 ⟹ x ∈ [ 2 E + t − 1 , 2 E + t + 1 ) ⟹ F ∈ { E , E + 1 } .
于是 x x x y y y 2 E 2^{E} 2 E F ≥ E F\ge E F ≥ E 2 E ∣ 2 F ∣ x 2^E\mid2^F\mid x 2 E ∣ 2 F ∣ x d = x − y ≥ 0 d=x-y\ge0 d = x − y ≥ 0 2 E 2^E 2 E
0 ≤ d ≤ 2 y − y = y < 2 t 2 E ⟹ d = k 2 E , 0 ≤ k ≤ 2 t − 1. 0\ \le\ d\ \le\ 2y-y=y\ \lt\ 2^{t}\,2^{E}\ \Longrightarrow\ d=k\,2^{E},\qquad 0\le k\le2^{t}-1 . 0 ≤ d ≤ 2 y − y = y < 2 t 2 E ⟹ d = k 2 E , 0 ≤ k ≤ 2 t − 1.
k k k t t t d d d t t t k ≥ 2 t − 1 k\ge2^{t-1} k ≥ 2 t − 1 d d d 1 ≤ k < 2 t − 1 1\le k\lt2^{t-1} 1 ≤ k < 2 t − 1 k k k k < 2 t − 1 k\lt2^{t-1} k < 2 t − 1 k = 0 k=0 k = 0 d = 0 d=0 d = 0 2 e min − t + 1 2^{e_{\min}-t+1} 2 e m i n − t + 1 本来就在网格上 的数当然不变：fl ( x − y ) = x − y \operatorname{fl}(x-y)=x-y fl ( x − y ) = x − y
∎ 推论 B.3.1 · 相消的完整账目
相近数相减时，减法本身零误差（本引理），结果的误差全部来自两个操作数在减法之前携带的误差 ：设 x , y x,y x , y ≤ u \le u ≤ u x − y x-y x − y ( ∣ x ∣ + ∣ y ∣ ) u ∣ x − y ∣ ≈ 2 ∣ x ∣ ∣ x − y ∣ u \dfrac{(\lvert x\rvert+\lvert y\rvert)u}{\lvert x-y\rvert}\approx\dfrac{2\lvert x\rvert}{\lvert x-y\rvert}u ∣ x − y ∣ (∣ x ∣ + ∣ y ∣) u ≈ ∣ x − y ∣ 2 ∣ x ∣ u
记号。 γ k : = k u 1 − k u \gamma_k:=\dfrac{ku}{1-ku} γ k := 1 − k u k u k u < 1 ku\lt1 k u < 1 k u ku k u
引理 B.4.1 · 乘积引理
若 ∣ δ j ∣ ≤ u \lvert\delta_j\rvert\le u ∣ δ j ∣ ≤ u ∏ j = 1 k ( 1 + δ j ) = 1 + θ \prod_{j=1}^{k}(1+\delta_j)=1+\theta ∏ j = 1 k ( 1 + δ j ) = 1 + θ ∣ θ ∣ ≤ γ k \lvert\theta\rvert\le\gamma_k ∣ θ ∣ ≤ γ k
证明（归纳）
k = 1 k=1 k = 1 ∣ θ ∣ = ∣ δ 1 ∣ ≤ u ≤ γ 1 \lvert\theta\rvert=\lvert\delta_1\rvert\le u\le\gamma_1 ∣ θ ∣ = ∣ δ 1 ∣ ≤ u ≤ γ 1 k − 1 k-1 k − 1 ∏ j = 1 k ( 1 + δ j ) = ( 1 + θ ′ ) ( 1 + δ k ) \prod_{j=1}^{k}(1+\delta_j)=(1+\theta')(1+\delta_k) ∏ j = 1 k ( 1 + δ j ) = ( 1 + θ ′ ) ( 1 + δ k ) ∣ θ ′ ∣ ≤ γ k − 1 \lvert\theta'\rvert\le\gamma_{k-1} ∣ θ ′ ∣ ≤ γ k − 1 1 + θ = 1 + θ ′ + δ k + θ ′ δ k 1+\theta=1+\theta'+\delta_k+\theta'\delta_k 1 + θ = 1 + θ ′ + δ k + θ ′ δ k
∣ θ ∣ ≤ γ k − 1 + u + u γ k − 1 = γ k − 1 + u ( 1 + γ k − 1 ) . \lvert\theta\rvert\le\gamma_{k-1}+u+u\gamma_{k-1}=\gamma_{k-1}+u(1+\gamma_{k-1}). ∣ θ ∣ ≤ γ k − 1 + u + u γ k − 1 = γ k − 1 + u ( 1 + γ k − 1 ) .
而直接展开 γ k − γ k − 1 = u ( 1 − k u ) ( 1 − ( k − 1 ) u ) = u ( 1 + γ k ) ( 1 + γ k − 1 ) ≥ u ( 1 + γ k − 1 ) \gamma_k-\gamma_{k-1}=\dfrac{u}{(1-ku)(1-(k-1)u)}=u(1+\gamma_k)(1+\gamma_{k-1})\ge u(1+\gamma_{k-1}) γ k − γ k − 1 = ( 1 − k u ) ( 1 − ( k − 1 ) u ) u = u ( 1 + γ k ) ( 1 + γ k − 1 ) ≥ u ( 1 + γ k − 1 ) ∣ θ ∣ ≤ γ k \lvert\theta\rvert\le\gamma_k ∣ θ ∣ ≤ γ k
∎ 定理 B.4.2 · 递归求和的误差界
S ^ n = fl ( ⋯ fl ( fl ( x 1 + x 2 ) + x 3 ) ⋯ + x n ) \hat S_n=\operatorname{fl}(\cdots\operatorname{fl}(\operatorname{fl}(x_1+x_2)+x_3)\cdots+x_n) S ^ n = fl ( ⋯ fl ( fl ( x 1 + x 2 ) + x 3 ) ⋯ + x n )
S ^ n = ∑ i = 1 n x i ( 1 + θ i ) , ∣ θ i ∣ ≤ γ n − 1 ⟹ ∣ S ^ n − ∑ i x i ∣ ≤ γ n − 1 ∑ i ∣ x i ∣ . \hat S_n=\sum_{i=1}^{n}x_i(1+\theta_i),\qquad\lvert\theta_i\rvert\le\gamma_{n-1}\qquad\Longrightarrow\qquad \Bigl\lvert\hat S_n-\sum_ix_i\Bigr\rvert\le\gamma_{n-1}\sum_i\lvert x_i\rvert . S ^ n = i = 1 ∑ n x i ( 1 + θ i ) , ∣ θ i ∣ ≤ γ n − 1 ⟹ S ^ n − i ∑ x i ≤ γ n − 1 i ∑ ∣ x i ∣ .
证明
S ^ k = ( S ^ k − 1 + x k ) ( 1 + δ k ) \hat S_k=(\hat S_{k-1}+x_k)(1+\delta_k) S ^ k = ( S ^ k − 1 + x k ) ( 1 + δ k ) ∣ δ k ∣ ≤ u \lvert\delta_k\rvert\le u ∣ δ k ∣ ≤ u k k k S ^ k = ∑ i ≤ k x i ∏ j = max ( i , 2 ) k ( 1 + δ j ) \hat S_k=\sum_{i\le k}x_i\prod_{j=\max(i,2)}^{k}(1+\delta_j) S ^ k = ∑ i ≤ k x i ∏ j = m a x ( i , 2 ) k ( 1 + δ j ) x k x_k x k k k k x i x_i x i ( 1 + δ k ) (1+\delta_k) ( 1 + δ k ) x i x_i x i n − 1 n-1 n − 1 ( 1 + δ j ) (1+\delta_j) ( 1 + δ j ) ∣ θ i ∣ ≤ γ n − 1 \lvert\theta_i\rvert\le\gamma_{n-1} ∣ θ i ∣ ≤ γ n − 1
∎ 三点读法。
界由 ∑ ∣ x i ∣ \sum\lvert x_i\rvert ∑ ∣ x i ∣ ∣ ∑ x i ∣ \lvert\sum x_i\rvert ∣ ∑ x i ∣ ：同号数据两者接近，界好用；正负相消（和远小于各项之和）时界可以远大于结果本身——这是「先小后大」「成对求和」之外，改写公式 （如 Welford）的动机。成对求和 把求和树深度降为 ⌈ log 2 n ⌉ \lceil\log_2n\rceil ⌈ log 2 n ⌉ x i x_i x i ⌈ log 2 n ⌉ \lceil\log_2n\rceil ⌈ log 2 n ⌉ ∣ \lvert ∣ ∣ ≤ γ ⌈ log 2 n ∣ ∑ ∣ x i ∣ \rvert\le\gamma_{\lceil\log_2n\rvert}\sum\lvert x_i\rvert ∣ ≤ γ ⌈ l o g 2 n ∣ ∑ ∣ x i ∣ n = 1 0 8 n=10^8 n = 1 0 8 γ n ≈ 1 0 − 8 \gamma_n\approx10^{-8} γ n ≈ 1 0 − 8 γ 27 ≈ 3 × 1 0 − 15 \gamma_{27}\approx3\times10^{-15} γ 27 ≈ 3 × 1 0 − 15 np.sum 内置此策略，NumPy 手写循环之差由此而来。随机符号 下误差的典型 大小是 O ( u n ) O(u\sqrt n,) O ( u n ) γ n \gamma_n γ n
算法。
s = x[0]; c = 0
for i = 1..n-1:
y = x[i] - c # 把上一轮丢失的低位补回来
t = s + y # 高位累加（可能丢失低位）
c = (t - s) - y # 记录本轮真正丢失的部分（Sterbenz：这步减法精确）
s = t
返回 s
机制。 t = s + y t=s+y t = s + y y y y c = ( t − s ) − y c=(t-s)-y c = ( t − s ) − y 精确地 等于这个丢失量（注意 t − s t-s t − s y y y y ′ = x i + 1 − c y'=x_{i+1}-c y ′ = x i + 1 − c 误差不许累积，当场记账、下轮还清 。严格误差界为 ∣ \lvert ∣ ∣ ≤ 2 u ∑ ∣ x i ∣ + O ( n u 2 ) \rvert\le2u\sum\lvert x_i\rvert+O(nu^2) ∣ ≤ 2 u ∑ ∣ x i ∣ + O ( n u 2 ) c c c u 2 u^2 u 2
逐位实例（四位十进制，加十个 1.000 1.000 1.000 1.000 × 1 0 4 1.000\times10^{4} 1.000 × 1 0 4 真值 10010 10010 10010 [ 10000 , 10000 ) [10000,10000) [ 10000 , 10000 ) { 10000 , 10010 , 10020 , … } {10000,10010,10020,\dots} { 10000 , 10010 , 10020 , … } 间距是 10 10 10 。
裸累加 ：10000 + 1 10000+1 10000 + 1 10000 10000 10000 10000 10000 10000 − 10 -10 − 10 Kahan （逐步跟踪 s , c s,c s , c
步 y = x − c y=x-c y = x − c t = s + y t=s+y t = s + y c = ( t − s ) − y c=(t-s)-y c = ( t − s ) − y s s s
1 1 − 0 = 1 1-0=1 1 − 0 = 1 10001 → 10000 10001\to10000 10001 → 10000 0 − 1 = − 1 0-1=-1 0 − 1 = − 1 10000 10000 10000 2 1 − ( − 1 ) = 2 1-(-1)=2 1 − ( − 1 ) = 2 10002 → 10000 10002\to10000 10002 → 10000 − 2 -2 − 2 10000 10000 10000 5 1 − ( − 4 ) = 5 1-(-4)=5 1 − ( − 4 ) = 5 10005 → 10000 10005\to10000 10005 → 10000 − 5 -5 − 5 10000 10000 10000 6 1 − ( − 5 ) = 6 1-(-5)=6 1 − ( − 5 ) = 6 10006 → 10010 10006\to10010 10006 → 10010 10 − 6 = + 4 10-6=+4 10 − 6 = + 4 10010 10010 10010 7 1 − 4 = − 3 1-4=-3 1 − 4 = − 3 10007 → 10010 10007\to10010 10007 → 10010 0 − ( − 3 ) = 3 0-(-3)=3 0 − ( − 3 ) = 3 10010 10010 10010 9 1 − 2 = − 1 1-2=-1 1 − 2 = − 1 10009 → 10010 10009\to10010 10009 → 10010 1 1 1 10010 10010 10010 10 1 − 1 = 0 1-1=0 1 − 1 = 0 10010 10010 10010 0 0 0 10010 10010 10010
前五步零头被「记账欠着」，第六步攒够一格一次跨过（c c c + 4 +4 + 4 多进了 ，随后逐步还账），第十步账面清零——返回 10010 10010 10010 。注意第 6 步的「先过头、再还账」：补偿求和允许中途反向修正，这正是它比「先小后大」这类纯排序技巧更强的地方（后者只是降低概率，前者给出保证）。
变体与现状。 Neumaier 变体（用 ∣ s ∣ \lvert s\rvert ∣ s ∣ ∣ y ∣ \lvert y\rvert ∣ y ∣ sum 内建 Neumaier 补偿、math.fsum 给出精确舍入的和——正文 1.2.5 思考题的实验差异全部由此解释。
（一）0.1 0.1 0.1 ( 0. 0011 ‾ ) 2 = 3 2 4 − 1 = 1 5 (0.\overline{0011})_2=\dfrac{3}{2^4-1}=\dfrac15 ( 0. 0011 ) 2 = 2 4 − 1 3 = 5 1
0.1 = ( 0.0 0011 ‾ ) 2 = 0.000110011001100110011 ⋯ 0.1=(0.0\overline{0011})_2=0.000110011001100110011\cdots 0.1 = ( 0.0 0011 ) 2 = 0.000110011001100110011 ⋯
——正如 1 / 3 1/3 1/3 1 / 10 1/10 1/10 存进去的那一刻就已失真 。
（二）双精度里到底存了什么。 尾数 53 53 53 0.1 × 2 55 = 3602879701896396.8 0.1\times2^{55}=3602879701896396.8 0.1 × 2 55 = 3602879701896396.8
fl ( 0.1 ) = 3602879701896397 × 2 − 55 = 0.1000000000000000055511151231257827 ⋯ \operatorname{fl}(0.1)=3602879701896397\times2^{-55}=0.1000000000000000055511151231257827\cdots fl ( 0.1 ) = 3602879701896397 × 2 − 55 = 0.1000000000000000055511151231257827 ⋯
0.2 0.2 0.2 2 2 2 fl ( 0.2 ) = 7205759403792794 × 2 − 55 \operatorname{fl}(0.2)=7205759403792794\times2^{-55} fl ( 0.2 ) = 7205759403792794 × 2 − 55
（三）相加恰好处在正中央。 精确和为
3602879701896397 × 2 − 55 + 7205759403792794 × 2 − 55 = 10808639105689191 × 2 − 55 = 5404319552844595.5 × 2 − 54 , 3602879701896397\times2^{-55}+7205759403792794\times2^{-55}=10808639105689191\times2^{-55}=\mathbf{5404319552844595.5}\times2^{-54}, 3602879701896397 × 2 − 55 + 7205759403792794 × 2 − 55 = 10808639105689191 × 2 − 55 = 5404319552844595.5 × 2 − 54 ,
而 [ 0.25 , 0.5 ) [0.25,0.5) [ 0.25 , 0.5 ) 2 − 54 2^{-54} 2 − 54 精确和恰好在两个相邻浮点数的正中央 ：
5404319552844595 × 2 − 54 ⏟ = 0.2999999999999999888 ⋯ < 0.3000000000000000166533 ⋯ < 5404319552844596 × 2 − 54 ⏟ = 0.3000000000000000444 ⋯ \underbrace{5404319552844595\times2^{-54}}_{=0.2999999999999999888\cdots}\quad\lt\ 0.3000000000000000166533\cdots\ \lt\quad\underbrace{5404319552844596\times2^{-54}}_{=0.3000000000000000444\cdots} = 0.2999999999999999888 ⋯ 5404319552844595 × 2 − 54 < 0.3000000000000000166533 ⋯ < = 0.3000000000000000444 ⋯ 5404319552844596 × 2 − 54
到两端的距离相等（都是 2 − 55 2^{-55} 2 − 55 平局，触发 round-half-to-even ：两个候选的尾数整数分别为 5404319552844595 5404319552844595 5404319552844595 5404319552844596 5404319552844596 5404319552844596
0.1 + 0.2 = 0.30000000000000004440892098500626 ⋯ ⟹ 0.1 + 0.2 == 0.3 为 False . 0.1+0.2=0.30000000000000004440892098500626\cdots\ \Longrightarrow\ \texttt{0.1 + 0.2 == 0.3}\ \text{为 False}. 0.1 + 0.2 = 0.30000000000000004440892098500626 ⋯ ⟹ 0.1 + 0.2 == 0.3 为 False .
顺带：fl ( 0.3 ) = 5404319552844595 × 2 − 54 = 0.29999999999999998889 ⋯ \operatorname{fl}(0.3)=5404319552844595\times2^{-54}=0.29999999999999998889\cdots fl ( 0.3 ) = 5404319552844595 × 2 − 54 = 0.29999999999999998889 ⋯ 0.3，存的也不是「数学的 0.3 0.3 0.3 1.2.5 三条后果之第三条的逐位版）。
1.2.7 实验二 (b) 的定量解释）（一）失效何时发生。 两遍法 s 2 = 1 n ∑ x i 2 − x ˉ 2 s^2=\dfrac1n\sum x_i^2-\bar x^2 s 2 = n 1 ∑ x i 2 − x ˉ 2 1.2.4 原则 (3) 的完美标本）：∑ x i 2 ≈ n ( x ˉ 2 + s 2 ) \sum x_i^2\approx n(\bar x^2+s^2) ∑ x i 2 ≈ n ( x ˉ 2 + s 2 ) u ⋅ n x ˉ 2 u\cdot n\bar x^2 u ⋅ n x ˉ 2 s 2 s^2 s 2
∣ 误差 ∣ s 2 ≈ u n x ˉ 2 n s 2 = u x ˉ 2 s 2 ⟹ 失效阈值： ∣ x ˉ ∣ s ≳ 1 u ≈ 1 0 8 . \frac{\lvert\text{误差}\rvert}{s^2}\ \approx\ \frac{u,n,\bar x^{2}}{n,s^{2}}=u,\frac{\bar x^{2}}{s^{2}}\qquad\Longrightarrow\qquad \textbf{失效阈值：}\ \frac{\lvert\bar x\rvert}{s}\ \gtrsim\ \frac{1}{\sqrt{u}}\approx10^{8}. s 2 ∣ 误差 ∣ ≈ n s 2 u n x ˉ 2 = u s 2 x ˉ 2 ⟹ 失效阈值： s ∣ x ˉ ∣ ≳ u 1 ≈ 1 0 8 .
正文实验的数据 1 0 9 + 1 , … , 1 0 9 + 8 10^9+1,\dots,10^9+8 1 0 9 + 1 , … , 1 0 9 + 8 x ˉ ≈ 1 0 9 \bar x\approx10^9 x ˉ ≈ 1 0 9 s ≈ 5.25 ≈ 2.29 s\approx\sqrt{5.25}\approx2.29 s ≈ 5.25 ≈ 2.29 ≈ 4.4 × 1 0 8 ≳ 1 0 8 \approx4.4\times10^{8}\gtrsim10^{8} ≈ 4.4 × 1 0 8 ≳ 1 0 8 0.0 0.0 0.0 8 × 1 0 18 8\times10^{18} 8 × 1 0 18 ∼ 1 0 2 \sim10^{2} ∼ 1 0 2 5.25 × 8 = 42 5.25\times8=42 5.25 × 8 = 42 后者是良性的 （相对误差不变），前者是致命的 （绝对误差不变而真实信号被平移压小）。
（二）Welford 单遍法的递推及其推导。 维护 ( μ k , M k ) (\mu_k,\ M_k) ( μ k , M k ) k k k
μ k = μ k − 1 + x k − μ k − 1 k , M k = M k − 1 + ( x k − μ k − 1 ) ( x k − μ k ) , s k 2 = M k k − 1 . \mu_k=\mu_{k-1}+\frac{x_k-\mu_{k-1}}{k},\qquad M_k=M_{k-1}+(x_k-\mu_{k-1})(x_k-\mu_k),\qquad s_k^2=\frac{M_k}{k-1}. μ k = μ k − 1 + k x k − μ k − 1 , M k = M k − 1 + ( x k − μ k − 1 ) ( x k − μ k ) , s k 2 = k − 1 M k .
推导。 离差的平移不变性（μ k = μ k − 1 + ( x k − μ k − 1 ) / k \mu_k=\mu_{k-1}+(x_k-\mu_{k-1})/k μ k = μ k − 1 + ( x k − μ k − 1 ) / k x i − μ k = ( x i − μ k − 1 ) − x k − μ k − 1 k x_i-\mu_k=(x_i-\mu_{k-1})-\frac{x_k-\mu_{k-1}}{k} x i − μ k = ( x i − μ k − 1 ) − k x k − μ k − 1
M k = ∑ i ≤ k ( x i − μ k ) 2 = ∑ i ≤ k ( x i − μ k − 1 ) 2 − k ( μ k − μ k − 1 ) 2 = M k − 1 + ( x k − μ k − 1 ) 2 − ( x k − μ k − 1 ) 2 k , M_k=\sum_{i\le k}(x_i-\mu_k)^2=\sum_{i\le k}(x_i-\mu_{k-1})^2-k(\mu_k-\mu_{k-1})^2=M_{k-1}+(x_k-\mu_{k-1})^2-\frac{(x_k-\mu_{k-1})^2}{k}, M k = i ≤ k ∑ ( x i − μ k ) 2 = i ≤ k ∑ ( x i − μ k − 1 ) 2 − k ( μ k − μ k − 1 ) 2 = M k − 1 + ( x k − μ k − 1 ) 2 − k ( x k − μ k − 1 ) 2 ,
最后一步对 i = k i=k i = k M k − 1 M_{k-1} M k − 1 k − 1 k-1 k − 1 x k − μ k = ( x k − μ k − 1 ) k − 1 k x_k-\mu_k=(x_k-\mu_{k-1})\frac{k-1}{k} x k − μ k = ( x k − μ k − 1 ) k k − 1 M k = M k − 1 + ( x k − μ k − 1 ) ( x k − μ k ) M_k=M_{k-1}+(x_k-\mu_{k-1})(x_k-\mu_k) M k = M k − 1 + ( x k − μ k − 1 ) ( x k − μ k )
为什么稳。 每步相减的两个量 ( x k − μ k − 1 ) (x_k-\mu_{k-1}) ( x k − μ k − 1 ) 大偏移 x ˉ \bar x x ˉ ，1 0 9 10^9 1 0 9 k k k O ( 1 ) O(1) O ( 1 )
本附录与正文的对应
B.1–B.2 对应 1.2.5 浮点数集与机器精度 ；B.3 对应 引理 1.2.1 （正文为证明思路）；B.4–B.5 对应 1.2.5 思考题 （sum vs fsum）；B.6 对应 1.2.5 三条后果之第三条；B.7 对应 实验二 (b) 的「报告要求」。系统论述见 Goldberg《What Every Computer Scientist Should Know about Floating-Point Arithmetic》与 Higham《Accuracy and Stability of Numerical Algorithms》Chapter 4。