language-agnostic - 典型 float 中无倒数的最小正整数是多少？

Question

一个常见的假设是 1/x * x == 1。在符合 IEEE 754 标准的常见硬件上打破这个的最小正整数是多少？

当乘法逆的假设失败时，写得不好的有理算术就会停止工作。由于包括 C 和 C++ 在内的许多语言在默认情况下使用舍入到零将 float 转换为整数，因此即使是一个小错误也可能导致积分结果相差一个。

快速测试程序会产生各种结果。

#include <iostream>

int main () {
    {
        double n;
        for ( n = 2; 1 / n * n == 1; ++ n ) ;
        std::cout << n << " (" << 1 - 1/n*n << ")\n";
        for ( ; (int) ( 1 / n * n ) == 1; ++ n ) ;
        std::cout << n << " (" << 1 - 1/n*n << ")\n";
    }
    {
        float n;
        for ( n = 2; 1 / n * n == 1; ++ n ) ;
        std::cout << n << " (" << 1 - 1/n*n << ")\n";
        for ( ; (int) ( 1 / n * n ) == 1; ++ n ) ;
        std::cout << n << " (" << 1 - 1/n*n << ")\n";
    }
}

关于 ideone.com使用 GCC 4.3.4 结果是

41 (5.42101e-20)
45 (5.42101e-20)
41 (5.42101e-20)
45 (5.42101e-20)

使用 GCC 4.5.1 会产生相同的结果，但报告的误差范围恰好为零。

在我的机器上(GCC 4.7.2 或 Clang 4.1)，结果是

49 (1.11022e-16)
49 (1.11022e-16)
41 (5.96046e-08)
41 (5.96046e-08)

这与 --fast-math 选项无关。使用 -mfpmath=387 令人惊讶地产生

41 (5.42101e-20)
41 (5.42101e-20)
41 (5.42101e-20)
41 (5.42101e-20)

值 5×10^-20 似乎暗示 epsilon 对应一个 64 位尾数，即使用 Intel 80 位扩展精度的内部计算。

这似乎高度依赖于 FPU 硬件。是否有适合测试的可靠值？

注意:我不关心什么语言标准或编译器保证 float 系统，尽管我不认为在任何常见的编程系统中有很多有意义的保证。我想知道数字与现实世界计算机之间的交互。

Answer 1

double :

1/41 = 0x1.8f9c18f9c18fap-6，41*0x1.8f9c18f9c18fap-6 = 0x1.000000000000028，四舍五入为 1。 1/45 = 0x1.6c16c16c16c17p-6，和 45*0x1.6c16c16c16c17p-6 = 0x1.00000000000002c，四舍五入为 1。

但是，

1/49 = 0x1.4e5e0a72f0539p-6，和 49*0x1.4e5e0a72f0539p-6 = 0x0.fffffffffffffa4，四舍五入为 0x0.fffffffffffff8 = 0x1.fffffffffffff0p-1

不过，49 确实有倒数!它是 0x1.4e5e0a72f053ap-6。

更一般地，如果 f 是 [1, 2) 中的 float ，则 f 具有倒数。在通常的四舍五入算法中，如果一个数字位于 [1 - 2^-54, 1 + 2^-53] 范围内，则该数字将四舍五入为 1。 请注意，最接近 1/f 的 double d 与 1/f 的距离小于 2^-54。如果 d > 1/f，那么我们就成功了； 1 < f*d < f*(1/f+2^-54) <= 1 + 2^-54 * f < 1 + 2^-53，因此 f*d 舍入为 1。如果 d < 1/f，则 f*d 可能舍入为 1 - 2^-53。如果是，则 f*d 位于 [1 - 2^-53, 1 - 2^-54)。如果取 e = 2^-53 + d，则 e*f > 1 且 e*f = d*f + 2^-53*f < 1 - 2< sup>-53 + 2^-52 = 1 + 2^-53，再次四舍五入为 1。

编辑:上述推理是错误的，因为连续两次 double 之间的步幅相差两倍。没有倒数的 double 的一个例子是 0x1.ffffffbfffffe。 0x1.0000002000001p-1 太小，但 0x1.0000002000002p-1 太大。没有倒数的最小整数示例是 237。1/237 大致为 0x1.1485f0e0acd3B68c6Bp-8，四舍五入为 0x1.1485f0e0acd58p-8。这个数字太小，而它后面的下一个双数太大。