15. 浮点数运算:问题和局限性

浮点数在计算机硬件中以二进制(base 2)分数表示。例如,十进制分数 0.625 的值为 6/10 + 2/100 + 5/1000,同样,二进制分数 0.101 的值为 1/2 + 0/4 + 1/8。这两个分数的值相同,唯一的区别是第一个用十进制分数表示法写成,第二个用二进制分数表示法写成。

不幸的是,大多数十进制分数不能精确地表示为二进制分数。结果是,通常情况下,您输入的十进制浮点数只是机器实际存储的二进制浮点数的近似值。

这个问题在十进制中更容易理解。考虑分数 1/3。您可以将其近似为十进制分数

0.3

或者,更好的是

0.33

或者,更好的是

0.333

等等。无论您愿意写下多少位数字,结果永远不会完全等于 1/3,而只会越来越接近 1/3。

同样,无论您愿意使用多少个二进制位,十进制值 0.1 都不能精确地表示为二进制分数。在二进制中,1/10 是一个无限循环的小数

0.0001100110011001100110011001100110011001100110011...

在任何有限的位数处停止,您都会得到一个近似值。在当今大多数机器上,浮点数使用二进制分数近似,分子使用从最高有效位开始的 53 位,分母为 2 的幂。在 1/10 的情况下,二进制分数为 3602879701896397 / 2 ** 55,它接近但并不完全等于 1/10 的真实值。

许多用户没有意识到这种近似,因为他们没有注意到值的显示方式。Python 只打印机器存储的二进制近似值的真实十进制值的十进制近似值。在大多数机器上,如果 Python 要打印存储的 0.1 的二进制近似值的真实十进制值,它将不得不显示

>>> 0.1
0.1000000000000000055511151231257827021181583404541015625

这比大多数人认为有用的位数要多,因此 Python 通过显示一个四舍五入的值来保持位数的可管理性

>>> 1 / 10
0.1

请记住,即使打印的结果看起来像 1/10 的精确值,实际存储的值是最近的可表示的二进制分数。

有趣的是,许多不同的十进制数字共享相同的最近近似二进制分数。例如,数字 0.10.100000000000000010.1000000000000000055511151231257827021181583404541015625 都由 3602879701896397 / 2 ** 55 近似。由于所有这些十进制值共享相同的近似值,因此它们中的任何一个都可以显示,同时仍然保留不变式 eval(repr(x)) == x

从历史上看,Python 提示符和内置的 repr() 函数将选择具有 17 位有效数字的数字,0.10000000000000001。从 Python 3.1 开始,Python(在大多数系统上)现在能够选择这些数字中最短的数字,并简单地显示 0.1

请注意,这是二进制浮点数的本质:这不是 Python 中的错误,也不是您代码中的错误。您将在所有支持硬件浮点数运算的语言中看到相同的情况(尽管某些语言可能不会显示默认情况下或在所有输出模式下的差异)。

为了获得更令人愉快的输出,您可能希望使用字符串格式化来生成有限数量的有效数字

>>> format(math.pi, '.12g')  # give 12 significant digits
'3.14159265359'

>>> format(math.pi, '.2f')   # give 2 digits after the point
'3.14'

>>> repr(math.pi)
'3.141592653589793'

重要的是要意识到,这在某种意义上是一种错觉:您只是对真实机器值的显示进行了四舍五入。

一个幻觉可能会导致另一个幻觉。例如,由于 0.1 不完全等于 1/10,因此将三个 0.1 值相加可能也不会完全等于 0.3。

>>> 0.1 + 0.1 + 0.1 == 0.3
False

此外,由于 0.1 无法更接近 1/10 的精确值,而 0.3 无法更接近 3/10 的精确值,因此使用 round() 函数进行预先舍入也无济于事。

>>> round(0.1, 1) + round(0.1, 1) + round(0.1, 1) == round(0.3, 1)
False

虽然这些数字无法更接近其预期的精确值,但 math.isclose() 函数可以用于比较不精确的值。

>>> math.isclose(0.1 + 0.1 + 0.1, 0.3)
True

或者,可以使用 round() 函数来比较粗略的近似值。

>>> round(math.pi, ndigits=2) == round(22 / 7, ndigits=2)
True

二进制浮点运算会带来很多这样的意外。 “0.1” 的问题将在下面的“表示误差”部分中详细解释。请参阅 浮点问题示例,以了解二进制浮点运算的工作原理以及实践中常见的各种问题的简要概述。还可以参阅 浮点数的风险,以更全面地了解其他常见的意外情况。

正如文章结尾所说,“没有简单的答案”。不过,不要过分担心浮点数!Python 浮点运算中的误差继承自浮点硬件,并且在大多数机器上,每次运算的误差不会超过 2**53 分之一。这对于大多数任务来说已经足够了,但您需要记住,它不是十进制运算,并且每次浮点运算都可能产生新的舍入误差。

虽然确实存在病态情况,但对于大多数浮点运算的非正式使用,如果您只是将最终结果的显示舍入到您期望的小数位数,那么您最终会看到您期望的结果。 str() 通常就足够了,而对于更精细的控制,请参阅 str.format() 方法在 格式字符串语法 中的格式说明符。

对于需要精确十进制表示的用例,请尝试使用 decimal 模块,该模块实现了适用于会计应用程序和高精度应用程序的十进制运算。

fractions 模块支持另一种形式的精确运算,该模块实现了基于有理数的运算(因此像 1/3 这样的数字可以精确表示)。

