15,921,989 members
Articles / Programming Languages / C++

Best Square Root Method - Algorithm - Function (Precision VS Speed)

Rate me:
15 Sep 2010CPOL4 min read 386K   2.5K   99   49
Square Root Methods Fast Algorithm Speed Precision computational Quake3 Fast Square Root Function Fast Gaming

Introduction

I enjoy Game Programming with Directx and I noticed that the most called method throughout most of my games is the standard sqrt method in the Math.h and this made me search for faster functions than the standard sqrt. And after some searching, I found lots of functions that were much much faster but it's always a compromise between speed and precision. The main purpose of this article is to help people choose the best square-root method that suits their program.

Background

In this article, I compare 14 different methods for computing the square root with the standard sqrt function as a reference, and for each method I show its precision and speed compared to the sqrt method.

1. Explaining how each method works
2. New ways to compute the square root

Using the Code

The code is simple, it basically contains:

1. main.cpp

Calls all the methods and for each one of them, it computes the speed and precision relative to the sqrt function.

2. SquareRootmethods.h

This Header contains the implementation of the functions, and the reference of where I got them from.

First I calculate the Speed and Precision of the sqrt method which will be my reference.

For computing the Speed, I measure the time it takes to call the sqrt function (M-1) times and I assign this value to `RefSpeed` which will be my reference.

And for computing the Precision, I add the current result to the previous result in `RefTotalPrecision` every time I call the method. `<code>RefTotalPrecisi`on will be my reference.

For measuring runtime duration (Speed) of the methods, I use the `CDuration` class found in this link, and I use `dur` as an instance of that class.

C++
```for(int j=0;j<AVG;j++)
{
dur.Start();

for(int i=1;i<M;i++)
RefTotalPrecision+=sqrt((float) i);

dur.Stop();

Temp+=dur.GetDuration();
}

RefTotalPrecision/=AVG;
Temp/=AVG;

RefSpeed=(float)(Temp)/CLOCKS_PER_SEC;      ```

And for the other methods I do the same calculations, but in the end, I reference them to the sqrt.

C++
```	for(int j=0;j<AVG;j++)
{
dur.Start();

for(int i=1;i<M;i++)
TotalPrecision+=sqrt1((float) i);

dur.Stop();

Temp+=dur.GetDuration();
}

TotalPrecision/=AVG;
Temp/=AVG;

Speed=(float)(Temp)/CLOCKS_PER_SEC;

cout<<"Precision = "
<<(double)(1-abs((TotalPrecision-RefTotalPrecision)/(RefTotalPrecision)))*100<<endl;         ```

NOTES:

1. I assume that the error in Precision whether larger or smaller than the reference is equal, that's why I use "abs".
2. The Speed is referenced as the actual percentage, while the Precision is referenced as a decrease percentage.

You can modify the value of `M` as you like, I initially assign it with 10000.

You can modify `AVG` as well, the higher it is, the more accurate the results.

C++
```#define M 10000
#define AVG 10   ```

Points of Interest

Precision wise, the sqrt standard method is the best. But the other functions can be much faster even 5 times faster. I would personally choose Method N# 14 as it has high precision and high speed, but I'll leave it for you to choose. :)

I took 5 samples and averaged them and here is the output:

According to the analysis the above Methods Performance Ranks (Speed x Precision) is:

NOTE: The performance of these methods depends highly on your processor and may change from one computer to another.

The METHODS

Sqrt1

Algorithm: Babylonian Method + some manipulations on IEEE 32 bit floating point representation

C++
```float sqrt1(const float x)
{
union
{
int i;
float x;
} u;
u.x = x;
u.i = (1<<29) + (u.i >> 1) - (1<<22);

// Two Babylonian Steps (simplified from:)
// u.x = 0.5f * (u.x + x/u.x);
// u.x = 0.5f * (u.x + x/u.x);
u.x =       u.x + x/u.x;
u.x = 0.25f*u.x + x/u.x;

return u.x;
}  ```

Sqrt2

Algorithm: The Magic Number (Quake 3)

C++
```#define SQRT_MAGIC_F 0x5f3759df
float  sqrt2(const float x)
{
const float xhalf = 0.5f*x;

union // get bits for floating value
{
float x;
int i;
} u;
u.x = x;
u.i = SQRT_MAGIC_F - (u.i >> 1);  // gives initial guess y0
return x*u.x*(1.5f - xhalf*u.x*u.x);// Newton step, repeating increases accuracy
}   ```

