我将问题简化为一个简单的测试用例：
#include <stdio.h>
int main(void) {
printf("%d bits: ", sizeof(void*) * 8);
volatile double d = 10.870191700000001;
volatile double x = 0.10090000000000002;
d += x * 30.07;
printf("%0.15f\n", d);
}
为避免编译器差异,我使用相同的编译器和交叉编译.在这种情况下,我在Core i5 CPU上使用Windows 7上的TDM-GCC 64位5.1.0.这是我的命令行：
gcc double_test.c -o double_test.exe -m32 -O0 && double_test.exe && gcc double_test.c -o double_test.exe -m64 -O0 && double_test.exe
输出是：
32 bits: 13.904254700000001 64 bits: 13.904254700000003
在这种情况下,错误是最小的,但在我的完整测试用例中,错误可以加起来并足以使我的输出加倍.
如何获得与32位输出匹配的精确位操作？
最接近我得到的东西是使用-ffloat-store,但是在这个片段中,它得到了像64位一样的32位执行,而我需要恰恰相反.但是,这对我的图书馆没有任何明显的影响.我还测试了-fexcess-precision = standard和-mfp-math选项无济于事.
既然你说你需要更精确的… 01结果,以及确定性,你遗憾的是你不能在你的32位构建中使用-msse2 -mfpmath = sse.寻找决定论的未来读者应该使用它.您可以使用-mfpmath = 387来请求gcc在64位模式下使用慢速/过时的x87数学运算,这不是默认值.调用约定在xmm寄存器中传递/返回FP args,因此这比在32位模式下更糟糕,有时需要额外的存储/重新加载.
peter@volta:/tmp$gcc -m64 -mfpmath=387 -O3 fp-prec.c -o fp-64-387 peter@volta:/tmp$./fp-64-387 64 bits: 13.904254700000001
我不确定当自动矢量化可能时gcc是否严格限制为x87.如果是这样,你会错过表现.
而BTW,在你的例子中,… 01是在将x + 30.07添加到d之前为x * 30.07保持80位临时值的额外精度的结果. (d是易失性的,但是d = stuff仍然相当于d = d的东西,因此x * 30.07不会首先四舍五入到64位双精度).
你可以使用long double,例如d = x *(long double)30.07强制80位临时存在.在x86-64 System V ABI Linux / OS X / * BSD /等中,long double是80位,但在x64 Windows上,它与64位double相同.所以这可能不适合你.
在这种情况下,您可以使用FMA获得相同的结果,该FMA在执行添加之前保持乘法的无限精度.在没有FMA支持的情况下,这在硬件上很慢,但是fma(d,30.07,x)将可靠地提供您想要的结果.
如果需要,请在需要精度的地方使用它.
如果在启用FMA的情况下进行编译,则可以内联到FMA指令. (例如-march =我的Skylake CPU上的原生)
即使不使用fma()math.h函数,gcc也会在优化时将mul添加表达式添加到FMA中. (与Clang不同,我认为默认情况下不使用-ffast-math进行FP_CONTRACT).请注意,我没有使用-march = 387
# your original source code, using an FMA instruction (native=skylake in my case) peter@volta:/tmp$gcc -m64 -march=native -O3 fp-prec.c -o fp-64-native peter@volta:/tmp$./fp-64-native 64 bits: 13.904254700000001
主要的相关部分是：
57e: c5 fb 10 44 24 08 vmovsd xmm0,QWORD PTR [rsp+0x8] # load x 584: c5 fb 10 0c 24 vmovsd xmm1,QWORD PTR [rsp] # load d 589: c4 e2 f1 99 05 d6 01 00 00 vfmadd132sd xmm0,xmm1,QWORD PTR [rip+0x1d6] # the 30.07 constant 592: c5 fb 11 04 24 vmovsd QWORD PTR [rsp],xmm0 # store d 597: c5 fb 10 04 24 vmovsd xmm0,QWORD PTR [rsp] # reload d 59c: e8 8f ff ff ff call 530 <printf@plt>
FP决定论一般很难.
另见https://randomascii.wordpress.com/2013/07/16/floating-point-determinism/和https://randomascii.wordpress.com/2012/03/21/intermediate-floating-point-precision/