- 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.