Sqrt3

C++
```float sqrt3(const float x)
{
union
{
int i;
float x;
} u;

u.x = x;
u.i = (1<<29) + (u.i >> 1) - (1<<22);
return u.x;
} ```

Sqrt4

Reference: I got it a long time a go from a forum and I forgot, please contact me if you know its reference.

Algorithm: Bakhsali Approximation

C++
```float sqrt4(const float m)
{
int i=0;
while( (i*i) <= m )
i++;
i--;
float d = m - i*i;
float p=d/(2*i);
float a=i+p;
return a-(p*p)/(2*a);
}  ```

Sqrt5

Reference: http://www.dreamincode.net/code/snippet244.htm

Algorithm: Babylonian Method

C++
```   float sqrt5(const float m)
{
float i=0;
float x1,x2;
while( (i*i) <= m )
i+=0.1f;
x1=i;
for(int j=0;j<10;j++)
{
x2=m;
x2/=x1;
x2+=x1;
x2/=2;
x1=x2;
}
return x2;
}   ```

Sqrt6

Algorithm: Dependant on IEEE representation and only works for 32 bits

C++
```double sqrt6 (double y)
{
double x, z, tempf;
unsigned long *tfptr = ((unsigned long *)&tempf) + 1;
tempf = y;
*tfptr = (0xbfcdd90a - *tfptr)>>1;
x =  tempf;
z =  y*0.5;
x = (1.5*x) - (x*x)*(x*z);    //The more you make replicates of this statement
//the higher the accuracy, here only 2 replicates are used
x = (1.5*x) - (x*x)*(x*z);
return x*y;
}  ```

Sqrt7

Algorithm: Dependant on IEEE representation and only works for 32 bits

C++
```float sqrt7(float x)
{
unsigned int i = *(unsigned int*) &x;
i  += 127 << 23;
// approximation of square root
i >>= 1;
return *(float*) &i;
}   ```

Sqrt8

Algorithm: Babylonian Method

C++
```double sqrt9( const double fg)
{
double n = fg / 2.0;
double lstX = 0.0;
while(n != lstX)
{
lstX = n;
n = (n + fg/n) / 2.0;
}
return n;
}  ```

Sqrt9

Algorithm: Babylonian Method

C++
``` double Abs(double Nbr)
{
if( Nbr >= 0 )
return Nbr;
else
return -Nbr;
}

double sqrt10(double Nbr)
{
double Number = Nbr / 2;
const double Tolerance = 1.0e-7;
do
{
Number = (Number + Nbr / Number) / 2;
}while( Abs(Number * Number - Nbr) > Tolerance);

return Number;
}   ```

Sqrt10

Algorithm: Newton's Approximation Method

C++
```double sqrt11(const double number)e
{
const double ACCURACY=0.001;
double lower, upper, guess;

if (number < 1)
{
lower = number;
upper = 1;
}
else
{
lower = 1;
upper = number;
}

while ((upper-lower) > ACCURACY)
{
guess = (lower + upper)/2;
if(guess*guess > number)
upper =guess;
else
lower = guess;
}
return (lower + upper)/2;

}  ```

Sqrt11

Algorithm: Newton's Approximation Method

C++
``` double sqrt12( unsigned long N )
{
double n, p, low, high;
if( 2 > N )
return( N );
low  = 0;
high = N;
while( high > low + 1 )
{
n = (high + low) / 2;
p = n * n;
if( N < p )
high = n;
else if( N > p )
low = n;
else
break;
}
return( N == p ? n : low );
}  ```

Sqrt12

Reference: http://cjjscript.q8ieng.com/?p=32

Algorithm: Babylonian Method

C++
```double sqrt13( int n )
{
// double a = (eventually the main method will plug values into a)
double a = (double) n;
double x = 1;

// For loop to get the square root value of the entered number.
for( int i = 0; i < n; i++)
{
x = 0.5 * ( x+a / x );
}

return x;
}    ```

Sqrt13

Reference: N/A

Algorithm: Assembly fsqrt

C++
```double sqrt13(double n)
{
__asm{
fld n
fsqrt
}
}
```

