启用优化后的不同浮点结果 - 编译器错误?

Different floating point result with optimization enabled - compiler bug?

提问人:Bear 提问时间:9/22/2011 最后编辑:Peter MortensenBear 更新时间:6/18/2020 访问量:21920

问:

下面的代码适用于 Visual Studio 2008,无论是否经过优化。但它只适用于没有优化 (O0) 的 g++。

#include <cstdlib>
#include <iostream>
#include <cmath>

double round(double v, double digit)
{
    double pow = std::pow(10.0, digit);
    double t = v * pow;
    //std::cout << "t:" << t << std::endl;
    double r = std::floor(t + 0.5);
    //std::cout << "r:" << r << std::endl;
    return r / pow;
}

int main(int argc, char *argv[])
{
    std::cout << round(4.45, 1) << std::endl;
    std::cout << round(4.55, 1) << std::endl;
}

输出应为:

4.5
4.6

但是带有优化 ( - ) 的 g++ 将输出:O1O3

4.5
4.5

如果我在 t 之前添加关键字,它会起作用,那么是否存在某种优化错误?volatile

在 g++ 4.1.2 和 4.4.4 上进行测试。

这是 ideone 的结果: http://ideone.com/Rz937

我在 g++ 上测试的选项很简单:

g++ -O2 round.cpp

更有趣的结果是,即使我在Visual Studio 2008上打开了选项,结果仍然是正确的。/fp:fast

进一步的问题:

我想知道,我应该总是打开这个选项吗?-ffloat-store

因为我测试的 g++ 版本是随 CentOS/Red Hat Linux 5 和 CentOS/Redhat 6 一起提供的

我在这些平台下编译了我的许多程序,我担心这会在我的程序中引起意想不到的错误。调查我所有的 C++ 代码和使用的库是否存在此类问题似乎有点困难。有什么建议吗?

有没有人对为什么打开Visual Studio 2008仍然有效感兴趣?似乎 Visual Studio 2008 在这个问题上比 g++ 更可靠?/fp:fast

C 优化 G C++-FAQ

评论

1赞 Steve Blackwell 9/22/2011
FWIW,我使用 MinGW 使用 g++ 4.5.0 获得了正确的输出。
1赞 Daniel A. White 9/22/2011
我想知道这是否与它有关:gcc.gnu.org/gcc-4.5/changes.html“GCC 已与 MPC 库集成。这使得 GCC 能够在编译时更准确地评估复杂的算术。
2赞 Daniel A. White 9/22/2011
ideone 使用 4.3.4 ideone.com/b8VXg
5赞 Jakub Wieczorek 9/23/2011
你应该记住,你的例程不太可能可靠地处理各种输出。与将双精度四舍五入为整数相比,这容易受到以下事实的影响:并非所有实数都可以表示,因此您应该期望得到更多这样的错误。
2赞 n. m. could be an AI 9/23/2011
对于那些无法重现错误的人:不要取消注释掉的调试 stmts,它们会影响结果。

答:

6赞 Puppy 9/22/2011 #1

不同的编译器具有不同的优化设置。根据 IEEE 754,其中一些更快的优化设置不会维护严格的浮点规则。Visual Studio 具有特定设置 , , , ,其中违反了有关可以执行的操作的标准。您可能会发现,此标志是控制此类设置中的优化的标志。您可能还会在 GCC 中找到类似的设置,该设置会更改行为。/fp:strict/fp:precise/fp:fast/fp:fast

如果是这种情况,那么编译器之间唯一的区别是,默认情况下,GCC 会在更高的优化中寻找最快的浮点行为,而 Visual Studio 不会在更高的优化级别上更改浮点行为。因此,它可能不一定是实际的错误,而是您不知道自己正在打开的选项的预期行为。

评论

