python - 公式在IEEE 754中失败的可能性

在我的电脑上，我可以检查一下

(0.1 + 0.2) + 0.3 == 0.1 + (0.2 + 0.3)

评估为False。

更一般而言，我可以使用以下模拟估算公式(a + b) + c == a + (b + c)在17%上统一且独立地选择a,b,c时的[0,1]大致失败:

import numpy as np
import numexpr

np.random.seed(0)
formula = '(a + b) + c == a + (b + c)'


def failure_expectation(formula=formula, N=10**6):
    a, b, c = np.random.rand(3, N)
    return 1.0 - numexpr.evaluate(formula).mean()
# e.g. 0.171744

我想知道是否有可能手动得出这个概率，例如使用浮点标准中的定义以及对均匀分布的一些假设。

给定以下答案，我认为至少到现在为止，原始问题的以下部分是遥不可及的。



(以下内容可能是numpy随机数生成的产物，但似乎仍然很有趣。)

奖励问题基于NPE的观察。我们可以使用以下代码来生成失败概率，以求在一系列范围[[-n,n] for n in range(100)]上的均匀分布:

import pandas as pd

def failures_in_symmetric_interval(n):
    a, b, c = (np.random.rand(3, 10**4) - 0.5) * n
    return 1.0 - numexpr.evaluate(formula).mean()

s = pd.Series({
    n: failures_in_symmetric_interval(n)
    for n in range(100)
})

情节看起来像这样:



特别地，当0是n的幂且似乎具有分形模式时，失败概率会下降到2。看起来每个“下降”的失败概率都等于先前的某个“峰值”的失败概率。对为什么会这样的任何解释都很好!

最佳答案

绝对有可能手工评估这些东西，但是我所知道的唯一方法是乏味的，并且涉及大量的个案枚举。

例如，对于确定(a + b) + c == a + (b + c)的概率的特定示例，该概率为53/64，在机器epsilon的几倍内。因此，不匹配的可能性为11/64，或大约17.19％，这与您从仿真中观察到的结果一致。

首先，请注意，在这种特殊情况下，存在一个主要的简化因素，那就是Python和NumPy的“[0，1]均匀”随机数的形式始终是n/2**53，对于range(2**53)中的某个整数n，并且在受底层梅森·Twister PRNG的约束，每个这样的数字都有可能发生。由于在2**62范围内大约有[0.0, 1.0] IEEE 754 binary64可表示的值，这意味着这些IEEE 754值中的绝大部分不是由random.random()(或np.random.rand())生成的。这个事实大大简化了分析，但也意味着它有点作弊。

这是一个不完整的草图，仅是提供了所涉及内容的概念。为了计算53/64的值，我不得不分为5种情况:

a + b

a + b = 1的情况。这里(a + b)+ c将是正确总和的正确舍入值，但a +(b + c)可能不是。根据a，b和c的最低有效位的奇偶性，我们可以进一步分割为子情况。让我们滥用术语，如果它的形式为带有奇数n的n/2**53，则叫“奇数”，如果它的形式为带有n偶数的n/2**53，则叫“偶数”，b和c也是类似的。如果b和c具有相同的奇偶校验(将发生一半的时间)，则将精确计算(b + c)，并且a +(b + c)必须再次匹配(a + b)+ c。对于其他情况，达成协议(protocol)的可能性在每种情况下均为1/2；细节非常相似，但是例如，在a为奇数，b为奇数且c为偶数的情况下，(a + b)+ c的计算准确，而在计算a +(b + c)的情况下，四舍五入误差，每个误差精确等于2**-53。如果这两个错误的方向相反，则将它们抵消，我们将达成协议(protocol)。如果没有，我们就不会。总体而言，在这种情况下达成协议(protocol)的可能性为3/4。

a + b> = 1且b + c

a + b> = 1且b + c> = 1，但a + b + c

情况a + b + c> =2。在这种情况下，由于我们将最终结果四舍五入为2**-53的四倍，因此有必要不仅查看a，b和c的奇偶校验，而且看一下最后两个有效位。我将为您省去一些细节，但达成协议(protocol)的可能性为13/16。

最后，我们可以将所有这些情况放在一起。为此，我们还需要知道三元组(a，b，c)在每种情况下着陆的概率。 a + b
因此我们的总计是1/3 * 1 + 1/6 * 3/4 + 1/6 * 3/4 + 1/6 * 21/32 + 1/6 * 13/16，如所声称的那样就是53/64。

最后一点:53/64几乎肯定不是一个正确的答案-要获得完美准确的答案，我们需要注意a + b，b + c或a + b + c的所有极端情况达到binade边界(1.0或2.0)。当然可以对上述方法进行改进，以精确计算出多少个2**109可能的三元组(a，b，c)满足(a + b)+ c == a +(b + c)，但在满足之前我上床 sleep 的时间。但是，极端情况应该构成案件总数的1/2**53的顺序，因此我们对53/64的估计至少应精确到15个小数位。

当然，上面遗漏了许多细节，但是我希望它能使人们对如何做到这一点有所了解。

关于python - 公式在IEEE 754中失败的可能性，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/57518436/