Sqrt14

Reference: N/A

Algorithm: Assembly fsqrt 2

C++
```double inline __declspec (naked) __fastcall sqrt14(double n)
{
_asm fld qword ptr [esp+4]
_asm fsqrt
_asm ret 8
} ```

1.3 (15 September 2010)

• Added Method N#14 (which is the best method till now)

1.2 (24 June 2010)

• Added the Methods Performance Rank

1.0 (31 March 2010)

• Initial release

I hope that this article would at least slightly help those who are interested in this issue.

Written By
Software Developer (Senior) Free Lancer
Egypt
BSc Computer Engineering, Cairo University, Egypt.

I love teaching and learning and I'm constantly changing
Feel free to discuss anything with me.

A Freelance UX designer, Code is my tool brush, I design for experience !

My Website

 First PrevNext
 Square Root Method for Integer Jens Grabner16-May-21 7:29 Jens Grabner 16-May-21 7:29
 how should version 14 look translated into AT&T? Member 139044695-Aug-18 18:01 Member 13904469 5-Aug-18 18:01
 Timings Member 1138475128-Apr-15 9:11 Member 11384751 28-Apr-15 9:11
 Re: Timings Member 80284076-Jun-16 9:16 Member 8028407 6-Jun-16 9:16
 great article - exactly what I need :-) Matth Moestl22-Apr-15 0:58 Matth Moestl 22-Apr-15 0:58
 Not sure methods are best Alexey KK6-May-14 1:17 Alexey KK 6-May-14 1:17
 GCC code for testing Member 1055101029-Apr-14 20:58 Member 10551010 29-Apr-14 20:58
 Re: GCC code for testing Member 1225871319-Jan-16 20:52 Member 12258713 19-Jan-16 20:52
 objective-c 10z28-Apr-14 8:04 10z 28-Apr-14 8:04
 another way to do it Marcos Lohmann20-Nov-13 4:40 Marcos Lohmann 20-Nov-13 4:40
 Re: another way to do it Marcos Lohmann20-Nov-13 5:41 Marcos Lohmann 20-Nov-13 5:41
 Further advice on precision calculation Conor Manning22-Feb-13 8:22 Conor Manning 22-Feb-13 8:22
 Good Read Bassam Abdul-Baki29-May-12 7:26 Bassam Abdul-Baki 29-May-12 7:26
 Re: Good Read Mahmoud Hesham El-Magdoub5-Jun-12 20:33 Mahmoud Hesham El-Magdoub 5-Jun-12 20:33
 actual time taken by the methods? utkarshs28-May-12 6:47 utkarshs 28-May-12 6:47
 Re: actual time taken by the methods? Mahmoud Hesham El-Magdoub28-May-12 7:41 Mahmoud Hesham El-Magdoub 28-May-12 7:41
 Re: actual time taken by the methods? utkarshs29-May-12 17:19 utkarshs 29-May-12 17:19
 Re: actual time taken by the methods? Mahmoud Hesham El-Magdoub5-Jun-12 20:40 Mahmoud Hesham El-Magdoub 5-Jun-12 20:40
 An Improvement to method #13 (x86) Lior Kogan5-Jul-10 9:58 Lior Kogan 5-Jul-10 9:58
 Re: An Improvement to method #13 (x86) Mahmoud Hesham El-Magdoub15-Sep-10 6:47 Mahmoud Hesham El-Magdoub 15-Sep-10 6:47
 Re: An Improvement to method #13 (x86) Member 1225871320-Jan-16 12:33 Member 12258713 20-Jan-16 12:33
 [My vote of 2] Precision calculation is broken Peter_in_278028-Jun-10 14:14 Peter_in_2780 28-Jun-10 14:14
 quick sqrt, precission 100 =asm=7-Jun-10 21:21 =asm= 7-Jun-10 21:21
 Re: quick sqrt, precission 100 Mahmoud Hesham El-Magdoub23-Jun-10 23:37 Mahmoud Hesham El-Magdoub 23-Jun-10 23:37
 Re: quick sqrt, precission 100 supercat924-Jun-10 5:11 supercat9 24-Jun-10 5:11
 Last Visit: 31-Dec-99 18:00     Last Update: 23-Jun-24 16:42 Refresh 12 Next ᐅ