解释了一种将双精度浮点数舍入为32位整数的快速方法

在阅读 Lua的源代码时,我注意到 Lua 使用一个宏将 double值舍入到32位 int值。宏在 Llimits.h头文件中定义如下:

union i_cast {double d; int i[2]};
#define double2int(i, d, t) \
{volatile union i_cast u; u.d = (d) + 6755399441055744.0; \
(i) = (t)u.i[ENDIANLOC];}

在这里,ENDIANLOC是根据 Endianness:0定义的,对于 little endian 是1,对于 big endian 体系结构是1; Lua 小心地处理 endianness。t参数被类似于 intunsigned int的整数类型所替代。

我做了一些研究,发现有一种更简单的宏格式,使用相同的技术:

#define double2int(i, d) \
{double t = ((d) + 6755399441055744.0); i = *((int *)(&t));}

或者,以 C + + 的风格:

inline int double2int(double d)
{
d += 6755399441055744.0;
return reinterpret_cast<int&>(d);
}

这个技巧可以在任何使用 IEEE 754的机器上使用(这意味着现在几乎所有的机器)。它既适用于正数,也适用于负数,四舍五入遵循 银行家守则。(这并不奇怪,因为它遵循 IEEE 754。)

我写了一个小程序来测试它:

int main()
{
double d = -12345678.9;
int i;
double2int(i, d)
printf("%d\n", i);
return 0;
}

正如预期的那样,它输出 -12345679

我想详细了解这个棘手的宏是如何工作的。神奇数字 6755399441055744.0实际上是251 + 252,或者1.5 × 252,而1.5的二进制数可以表示为1.1。当任何32位整数被加到这个神奇的数字ー

我迷路了

更新

  1. 正如@Mystilic 所指出的,这种方法并不局限于32位的 int,它也可以扩展为64位的 int,只要数字在252的范围内。(尽管宏需要进行一些修改。)

  2. 有些材料说,这种方法不能用于 直接3D

  3. 在使用 Microsoft x86汇编程序时,有一个用汇编代码编写的更快的宏(以下内容也是从 Lua source 中提取的) :

     #define double2int(i,n)  __asm {__asm fld n   __asm fistp i}
    
  4. 对于单精度数,有一个类似的幻数: 1.5 × 223

17046 次浏览

double浮点类型的值表示如下:

double representation

它可以看作是两个32位的整数; 现在,所有版本的代码中的 int(假设它是一个32位的 int)就是图中右边的那个,所以你最后要做的就是取最低的32位尾数。


现在,来看看这个神奇的数字: 6755399441055744是251 + 252; 加上这个数字,double就会进入252和253之间的“甜蜜范围”,就像 由维基百科解释一样,它有一个有趣的属性:

在252 = 4,503,599,627,370,496和253 = 9,007,199,254,740,992之间,可表示的数正好是整数。

这是因为尾数是52位宽。

另一个关于加入251 + 252的有趣事实是,它只影响尾数的最高位ーー这两个位无论如何都会被丢弃,因为我们只取其最低的32位。


最后但并非最不重要的: 标志。

IEEE 754浮点数使用数量级和符号表示,而“正常”计算机上的整数使用2的补数算法; 这里是如何处理的?

我们只讨论了正整数; 现在假设我们处理的是一个32位 int表示的范围内的负数,所以小于(- 231 + 1) ; 称之为-a。这样一个数字显然是通过加上这个神奇的数字而得到的,得到的结果是252 + 251 + (- a)。

现在,如果我们把尾数解释为2的补数表示,我们会得到什么?它必须是2的补数和(252 + 251)和(- a)的结果。同样,第一项只影响上两位,位0-50中剩下的是(- a)的2的补数表示(同样,减去上两位)。

因为将2的补数减少到一个更小的宽度仅仅是通过去掉左边的多余位来完成的,所以取较低的32位就可以得到32位的补数算法中的正确(- a)。

下面是上述 Lua 技巧的一个更简单的实现:

/**
* Round to the nearest integer.
* for tie-breaks: round half to even (bankers' rounding)
* Only works for inputs in the range: [-2^51, 2^51]
*/
inline double rint(double d)
{
double x = 6755399441055744.0;  // 2^51 + 2^52
return d + x - x;
}

这个技巧适用于绝对值 < 2 ^ 51的数字。

这是一个测试它的小程序: 网址: ideone.com

#include <cstdio>


int main()
{
// round to nearest integer
printf("%.1f, %.1f\n", rint(-12345678.3), rint(-12345678.9));


// test tie-breaking rule
printf("%.1f, %.1f, %.1f, %.1f\n", rint(-24.5), rint(-23.5), rint(23.5), rint(24.5));
return 0;
}


// output:
// -12345678.0, -12345679.0
// -24.0, -24.0, 24.0, 24.0

这种“技巧”来自较老的 x86处理器,它使用8087指令/接口作为浮点数。在这些机器上,有一个将浮点数转换为整数“拳头”的指令,但它使用的是当前的 fp 舍入模式。遗憾的是,C 规范要求 fp-> int 转换截断为零,而所有其他 fp 操作都是四舍五入到最接近的值,因此
Fp-> int 转换需要首先改变 fp 舍入模式,然后做一个拳头,然后恢复 fp 舍入模式。

现在,在最初的8086/8087上,这还不算太糟糕,但是在后来的处理器上,开始变得超标量和乱序执行,改变 fp 舍入模式通常会序列化 CPU 核心,并且非常昂贵。所以在像 Pentium-III 或者 Pentium-IV 这样的 CPU 上,这个总成本是相当高的——一个普通的 fp-> int 转换要比这个 add + store + load 技巧贵10倍甚至更多。

但是,在 x86-64上,浮点运算是使用 xmm 指令完成的,并且转换的成本也是如此
Fp-> int 非常小,所以这个“优化”可能比正常的转换慢。

如果这有助于形象化,那 Lua 的神奇价值

  (2^52+2^51, or base2 of 110 then [50 zeros]

巫术

  0x  0018 0000 0000 0000 (18e12)

八进制

  0 300 00000 00000 00000 ( 3e17)