这个 sqrt 近似内联汇编函数是如何工作的？

Question

通读 3D 游戏编程大师的技巧，我发现了这个用内联汇编编写的排序函数：

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);
}

Answer 1

float中的0x3F800000为1。这是因为float的存储方式。您可以在 https://gregstoll.dyndns.org/~gregstoll/floattohex/.

处看到直观表示

这是一个很好的近似值，我相信 sqrt。其起源来自游戏 Quake for inverse sqrt (https://en.wikipedia.org/wiki/Fast_inverse_square_root#Aliasing_from_floating_point_to_integer_and_back).

Answer 2

下面是一个实际运作机制的示例：

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

Answer 3

很多人指出0x3f800000是1.0的表示。虽然这是事实，但它与计算的工作方式无关。要理解它，您需要知道 non-negative 浮点数是如何存储的。 f = (1+m)*2^x，0 <= m < 1 和 m 是尾数，x 是指数。另请注意，x 存储时带有偏差，因此二进制文件中的实际内容是 x+127。 32 位值由符号位（在我们的例子中为零）和随后的 8 位指数存储 x+127 以及最后的 23 位尾数 m 组成。（参见wikipedia article）。

应用一些基础数学，

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)。如果 m 几乎 1 意味着 f 接近，但小于 2 的奇次方，则此最大相对误差为 (1.5 - sqrt(2))/sqrt(2) ~ 6%。以 f=7.99 为例。该公式为我们提供了 2.998 而不是 2.827，它确实有 6%.

的错误

现在，如果指数是奇数，那么最低有效位将是 1 并且当移入尾数时将导致增加一半。因此，我们得到 sqrt(f) = (1.5+m/2)*2^((x-1)/2)。这个的最大错误实际上是在 m=0 时，那将是 (1.5/sqrt(2)-sqrt(1))/sqrt(1)，这又是在 6% 附近。对于从上面看接近 2 的奇次方的数字，会出现这种情况。

如果输入值恰好接近 2 的奇次方，这两种情况加在一起意味着最严重的不准确度约为 6%。对于 2 的偶次方，结果是准确的。

Answer 4

浮点数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 + rx + (r(r - 1)/2)*(x^2) + ....

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

但是，尾数和指数合并为一个数，右移除指数和尾数BOTH。

要评估误差，您可以考虑二项式级数的第二项 - (m^2)/8。

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

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

这个 sqrt 近似内联汇编函数是如何工作的？

How does this sqrt approximation inline assembly function work?

optimization

x86

assembly

sar

square-root