有限差分中步长h的选取方法及可重复应用的稳定数值微分算法咨询

一、针对指定精度N的步长h确定 #

数值差分的误差来自两部分：截断误差（由差分公式的近似性导致）和舍入误差（由浮点数计算的精度限制导致），两者随h的变化趋势相反：

• h越大，截断误差越高（前向差分(f(x+h)-f(x))/h的截断误差为O(h)）；
• h越小，f(x+h)与f(x)的差值越接近，浮点数相减会损失有效数字，舍入误差急剧上升。

1. 理论估算最优步长 #

假设机器浮点精度为ε（双精度环境下约为1e-16），函数f在x处的二阶导数绝对值上限为M₂，前向差分的总误差可近似为：
总误差 ≈ (M₂·h)/2 + (ε·|f(x)|)/h
对h求导并令导数为0，可得到最优步长的估算值：
h_opt ≈ sqrt(2·ε·|f(x)| / M₂)
如果无法提前获取M₂，可以先通过少量试算估计f的二阶导数范围，再代入计算。

2. 实验验证法 #

直接测试不同h对应的差分结果：

• 从较大的h开始逐步减小，观察结果的变化：当结果开始稳定在某个区间（波动小于指定精度要求的10^-N），此时的h处于截断误差主导的合理区间；
• 继续减小h，当结果开始出现无规律的波动（舍入误差主导），停止测试，取稳定区间内的中间值作为最终h。

3. 收敛性与稳定性验证 #

• 收敛性验证：用h和h/2分别计算差分结果，若为前向差分（一阶收敛），两次结果的误差应大致减半；若使用中心差分(f(x+h)-f(x-h))/(2h)（二阶收敛），误差应大致减至1/4。若符合该规律，说明方法收敛。
• 稳定性验证：在合理h范围内取多个值计算，若结果的波动幅度小于指定精度N的要求，说明步长选择稳定；若波动过大，说明h过小（舍入误差主导）或过大（截断误差超标）。

二、可重复1000次的稳定数值微分方法 #

普通有限差分重复多次微分后，误差会快速累积放大，无法满足1000次重复的需求，推荐以下几种稳定方法：

1. Chebyshev微分矩阵法 #

基于Chebyshev正交多项式构造微分矩阵，将微分操作转化为矩阵乘法。多次微分只需对矩阵进行幂运算，Chebyshev多项式的正交性能有效抑制误差累积，即使重复数百次甚至上千次微分，仍能保持较高的稳定性。

2. 样条插值微分法 #

先对目标函数进行高精度的样条拟合（比如B样条或自然样条），样条函数本身是分段多项式，其导数可以通过解析方式精确计算。重复求导时，只需对样条的导数继续求导，整个过程仅涉及多项式的精确运算，误差仅来自初始拟合阶段，能稳定支持上千次重复微分。

3. 自动微分（AD） #

通过链式法则将导数计算拆解为一系列基本运算的导数组合，没有有限差分的截断误差，仅存在浮点数舍入误差。重复多次微分时，只需逐步展开链式法则，精度与原函数计算一致，稳定性远优于传统有限差分方法。

内容的提问来源于stack exchange，提问作者ShoutOutAndCalculate