提问人:Bear 提问时间:9/22/2011 最后编辑:Peter MortensenBear 更新时间:6/18/2020 访问量:21920
启用优化后的不同浮点结果 - 编译器错误?
Different floating point result with optimization enabled - compiler bug?
问:
下面的代码适用于 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++ 将输出:O1
O3
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
答:
不同的编译器具有不同的优化设置。根据 IEEE 754,其中一些更快的优化设置不会维护严格的浮点规则。Visual Studio 具有特定设置 , , , ,其中违反了有关可以执行的操作的标准。您可能会发现,此标志是控制此类设置中的优化的标志。您可能还会在 GCC 中找到类似的设置,该设置会更改行为。/fp:strict
/fp:precise
/fp:fast
/fp:fast
如果是这种情况,那么编译器之间唯一的区别是,默认情况下,GCC 会在更高的优化中寻找最快的浮点行为,而 Visual Studio 不会在更高的优化级别上更改浮点行为。因此,它可能不一定是实际的错误,而是您不知道自己正在打开的选项的预期行为。
评论
-ffast-math
-O
-ffast-math
g++ 4.4.3
-ffast-math
4.5
0
4.5
-O1
-O2
-O0
-O3
-O1,2,3
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 寄存器,因此不使用扩展精度,并且不会出现此问题。float
double
gcc
编译器选项 -mfpmath
控制这一点。
评论
inf
对于那些无法重现错误的人:不要取消注释掉的调试 stmts,它们会影响结果。
这意味着问题与调试语句有关。在输出语句期间将值加载到寄存器中似乎会导致舍入错误,这就是为什么其他人发现您可以使用以下方法解决此问题的原因-ffloat-store
进一步的问题:
我想知道,我应该总是打开选项吗?
-ffloat-store
要轻率,必须有原因让某些程序员不开机,否则该选项将不存在(同样,某些程序员确实开机也一定有原因)。我不建议总是打开它或总是关闭它。打开它会阻止某些优化,但关闭它允许你得到的那种行为。-ffloat-store
-ffloat-store
但是,一般来说,二进制浮点数(如计算机使用)和十进制浮点数(人们熟悉)之间存在一些不匹配,并且这种不匹配会导致与您得到的行为相似的行为(需要明确的是,您得到的行为不是由这种不匹配引起的,但类似的行为可以是)。问题是,由于您在处理浮点数时已经有一些模糊之处,因此我不能说这会使它变得更好或更糟。-ffloat-store
相反,你可能想研究你试图解决的问题的其他解决方案(不幸的是,Koenig 没有指出实际的论文,而且我无法真正找到一个明显的“规范”位置,所以我必须把你送到谷歌)。
如果您不是为了输出目的而四舍五入,我可能会查看 (in ) 和 (in )。考虑原始函数,我相信用对这个函数的调用替换对的调用会更干净:std::modf()
cmath
std::numeric_limits<double>::epsilon()
limits
round()
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
.5
std::numeric_limits<double>::epsilon()
std::numeric_limits<double>::epsilon()
std::numeric_limits<double>::epsilon()
如今,您应该考虑 std::nearbyint()。
评论
x - nextafter(x, INFINITY)
epsilon()
-ffloat-store
-mfpmath=sse -msse2
double
float
输出应为:4.5 4.6 如果您具有无限精度,或者您使用的是基于十进制而不是基于二进制的浮点表示形式的设备,则输出就是这样。但是,你不是。大多数计算机使用二进制 IEEE 浮点标准。
正如马克西姆·叶戈鲁什金(Maxim Yegorushkin)在他的回答中已经指出的那样,部分问题在于您的计算机内部使用了80位浮点表示。不过,这只是问题的一部分。问题的基础是 n.nn5 形式的任何数字都没有精确的二进制浮点表示。这些极端情况总是不准确的数字。
如果你真的希望你的舍入能够可靠地舍入这些极端情况,你需要一个舍入算法来解决 n.n5、n.nn5 或 n.nnn5 等(但不是 n.5)总是不准确的事实。查找确定某个输入值是向上舍入还是向下舍入的极端情况,并根据与此极端情况的比较返回向上舍入或向下舍入的值。而且,您确实需要注意,优化编译器不会将发现的极端情况放在扩展精度寄存器中。
或者你可以忍受这样一个事实,即极端情况有时会错误地四舍五入。
就我个人而言,我遇到了同样的问题——从 gcc 到 VS。在大多数情况下,我认为最好避免优化。唯一值得使用的时候是处理涉及大型浮点数据数组的数值方法。即使在拆解之后,我也经常对编译器的选择感到不知所措。很多时候,使用编译器内部函数或自己编写程序集会更容易。
如果要编译为不包含 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 的其他性能优势。
评论
-msse2
-msse
我深入研究了这个问题,我可以带来更多的精度。首先,根据 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!double
long 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/128
libdfp
评论