
Integer Square Roots - chanux
http://www.embedded.com/98/9802fe2.htm
======
chanux
If you feel bored at the start of the article, Scroll down to "The friden
algorithm".

Flashback: I found the Friden algo part of this article printed and abandoned
(wasted) in a common room in uni. Took home one copy (There were few for some
reason) and loved it. Today I found embedded.com article on HN & looked for
this article. And thought HNers would love.

~~~
imp
Thanks for the tip. I skipped down to the friden algorithm and it was an
interesting story.

------
silentbicycle
> "When I'm starting on a new algorithm, I sometimes find it's both fun and
> instructional to begin by trying to see just how simple-minded an algorithm
> I can write. I think I do this just to prove to myself that the problem can
> at least be solved. From there, anything else is mere refinement."

This is excellent advice.

------
chmike
Yeah. Solved it 16 years ago as an intellectual game. comp.lang.c was my
Hacker News at the time ;)

See
[http://groups.google.fr/group/comp.lang.c/browse_thread/thre...](http://groups.google.fr/group/comp.lang.c/browse_thread/thread/229c0a6f0d0c031c/5aacf5997b615c37)

The method is very simple. We know how to compute the square of a value. So
start by finding the biggest integer value with a single bit set whose square
is smaller than the input value. Then set the bit just below to one, compute
the square of the value and compare it with the input value. If it is bigger,
the bit should be 0 otherwise it should be left to 1. And so on.

The thing is to avoid the multiplication to compute the square. So if you have
a value n and you want to set a bit to one. We compute in fact the square of
(N+B) where B is an integer with the bit to set to one, set to one. We then
have (N+B)² = N² + 2NB + B² The last term is easy to compute from the previous
B with the previous bit set. We simply divide by 4. 2NB is computed by a
simple shift and N² is the result of the previous computation.

For fixed point, we benefit from the fact that while we progress in computing
the square root, we have an increasing number of bits set to 0 in the most
significant bits of the values. I then simply shift up to get usable bits in
the less significant bits.

The only limitation is that we need the sign bit of the fixed point. So the
computation won't work with big unsigned fixed points. This is solved by
unwinding the first iteration.

It took me some time to understand all this. All I had was a short code to
compute integer square root and I had no idea how it worked.

------
lkozma
I believe after the first formula in the text:

 _In the vicinity of the correct root, the method converges exponentially,
meaning that it gives twice as many digits of accuracy at each step_

should instead be _converges quadratically_

~~~
dreish
No, it is exponential in the number of digits of accuracy as a function of the
number of iterations performed.

~~~
lkozma
Yes, but I think this is what is called quadratic convergence if the exponent
is 2, as it is here. See for ex.
[http://sepwww.stanford.edu/public/docs/sep97/paul1/paper_htm...](http://sepwww.stanford.edu/public/docs/sep97/paul1/paper_html/node5.html)

~~~
jonsen
Don't cross purpose each other.

In fact it's somewhat interesting, that quadratic error decay equals
exponential precision growth:

    
    
      0.1^2 = 0.01
      0.01^2 = 0.0001
      0.0001^2 = 0.00000001
      ...
    

Value equals error; number of zeros equals precision _measured as number of
digits_.

------
Locke1689
I wont lie, I went straight for the Taylor approximation. According to
wikipedia:

"As an iterative method, the order of convergence is equal to the number of
terms used. With 2 terms, it is identical to the Babylonian method; With 3
terms, each iteration takes almost as many operations as the Bakhshali
approximation, but converges more slowly. Therefore, this is not a
particularly efficient way of calculation."

Heh, woops. My second inkling was the way slide rules do it: exponentials and
logarithms ( sqrt(n)=exp(.5*ln(n)) ) Turns out that's more on target:

"Pocket calculators typically implement good routines to compute the
exponential function and the natural logarithm, and then compute the square
root of S"

~~~
jonsen
_Pocket calculators typically implement good routines_

Yeah, they are accurate, but don't need to be fast. The blink of an eye is a
long time.

