我们先从一个简单的问题入手。

在二维平面中有一个向量，它代表了我们生活在二维空间中的主角当下运动的状态。向量的方向即运动方向，向量的大小即速度大小。如果我们不在于主角的移动速度，只是想得到他的运动方向，也就是我们常说的方向向量，我们该如何计算？

很简单，我们利用勾股定理 a^2 + b^2 = c^2 即可求得向量长度的二次方，随后我们再对结果做平方根运算，得到的值为向量的大小，最后我们再用原向量除以向量大小，得到的结果即为方向向量，这个计算方法通常被叫做归一化（normalize），得到的向量我们也可以叫做单位向量。

方向向量不仅是在游戏中的物理模拟非常重要，对于渲染来说也是一样的。它能用来指明光线与三角面的方向，而方向则决定了光照程度、材质质感等关键效果的计算结果。

我们用简单的式子来描述一下这个计算过程：V / ||V||，其中 V 为任意向量，||V|| 为向量大小。也就是说，我们只要得到一个向量大小的倒数，即 1 / ||V||，就可以快速计算出方向向量。而 ||V|| 又涉及到平方根操作，所以最终我们要计算的其实是平方根的倒数。

当然我们可以把倒数与平方根运算分开操作，但是在游戏开发的蛮荒时代，平方根运算并不容易，倒数的操作也不是非常快速的简单运算。被机能困在旧时代的开发者们渴望得到一个快速算法，能够解决游戏开发中涉及到的大量的平方根倒数运算。

那么这个问题该怎么解决呢？

平方根该怎么计算？当然你不能使用计算器。

举个例子来说，我们就求 √2 吧。

我们可以做如下推导：

x = √2

=> x^2 = 2

=> 0 = x^2 - 2

所以我们可知，函数 f(x) = x^2 - 2 与 x 轴的交点即为 √2 的值，函数图像如下：

观察图像，我们可以推测这个值介于 1 与 1.5 之间。

上面的函数与方程只能让我们得到 √2 这种形式的数值表达形式，那么我们怎么求得 √2 具体值呢？

如上图所示，我们画一条线，这条线与函数 f(x) 的曲线相切于点（2,2），并且这条线与 x 轴的交点十分的接近我们要求得的那个点，这是一个好现象。也就是说，我们只要找出这条线的函数形式，并求得这条线与 x 轴的交点，我们就能的到一个 √2 的近似值。

那么我们怎么求得这条线的函数呢？

首先我们知道这条线经过点（2,2），那么如果我们还能得到线的斜率，通过点斜式即可求得函数。

因为这条线本质上是函数 f(x) 在点（2,2）切线，所以它的斜率就是 f(x) 在这一点的导数。

所谓的函数的导数，它表示这个函数在某一点上变化的趋势与程度。这一概念其实与经过切点的切线是保持一致的。切线的倾斜程度（即斜率）不就是代表了这条线上的值的变化程度吗？

我们用 f '(x) 来表示函数 f(x) 的导函数，对于 f(x) = x^2 - 2 来说，f '(x) = 2x。

其实对于类似 f(x) = ax^n + b 这种形式的函数来说，导函数的形式非常简单，就是 x 的指数减一，得到的值与 x 的系数相乘，偏移量 b 因为是常数，所以会被省略。按照这个规则，f(x) = ax^n + b 的导数就是 f ‘(x) = an * x^(n-1)，所以 f(x) = x^2 - 2 导数即为 f '(x) = 2x。

我们将点（2,2）带入函数 f '(x) = 2x 求得此处的导数为 4，即经过此处的切线的斜率为 4。

现在切线的斜率与切线经过的点都有了，带入点斜式我们可得：

y - b = k(x - a)；a，b 为函数图像经过的点

y - 2 = 4(x - 2);

y = 4x - 6;

即切线函数为 g(x) = 4x - 6，并计算 0 = 4x - 6，所以切线与 x 轴的交点为（1.5, 0）。那么我们要求的 √2 近似值即为1.5。

