2. 数值微分基础

数值微分要解决的问题很朴素：给定一个只能通过函数值访问的 f(x)，如何估计它在某点的导数。解析求导需要函数表达式，而工程中我们常常只有采样数据、黑箱程序或代价高昂的仿真器。此时导数只能由若干函数值的组合来近似，这类近似统称为差商。本节从差商定义出发，推导前向、后向与中心差分公式，并重点讨论一个贯穿全课程的核心矛盾：步长 h 越小，截断误差越小，但舍入误差越大。

差商的定义与三种基本形式

设函数 f 在点 x 附近充分光滑，取步长 h > 0。导数的定义是差商的极限：

f'(x) = lim(h→0) [f(x+h) − f(x)] / h

把极限去掉，就得到三种最常用的近似。前向差分取 x 与 x+h 两点：

D₊f(x) = [f(x+h) − f(x)] / h

后向差分取 x−h 与 x 两点：

D₋f(x) = [f(x) − f(x−h)] / h

中心差分取 x−h 与 x+h 两点，把 x 放在中间：

D₀f(x) = [f(x+h) − f(x−h)] / (2h)

三者都写成“函数值之差除以自变量之差”的形式，区别只在于采样点的位置与权重。前向与后向是单侧差分，中心差分是双侧差分。这个结构差异直接决定了它们的精度阶数。

用泰勒展开推导截断误差

截断误差来自“用差商代替导数”这一数学近似本身，与计算机浮点无关。把 f(x±h) 在 x 处展开：

f(x+h) = f(x) + h f'(x) + h²/2 f''(x) + h³/6 f'''(x) + O(h⁴)
f(x−h) = f(x) − h f'(x) + h²/2 f''(x) − h³/6 f'''(x) + O(h⁴)

对前向差分，代入并整理得：

D₊f(x) = f'(x) + h/2 f''(x) + O(h²)

所以前向差分的截断误差主项是 (h/2) f''(x)，为 一阶精度 O(h)。后向差分同理：

D₋f(x) = f'(x) − h/2 f''(x) + O(h²)

也是一阶精度，只是误差主项符号相反。这说明前向与后向的误差方向相反，二者平均可以消去一阶项——这正是中心差分的来源。把两式相减：

f(x+h) − f(x−h) = 2h f'(x) + h³/3 f'''(x) + O(h⁵)
D₀f(x) = f'(x) + h²/6 f'''(x) + O(h⁴)

中心差分的截断误差主项是 (h²/6) f'''(x)，为 二阶精度 O(h²)。直观理解：中心差分对 x 两侧对称采样，一阶误差项因对称而抵消，只剩下三阶导数项。这就是为什么在同样步长下，中心差分通常比单侧差分精确得多。

舍入误差：小步长的另一面

截断误差随 h 减小而减小，看起来 h 越小越好。但计算机用有限位浮点数存储函数值，f(x+h) 与 f(x) 各自带有约 ε·|f| 量级的舍入误差，其中 ε 是机器精度（双精度约 2.2×10⁻¹⁶）。差商中两个相近的数相减，会产生灾难性抵消：有效数字大量丢失，误差被 1/h 放大。

以中心差分为例，设函数值误差上界为 ε_f，则舍入误差约为 ε_f / (2h)。总误差是两者之和：

E(h) ≈ (h²/6)|f'''(x)| + ε_f / (2h)

第一项随 h 减小而下降，第二项随 h 减小而上升，总误差存在一个最优步长。对 h 求导令其为零，可得最优步长量级：

h* ≈ (3ε_f / |f'''(x)|)^(1/3)

对前向差分，总误差为 (h/2)|f''| + 2ε_f/h，最优步长量级为 h* ≈ 2√(ε_f/|f''|)。注意指数不同：中心差分最优步长约按 ε^(1/3) 缩放，前向差分约按 ε^(1/2) 缩放。以双精度 ε≈10⁻¹⁶ 估算，中心差分的最优步长在 10⁻⁵ 到 10⁻⁶ 量级，前向差分则在 10⁻⁸ 量级，且中心差分能达到的最小误差更小。这解释了为什么实际代码中不应盲目把 h 设成 1e-15。

一个可复现的示例

取 f(x)=sin(x)，在 x=1 处求导，解析值为 cos(1)≈0.540302306。用双精度分别计算中心差分在不同步长下的结果（以下为按公式推算的示例值，用于说明趋势）：

可以看到误差先随 h 减小而下降，在 h≈1e-5 附近达到谷底，之后因舍入误差回升。这条“U 形”误差曲线是数值微分的典型特征，也是选择步长的直接依据。

工程实现要点

实际编码时，优先使用中心差分，并把步长取在 sqrt(ε) 到 cbrt(ε) 之间再结合问题尺度缩放。若函数值量级为 |f|，可用 h = max(|x|,1)·ε^(1/3) 作为起点，再通过试算几个步长观察误差是否稳定来微调。对于边界点无法取双侧样本时，退化为单侧差分，此时应意识到精度降为一阶。若需要更高精度，可进一步使用五点差分公式（四阶精度），或采用 Richardson 外推把低阶误差项消去。此外，若函数值本身来自带噪测量，噪声会被 1/h 放大，此时应先用平滑或正则化手段抑制噪声，再讨论差分格式——这正是后续噪声抑制内容的出发点。