c# - 浮点运算以及x86和x64上下文

Question

我们正在VisualStudio进程上下文（x86上下文）和VisualStudio上下文（x64上下文）之外运行一些代码。我注意到以下代码在两种情况下都提供了不同的结果（x86中为100000000000，x64中为99999997952）

float val = 1000f;
val = val * val;
return (ulong)(val * 100000.0f);

我们需要以可靠的方式从浮点值中获得ulong值，无论上下文和ulong值如何，它仅用于哈希目的。我在x64和x86上下文中测试了此代码，并确实获得了相同的结果，它看起来很可靠：

float operandFloat = (float)obj;
byte[] bytes = BitConverter.GetBytes(operandFloat);
Debug.Assert(bytes.Length == 4);
uint @uint = BitConverter.ToUInt32(bytes, 0);
return (ulong)@uint;

该代码可靠吗？

Answer 1

正如其他人在评论中推测的那样，您观察到的差异是执行浮点运算时差分精度的结果，这是由于32位和64位内部版本执行这些操作的方式之间的差异引起的。

您的代码由32位（x86）JIT编译器转换为以下目标代码：

fld   qword ptr ds:[0E63308h]  ; Load constant 1.0e+11 onto top of FPU stack.
sub   esp, 8                   ; Allocate 8 bytes of stack space.
fstp  qword ptr [esp]          ; Pop top of FPU stack, putting 1.0e+11 into
                               ;  the allocated stack space at [esp].
call  73792C70                 ; Call internal helper method that converts the
                               ;  double-precision floating-point value stored at [esp]
                               ;  into a 64-bit integer, and returns it in edx:eax.
                               ; At this point, edx:eax == 100000000000.

请注意，优化器已将算术运算（(1000f * 1000f) * 100000f）折叠为常数1.0e + 11。它已将该常量存储在二进制文件的数据段中，并将其加载到x87浮点堆栈的顶部（fld指令）。然后，代码通过sub压缩堆栈指针（esp）分配8字节的堆栈空间（足以容纳64位双精度浮点值）。 fstp指令将值从x87浮点堆栈的顶部弹出，并将其存储在其内存操作数中。在这种情况下，它将其存储到我们刚刚分配给堆栈的8个字节中。所有这些改组是毫无意义的：它可能只是将浮点常量1.0e + 11直接加载到内存中，而不是通过x87 FPU来回进行行程，但是JIT优化器并不完美。最后，JIT发出代码以调用内部帮助器函数，该函数将存储在内存（1.0e + 11）中的双精度浮点值转换为64位整数。按照32位Windows调用约定的惯例，在寄存器对edx:eax中返回64位整数结果。这段代码完成后，edx:eax包含64位整数值100000000000（即1.0e + 11），完全符合您的期望。

（希望这里的术语不太混乱。请注意，有两个不同的“堆栈”。x87FPU有一系列寄存器，它们像堆栈一样被访问。我将其称为FPU堆栈。然后，您可能熟悉的堆栈，它存储在主存储器中，并可以通过堆栈指针esp访问。）


但是，64位（x86-64）JIT编译器的处理方式略有不同。此处的最大区别是64位目标始终使用SSE2指令进行浮点运算，因为所有支持AMD64的芯片也都支持SSE2，并且SSE2比旧的x87 FPU更高效，更灵活。具体来说，64位JIT将您的代码转换为以下内容：

movsd  xmm0, mmword ptr [7FFF7B1A44D8h]  ; Load constant into XMM0 register.
call   00007FFFDAC253B0                  ; Call internal helper method that converts the
                                         ;  floating-point value in XMM0 into a 64-bit int
                                         ;  that is returned in RAX.

事情马上就出了问题，因为第一条指令加载的常数值为0x42374876E0000000，这是99999997952.0的二进制浮点表示形式。问题不在于正在转换为64位整数的辅助函数。相反，它是JIT编译器本身，特别是优化器例程，用于预先计算常量。

为了深入了解问题的出处，我们将关闭JIT优化，然后查看代码如下：

movss    xmm0, dword ptr [7FFF7B1A4500h]  
movss    dword ptr [rbp-4], xmm0  
movss    xmm0, dword ptr [rbp-4]  
movss    xmm1, dword ptr [rbp-4]  
mulss    xmm0, xmm1  
mulss    xmm0, dword ptr [7FFF7B1A4504h]  
cvtss2sd xmm0, xmm0  
call     00007FFFDAC253B0 

第一条movss指令将一个单精度浮点常量从内存加载到xmm0寄存器中。但是，这次，该常数为0x447A0000，它是1000的精确二进制表示形式-代码中的初始float值。

第二条movss指令转过来，并将该值从xmm0寄存器存储到内存中，而第三条movss指令将刚存储的值从存储器中重新加载回xmm0寄存器中。 （告诉您这是未经优化的代码！）它还将来自内存的相同值的第二个副本加载到xmm1寄存器中，然后将mulss和xmm0中的两个单精度值相乘（xmm1）。一起。这是您的val = val * val代码的字面翻译。精确地，此操作的结果（以xmm0结尾）为0x49742400或1.0e + 6。

第二个mulss指令执行val * 100000.0f操作。它隐式加载单精度浮点常量1.0e + 5，并将其与xmm0中的值相乘（回想一下，该值为1.0e + 6）。不幸的是，此操作的结果与您期望的不一样。实际上是9.9999998e + 10，而不是1.0e + 11。为什么？因为不能将1.0e + 11精确地表示为单精度浮点值。最接近的表示形式是0x51BA43B7或9.9999998e + 10。

最后，cvtss2sd指令将xmm0中的（错误！）标量单精度浮点值执行就地转换为标量双精度浮点值。 Neitsa在对问题的评论中建议，这可能是问题的根源。实际上，正如我们所看到的，问题的根源是上一条指令，该指令执行乘法。 cvtss2sd只是将已经不精确的单精度浮点表示形式（0x51BA43B7）转换为不精确的双精度浮点表示形式：0x42374876E0000000或99999997952.0。

这正是JIT编译器执行的一系列操作，以生成初始的双精度浮点常量，该常量以优化的代码加载到xmm0寄存器中。

尽管我在整个答案中一直暗示要归咎于JIT编译器，但事实并非如此！如果您在针对SSE2指令集时用C或C ++编译了相同的代码，则将获得完全相同的不精确结果：99999997952.0。 JIT编译器的性能与预期的一样—如果将人们的期望正确地校准为不精确的浮点运算！


那么，这个故事的寓意是什么？有两个。首先，浮点运算很棘手，there is a lot to know about them。其次，鉴于此，在进行浮点运算时，请始终使用现有的最高精度！

32位代码产生正确的结果，因为它使用双精度浮点值进行操作。使用64位，可以精确表示1.0e + 11。

64位代码产生不正确的结果，因为它使用的是单精度浮点值。仅可使用32位，因此无法精确表示1.0e + 11。

如果使用double类型开头，则不会出现此问题：

double val = 1000.0;
val = val * val;
return (ulong)(val * 100000.0);

这样可以确保在所有体系结构上获得正确的结果，而无需像问题中所建议的那样丑陋，不可移植的位操作hack。 （由于它不能解决问题的根源，即您不能用32位单精度float直接表示您想要的结果，因此仍然不能保证正确的结果。）

即使必须将输入作为单精度float，也要立即将其转换为double，然后在双精度空间中进行所有后续的算术运算。那仍然可以解决这个问题，因为初始值1000可以精确地表示为float。