1.5 这个值明显并不精确，它有点过于“近似”了。从图上我们也能观察到它与真正的 √2 的值还是有一定差距的。但是我们已经掌握了求近似值的方法，那么我们是否可以重复上述过程，来进一步接近正确值呢？

答案当然是肯定的。我们已经得到了一个接近 √2 的值 1.5，那我们可以尝试以这个值为基础，再做一次上面的运算。

首先我们将 1.5 带入 f(x)，求得一个新的切点：

f(x) = x^2 - 2 = 1.5^2 - 2 = 0.25

新切点即为（1.5, 0.25）。

再导入导函数求切线斜率：

f '(x) = 2x = 2 * 1.5 = 3

最后带入点斜式，得到切线函数：

y - b = k(x - a)；

y - 0.25 = 3(x - 1.5);

y = 3x - 4.25;

3x = 4.25;

x ≈ 1.41666

我们也可以从图像中观察到这个结果：

我们将这几条线与 x 轴相交的部分放大观察：

正如我们设想的那样，运算的得到的结果作为下一次运算的初始值，递归运算下去，我们就能不断的接近那个我们想要的值。

我们再从全局的角度观察一下我们的计算的流程。我们有两个关键函数，f(x) 和它的导函数 f '(x)，用它们我们可以求得切点与斜率，再带入点斜式求得切线函数，最后得到切线与 x 轴的交点，即为我们想要的值。从这个过程来看，我们好像可以总结出一个更简单的公式，来表示这个计算过程，我们来尝试一下。

点斜式 y - b = k(x - a)，因为我们最后要求得的值为函数与 x 轴交点中的横坐标值，即点斜式中的 a 值；又因为 x 轴交点的纵坐标必为 0，所以我们将点斜式变化为：

a = x - y/k

其中 x 与 y 组成我们当前已知的切点（例如上述第一次计算过程中的点 2,2），y 由 f(x) 计算得到，所以 y = f(x)；k 为斜率，由导函数 f ’(x) 计算的到，所以 k = f ’(x)。

我们还知道 a 为我们要求得的最新的近似值，我们设它为Xn+1；x 为我们当前已经求得的近似值，我们设它为Xn。

将上述各个值的转换导入式子，我们可得：

Xn+1 = Xn - f(Xn) / f '(Xn)；

带入具体的实例测试一下？

第一个切点我们选（2,2），并且我们有 f(x) = x^2 - 2，f '(x) = 2x，带入公式即为：

Xn+1 = 2 - (2^2 - 2)/(2*2) = 2 - 1/2 = 1.5。

恭喜！我们又一次发明了牛顿的咒语：牛顿迭代法。

在这个过程中，我们不难发现牛顿灵活的运用了近似值这个概念，以近似值的方式，不断逼近精确值。与其说是算法， 不如说是一种策略。

对于这条所谓的咒语来说，它不仅简单，而且还特别适合用代码实现，我们可以用 for 或 while 为循环计算提供动力，直到计算出我们想要的精度，一次平方根计算就完成了。

但是，我们开头也提到过，这是一个蛮荒的时代，计算机的算力有限，而我们要完成的却有很多，我们很难使用一个依赖循环次数得到确切值的算法。循环就像一个吞噬算力的猛兽，它很快就会让机器的计算能力变得捉襟见肘。

如果只需要念一次咒语就好了。

没错，只进行一次牛顿迭代法。不是不可能，只不过需要我们首次提供的 Xn 就足够精确罢了，足够精确的 Xn 得到更精确的 Xn+1。

这好像是个死循环。

如果卡马克的魔法没有诞生的话。

假设我们有一个数字，30000，用科学计数法表示它即为，3 * 10^4。不难发现 10 的指数其实就是数值中 0 的个数。

十进制的规则深入我们的骨髓，其中的道理也无需多言。

那么二进制呢？

