8

考虑以下代码:

// Filename fputest.cpp

#include <cmath>
#include <cstdio>

int main()
{
    double x;
    *(__int64 *) &x = 0xc01448ec3aaa278di64; // -5.0712136427263319
    double sine1 = sin(x);
    printf("%016llX\n", sine1);
    double sine2;
    __asm {
    fld x
    fsin
    fstp sine2
    }
    printf("%016llX\n", sine2);
    return 0;
}

使用 Visual C++ 2012 ( cl fputest.cpp) 编译并执行程序时,输出如下:

3FEDF640D8D36174
3FEDF640D8D36175

问题:

  • 为什么这两个值不同?
  • 是否可以发布一些编译器选项,以便计算的正弦值完全相同?
4

4 回答 4

9

这个问题不是由 long double 到 double 的转换引起的。这可能是由于sin数学库中的例程不准确。

fsin指令指定为在其范围内的操作数(根据英特尔 64 和 IA-32 架构软件开发人员手册,2011 年 10 月,第 1 卷,8.3.10)在 1 ULP(长双精度格式)内产生结果到最近模式。在英特尔酷睿 i7 上,fsin提问者的值,-5.07121364272633190495298549649305641651153564453125 或 -0x1.448ec3aaa278dp+2,产生 0xe.fb206c69b0ba402p-4。我们可以很容易地从这个十六进制看到最后 11 位是 100 0000 0010。这些是从 long double 转换时将被舍入的位。如果它们大于 100 0000 0000,则数字将向上取整。他们更大。因此,将这个 long double 值转换为 double 的结果是 0xe.fb206c69b0ba8p-4,等于 0x1.df640d8d36175p-1 和 0.93631021832247418590355891865328885614871978759765625。另请注意,即使结果低一个 ULP,最后 11 位仍将大于 100 0000 0000 并且仍会向上取整。因此,此结果在符合上述文档的 Intel CPU 上应该不会有所不同。

sin将此与使用产生正确舍入结果的理想例程直接计算双精度正弦进行比较。的值的正弦近似0.93631021832247413051857150785044253634581268961333520518023697738674775240815140702992025520721336793516756640679315765619707343171517531053811196321335899848286682535203710849065933755262347468763562(的Maple 10计算的)。最接近此的双精度为 0x1.df640d8d36175p-1。这与我们通过将fsin结果转换为 double 获得的值相同。

因此,差异不是由 long double 转换为 double 引起的;将 long double结果转换为 double 会产生与理想双精度例程fsin完全相同的结果。sin

sin对于提问者的 Visual Studio 包所使用的例程的准确性,我们没有规范。在商业图书馆中,允许 1 个 ULP 或多个 ULP 的错误是很常见的。观察正弦与双精度值舍入的点有多接近:距双精度值的距离为 0.498864 ULP(双精度 ULP),因此距舍入更改的点为 0.001136 ULP。因此,即使例程中非常轻微的不准确sin也会导致它返回 0x1.df640d8d36174p-1 而不是更接近的 0x1.df640d8d36175p-1。

因此,我推测差异的来源是sin例程中的一个非常小的不准确。

于 2012-09-01T12:23:22.073 回答
1

(注意:正如评论中提到的,这在 VC2012 上不起作用。我把它留在这里作为一般信息。我不建议依赖任何依赖于优化级别的东西!)

我没有 VS2012,但是在 VS2010 编译器上你可以/fp:fast在命令行上指定然后我得到相同的结果。这会导致编译器生成“快速”代码,这些代码不一定完全遵循 C++ 中所需的舍入和规则,但与您的汇编语言计算相匹配。

我不能在 VS2012 中尝试这个,但我想它有相同的选项。

这似乎也只能在优化的构建中/Ox作为一个选项起作用。

于 2012-09-01T07:24:07.160 回答
1

请参阅为什么即使 x == y 也是 cos(x) != cos(y)?

正如大卫在评论中提到的那样,差异来自将 FP 寄存器中的数据移动到不同大小的内存位置(寄存器/ram)。而且也不总是分配;即使是另一个附近的浮点操作也足以刷新 FP 寄存器,从而使任何试图保证特定值的尝试都徒劳无功。如果您需要进行比较,您可以通过将所有结果强制到内存位置来缓解其中的一些问题,如下所示:

float F1 = sin(a); float F2 = sin(b); if (F1 == F2)

但是,即使这样也可能行不通。最好的方法是接受任何浮点运算只会“大部分准确”,并且从程序员的角度来看,即使重复执行相同的运算,这个错误实际上也是不可预测的和随机的。代替

if (F1 == F2)

你应该使用一些东西

if (isArbitrarilyClose(F1, F2))

或者

if (absf(F1 - F2) <= n)

哪里n是一个很小的数字。

于 2012-09-01T07:56:52.590 回答
0

VS2012 中的代码生成器进行了重大更改以支持自动矢量化。部分更改是 x86 浮点数学现在在 SSE2 中完成,不再使用 FPU,这是必要的,因为 FPU 代码无法向量化。SSE2 使用 64 位精度而不是 80 位精度进行计算,因此结果很有可能由于舍入而偏离一位。也是@J99在VS2010中可以与/fp:fast得到一致结果的原因,其编译器仍然使用FPU,而/fp:fast直接使用FSIN结果。

此功能有很多好处,请查看链接网址上的 Jim Hogg 的视频,了解如何利用它。

于 2012-09-01T12:49:43.450 回答