如果您是浮点运算的重度用户,您应该查看 NumPy 包以及 SciPy 项目提供的用于数学和统计运算的许多其他包。请参阅 <https://scipy.org.cn>。

在您确实需要知道浮点数的精确值时,Python 提供了一些工具来帮助您应对这些罕见的情况。 float.as_integer_ratio() 方法将浮点数的值表示为分数。

>>> x = 3.14159
>>> x.as_integer_ratio()
(3537115888337719, 1125899906842624)

由于该比率是精确的,因此它可以用来无损地重新创建原始值。

>>> x == 3537115888337719 / 1125899906842624
True

float.hex() 方法以十六进制(16 进制)表示浮点数,同样给出计算机存储的精确值。

>>> x.hex()
'0x1.921f9f01b866ep+1'

此精确的十六进制表示可用于精确地重建浮点值。

>>> x == float.fromhex('0x1.921f9f01b866ep+1')
True

由于表示是精确的,因此它对于跨不同版本的 Python(平台独立性)可靠地移植值以及与支持相同格式的其他语言(如 Java 和 C99)交换数据非常有用。

另一个有用的工具是 sum() 函数,它有助于减轻求和过程中的精度损失。它在将值添加到运行总计的过程中使用扩展精度进行中间舍入步骤。这可能会对整体精度产生影响,从而使误差不会累积到影响最终总计的程度。

>>> 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 == 1.0
False
>>> sum([0.1] * 10) == 1.0
True

math.fsum() 更进一步,它在将值添加到运行总计的过程中跟踪所有“丢失的数字”,以便结果只有一个舍入。这比 sum() 慢,但在极少数情况下,当大量级输入大部分相互抵消,留下接近零的最终总计时,它会更准确。

>>> arr = [-0.10430216751806065, -266310978.67179024, 143401161448607.16,
...        -143401161400469.7, 266262841.31058735, -0.003244936839808227]
>>> float(sum(map(Fraction, arr)))   # Exact summation with single rounding
8.042173697819788e-13
>>> math.fsum(arr)                   # Single rounding
8.042173697819788e-13
>>> sum(arr)                         # Multiple roundings in extended precision
8.042178034628478e-13
>>> total = 0.0
>>> for x in arr:
...     total += x                   # Multiple roundings in standard precision
...
>>> total                            # Straight addition has no correct digits!
-0.0051575902860057365

15.1. 表示误差

本节将详细解释“0.1”示例,并展示如何自行对类似情况进行精确分析。假设您已具备二进制浮点数表示的基本知识。

表示误差指的是某些(实际上是大多数)十进制小数无法精确地表示为二进制(以 2 为基数)小数。这是 Python(或 Perl、C、C++、Java、Fortran 以及许多其他语言)经常无法显示您期望的精确十进制数字的主要原因。

为什么会出现这种情况?1/10 无法精确地表示为二进制小数。自 2000 年以来,几乎所有机器都使用 IEEE 754 二进制浮点运算,并且几乎所有平台都将 Python 浮点数映射到 IEEE 754 binary64 “双精度”值。IEEE 754 binary64 值包含 53 位精度,因此在输入时,计算机会努力将 0.1 转换为它所能得到的J/2**N形式的最接近的分数,其中J 是一个包含正好 53 位的整数。重写

1 / 10 ~= J / (2**N)

J ~= 2**N / 10

并回想J 恰好包含 53 位(>= 2**52< 2**53),N 的最佳值为 56

>>> 2**52 <=  2**56 // 10  < 2**53
True

也就是说,56 是N 的唯一值,它使J 恰好包含 53 位。然后,J 的最佳可能值为该商四舍五入后的值

>>> q, r = divmod(2**56, 10)
>>> r
6

由于余数大于 10 的一半,因此通过向上取整可以获得最佳近似值

>>> q+1
7205759403792794

因此,IEEE 754 双精度中 1/10 的最佳可能近似值为

7205759403792794 / 2 ** 56

将分子和分母都除以 2 将分数简化为

3602879701896397 / 2 ** 55

请注意,由于我们向上取整,这实际上略大于 1/10;如果我们没有向上取整,商将略小于 1/10。但在任何情况下,它都不能完全等于 1/10!

因此,计算机永远不会“看到” 1/10:它看到的是上面给出的精确分数,这是它所能得到的最佳 IEEE 754 双精度近似值

>>> 0.1 * 2 ** 55
3602879701896397.0

如果我们将该分数乘以 10**55,我们可以看到小数点后 55 位的值

>>> 3602879701896397 * 10 ** 55 // 2 ** 55
1000000000000000055511151231257827021181583404541015625

这意味着存储在计算机中的精确数字等于十进制值 0.1000000000000000055511151231257827021181583404541015625。许多语言(包括旧版本的 Python)不会显示完整的十进制值,而是将结果四舍五入到 17 位有效数字

>>> format(0.1, '.17f')
'0.10000000000000001'

fractionsdecimal 模块使这些计算变得容易

>>> from decimal import Decimal
>>> from fractions import Fraction

>>> Fraction.from_float(0.1)
Fraction(3602879701896397, 36028797018963968)

>>> (0.1).as_integer_ratio()
(3602879701896397, 36028797018963968)

>>> Decimal.from_float(0.1)
Decimal('0.1000000000000000055511151231257827021181583404541015625')

>>> format(Decimal.from_float(0.1), '.17')
'0.10000000000000001'