为什么在逼近导数时不同的步长会产生很大的数值误差？答案

【问题标题】：Why do different step sizes produce large numerical error while approximating derivatives?为什么在逼近导数时不同的步长会产生很大的数值误差？
【发布时间】：2021-04-15 07:36:31
【问题描述】：

我正在尝试一些代码来帮助我近似任意函数的导数。我看到另一个post上有四个选项：

有限的差异
自动导数
符号微分
手动计算导数

我发现我的方法最适合第一个选项，它有注释“容易出现数字错误”。所以我知道这种方法并不准确，这很好。

话虽如此，我对不同数据类型可以存储的大小数字进行了一些研究，并在post 中发现它可以非常小，（大约 10^{–308 sup>) 并且“在正常范围内，基本运算的结果将在格式的正常精度范围内准确”。}

话虽如此，对于以下代码 sn-p，我似乎得到了非常糟糕的结果，我在其中探索了不同大小的间隔；最小差异不应小于10^–27（10^–9，立方），远大于极限值。我会很感激也许更具体的回应？

epsilon = 0.01 # is "small" w.r.t. to 3

def approx_derivative(func): # rough derivative factory function
  return lambda x : (func(x + epsilon) - func(x)) / epsilon

while epsilon > 10**-9:
  nth_deriv = lambda x : x ** 3 # 0th derivative
  for i in range(5): # should read about 27, 27, 18, 6, 0
    print(nth_deriv(3), end=', ')
    nth_deriv = approx_derivative(nth_deriv) # take derivative
  print('\n')
  epsilon *= 0.1

输出是：

27, 27.090099999999495, 18.0599999999842, 6.000000002615025, -3.552713678800501e-07, 

27, 27.009000999996147, 18.00600000123609, 6.000000496442226, -0.007105427357601002, 

27, 27.00090001006572, 18.000599766310188, 6.004086117172847, -71.05427357601002, 

27, 27.000090000228735, 18.000072543600254, 3.5527136788005005, 355271.36788005003, 

27, 27.000009005462285, 17.998047496803334, 0.0, 3552713678.8005, 

27, 27.000000848431675, 18.11883976188255, 0.0, -35527136788004.99, 

27, 27.0000001023618, 0.0, 0.0, 3.552713678800497e+17, 

27, 27.000002233990003, 0.0, 0.0, 0.0,

正如我们在前几个示例中看到的那样，结果并不准确，但非常好。但是，对于某些间隔大小，某些值被放大了；其他人为0；有些是完全错误的，比如给出一半的值，尽管直觉上它们应该对较小的 epsilon 变得更准确。我可以将此错误归因于哪些主要因素？我应该注意/注意什么？是否有我应该担心用块捕获的错误（例如除以 0）？

对于使用浮点数进行计算，通常认为“最佳”的 epsilon 值是否存在？或者是否有根据您的输入选择合适大小的 epsilon 的“经验法则”？与我实现的衍生品相比，是否有首选的衍生品定义？

【问题讨论】：

标签： python precision numeric ieee-754 derivative

【解决方案1】：

首先，浮点类型可表示的最小值在此无关紧要。精度是一个问题。 Python 规范没有具体说明 Python 实现使用什么浮点格式，但许多使用 IEEE-754 binary64 格式，其有效位有 53 位，这意味着值的最小可表示变化约为 2⁻⁵² 乘以所表示的数字的大小。我们将研究这对您的计算有何影响。

计算导数 f' 的近似值，计算 (x+e)³−x³。结果是 3ex² + 3e²x + e³，大约是 3ex²，然后除以 e 得到大约 3 x²。请注意，涉及到大小约为 x³ 的两个项，它们相差大约 3ex²。由于您在 x = 3 处取导数，因此 3ex² 恰好等于 ex³，所以它比 x³ 低 e 倍。只要 e 是适度的，x³ 附近的两个二进制 64 数字很容易相差 ex³ 不失准确性。

计算二阶导数 f'' 的近似值，计算 ((x+2e)³-( x+e)³)-((x+e)^{3-x³) = (x+2e)³ −2(x+e)³+x³ = 6e²x + 6e³。一路上，这被 e² 划分，但我们现在可以忽略它进行错误分析。原始项的大小约为 x³，但我们减法的结果约为 6e²x。所以这取决于数字能够相差大约 e²。由于 2^-52 大约是 2•10^-16，所以当 e 大约是 10^-8 时，这会大大失败，在此之前错误较少。当 e 在这个大小附近时，binary64 格式根本无法准确地记录数字之间的差异，以在以后减去它们时保持准确性。}

对于三阶导数，您需要表示 e³ 周围相对于主要幅度的差异，因此它在 e = 周围相当失败10⁻⁵，此后e³接近2•10⁻¹⁶的极限。 p>

【讨论】：