Tech Support > Computers & Technology > Programming > square root algorithm
square root algorithm
Posted by Gustavo Rios on August 28th, 2003


Hi folks.

I am in need for an algorithm to evaluates the square root for a given
number (integral). the specification is the following:

my_sqrt(x) ^ 2 <= x < (my_sqrt(x) + 1) ^ 2

Since, my requirement is for a very fast algorithm i seek for your
help. Is there a pretty fast method to evaluate this ?

Posted by Ricardo Gibert on August 28th, 2003


Aside from being integral, what other constraints are there on x? Lower limit on the value of x is assumed to be zero of course, but
upper limit on the value of x?

BTW, if you enter: "integer square root" into www.google.com (including the quotation marks), you should get lots of useful hits.


"Gustavo Rios" <gustavo.rios@terra.com.br> wrote in message news:c37fd0bd.0308280919.3518a3cc@posting.google.c om...


Posted by Corey Murtagh on August 28th, 2003


Gustavo Rios wrote:

A possible bit-shifting solution:

int my_sqrt(int x)
{
if (x < 0)
return 0;
unsigned int bit = 1;
for (; bit < x; bit <<= 1);

int res = 0;
for (bit >>= 1; bit; bit >>= 1)
{
int next = res | bit;
if (next * next <= x)
res = next;
}
return res;
}

Untested, but should give a reasonable approximation. You can skip the
initial calculation of bit and just set it to 0x80000000 (adjust for
your integer width) if you like.

Would be interested to hear if you find an integer sqrt that outperforms
the floating point version on modern systems.

--
Corey Murtagh
The Electric Monk
"Quidquid latine dictum sit, altum viditur!"


Posted by pete on August 28th, 2003


Gustavo Rios wrote:
The algorithm in isqrt() is fast on systems that don't
have multiplication operations, like Microchip PIC's.
I don't know a name for that algorithm.
irt() uses Newton's or Hero's method.



/* BEGIN, isqrt.c */

#include <stdio.h>

unsigned isqrt(unsigned, unsigned *);
unsigned irt(unsigned);

int main(void)
{
unsigned n, *rem, root;

for (n = 0; n < 17; ++n) {
root = irt(n);
printf("The square root of %u is %u.\n", n, root);
}
putchar('\n');
for (n = 15; n < 17; ++n) {
root = isqrt(n, rem);
if (*rem) {
printf("%u is not a perfect square.\n", n);
} else {
printf("%u is a perfect square.\n", n);
printf("The square root of %u is %u.\n", n, root);
}
}
return 0;
}

unsigned isqrt(unsigned n, unsigned *remainder)
{
unsigned power_of_4, root, sum;

power_of_4 = ((unsigned)-1 >> 2 - (unsigned)-1 % 3) + 1;
root = 0;
do {
sum = root + power_of_4;
if (n >= sum) {
n -= sum;
root = sum + power_of_4;
}
root >>= 1;
power_of_4 >>= 2;
} while (power_of_4);
if (remainder) {
*remainder = n;
}
return root;
}

unsigned irt(unsigned n)
{
unsigned a, b;

a = n;
for (b = (1 == a) + a >> 1; a > b; b = n / a + a >> 1) {
a = b;
}
return a;
}

/* END, isqrt.c */

--
pete

Posted by Thad Smith on August 29th, 2003


pete wrote:

Great piece of code. Thanks!


You need to add:

unsigned remainder;
rem = &remainder;

in order for rem to be defined. Otherwise rem might be NULL, in which
case no remainder is returned, or might be any address, in which case
anything could be clobbered.

Thad

Posted by pete on August 29th, 2003


Thad Smith wrote:
The function is designed to take NULL
as an argument for remainder, for times when you
don't care whether the number is a perfect square.

If you put in just any address, then you get what you deserve.

--
pete

Posted by Gustavo Rios on September 9th, 2003


gustavo.rios@terra.com.br (Gustavo Rios) wrote in message news:<c37fd0bd.0308280919.3518a3cc@posting.google. com>...
After some research, i could write a function to perform return the
desired value.

I have posted it: http://www.rootshell.be/~grios

It is in adk/msc/nsqrt.c

and be warned, before compiling anything, adjust adk/include/__conf.h
to suite your environment. And if after all you need futher help, feel
free to get in touch.


Similar Posts