我在 Windows 下使用 Visual Studio 和 Intel 编译器时遇到双精度问题。
#include <stdio.h>
#include <math.h>
int main()
{
double result = 42;
int rounds = 14;
for(int a=1; a<=rounds; a++)
{
result = sqrt(sqrt(result));
printf("Round %2.i is %.17g\n", a, result);
}
for(int a=1; a<=rounds; a++)
{
result = pow(result, 4);
printf("Round %2.i is %.17g\n", a, result);
}
return 0;
}
使用此代码,我在 Linux 下使用这些编译器设置生成:
Linux下Intel编译器:
icc -O2 -fp-model double test.cpp -o test
linux下的GCC编译器:
g++ -O2 test.cpp -o test
以及 Windows 下的 Matlab 和 Scilab 正确结果:
Round 1 is 2.5457298950218306
Round 2 is 1.2631446315995756
Round 3 is 1.0601401197016873
Round 4 is 1.0147073765340913
Round 5 is 1.0036567375971333
Round 6 is 1.0009129334669549
Round 7 is 1.0002281552726189
Round 8 is 1.0000570339386641
Round 9 is 1.0000142581797196
Round 10 is 1.0000035645258711
Round 11 is 1.0000008911302767
Round 12 is 1.0000002227824947
Round 13 is 1.000000055695619
Round 14 is 1.0000000139239045
Round 1 is 1.000000055695619
Round 2 is 1.0000002227824947
Round 3 is 1.0000008911302765
Round 4 is 1.0000035645258705
Round 5 is 1.0000142581797171
Round 6 is 1.0000570339386543
Round 7 is 1.0002281552725802
Round 8 is 1.0009129334668005
Round 9 is 1.0036567375965142
Round 10 is 1.0147073765315875
Round 11 is 1.0601401196912235
Round 12 is 1.263144631549705
Round 13 is 2.545729894619797
Round 14 is 41.999999973468661
Linux 下使用此 Intel 编译器选项
icc -O2 test.cpp -o test
还有这个linux下的gcc编译器选项
g++ -O2 -ffast-math test.cpp -o test
和所有 Visual Studio 2012 编译器设置我得到这个错误的结果:
Round 1 is 2.5457298950218306
Round 2 is 1.2631446315995756
Round 3 is 1.0601401197016873
Round 4 is 1.0147073765340913
Round 5 is 1.0036567375971333
Round 6 is 1.0009129334669549
Round 7 is 1.0002281552726189
Round 8 is 1.0000570339386641
Round 9 is 1.0000142581797196
Round 10 is 1.0000035645258711
Round 11 is 1.0000008911302767
Round 12 is 1.0000002227824947
Round 13 is 1.000000055695619
Round 14 is 1.0000000139239045
Round 1 is 1.000000055695619
Round 2 is 1.0000002227824947
Round 3 is 1.0000008911302767
Round 4 is 1.0000035645258711
Round 5 is 1.0000142581797198
Round 6 is 1.0000570339386647
Round 7 is 1.0002281552726222
Round 8 is 1.0009129334669684
Round 9 is 1.0036567375971872
Round 10 is 1.0147073765343093
Round 11 is 1.0601401197025984
Round 12 is 1.2631446316039172
Round 13 is 2.5457298950568319
Round 14 is 42.000000002309839
在 Windows 下,以下设置将被忽略,我总是得到错误的结果。
Windows 下的 Intel 编译器选项:
icl /O2 /fp:double test.cpp -o test.exe
Windows 下的 Visual Studio 编译器选项:
/fp:precise
/fp:strict
/fp:fast
产生了所有相同的“错误数字”...
为什么会出现此错误?如何修复它,使我在所有平台上都获得相同的工作双精度?
非常感谢和问候
娜塔莉
编辑:
为了进行基准测试,我循环了代码:
#include <stdio.h>
#include <math.h>
int main()
{
double x = 1;
double result = 0;
int rounds = 12;
while (x!=result)
{
x++;
result=x;
for(int a=1; a<=rounds; a++)
{
result = sqrt(sqrt(result));
}
for(int a=1; a<=rounds; a++)
{
result = pow(result, 4);
}
}
printf("The next possible with %2.i rounds is %.0f\n", rounds, result);
return 0;
}
使用可用的 Linux 编译器编译: 英特尔编译器:
icc -O2 -fp-model double speedtest3.cpp -o speedtest3
结果:
The next possible with 12 rounds is 3671078
real 0m1.950s
user 0m1.947s
sys 0m0.003s
GCC 编译器:
g++ -O2 speedtest3.cpp -o speedtest3
结果:
The next possible with 12 rounds is 3671078
real 0m3.445s
user 0m3.442s
sys 0m0.004s
带 pow 的代码:
#include <stdio.h>
#include <math.h>
int main()
{
double x = 1;
double result = 0;
int rounds = 12;
while (x!=result)
{
x++;
result=x;
for(int a=1; a<=rounds; a++)
{
result = pow(result, 0.25);
}
for(int a=1; a<=rounds; a++)
{
result = pow(result, 4);
}
}
printf("The next possible with %2.i rounds is %.0f\n", rounds, result);
return 0;
}
使用相同选项编译: 英特尔编译器结果:
The next possible with 12 rounds is 3671078
real 0m2.887s
user 0m2.885s
sys 0m0.004s
GCC 编译器结果:
The next possible with 12 rounds is 3671078
real 0m5.905s
user 0m5.905s
sys 0m0.008s
结果:
sqrt 比 pow 快近 2 倍 intel 编译器比 gcc 编译器快近 2 倍
他们都是正确的。
在处理任何精度的浮点数学时,您不能要求结果是精确到最后一位的任何数字。期望准确地获得 41.999999973468661 与期望准确地获得 42.000000000000000 一样被误导。
以下是如何获得不同结果的示例:C++ 中的双精度浮点运算是(通常取决于编译器和平台)64 位。然而,Intel CPU 的浮点硬件使用80 位在内部计算数字。在某个时刻,这些值会从 FPU 中提取并修剪为 64 位。 “在某个时刻”的含义取决于编译器及其优化器认为合适或最佳的内容,因此,不同的编译器/设置将导致计算的不同部分进行舍入,从而影响结果的确切数字。
即使在不同的计算机上使用相同的二进制程序也可能合法地改变结果,并且仍然是“正确的”。例如,曾经有一些英特尔处理器根本没有硬件浮点单元(很久以前,最后一个是486SX)。在这些情况下,您的程序将在一个单独的纯软件库上运行计算,该库在 80 位精度下根本不执行任何操作。
附录:
我想正确看待事情。 Wolfram Alpha 表示纽约和旧金山之间的距离约为 2577 英里。在这个比例下,您如此担心两个数字之间的差异就像说您自己对距离的估计“错误了 1/10th 英寸”(*)。如果这种差异对您的问题很重要,那么浮点算术是不适合该工作的工具。
* 相对误差为 (42.000000002309839 - 41.999999973468661)/41.999999973468661 = 6.8669471471949833966026905617771e-10。 2577 英里 = 163278720 英寸。该测量的等效散度为 (163278720 英寸 * 6.8669471471949833966026905617771e-10) = 0.1121226... 英寸。所有计算均在 Windows 7 上的 Windows 计算器中完成,该计算器使用定点运算。
查看 pow(double,double) 的一种实现方式,这里。来自不同供应商的不同编译器必然包含不同的实现,从而产生差异。
正如所建议的,更好的选择是使用
double h = x*x;
x = h*h;
你所说的“不正确”更准确,比较一下:
delta is 2.3098394308362913e-09 42.00000000230983 - 42
delta is 2.6531338903623691e-08 42 - 41.999999973468661
稍后你想证明或展示什么?重复执行一些数值计算通常比直接计算结果更糟糕;因此 pow( 42, 1.0/(1024*1024*64)) 理论上比在循环中求 4 次根更可取,相反的过程也是如此: pow( result, 1024*1024*64 ),这可能使 pow 能够使用最少次数的乘法。
无论如何,42 的无数次方无法表示为机器数,因此它的无数次方必然会偏离 42。您可以使用初等数学来估计误差。