使用 numpy.log 时如何理解不准确的结果？

Question

我想计算一个非常接近 1 但不完全是 1 的值的自然对数。

例如，np.log(1 + 1e-22) 是 0 而不是某个非零值。但是，np.log(1 + 1e-13)不为零，计算为9.992007221625909e-14.

在使用 numpy 函数与将 dtype 定义为 np.longdouble 的 numpy 数组时，我如何理解这种精度权衡？

我使用的系统上numpy(v1.22.2)的浮点精度信息：

>>> np.finfo(np.longdouble)
finfo(resolution=1e-15, min=-1.7976931348623157e+308, max=1.7976931348623157e+308, dtype=float64)

>>> 1 + np.finfo(np.longdouble).eps
1.0000000000000002

Answer 1

对于实际使用，请查看 np.log1p(x)，它计算 log(1 + x) 时没有来自 1 + x 的舍入误差。来自文档：

For real-valued input, log1p is accurate also for x so small that 1 + x == 1 in floating-point accuracy.

即使 log(1 + 1e-13) 的看似有效的示例也不同于 64 位浮点数小数点后第三位的真实值，以及 128 位浮点数的第六位（真实值来自 WolframAlpha):

>>> (1 + 1e-13) - 1
9.992007221626409e-14
>>> np.log(1 + 1e-13)
9.992007221625909e-14
>>> np.log(1 + np.array(1e-13, dtype=np.float128))
9.999997791637127032e-14
>>> np.log1p(1e-13)
9.9999999999995e-14
>>> 9.999999999999500000000000033333333333330833333333333533333*10**-14
9.9999999999995e-14

Answer 2

使用Numpy完成@yut23的好解。如果您需要处理不适合 Numpy 定义的本机类型的非常小的浮点数或接近 1 且精度超过 ~10 位的数字，那么您可以使用 decimal 包 .它 比原生 float 慢，但它可以为您提供 任意精度 。问题是它不直接支持自然对数（即 log）函数，因为它基于超越欧拉数（即 e），很难以任意精度计算（在至少在精度很高的时候不会）。幸运的是，您可以根据基于 10 的对数和基于现有数字数据库（如 this one）预先计算的欧拉数计算自然对数（我想 10000 位数字应该足以满足您的需要；））。这是一个例子：

import decimal
from decimal import Decimal

decimal.setcontext(decimal.Context(prec=100)) # 100 digits of precision
e = Decimal('2.71828182845904523536028747135266249775724709369995957496696762772407663035354759457138217852516643')
result = (Decimal("1")+Decimal("1e-13")).log10() / e.log10()

# result = 9.999999999999500000000000033333333333330833333333333533333333333316666666666668095238095237970238087E-14

result精度为99位（只有最后一位不正确）

使用 numpy.log 时如何理解不准确的结果？

How to understand inaccurate results while using numpy.log?

python

numpy