Libc hypot 函数似乎为 double 类型返回不正确的结果......为什么？

Question

#include <tgmath.h>
#include <iostream>
int main(int argc, char** argv) {

        #define NUM1 -0.031679909079365576
        #define NUM2 -0.11491794452567111

        std::cout << "double precision :"<< std::endl;
        typedef std::numeric_limits< double > dbl;
        std::cout.precision(dbl::max_digits10);
        std::cout << std::hypot((double)NUM1, (double)NUM2);
        std::cout << " VS sqrt :" << sqrt((double )NUM1*(double )NUM1 
                                  + (double )NUM2*(double )NUM2) << std::endl;

        std::cout << "long double precision :"<< std::endl;
        typedef std::numeric_limits<long double > ldbl;
        std::cout.precision(ldbl::max_digits10);
        std::cout << std::hypot((long double)NUM1, (long double)NUM2);
        std::cout << " VS sqrt :" << sqrt((long double )NUM1*(long double )NUM1 + (long double )NUM2*(long double )NUM2);
}

在 Linux 下返回（Ubuntu 18.04 clang 或 gcc，无论优化，glic 2.25）：

双精度： 0.1192046585217293 VS sqrt：0.11920465852172932

long 双精度： 0.119204658521729311251 VS sqrt：0.119204658521729311251

根据 cppreference ：

实现通常保证小于 1 ulp（最后一个单位）的精度：GNU、BSD、Open64 std::hypot(x, y) 等价于 std::abs(std::complex(x,y)) POSIX 规定只有当两个参数都低于正常且正确的结果也低于正常时才会发生下溢（这禁止幼稚的实现）

所以，hypot((double)NUM1, (double)NUM2) 应该返回 0.11920465852172932，我想（作为天真的 sqrt 实现）。 在 Windows 上，使用 MSVC 64 位，就是这种情况。

为什么我们使用 glibc 会看到这种差异？如何解决这种不一致？

Answer 1

• 0.11920465852172932 由 0x1.e84324de1b576p-4 表示（作为双精度）
• 0.11920465852172930 由 0x1.e84324de1b575p-4 表示（作为双精度）
• 0.119204658521729311251 是 long-double 结果，我们可以假设它在小数点后几位是正确的。即确切的结果更接近四舍五入的结果。

那些 FP 位模式仅在尾数的低位（即有效位）上有所不同，确切的结果在它们之间。 因此它们每个都有小于 1 ulp 的舍入误差，实现了典型实现（包括 glibc）的目标。

与 IEEE-754 “基本”操作 (add/sub/mul/div/sqrt) 不同，hypot 不需要“正确舍入”。这意味着

碰巧在这种情况下，简单的计算方法产生了正确舍入的结果，而 glibc 的 std::hypot 的“安全”实现（在相加前对小数进行平方时必须避免下溢）产生的结果 >0.5 但是

---

您没有指定是否在 32 位模式下使用 MSVC。

大概 32 位模式将使用 x87 进行 FP 数学运算，从而提供额外的临时精度。尽管某些 MSVC 版本的 CRT 代码在每次操作后将 x87 FPU 的内部精度设置为舍入到 53 位尾数，因此它的行为类似于使用实际 double 的 SSE2，只是指数范围更广。见Bruce Dawson's blog post。

所以我不知道除了运气之外，MSVC 的std::hypot 是否有任何理由得到正确舍入的结果。

注意MSVC中的long double与64位double的类型相同；该 C++ 实现不公开 x86 / x86-64 的 80 位硬件扩展精度类型。 （64 位尾数）。