第二周记

前言：

10.20

今天体测跑了一公里，拿下了三分三十的成绩，有时感觉自己当时走体育会不会比现在混得好（）

下午心理课颓废整场，其实越学就越想学（至少对于我是这样），主要我现在真不知道怎么学，没有老师带领，做题像大海捞针，做题也不发答案，也不知道是我的问题还是学校的问题，tnd高中也是被折磨成斯德哥尔摩了。

听说导员当年四级裸考过线（高考英语118），感觉又有点放松了，但是心里还是毛毛的。

不会高数。不会高数。不会高数。不会高数。不会高数。不会高数。不会高数。

晚上终于更新了“快速平方根取倒数算法”。还是很震撼的，这个算法刷新了我对“卡常”的认知：之前总是把常数优化和算法优化区别来看，觉得“卡常”是相当不优雅的行为，实际上充分调用“计算机”的所有资源本身就是一种艺术，内核级的优化绝对不是几个常规操作就能解决的。

再说一点：虽然上文为“卡常”“正名”，但是在高级语言中不要老是觉得“卡常”有多万能，别想着直接加个循环展开就能 \(n^2\) 过百万，写编译器的一定比你更懂优化。

有用的“卡常”也就是快读快输，递归改迭代，数组模拟这些吧，当然还有指令集之后可能会学一下？比如这位大佬就很厉害。

淦，键盘没电了，睡觉睡觉。

10.21

学习笔记：

一些杂项

快速平方根取倒数算法

很久以前就想写一下这个快速平方根取倒数算法了，这个好像还很有历史“渊源”来着，据说给当时“雷神之锤”的开发带来了巨大优化。

首先我们思考一下，平方根应该用哪种算法求，显然这是一个只能逼近的值，因此我们想到了“牛顿迭代”算法：

\[x_n=x_{n-1}-\frac{f(x_{n-1})}{f'(x_{n-1})} \]

考虑函数 \(f(x)=\frac{1}{x^2}-a\) ，先选定一个初始值 \(x_0\) 然后：

\[\begin{aligned} x_n&=x_{n-1}-\frac{ax_{n-1}^3-x_{n-1}}{2}\\ &=\frac{3x_{n-1}-ax_{n-1}^3}{2} \end{aligned} \]

迭代求解。

显然复杂度是与求解的精度相关的，迭代次数和答案的精度正相关。那么，能否选定某个恰当的 \(x_0\)，使得迭代次数得到大幅的减少呢？答案是肯定的。

早期游戏工程师为了简化复杂度（众所周知早期算法工程师常从底层最大限度的优化算法，极致压榨 \(PC\) 机效率），从浮点数的性质入手：

\(32\) 位浮点数的 IEEE754 存储形式是这样的:

• \(S\) 是符号位，在开平方根是不必理会，因为计算机系统运算是基于实数的，在实数中对负数开平方是没有意义的
• \(E\) 表示阶码，但是有 \(127\) 的平移用来表示正负，故 \(E−127\) 是阶码表示的实际值
• \(M\) 表示尾码，表示原二进制浮点数第一个 \(1\) 之后的所有数位，基本可以认为隐藏了一位 \(1\) ,实际值为 \(1.xxxxxxx\)。

那么类似于十进制的科学计数法，我们可以将浮点数表示为：

\[y=(1+ \frac{M}{2^{23}})⋅2^{E−127} \]

对于浮点数 \(x=\frac{1}{\sqrt{a}}\)，我们不妨取对数；

\[\begin{aligned} -\frac{1}{2}\log_2a&=-\frac{1}{2}[\log_2(1+\frac{M_a}{2^{23}})+E_a-127] \end{aligned} \]

我们观察到 \(\frac{M_a}{2^{23}}\in(0,1)\)，故 \(\log_a(1+x)\) 可近似为 \(x+\epsilon\) （\(\epsilon\) 为修正值），即：

\[\begin{aligned} \log_2\frac{1}{\sqrt{a}}&=-\frac{1}{2}(\frac{M_a}{2^{23}}+E_a+\epsilon-127)\\ \log_2\frac{1}{\sqrt{a}}&=-\frac{1}{2}(\frac{M_a+2^{23}E_a}{2^{23}}+\epsilon-127)\\ \end{aligned} \]

观察这个部分 \(M_a+2^{23}E_a\)，非常巧妙地，跟据浮点数的存储规则，我们发现这个部分就是浮点数 \(a\) 强制转换为整数型后的表示，我们把这个部分称为 \(X_a\)。接下来，我们将 \(x\) 和 \(a\) 分别表示并联立：

\[\left\{ \begin{aligned} \log_2\frac{1}{\sqrt{a}}&=-\frac{1}{2}(\frac{X_a}{2^{23}}+\epsilon-127)\\ \log_2x&=\frac{X_x}{2^{23}}+\epsilon-127 \end{aligned} \right. \]

化简得；

\[\begin{aligned} \frac{X_x}{2^{23}}+\epsilon-127&=-\frac{1}{2}(\frac{X_a}{2^{23}}+\epsilon-127)\\ X_x&=(381+\epsilon)2^{23}-\frac{1}{2}X_a \end{aligned} \]

接下来，为了使修正值更加正确，则应该使得 \(\log2(1+x)\) 与 \(x+\epsilon\) 在区间 \([0,1]\) 上平均差值为 \(0\)，即：

\[\begin{aligned} \int_0^1\log_2(1+x)dx-(x+\epsilon)dx&=0\\ \int_0^1\log_2(1+x)dx-\int_0^1(x+\epsilon)dx&=0\\ \int_0^1\log_2(1+x)dx&=\int_0^1(x+\epsilon)dx\\ \frac{1}{\ln2}\int_1^2\ln{u}du&=\frac{1}{2}+\epsilon\\ \frac{1}{\ln2}(2ln{2}-1)&=\frac{1}{2}+\epsilon\\ 2-\frac{1}{ln{2}}&=\frac{1}{2}+\epsilon\\ \epsilon&=\frac{3}{2}-\frac{1}{\ln2}\\ \end{aligned} \]

最后将\(2^{22}\times(381-\epsilon)\) 十六进制硬编码为 0x5f3c551d，而 0x5f3759df 这个经典算法中的”魔数“怎么来的希腊奶（事实上确实是经典算法中的“魔数”能带来更加精确的答案。我猜是暴力二分试出来的，毕竟是程序师不是数学家，费这劲干啥）。

综上，我们就完整的得到了精度非常高的初始值 \(x_0\)，接着将这种初始值代入牛顿法仅需 \(1\) 次迭代就能达到相当高的精度，大大降低了复杂度。

上个原版代码；

float Q_rsqrt( float number )
{
    long i;
    float x2, y;
    const float threehalfs = 1.5F;

    x2 = number * 0.5F;
    y  = number;
    i  = * ( long * ) &y;                       // evil floating point bit level hacking
    i  = 0x5f3759df - ( i >> 1 );               // what the fuck?
    y  = * ( float * ) &i;
    y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
//    y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

    return y;
}