5

阅读The Tricks of the 3D Game Programming Gurus时,我遇到了这个用内联汇编编写的排序函数:

inline float FastSqrt(float Value)
{
    float Result;

    _asm
    {
        mov eax, Value
        sub eax, 0x3F800000
        sar eax, 1
        add eax, 0x3F800000
        mov Result, eax
    }

    return(Result);
}

这是实际平方根的近似值,但准确度足以满足我的需要。

这实际上是如何工作的?这是什么神奇的0x3F800000价值?我们如何通过减法、旋转和加法来获得平方根?

下面是它在 C/C++ 代码中的样子:

inline float FastSqrt_C(float Value)
{
    float Result;

    long Magic = *((long *)&Value);
    Magic -= 0x3F800000;
    Magic >>= 1;
    Magic += 0x3F800000;
    Result = *((float *)&Magic);

    return(Result);
}
4

4 回答 4

10

很多人都指出,0x3f800000是 的表示1.0。虽然这是真的,但它与计算的方式无关。要理解它,您需要知道非负浮点数是如何存储的。f = (1+m)*2^x,0 <= m < 1并且m是尾数,x指数。另请注意,它x是带有偏差存储的,因此二进制文件中的实际内容是x+127. 32 位值由符号位(在我们的例子中为零)、8 位指数存储x+127和最后 23 位尾数组成m。(参见维基百科文章)。

应用一些基本的数学,

sqrt(f) = sqrt((1+m)*2^x)
        = sqrt(1+m)*sqrt(2^x)
        = sqrt(1+m)*2^(x/2)

所以,作为一个粗略的近似,我们需要将指数减半,但由于偏差,我们不能只做x/2我们需要(x-127)/2 + 127的。这种127移位到合适的位位置是神奇的0x3f800000

除以 2 是通过右移一位来实现的。由于这对整个浮点数起作用,因此它对尾数也有副作用。

首先,假设原始指数是偶数。然后移出的最低有效位为零。因此,尾数也减半,所以我们最终得到的是:sqrt(f) = (1+m/2)*2^(x/2). 我们得到了正确的指数,但尾数(1+m/2)不是sqrt(1+m). 最大相对误差是(1.5 - sqrt(2))/sqrt(2) ~ 6%m几乎1意义f接近但小于 的奇数幂时发生的2。举个例子f=7.99。该公式为我们提供了 about2.998而不是2.827which 确实有 的错误6%

现在,如果指数是奇数,那么最低有效位将是1,当移入尾数时将导致增加一半。因此,我们得到sqrt(f) = (1.5+m/2)*2^((x-1)/2)。最大的错误实际上是 when m=0,那将是(1.5/sqrt(2)-sqrt(1))/sqrt(1)which 再次出现6%。这发生在从上方接近 2 的奇数次幂的数字上。

如果输入值恰好接近 2 的奇数幂,则这两种情况相结合意味着最差的不准确度约为 6%。对于偶数次方,结果是准确的。

于 2017-01-22T00:39:35.127 回答
0

浮点数 f = (1 + m)* [2^(e+127)],其中 m 是尾数部分,e 是指数部分。

因此: sqrt(f) = (f)^(1/2) = ((1 + m)* [2^(e+127)] )^(1/2)

-> ((1 + m)* [2^(e+127)] )^(1/2) = (1 + m)^(1/2) * 2^((e + 127)/2)

在指数部分,2^((e + 127)/2):

2^((e + 127)/2) = 2^( (e-127/2) + 127)

因此,在浮动表示中,它是 (e - 0x3F800000) /2 + 0x3F800000

在尾数部分,(1 + m)^(1/2):

从二项式级数公式, (1 + x)^r = 1 + r x + (r (r - 1)/2)*(x^2) +....

因此,(1 + m)^(1/2) 等于 (1 + m/2 - (m^2)/8 + ...) 它大约等于 1 + m/2(一阶的典型近似值)因此, 尾数部分应除以 2。

但是,尾数和指数组合为 A 数,右移将指数和尾数 BOTH 分开。

要评估错误,您可以考虑二项式系列的第二项 - (m^2)/8。

因为 m 总是小于 1,我将 m 替换为 0.9999 (0.5 + 0.25 + 0.125 + ...)

(m^2)/8 = 0.12497500125,这是最坏的情况。

于 2018-04-20T19:34:27.233 回答
0

浮点数中的 0x3F800000 为 1。这是因为浮点数的存储方式。你可以在https://gregstoll.dyndns.org/~gregstoll/floattohex/看到一个视觉表示。

这是一个很好的近似值,我相信 sqrt。这起源于反 sqrt 游戏 Quake ( https://en.wikipedia.org/wiki/Fast_inverse_square_root#Aliasing_from_floating_point_to_integer_and_back )。

于 2017-01-21T23:09:41.670 回答
0

下面是这个机制的一个例子:

FastSqrt(4.0) == 2.0

4.0 to hex -> 0x40800000
0x40800000 - 0x3f800000 = 0x1000000
0x1000000 to binary -> 00000001 00000000 00000000 00000000
shift toward the lsb (sar) -> 00000000 10000000 00000000 00000000
00000000 10000000 00000000 00000000 back to hex -> 0x00800000
0x00800000 + 0x3f800000 = 0x40000000
0x40000000 to dec -> 2.0
于 2017-01-21T23:25:38.287 回答