Root-finding algorithm: Difference between revisions

Content deleted Content added
m remove extlink again
No edit summary
Line 29:
* [[Broyden's method]]
* [[Fixed-point iteration]]
 
===A fast implementation of a square root algorithm using integer shift===
 
By Ricardo Pereira
 
I have just figured out that the square root of a number
is very close to the same number less half the number of bits on it
 
Looking in the books I saw that a square root has always half the number of
digits , then just making a few tests I have observed the truth :
 
The square root of a number is very close to the same number less
half the bits on it
 
With this information I have constructed this algorithm , here the
approximation initiate not at an arbitrary value but exactly at the
approximation based on the half bits definition
 
The second part of the algorithm is the same Babylonian square root method
, what make it different is that it initiate already very close to the final
square root value , giving a good improvement in speed , and the interations
in the approximation loop never will be more than 7 , limiting it to a
fixed value no matter what is the value of the number
 
Without this approximation the number of interations is undefined and
varies depending on the number
 
This new method gives at least 50 percent speed improvement over the
original method
 
Comparing this method with the fsqrt function that is part of the x86
processor this implementation is only 5 times slower than the hardware
implementation, and yet I think that there is too much room for
improvements on this square root algorithm
 
//////////////////////////////////////////////////
 
#include <windows.h>
#include <stdio.h>
#include <math.h>
#include <unistd.h>
#include <limits.h>
#include <time.h>
 
#undef NDEBUG
#include <assert.h>
#include <string.h>
#include <mmsystem.h>
#include <float.h>
#include <stdlib.h>
 
 
double
sqrt2 (double val)
{
 
double temp;
double temp2;
double temp3;
 
/*variable to see whether we have a good appproximation */
double lastval;
 
 
if (1 == val)
{
 
return 1;
}
 
/*initiate at 2 , later it will be adjusted to a more close square root */
temp3 = 2;
 
 
/*be sure that lastval is different than val when entering the approximation function */
lastval = -val;
 
 
/*no reason to execute the approximation in values less than 4 */
if (4 < val)
{
 
unsigned int p;
unsigned int p2;
int counter = 0;
 
p = (int) val;
 
p2 = p;
 
while (p != 0)
{
 
/*count the number of bits in the integer portions of the number */
 
p = p >> 1;
 
counter++;
 
}
 
 
/*split the number */
counter = ((counter) / 2);
 
/*remove half the bits on it */
p2 = p2 >> counter;
 
 
/*here we have the approximation , the second part of
the code will start already very close to the final square root */
temp3 = p2;
 
}
 
 
while (1)
{
/*aproximation based on the code t=(t/n)/2 */
temp2 = val / temp3;
 
temp3 = (temp2 + temp3) * 0.5;
 
if (lastval == temp3)
{
/*if the last approximation is the same then we have found the
best possible approximation for the double type precision */
goto achou;
 
}
 
lastval = temp3;
 
}
 
achou:
 
return temp2;
 
}
 
int
main ()
{
 
double val = 2;
 
val = sqrt2 (val);
 
printf ("%.16f \n", val);
 
val = 100000000;
 
val = sqrt2 (val);
 
printf ("%.16f \n", val);
 
val = 555555555;
 
val = sqrt2 (val);
 
printf ("%.16f \n", val);
 
val = 999999999;
 
val = sqrt2 (val);
 
printf ("%.16f \n", val);
 
return 0;
 
}
/////////////////////////////////////////////////
 
For more information about this algorithm , visit:
http://rspsoftware.esparta.8x.com.br/squareroot.htm
maquisistem@wln.com.br
 
==Finding roots of polynomials==