假设我们有一个二进制数 10000（十进制值为 16），我们类比上面的十进制规则，那么二进制 10000 的科学计数法表示形式即为 1 * 2^4。此处与十进制不同的是，二进制是以 2 为底的指数形式。你可以仔细品味一下乘法的本质，再从十进制类比到二进制。

这一切都说得通，而且非常自然。

我们再进一步。

十进制数 31000，用科学计数法表示即为 3.1 * 10^4，观察这个表达形式，我们不难发现，10 的指数本质其实是在表明小数点的移动位数。

那么类比到二进制中，二进制数 11000 （十进制值为 24）的科学计数法表达形式即为 1.1 * 2^4，符合我们上面观察到的规律。

这就是二进制的位移规则。

如上图所示，浮点数在计算机中的二进制存储模式被分为三个部分，首位是符号位，用来标识数字的正负号；后八位是偏执指数位，用来表示在科学计数法下小数点偏移的位数；最后的二十三为则表示在科学计数法下，小数点后面的数字排列，这部分融合了原数字的整数与小数两个部分。

我们以一个具体的数字来举例，二进制 111 的科学计数法表达形式为：

1.11 * 2^2

以浮点数编码规则编码后，在计算机中的二进制形式为：

其中指数位的计算规则是，小数点偏移量 + 127。以 111 为例，科学计数法表达后小数点向左偏移两位，所以指数位为 2 + 127 = 129，二进制表达即为 10000001。

至于尾数位，因为二进制科学计数法的规则就是小数点会停留在首个不为 0 的位数后面，所以我们只需要去掉科学计数法中的“1.”，剩下的部分就是尾数位的本体。

以上就是我们对浮点数编码规则的理解。

关于浮点数，我们前面所讨论的一直都是科学计数法到浮点数二进制形式的转换，那我们就忍不住好奇，将这个过程反过来呢？浮点数二进制形式如何转化为科学计数法呢？

很简单，对于偏置指数位来说，我们只需要用偏执指数位 - 127，结果就是科学计数法中的小数点偏移量，就上述中的例子来说：

100000001 - 01111111

= 129 - 127

= 2

对于尾数位来说，反向推到科学计数法的表达形式时，需要在尾数位前面加上“1.”，就上述中的例子来说，也就是在 “11000……” 前面加上 “1.”，即 1.11，0 被省略。

这还不够，根据上面的文字描述，我们可以试着总结出一条公式。

我们设尾数位为 M，因为尾数位共 23 位，尾数位反推至科学计数法时，小数点要出现在尾数位形式之前，那么根据我们前面提到的二进制位移规则，尾数位小数点向左移动 23 位，即 M / 2^23，最后别忘了，我们还要加上一个“1.”，所以最终公式应该为：1 + M / 2^23；

接着我们设偏置指数位为E，则科学计数法中小数点偏移量的计算为：E - 127。

结合尾数位部分的公式，我们可以得到（这里我们省略的符号位，因为加入符号位很简单，而且后续我们的讨论也不需要符号位）：

科学计数法表达形式 = (1 + M / 2^23) * 2^(E - 127)

对比一下例子中的科学计数法表达：1.11 * 2^2

那么这个公式有什么用呢？我们先来随意摆弄一下它看看会发生什么。

因为公式很明显被乘法分为两个部分，并且其中一部分还是 2^(E - 127) 这种 2 的指数形式，那么我们不防试试对数？

对数的运算法则很有意思：

它能将乘法的两部分，变为加法（第一条）；还能将指数转变为乘法（第四条）。这两条都存在于我们的公式之中，那我们就来试试。

设 a = (1 + M / 2^23) * 2^(E - 127)，我们对两边都进行以 2 为底的对数运算，可得：

㏒₂a = ㏒₂((1 + M / 2^23) * 2^(E - 127))

㏒₂a = ㏒₂(1 + M / 2^23) + ㏒₂2^(E - 127)

㏒₂a = ㏒₂(1 + M / 2^23) + ㏒₂2 * E - 127

