From 7025922124b84c231355592b0cbb1bb3be33ad70 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?=E1=B5=87=E1=B5=92?= <40080007@qq.com>
Date: Tue, 23 Apr 2024 13:45:27 +0800
Subject: [PATCH] feat: fast inverse square root

---
 _posts/2024-04-23-Fast-inverse-square-root.md | 145 ++++++++++++++++++
 1 file changed, 145 insertions(+)
 create mode 100644 _posts/2024-04-23-Fast-inverse-square-root.md

diff --git a/_posts/2024-04-23-Fast-inverse-square-root.md b/_posts/2024-04-23-Fast-inverse-square-root.md
new file mode 100644
index 0000000..d5da237
--- /dev/null
+++ b/_posts/2024-04-23-Fast-inverse-square-root.md
@@ -0,0 +1,145 @@
+---
+layout: post
+title: 快速平方根倒数
+date: 2024-04-23 13:00:00 +0800
+categories: algorithm
+tags: bit-hack
+published: true
+---
+
+* content
+{:toc}
+
+## 雷神之锤3源码
+
+```C
+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;
+}
+```
+
+## 牛顿迭代
+
+
+
+## IEEE 754
+
+科学计数法
+
+normalised numbers
+denormalised numbers
+NaN
+0 & -0
+
+由 IEEE 754 可知 float 的 long representation 为 $$(E*2^{23} + M)$$
+
+## 推导
+
+x 为 float，令 y 为 x 的平方根倒数：
+
+$$
+\begin{aligned}
+& \quad\,y = \tfrac{1}{\sqrt{x}} = {x^{-\tfrac{1}{2}}}\\
+&\Rightarrow\log_2 (y)=-\tfrac{1}{2}\log_2 (x)
+\end{aligned}
+$$
+
+令 $$E_y$$ $$M_y$$  $$E_x$$ 和 $$M_x$$ 分别为 y x 的指数和小数部分：
+
+$$
+\begin{aligned}
+&\Rightarrow \log_2 (2^{(E_y-127)}*(1.0+\tfrac{M_y}{2^{23}})) = -\tfrac{1}{2}\log_2 (2^{(E_x-127)}*(1.0+\tfrac{M_x}{2^{23}}))\\
+&\Rightarrow {(E_y-127)}+\log_2 (1.0+\tfrac{M_y}{2^{23}}) = -\tfrac{1}{2}(E_y-127)-\tfrac{1}{2}\log_2 (1.0+\tfrac{M_x}{2^{23}})
+\end{aligned}
+$$
+
+代入 $$ \log(1+x) \approx x+\mu, x \in [0,1] $$：
+
+$$
+\begin{aligned}
+&\Rightarrow {(E_y-127)} + \tfrac{M_y}{2^{23}} + \mu = -\tfrac{1}{2}(E_y-127) - \tfrac{1}{2}(\tfrac{M_x}{2^{23}} + \mu)\\
+&\Rightarrow E_y*{2^{23}} + M_y = \tfrac{3}{2}(127 - \mu)2^{23} - \tfrac{1}{2}(E_x*2^{23} + M_x)
+\end{aligned}
+$$
+
+记 $$y^{\prime}$$ $$x^{\prime}$$ 分别为 float 数 y 和 x 的 long representation：
+
+$$
+\begin{aligned}
+&\Rightarrow y^{\prime} = \tfrac{3}{2}(127 - \mu)2^{23} - \tfrac{1}{2}x^{\prime}\\
+&\Rightarrow y^{\prime} = \tfrac{3}{2}(127 - \mu)2^{23} - (x^{\prime} >> 1)
+\end{aligned}
+$$
+
+$$ \dot y $$
+
+$$ \dot{y} $$
+
+$$ \ddot y $$
+
+$$ \ddot{y} $$
+
+$\log_2 10$
+
+$$\lg_2 10$$
+
+$$\lg 10$$
+
+$$lg 10$$
+
+$$\ln 10$$
+
+$$ln 10$$
+
+$$\log_2 10$$
+
+$$log_2 10$$
+
+$\sqrt{3x-1}+(1+x)^2$
+
+$$\sqrt3$$
+
+$$\sqrt{3}$$
+
+$$\sqrt{3x-1}$$
+
+$$\sqrt{3x-1}+(1+x)^2$$
+
+$10^{10}$
+
+$\frac{4n^3-6n^2-10n-2k^3-3k^2+8k+12kn+c}{12}$
+
+$$\frac{4n^3-6n^2-10n-2k^3-3k^2+8k+12kn+c}{12}$$
+
+$$\frac{4}{12}$$
+
+$$\tfrac{4}{12}$$
+
+$$\dfrac{4}{12}$$
+
+$$\cfrac{4}{12}$$
+
+## 总结
+
+现在已经有了平方根运算器，已经不需要这么写了。
+
+<!-- https://www.youtube.com/watch?v=p8u_k2LIZyo -->
+<!-- https://www.zhihu.com/question/26287650 -->
+<!-- https://en.wikipedia.org/wiki/Fast_inverse_square_root -->
+<!-- https://zh.wikipedia.org/wiki/%E5%B9%B3%E6%96%B9%E6%A0%B9%E5%80%92%E6%95%B0%E9%80%9F%E7%AE%97%E6%B3%95 -->
+<!-- https://zhuanlan.zhihu.com/p/400064205 -->
+<!-- https://brilliant.org/wiki/newton-raphson-method/ -->
+<!-- https://www.cnblogs.com/ywsun/p/14271547.html -->
\ No newline at end of file