4赞 Mat 9/22/2011
GCC 有一个开关,它没有被任何优化级别打开,因为引用:“它可能导致依赖于 IEEE 或 ISO 数学函数规则/规范的精确实现的程序输出不正确。-ffast-math-O
0赞 NPE 9/22/2011
@Mat:我已经尝试了其他一些东西,但我仍然无法重现这个问题。-ffast-mathg++ 4.4.3
0赞 Kerrek SB 9/22/2011
不错:在这两种情况下,我都得到了大于 的优化水平。-ffast-math4.50
0赞 Kerrek SB 9/23/2011
(更正:我在 GCC 4.4.3 中得到了 和,但不是在 和中,而是在 GCC 4.6.1 中。4.5-O1-O2-O0-O3-O1,2,3
99赞 Maxim Egorushkin 9/23/2011 #2

Intel x86 处理器内部使用 80 位扩展精度,而通常为 64 位宽。不同的优化级别会影响 CPU 浮点值保存到内存中的频率,从而从 80 位精度四舍五入到 64 位精度。double

使用 gcc 选项可获得具有不同优化级别的相同浮点结果。-ffloat-store

或者,使用该类型,该类型在 gcc 上通常为 80 位宽,以避免从 80 位精度舍入到 64 位精度。long double

man gcc说明了一切:

   -ffloat-store
       Do not store floating point variables in registers, and inhibit
       other options that might change whether a floating point value is
       taken from a register or memory.

       This option prevents undesirable excess precision on machines such
       as the 68000 where the floating registers (of the 68881) keep more
       precision than a "double" is supposed to have.  Similarly for the
       x86 architecture.  For most programs, the excess precision does
       only good, but a few programs rely on the precise definition of
       IEEE floating point.  Use -ffloat-store for such programs, after
       modifying them to store all pertinent intermediate computations
       into variables.

在 x86_64 版本中,编译器默认使用 SSE 寄存器,因此不使用扩展精度,并且不会出现此问题。floatdouble

gcc 编译器选项 -mfpmath 控制这一点。

评论

20赞 Mark Ransom 9/23/2011
我认为这就是答案。常数 4.55 转换为 4.549999999999999,这是 64 位中最接近的二进制表示;乘以 10 并再次四舍五入到 64 位,得到 45.5。如果通过将舍入步骤保留在 80 位寄存器中来跳过舍入步骤,则最终会得到 45.49999999999999。
0赞 Bear 9/23/2011
谢谢,我什至不知道这个选项。但是我想知道,我应该总是打开-ffloat-store选项吗?因为我测试的 g++ 版本是随 CentOS/Redhat 5 和 CentOS/Redhat 6 一起提供的。我在这些平台下编译了许多程序,我担心这会导致我的程序中出现意想不到的错误。
5赞 Mark Ransom 9/23/2011
@Bear,debug 语句可能会导致该值从寄存器刷新到内存中。
2赞 Maxim Egorushkin 9/23/2011
@Bear,通常情况下,您的应用程序应该受益于扩展的精度,除非当 64 位浮点数预计会溢出或溢出并产生 .没有好的经验法则,单元测试可以给你一个明确的答案。inf
2赞 plugwash 11/26/2015
@bear 作为一般规则,如果你需要完全可预测的结果和/或确切的人类在纸上做总和会得到的结果,那么你应该避免浮点。-ffloat-store 消除了一个不可预测性的来源,但它不是灵丹妙药。
4赞 Max Lybbert 9/23/2011 #3

对于那些无法重现错误的人:不要取消注释掉的调试 stmts,它们会影响结果。

这意味着问题与调试语句有关。在输出语句期间将值加载到寄存器中似乎会导致舍入错误,这就是为什么其他人发现您可以使用以下方法解决此问题的原因-ffloat-store

进一步的问题:

我想知道,我应该总是打开选项吗?-ffloat-store

要轻率,必须有原因让某些程序员不开机,否则该选项将不存在(同样,某些程序员确实开机也一定有原因)。我不建议总是打开它或总是关闭它。打开它会阻止某些优化,但关闭它允许你得到的那种行为。-ffloat-store-ffloat-store

但是,一般来说,二进制浮点数(如计算机使用)和十进制浮点数(人们熟悉)之间存在一些不匹配,并且这种不匹配会导致与您得到的行为相似的行为(需要明确的是,您得到的行为不是由这种不匹配引起的,但类似的行为可以是)。问题是,由于您在处理浮点数时已经有一些模糊之处,因此我不能说这会使它变得更好或更糟。-ffloat-store

相反,你可能想研究你试图解决的问题的其他解决方案(不幸的是,Koenig 没有指出实际的论文,而且我无法真正找到一个明显的“规范”位置,所以我必须把你送到谷歌)。


如果您不是为了输出目的而四舍五入,我可能会查看 (in ) 和 (in )。考虑原始函数,我相信用对这个函数的调用替换对的调用会更干净:std::modf()cmathstd::numeric_limits<double>::epsilon()limitsround()std::floor(d + .5)

// this still has the same problems as the original rounding function
int round_up(double d)
{
    // return value will be coerced to int, and truncated as expected
    // you can then assign the int to a double, if desired
    return d + 0.5;
}

我认为这暗示了以下改进:

// this won't work for negative d ...
// this may still round some numbers up when they should be rounded down
int round_up(double d)
{
    double floor;
    d = std::modf(d, &floor);
    return floor + (d + .5 + std::numeric_limits<double>::epsilon());
}

一个简单的注释:被定义为“与 1 相加的最小数字,产生一个不等于 1 的数字”。您通常需要使用相对 epsilon(即,以某种方式缩放 epsilon 以解释您正在处理“1”以外的数字这一事实)。和 的总和应该接近 1,因此对加法进行分组意味着它的大小与我们正在做的事情大致相同。如果有的话,会太大(当所有三个的总和都小于 1 时),并且可能会导致我们在不应该的时候四舍五入一些数字。std::numeric_limits<T>::epsilon()d.5std::numeric_limits<double>::epsilon()std::numeric_limits<double>::epsilon()std::numeric_limits<double>::epsilon()


如今,您应该考虑 std::nearbyint()。

评论

0赞 Peter Cordes 11/1/2016
“相对 epsilon”称为 1 ulp(最后为 1 个单位)。 与 X 的 1 个 ULP 相关(但不要使用它;我敢肯定有极端情况,我只是编造的)。的 cppreference 示例有一个缩放它以获得基于 ULP 的相对误差的示例x - nextafter(x, INFINITY)epsilon()
2赞 Peter Cordes 11/1/2016
顺便说一句,2016 年的答案是:首先不要使用 x87。使用 SSE2 数学(64 位二进制文件,或用于制作粗糙的旧 32 位二进制文件),因为 SSE/SSE2 具有没有额外精度的临时文件。 XMM 寄存器中的 var 实际上是 IEEE 64 位或 32 位格式。(与 x87 不同,x87 的寄存器始终为 80 位,存储到内存会四舍五入到 32 位或 64 位。-ffloat-store-mfpmath=sse -msse2doublefloat
10赞 David Hammen 9/23/2011 #4

输出应为:4.5 4.6 如果您具有无限精度,或者您使用的是基于十进制而不是基于二进制的浮点表示形式的设备,则输出就是这样。但是,你不是。大多数计算机使用二进制 IEEE 浮点标准。

正如马克西姆·叶戈鲁什金(Maxim Yegorushkin)在他的回答中已经指出的那样,部分问题在于您的计算机内部使用了80位浮点表示。不过,这只是问题的一部分。问题的基础是 n.nn5 形式的任何数字都没有精确的二进制浮点表示。这些极端情况总是不准确的数字。

如果你真的希望你的舍入能够可靠地舍入这些极端情况,你需要一个舍入算法来解决 n.n5、n.nn5 或 n.nnn5 等(但不是 n.5)总是不准确的事实。查找确定某个输入值是向上舍入还是向下舍入的极端情况,并根据与此极端情况的比较返回向上舍入或向下舍入的值。而且,您确实需要注意,优化编译器不会将发现的极端情况放在扩展精度寄存器中。

请参阅Excel如何成功地舍入浮点数,即使它们不精确?

或者你可以忍受这样一个事实,即极端情况有时会错误地四舍五入。

-1赞 cdcdcd 11/1/2016 #5

就我个人而言,我遇到了同样的问题——从 gcc 到 VS。在大多数情况下,我认为最好避免优化。唯一值得使用的时候是处理涉及大型浮点数据数组的数值方法。即使在拆解之后,我也经常对编译器的选择感到不知所措。很多时候,使用编译器内部函数或自己编写程序集会更容易。

4赞 tmandry 3/15/2018 #6

如果要编译为不包含 SSE2 的 x86 目标,则接受的答案是正确的。所有现代 x86 处理器都支持 SSE2,因此,如果可以利用它,则应:

-mfpmath=sse -msse2 -ffp-contract=off

让我们来分析一下。

-mfpmath=sse -msse2.这通过使用 SSE2 寄存器执行舍入,这比将每个中间结果存储到内存中要快得多。请注意,这已经是 gcc 上 x86-64 的默认值。来自 GCC wiki

在支持 SSE2 的更现代的 x86 处理器上,指定编译器选项可确保所有浮点和双精度运算都在 SSE 寄存器中执行并正确舍入。这些选项不会影响 ABI,因此应尽可能使用这些选项以获得可预测的数值结果。-mfpmath=sse -msse2

-ffp-contract=off.但是,控制舍入不足以实现完全匹配。FMA(融合乘加)指令可以改变舍入行为与非融合对应指令相比,因此我们需要禁用它。这是 Clang 上的默认值,而不是 GCC。正如这个答案所解释的:

FMA 只有一个舍入(它有效地为内部临时乘法结果保持无限精度),而 ADD + MUL 有两个舍入。

通过禁用 FMA,我们得到的结果在调试和发布时完全匹配,但代价是一些性能(和准确性)。我们仍然可以利用 SSE 和 AVX 的其他性能优势。

评论

0赞 ManuelAtWork 12/3/2021
+1 将 FMA 指令作为潜在的麻烦来源。此外,SSE 指令在某些平台上可能存在错误,因此请尽可能优先使用。-msse2-msse
1赞 calandoa 6/16/2018 #7

我深入研究了这个问题,我可以带来更多的精度。首先,根据 gcc 在 x84_64 上,4.45 和 4.55 的精确表示如下(使用 libquadmath 打印最后一个精度):

float 32:   4.44999980926513671875
double 64:  4.45000000000000017763568394002504646778106689453125
doublex 80: 4.449999999999999999826527652402319290558807551860809326171875
quad 128:   4.45000000000000000000000000000000015407439555097886824447823540679418548304813185723105561919510364532470703125

float 32:   4.55000019073486328125
double 64:  4.54999999999999982236431605997495353221893310546875
doublex 80: 4.550000000000000000173472347597680709441192448139190673828125
quad 128:   4.54999999999999999999999999999999984592560444902113175552176459320581451695186814276894438080489635467529296875

所述,该问题出在FPU寄存器的80位大小上。

但是为什么这个问题在 Windows 上从未发生过呢?在 IA-32 上,x87 FPU 配置为使用 53 位尾数的内部精度(相当于 64 位的总大小:)。对于 Linux 和 Mac OS,使用默认精度 64 位(相当于 80 位的总大小:)。因此,通过更改 FPU 的控制字(假设指令序列会触发错误),在这些不同的平台上,问题应该是可能的,也可能不是。该问题已作为错误 323 报告给 gcc(至少阅读评论 92!doublelong double

要在 Windows 上显示尾数精度,您可以使用 VC++ 以 32 位编译:

#include "stdafx.h"
#include <stdio.h>  
#include <float.h>  

int main(void)
{
    char t[] = { 64, 53, 24, -1 };
    unsigned int cw = _control87(0, 0);
    printf("mantissa is %d bits\n", t[(cw >> 16) & 3]);
}

在 Linux/Cygwin 上:

#include <stdio.h>

int main(int argc, char **argv)
{
    char t[] = { 24, -1, 53, 64 };
    unsigned int cw = 0;
    __asm__ __volatile__ ("fnstcw %0" : "=m" (*&cw));
    printf("mantissa is %d bits\n", t[(cw >> 8) & 3]);
}

请注意,使用 gcc 时,您可以使用 设置 FPU 精度,尽管在 Cygwin 中会忽略它。但请记住,它会修改尾数的大小,但不会改变指数的大小,让大门向其他类型的不同行为敞开。-mpc32/64/80

在x86_64架构上,SSE 是按照 tmandry 所说的使用,因此除非您强制使用旧的 x87 FPU 进行 FP 计算,或者除非您以 32 位模式编译(您将需要 multilib 包),否则不会出现问题。我可以在 Linux 上使用不同的标志和 gcc 版本组合重现该问题:-mfpmath=387-m32

g++-5 -m32 floating.cpp -O1
g++-8 -mfpmath=387 floating.cpp -O1

我在 Windows 或 Cygwin 上使用 VC++/gcc/tcc 尝试了一些组合,但该错误从未出现。我想生成的指令顺序是不一样的。

最后,请注意,使用 4.45 或 4.55 防止此问题的一种奇特方法是使用 ,但支持确实很少......我花了很多时间只是为了能够用 !_Decimal32/64/128libdfp