㏒₂a = ㏒₂(1 + M / 2^23) + E - 127

唉？有趣，2 的指数形式消失了，只剩下 E - 127，从某种角度来说，整体运算变得简单了。那 ㏒₂(1 + M / 2^23) 部分呢？

我们知道 M 作为尾数位的二进制形式，具体二进制取决于我们浮点数的具体数值，所以这里可以把它看做是一个变量，M / 2^23 自然也是随着 M 的变化而变化的，我们不防暂时先将他们看做一个整体 x。于是 ㏒₂(1 + M / 2^23) 就变成了 ㏒₂(1 + x)。这就又简化了我们的式子，因为 ㏒₂x 本身就是一个很常见的对数函数，图像如下：

而 ㏒₂(1 + x) 只需对图中图像稍作位移：

这还没完，因为我们知道 M / 2^23 必然大于 0 小于 1，所以 ㏒₂(1 + x) 中，x 的取值范围也在 0 至 1 之间，那么我们放大这个图像的细节到 x 的 0 至 1 区间：

从图像中可以看到，这条曲线在 0 到 1 之间的部分，其实近乎于直线，再回想一下我们从牛顿的咒语中学到了什么？

没错，近似值。在我们实际的工程应用中，很多数值都不是十分精确的，或者说对数值的精确度要求是变化的。所以既然它接近于一条支线，那我们为什么不就把它当做一条支线呢？

因为曲线经过点原点与点（1，1），所以我们可以把它看做 y = x：

从图上可以观察到，y = x 与 ㏒₂(1 + x) 在 0 至 1 区间的值，真的非常相近。

那么，因为：

㏒₂(1 + M / 2^23) ≈ M / 2^23

所以：

㏒₂a = ㏒₂(1 + M / 2^23) + E - 127

≈ M / 2^23 + E - 127

到这里还不够，我们再稍作变形，将 M / 2^23 + E - 127 中 M / 2^23 + E 提取 1 / 2^23 可得：

1 / 2^23 * (M + E * 2^23) - 127

最后我们终于来到了这里。

观察式子中的 M + E * 2^23，你看它像什么？

还是 1.11 * 2^2 的二进制形式为例，它的二进制形式为：

其中，M = 11000000000000000000000；

E = 10000001；

那么，E * 2^23 根据二进制位移规则，就是将小数点向右移动 23 位，就是说要在 E 后面添加 23 个 0：

E * 2^23 = 10000001 00000000000000000000000；

于是可得：

M + E * 2^23

=

11000000000000000000000

+

10000001 00000000000000000000000

=

10000001 11000000000000000000000

原来，M + E * 2^23 就是这个 M 和 E 所代表的那个数字的二进制形式。

那么根据我们的推导过程：

㏒₂a = ㏒₂(1 + M / 2^23) + E - 127

≈ M / 2^23 + E - 127

= 1 / 2^23 * (M + E * 2^23) - 127

我们可以得到一个结论：

对于数字 a，设 A 为浮点数 a 二进制形式，那么：

㏒₂a ≈ 1 / 2^23 * A - 127

也就是说我们得到了一个快速计算对数的算法，只要利用浮点数二进制就可以快速计算。

还记得我们最开始的目的吗？我们是为了计算平方根倒数，那这个快速计算对数的方法和平方根倒数计算又有什么关系呢？别着急，我们接着来。

我们设 a 的平方根倒数，即 1/√a ，为 y，那么我们有：

㏒₂y = ㏒₂a^(-1/2)

= -1/2 * ㏒₂a = -1/2 * ((1 / 2^23) * A - 127)

也就是说：

㏒₂y = -1/2 * ((1 / 2^23) * A - 127)

y 虽然是 1/√a 的符号表现，但是，y 最终一定有一个值对不对？它既然是一个值，那就必然可以表达为浮点数的二进制形式，我们就设浮点数 y 的二进制形式为 Y。这个时候我们就可以 ㏒₂y 也转化为二进制快速计算对数的形式：

