std::pow 在 32 位和 64 位应用程序中产生不同的结果

Question

我发现一些复杂计算的结果不匹配。当我彻底观察中间结果时，它的 std::pow 函数会造成这种不匹配。 以下是输入/输出。

long double dvalue = 2.7182818284589998;
long double dexp = -0.21074699576017999;
long double result = std::powl( dvalue, dexp); 

64 位 -> 结果 = 0.80997896907296496 和 32 位 -> 结果 = 0.80997896907296507

我正在使用 VS2008。 我尝试过使用 long double 并返回 long double 的 pow 函数的其他变体，但仍然看到相同的差异。

double pow( double base, double exponent );

long double powl( long double base, long double exponent );

我已经阅读了一些关于此的信息：

Intel x86 处理器在内部使用 80 位扩展精度，而 double 通常是 64 位宽。不同的优化级别会影响 来自 CPU 的浮点值多久保存到内存中，以及 因此从 80 位精度四舍五入到 64 位精度。或者， 使用 long double 类型，在 gcc 上通常为 80 位宽 避免从 80 位精度舍入到 64 位精度。

谁能让我清楚地理解差异以及克服这种差异的方法。

Answer 1

可能发生的情况是 32 位构建使用 80 位 FPU 寄存器进行计算，而 64 位构建使用使用 64 位值的 SIMD 操作，导致略有差异。请注意，两个答案都同意小数点后 14 位，这大约是 64 位浮点值所希望的最佳值。

Visual C++ 提供compiler options，可让您说出在浮点运算方面您更喜欢速度、一致性还是精度。使用这些选项（例如，/fp:strict），如果这对您很重要，您可能会在两个构建之间获得一致的值。

另请注意，VC++2008 相当旧。较新的版本修复了许多错误，包括一些与浮点相关的错误。 （自 2008 年以来，开源软件中strtod 的流行实现已经检测到并修复了错误。）除了 80 位和 64 位操作之间的精度差异之外，您还可能遇到解析和显示错误。尽管如此，浮点还是很难的，bugs persist.

Answer 2

了解浮点计算最重要的一点是它们（几乎总是）不精确。大多数数字不能精确地表示为浮点数。并且即使计算结果可以准确表示，实际计算的结果可能仍然不完全正确。

解决这个问题的方法是编写不依赖于获得精确结果的代码。例如，您应该几乎从不测试浮点数是否相等。或者，如果您需要测试一个数字是否为正，您的程序可能需要拒绝极小的正数（它们近似为负）或接受极小的负数（它们近似为正）。

同样，你应该尽量避免数值不稳定的算法，因为这些小错误会很快爆发；相反，你应该尝试使用数值稳定的算法，因为它们是容错的。

如何做好数值计算是一门完整的研究领域！

Answer 3

您使用了 double 类型的文字，而不是 long double（您忘记了后缀）。这意味着当您编写 2.7182818284589998（double 的不可能值）时，编译器必须在 2.718281828458999793696193592040799558162689208984375 和 2.71828182845899934960698374197818338871002197265625 之间进行选择，

当您编写 -0.21074699576017999（double 的另一个不可能的值）时，编译器必须在 -0.2107469957601799948054832611887832172214984893798828125 和 -0.210746995760179967049907645559869706630706787109375 之间进行选择。

在默认舍入到最接近的情况下，您存储在dvalue 和dexp 中的值是2.718281828458999793696193592040799558162689208984375 和-0.2107469957601799948054832611887832172214984893798828125（在长双精度中存储双精度不会改变它的值）

pow 的结果应该接近0.8099789690729650165287354526069381795064774873497553965297999359066924950079080502973738475702702999114990234375，然后必须将其放在返回类型中，在您的情况下应该是long double（除了MSVC 没有将它们与double 区分开来，到目前为止我记得和你的结果显示）

将结果放在 64 位 double 中，我们必须在 0.80997896907296496049610823320108465850353240966796875 和 0.80997896907296507151841069571673870086669921875 之间进行选择。

正确答案（四舍五入）是0.80997896907296507151841069571673870086669921875，这正是您在“32 位结果”中得到的结果，截断为0.80997896907296507。

您的“64 位结果”似乎正好是另一个 64 位 double 值，从正确的结果以错误的方式四舍五入（并截断为 0.80997896907296496）。我认为 QoI 错误：gcc、clang、intel 和 oracle 都给出了唯一一个正确的结果（即使它们不必：IEEE 对 pow 的精度要求允许超过 0.5 ulp 的错误）

顺便说一句，如果您的 pow 返回一个 Intel 80 位长的 double，它必须介于 0.8099789690729650164951504420773886749884695746004581451416015625 和 0.809978969072965016549360550701663896688842214643955230712890625 之间，后者是最接近的。

Answer 4

来自wikipedia page on long double

在 x86 架构上，大多数 C 编译器都将 long double 实现为 x86 硬件支持的 80 位扩展精度类型 [...]。 Microsoft Visual C++ for x86 是个例外， 这使得 long double 成为 double 的同义词。

所以当你在 32 位 long double = double 上编译时，但在 x64 上 long double 实际上是一个 80 位浮点，所以结果是不同的。