引子:浮点数
计算机中的浮点数被存放为三个部分,分别放置在一个32位存储的三个段中,如图:
| S | E | M |
|---|---|---|
| 0 | 0000 0000 | 23 个 0 |
并且这三段内存表达的数为:S,E,M.
计算机计算对应的浮点数的方法为:
\[(-1)^S \times {1 + M \over 2^{23}}\times 2^{E-127}\]我们可以举个例子, 3.14在计算机中存储为:
其中S=0,E=128,M=4781507.
代入上公式计算:
\[{(-1)}^{0} \times {1 + {4781507 \over 2^{23}}} \times 2^{128-127} = 3.1400001049041748046875\]更快的平方根倒数计算
计算机图形学中，常常会使用一些数学方法，在可接受的精度范围内，加速一些经常进行的数学运算。比如对三维空间进行二维投影时，这就需要对这两个几何体面上的大量点，进行距离\(\sqrt{a^2 + b^2 + c^2}\)倒数的运算投影到摄影机空间，使得物体在二维平面也可以符合近大远小的透视关系：
\(\alpha = {a \over \sqrt{a^2+b^2+c^2}}(1)\) \(\beta = {b \over \sqrt{a^2+b^2+c^2}}(2)\) \(\zeta = {c \over \sqrt{a^2+b^2+c^2}}(3)\)
在上面的过程中，形如下式的函数将被大量计算：
\[y = {1 \over \sqrt{a^2+b^2+c^2}}(4)\]而且这些计算都是浮点数计算，每秒成千上万次类似的浮点运算，会导致对这个过程进行优化将获得很大的性能提升。
首先，我们对式1进行化简：
\(a^2+b^2+c^2=x\) \(d = {1 \over \sqrt{x}}\)
两端同时取对数：
\[\log_{2}{d} = - {1 \over 2 } \log_2{x}\]让我们回顾一下我们的目标,我们现在拥有的是一个浮点表示的\(I_x\)(\(S_x\),\(E_x\),\(M_x\)),我们需要拿到一个浮点表示的\(I_y\).
另，因为知道浮点数的近似表达公式为：
\[(1+m)\times 2^e\]代入上式：
\[\log_2({(1+m_y)\times 2^{e_y}}) = - {1 \over 2 } \log_2({(1+m_x)\times 2^{e_x}})\]对\(\log_2({1+m_y})\)进行线性化：
\[\log_2{(1+m_y)} \approx m+\alpha\]代入原式： \(m_y+\alpha + e_y = - {1 \over 2 } (m_x+\alpha) + e_x\) 因为\(m = {M \over 2^{23}}\),\(e = {E-127}\)
得: \({M_y \over 2^{23}}+\alpha + {E_y-127} = - {1 \over 2 } ({M_x \over 2^{23}}+\alpha) + {E_x-127}\)
再把整个表达式乘以\(2^{23}\) ,得到：
\[M_y+E_y 2^{23} = -{1\over2}(M_x+E_x 2^{23})-{3\over2}(\alpha-127)2^{23}\]左侧就是浮点表示对应的整形数的值\(I_y\),右侧为浮点表示对应的整形数的值\(I_x\)
最终段简化公式为: \(I_y \approx R - {1\over 2} I_x\)
基于这种算法，早在1999年就在游戏《雷神之锤》被应用，源码函数如下：
float Q_rsqrt( float number )
{
long i;
float x2, y;
const float threehalfs = 1.5F;
x2 = number * 0.5F;
y = number;
i = * ( long * ) &y; // 浮点表示法对的整型值I_x
i = 0x5f3759df - ( i >> 1 ); // 获得I_y
y = * ( float * ) &i; // 获得浮点y
y = y * ( threehalfs - ( x2 * y * y ) ); // 1st Newton iteration
// y = y * ( threehalfs - ( x2 * y * y ) ); // 2nd Newton iteration, this can be removed
return y;
}
最后，在自己的PC上实验一下，感受一下知识的力量吧！
运行一次运算的速度从73个tick,提升到了5tick