Fast Inversed Squareroot Function - Quake 3

Fast Inversed Squareroot Function – Quake 3

This technique gained popularity in the early days of computer graphics when processing power was limited. It’s an approximation method used to calculate the reciprocal of a square root (1 / sqrt(x)) faster than the standard division and square root functions. This was implemented in Quake 3 game.

The source code

				
					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;
}

Let’s dig into the source code

The IEEE 754 Standard

💡 Wikipedia Single-precision floating-point format

According to this standard, a float number is represented as 3 parts of 32 bits: Sign bit (x1), Exponential bits (x8), Mantissa bits (x23).
Then, a positive float number $x$ is calculated by: $x = (1 + \frac{M}{2^{23}}) * 2^{E - 127} (*)$

And, bit representation of $x$ would be: $M + 2^{23} * E (* *)$

Extract $x$ binary information from the $\log_{2} x $

For $x$ small $x \in [0, 1)$ , we could approximate $x$ by: $x \approx \log_{2} (1 + x)$ , then we can claim that:

$\log_{2} (1 + x) = x + μ; (μ is empirically chosen at 0.043)$

Looking at an arbitrary value $x$ , applying (*), we have:

$\log_{2} x = \log_{2} ((1 + \frac{M}{2^{23}}) * 2^{E - 127}) = \log_{2} (1 + \frac{M}{2^{23}}) + E - 127$

Since $\frac{M}{2^{23}} \in [0, 1)$ , we now apply the approximation rule:

$\log_{2} (1 + \frac{M}{2^{23}}) = \frac{M}{2^{23}} + μ$

Then:

$\begin{aligned} \log_{2} x & = \frac{M}{2^{23}} + μ + E - 127 \\ = \frac{M}{2^{23}} + \frac{2^{23} * E}{2^{23}} + μ - 127 \\ = \frac{1}{2^{23}} (\underset{bit representation (**)}{\underset{⏟}{M + 2^{23} * E}}) + μ - 127 \end{aligned}$

Calculate $\frac{1}{\sqrt{x}}$ using $\log_{2} $

The meaning of the line of code: -(i >>1)

$\log_{2} (\frac{1}{\sqrt{x}}) = \log_{2} (x^{- \frac{1}{2}}) = - \frac{1}{2} \log_{2} (x)$

Let $Γ = \frac{1}{\sqrt{x}} ⟹ \log_{2} (Γ) = \log_{2} (\frac{1}{\sqrt{x}}) = - \frac{1}{2} \log_{2} (x) $

We replace the logarithm with the bit representation for both sides:

$\begin{aligned} \frac{1}{2^{23}} (M_{Γ} + 2^{23} * E_{Γ}) + μ - 127 & = - \frac{1}{2} (\frac{1}{2^{23}} (M_{x} + 2^{23} * E_{x}) + μ - 127) \\ ⟺ M_{Γ} + 2^{23} * E_{Γ} & = \underset{0x5f3759df}{\underset{⏟}{\frac{3}{2} 2^{23} (127 - μ)}} - \underset{i >> 1}{\underset{⏟}{\frac{1}{2} (\underset{x}{\underset{⏟}{M_{x} + 2^{23} * E_{x})}}}} \end{aligned}$

Here reveals the secret of this line of code:

				
					i  = 0x5f3759df - ( i >> 1 );  // what the fuck?

Newton Iteration Method for Finding Approximation of Root

$x_{n e w} = x_{o l d} - \frac{f (x)}{f^{'} (x)}$

Let $f (y) = \frac{1}{y^{2}} - x$ , then: $f (y) = 0 ⟺ y = \frac{1}{\sqrt{x}}$ . We try to solve $f (y) = 0$ using Newton Method. First we compute the derivative of $f (y)$ :

$f^{'} (y) = \frac{y - x y^{3}}{- 2}$

Then

$y_{n e w} = y - \frac{f (y)}{f^{'} (y)} = y - \frac{y - x y^{3}}{- 2} = y (\frac{3}{2} - \frac{x}{2} y^{2})$

This is the meaning of the last line of code:

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