㏒₂y = (1 / 2^23) * Y - 127

而 Y 是浮点数 y 的二进制形式，也就代表了 y 的值，即 1/√a，也就是我们想要求得的 a 的平方根倒数。

终于，一切都说的通了：

(1 / 2^23) * Y - 127 = -1/2 * ((1 / 2^23) * A - 127)

Y = 381 * 2^22 - 1/2 * A

Y = 1598029824 - 1/2 * A

如果我们用符号 “>>” 表示二进制向右的移动，那么考虑到我们前面所提到的二进制位移规则， 1/2 * A 就等同于 A >> 1，它表示将二进制 A 向右移动一位；如果我们将 1598029824 转化为十六进制，那么最终的式子就变为：

Y = 5F400000 - A >> 1

根据最终我们得到得的这个式子，我们得出一个结论：

将浮点数 a 的二进制形式所代表的值乘以0.5，再与数值 1598029824 作差，结果的绝对值即为 A 的平方根倒数。

这就是卡马克在雷神之锤 3 中所写下的快速平方根倒数算法。它巧妙的将浮点数的编码规则与平方根倒数的计算相结合，一切都显得如此自然又仿佛是巧合。

我们看一下具体的代码实现：

其中第七行代码，先将 y 值转化为 long 型，因为整型的编码就是简单的二进制与十进制的转换，不涉及类似于浮点数编码规则对二进制数的特别解读，所以转化为 long 型，就相当于是原封不动的保存了浮点数 y 的二进制。

第八行代码与我们上面推导出的 Y = 5F400000 - A >> 1 相似，不同的点在于减数，代码中的减数 0x5f3759df 看起来比我们的 5F400000 要更加精确。这也是事实，还记得我们上面提到过的 y = x 与 ㏒₂(1 + x) 在 0 至 1 区间的图像对比吗？

我们能看到，虽然两个函数的图像在 0 至 1 区间内相似，但是归根结底还是不同的，然而我们的所有推导，都是建立在他们非常相似的基础上，它们越相似，我们的误差就越小，算法的精度就越高。如果我们将 y = x 变为 y = x + 0.05，则图像如下：

这个 0.05 我们可以叫做修正值。这个修正值在我们的计算过程中，最终会以常数的形式与 5F400000 相加，实际代码中的 0x5f3759df 就是这么来的。这个修正值可以按照一定的策略来做尝试，比如我们可以使用二分法，快速缩减尝试范围，最红确定一个相对精确的修正值。

最后第十一行代码，使用了牛顿迭代法进一精确我们得到的数值。具体代码的逻辑为：

Y = 1/√x

=> Y^2 = 1/x

=> 1/Y^2 - x = 0

即，我们有函数 f(Y) = 1/Y^2 - x，此函数当 x 为 0 时的图像如下：

由图可知，当 x = 0 时，函数图像无限接近 x 轴，但并不与 x 轴相交；如果 x >= 0，我们可以想象，函数图像会向下位移，之后与 x 轴相交处，即是我们要求得的 x 的平方根倒数。

因为 x 是我们想要求得平方根倒数的目标数值，所以对于每次计算来说，这都是一个常数，那么根据牛顿迭代法可得：

已知，f '(Y) = -2/Y^3

所以，

Yn+1 = Yn - f(Yn) / f '(Yn)

=> Yn - (1/Yn^2 - x) / (-2/Y^3)

=> Yn * (3/2 - x/2 * Yn^2)

此式子就对应代码中：

y = y * ( threeHalfs - ( x2 * y * y ) );

这与我们在牛顿的魔法一节最后提到的，以一个相对精确值得到更精确值的目标相对应。

最后值得一说的是，其实这个快速计算平方根倒数的算法并非出自卡马克之手，不然卡马克的注释也不会是 WTF？了。

对于这个算法最初的发明者，我们可以参考 wiki 上的说法：