🚀 go-pugleaf

> I sent him an integer square root code using Newton's method.  I've done
> this problem before and I'm sure there are better methods out there.  Does
> somebody want to post/email something better than Newton's?

B. Reese's description of Heron's algorithm is interesting, although the
obvious consequence of extension beyond 16 bits is increasing the number of
nested if's (or case labels), and if you ever go to 64-bit integers that
could really get out of control...  Here are two algorithms which I discovered
long ago, which may also be worthy of consideration.

The first old mentod uses base 2 arithmetic instead of base 10, but basically
just follows the method we all learned in high-school algebra class; I have
seen it mentioned in some old computer architecture/algorithm books, and saw a
particularly good implementation done by Eduardo D. Lara, at Link (the flight
simulation division of Singer) in Sunnyvale, California.  I later implemented
the same algorithm, for 32-bit integers, on an early Intel 8-bit chip.

Just after that I attended a talk by Tien Chi Chen at UC Berkeley, where he
described an algorithm that was much faster -- the complexity is comparable
to a multiplication algorithm, where the number of add-like operations is a
function of the number of 1's in the operand (here it's not the original bit
count, because the source operand gets subtractively modified during succes-
sive algorithm steps, but statistically it's still like a small constant
multiple of straightforward multiplication time).  In addition to square
root, he does division, log,* and exponential functions -- all based on
shift-and-add/subtract operations with execution time/complexity comparable
to integer multiplication.  (He operates on the mantissa of floating-point
numbers, and the same approach applies to integers.)

Chen published a paper describing this algorithm,
	T. C. Chen [1972] The Automatic Generation of Exponentials,
		Logarithms, Ratios and Square Roots.
		IBM Journal of Research and Development 16: 380-388
		(July 1972)
His paper has been cited by several other authors, and if you have heard of
the "CORDIC" technique you will recognize a lot of resemblance to Chen's
algorithm.  His major "value added" was hardware to locate the leftmost 1
in a binary integer representation; I know that sort of thing existed in
several minicomputers in the 70's and 80's, and hope it's available on most
machines now -- without it you'll just have to do that in software.

Because it depends on locating the leftmost bit, etc., I doubt that it could
be implemented in C (certainly not in Fortran), raising another issue about
Hardware Support for Numeric Algorithms, which I note has mostly become "how
high-level languages need to be expanded to comply with mathematical needs" --
but even in assembly language, with or without a "locate leftmost bit" opcode,
it's quite easy to write this (Chen modeled it in APL for verification
purposes, but I believe he relied on assembler for efficient execution).

Good luck...
DW

  * about Chen's log -- all the other algorithms, given that there isn't an
    exact integer result for an integer operand, have relative errors which
    are less than the least significant bit of the operand.  However, Chen's
    Log has an _absolute_ error near the least significant bit, meaning that
    a log of a floating-point mantissa which is very near to 1.0 will have a
    surprisingly large relative error because log(1+epsilon) is close to
    epsilon (if epsilon is very small) and the log error is also near epsilon.