如何求一个正数$n$的平方根倒数$\frac{1}{\sqrt n}$

牛顿迭代法它是一种在实数域和复数域上近似求解方程的方法。方法使用函数$f(x)$的泰勒级数前面几项来寻找方程$f(x)=0$的根。

首先,选择一个接近函数$f(x)$零点的$x_{0}$,计算相应的$f(x_{0})$和切线斜率$f'(x_{0})$(这里$f'$表示函数$f$的导数)。然后我们计算穿过点$(x_{0},f(x_{0}))$并且斜率为$f'(x_{0})$的直线和$x$轴的交点的$x$坐标,也就是方程$0=(x-x_{0})\cdot f'(x_{0})+f(x_{0})$的解。

我们将新求得的点的$x$坐标命名为$x_{1}$,通常$x_{1}$会比$x_{0}$更接近方程$f(x)=0$的解。因此我们现在可以利用$x_{1}$开始下一轮迭代。迭代公式可化简为$x_{n+1}=x_{n}-{\frac {f(x_{n})}{f'(x_{n})}}$

若$\frac{1}{\sqrt n}$为函数$f(x)$的根,$x=\frac{1}{\sqrt n} \Rightarrow n = \frac{1}{x^2} \Rightarrow f(x) = \frac{1}{x^2}-N$

通过已知的$f(x)$得到具体迭代公式为$x_{t+1} = x_t - \frac{x^{-2}-n}{-2x^{-3}} = \frac{3}{2}x_{t}-\frac{n}{2}x_t^3$

为了减少迭代次数,初始时找到选择一个接近函数$f(x)$零点的$x_{0}$尤为重要。

一个浮点数 (Value) 的表示其实可以这样表示:

$Value =sign \times exponent \times fraction$

也就是浮点数的实际值,等于符号位(sign bit)乘以指数偏移值(exponent bias)再乘以分数值(fraction)。

IEEE754浮点数标准中,32位浮点数的组成:最高位$S$为符号位,接下来8位数为指数$E$,最后低位23为分数值$F$。($S,E,F$均以无符号型二进制数表示)。

Untitled

浮点数的数值表示为$(-1)^{S+1}(1+\frac{F}{2^{23}})2^{E-127}$

我们只考虑正数,所以符号位$S=0$。

对于一个32位的正浮点数$T$,若将二进制存储视为32位整型数$T_b$,可知$T_b = 2^{23}E+F$

我们尝试寻找$T$和$T_b$的关系。二者取对数后仍相等,$log_2(T) = log_2(1+\frac{F}{2^{23}})+E-127$

由于$\frac{F}{2^{23}} \in [0, 1)$,通过观察$log_2(x+1)$的函数图像

Untitled

Untitled

发现,$0\le x < 1$时,$x\le log_2(x+1) \le x-log_2(ln(2))-\frac{1}{ln(2)}+1 \approx x+0.0860713320559$ 我们可以近似认为$log_2(x+1) \approx x+u$, $u$是一个附加常数。

$log_2(T) = log_2(1+\frac{F}{2^{23}})+E-127 \approx \frac{F}{2^{23}} + u + E - 127 = \frac{1}{2^{23}}(F+2^{23}E) + u - 127$

32位浮点数与其二进制之间存在一个近似的关系$log_2(T) \approx \frac{T_b}{2^{23}} + u - 127$

类似的,64位浮点数与其二进制之间存在一个近似的关系$log_2(T) \approx \frac{T_b}{2^{52}} + u - 1023$

利用这个关系来求$\frac{1}{\sqrt n}$的近似值,设$y = \frac{1}{\sqrt n}$

$log_2(y) = log_2 n^{-\frac{1}{2}} \Rightarrow \frac{y_b}{2^{23}} + u - 127 = -\frac{1}{2}(\frac{n_b}{2^{23}} + u - 127) \Rightarrow y_b = 3(127-u)2^{22} - \frac{n_b}{2}$

在c语言中,float型浮点数n在取地址后,可以用任意类型指针变量存储。一般情况下我们用float* 指针存储,通过指针访问所指内容就是float型变量,但是如果用long* 存储float型地址,那么通过指针访问会认为所指内容是long型,也就得到了浮点数的存储的二进制结构。

通过指针得到$n$的二进制$n_b$,通过神奇的运算$3(127-u)2^{22} - \frac{n_b}{2}$得到平方根倒数的二进制$y_b$,再通过指针转为浮点数$y$。$y \approx \frac{1}{\sqrt n}$