Calculating sqrt of integers
Marc Verbruggen
--
-----------------------------------------------------------
Marc Verbruggen tel : 03.240.94.19
Alcatel Bell Telephone fax : 03.240.99.50
F. Wellesplein 1 email : mv...@se.alcbel.be
B-2018 Antwerp
Belgium
Newsgroups: sci.math
Sender: mv...@btma56.nohost.nodomain (Marc Verbruggen)
From: mv...@se.alcbel.be (Marc Verbruggen)
Path: se.alcbel.be!mvbr
Distribution: world
Reply-To: mv...@se.alcbel.be
Organization: Alcatel Bell Telephone
Subject:Calculatring sqrt of an integer
--
Recently, someone explained me a method to calculate the square root of any
integer up to any required precision. In the following explanation I use
the following notation:
x="a"+"b" means x is a number created by concatenating a and b.
e.g. z="3"+"4" yields z = 34
The result will be notated as r(1)r(2)r(3)...
e.g sqrt(2)=1.414..., so r(1)=1, r(2)=4, r(3)=1, ...
Take the integer I and split in groups of 2 digits, starting from the right,
resulting in parts S(M), S(M-1) ... S(1), where S(M) is 1 or 2 digits
depending on the number of digits in I. e.g 123 yields S(1)=23 and S(2)=1.
1st_step : take largest number r(1) (<10) for which r(1)**2 < S(M)
put z = r(1)
x = S(M)
i = M
j = 1
root = ""
next_step : root = "root"+"r(j)"
d = x - z*r(j)
j = j + 1
x = "d" + "S(i-1)" note: if i-1 = 0 then x = "d" + "00"
and subsequent r(j) form the decimal part
y = 2*root
find largest r(j) (<10) such that z="y"+"r(j)"
and z*r(j) < x
goto next_step
In this way you can calculate sqrt with paper and pencil.
Now my questions :
(1) What are the basic mathematical reasons why this works as such
(2) Does it also works in a similar way for I**(1/N) with N>2 ? E.g.
splitting in groups of three digits, etc...
-----------------------------------------------------------
Marc Verbruggen tel : 03.240.94.19
Alcatel Bell Telephone fax : 03.240.99.50
F. Wellesplein 1 email : mv...@se.alcbel.be
B-2018 Antwerp
Belgium
USENET News System
>
> Recently, someone explained me a method to calculate the square root of any
> integer up to any required precision. In the following explanation I use
> the following notation:
>
>
> (1) What are the basic mathematical reasons why this works as such
> (2) Does it also works in a similar way for I**(1/N) with N>2 ? E.g.
> splitting in groups of three digits, etc...
>
and b is the rest.
2. It can be done for cube roots but is rather messy.
3. Neither is very efficient, Newton-Raphson gives a better
algorithm. For the square root of A this amounts to
next approx = (prev approx + A/(prev approx))/2
For nth roots, it amounts to new = (n-1)old/n + A/n / old^(n-1)
(I don't claim these are the very best).
Terry Moore
Michael L. Brundage
[procedure omitted]
>In this way you can calculate sqrt with paper and pencil.
>Now my questions :
>(1) What are the basic mathematical reasons why this works as such
>(2) Does it also works in a similar way for I**(1/N) with N>2 ? E.g.
>splitting in groups of three digits, etc...
I'd like to know (2) also ... I know there is a procedure for computing cube
roots on abacus, but I found that procedure a lot more complicated that for
taking square roots. Is there a simpler procedure for taking them by hand?
(N=3 in particular, but N in general would be cool too)
>-----------------------------------------------------------
>Marc Verbruggen tel : 03.240.94.19
>Alcatel Bell Telephone fax : 03.240.99.50
>F. Wellesplein 1 email : mv...@se.alcbel.be
>B-2018 Antwerp
>Belgium
>
--
mo...@cco.caltech.edu | "Let not your heart be troubled,
or | believe in God, believe also in Me."
mich...@iago.caltech.edu | John 14:1 (NASVB)
Steve Monson
OK. There have been hints enough posted. I remember when I was in the 9th
grade I asked my math teacher if he had any idea how to do cube roots by hand,
and he did not. I looked up the articles on math in the Encyclopaedia
Britannica, and found an example worked out, but with apparently no explana-
tion. It just said it was based on the expansion of (a+b)^n
and that was that. After an hour or two of trying to decipher the example, I
finally realized what they were talking about, and it's really cinchy (but
gets messy in a hurry) to extract roots to any order. It's just generalized
long division, in a way.
When we divide, we make a guess at first, and do subtractions to reduce the
rest of the number, till we get down to where the dividend is less than the
divisor. Then we are done. So. On to square roots.
2 2 2 2
(a + b) = a + 2ab + b = a + (2a+b)b
In our example, we guess our first a, and go on from there. Say we want to
extract the square root of 472843. Well, our first guess is 600, since 600^2
is less than 472843, and 700^2 is greater. Then we want to refine it. 600
is our "a" in the polynomial; we need to guess successive b's till we get as
close to the value as we desire.
_6_8_7
\/472843
a=60 36 (discard trailing zereos)
--
2a: 120 1128
b: 8
2a+b: 128 1024 1024 = 8*128 = b(2a+b)
----
a=680 10443
2a: 1360
b: 7
2a+b: 1367 9569
-----
87400
and so on to further digits. Now, for cube roots, we have
3 3 2 2 3 3 2 2
(a+b) = a + 3a b + 3ab + b = a + (3a + 3ab + b )b
So, to extract cube roots, we just do what we did with square roots, but
with a bit more arithmetic:
3 _3__3__1__5
2 \/36429280875
a=30, 3a : 2700 |27 subtract our a^3
3ab: 270 |-- guess for next b is 3
b^2: 9 | 9429
----- |
2979 | 8937 (8937 = 2979 * 3)
----
a=330 |492280
3a^2: 326700 |
3ab: 990 | next b is 1
b^2 1 |
------ |
327691 | 327691
------
a=3310 |164589875
3a^2: 32868300 |
3ab: 49650 | next b is 5
b^2 25 |
-------- |
32917975 164589875
---------
0 (exact cube root)
Note that throughout the root extraction, we still have to guess the next "b"
just as we do in long division. But that is how it's done.
Steve
-- mon...@diablo.amd.com -- (512) 462-4013
__ | signature designed by | You have to be able to tell a tracheotomist
(_ | (and ripped off from) | from a cutthroat.
__)teve | Stephen Wayne Miller | Gregory Benford, _Artifact_
Torben AEgidius Mogensen
>In article <9...@se.alcbel.be>, mv...@se.alcbel.be (Marc Verbruggen) writes:
>>
>> Recently, someone explained me a method to calculate the square root of any
>> integer up to any required precision. In the following explanation I use
>> the following notation:
>>
>> ...method deleted...
>> Now my questions :
>>
>> (1) What are the basic mathematical reasons why this works as such
>> (2) Does it also works in a similar way for I**(1/N) with N>2 ? E.g.
>> splitting in groups of three digits, etc...
>>
>1. It uses the formula for sqr(a+b), a is the part already calculated
> and b is the rest.
Below is shown a derivation of the algorithm from the sum formula (for
cube roots, though).
>2. It can be done for cube roots but is rather messy.
Below is shown a cube root version of the algorithm.
>3. Neither is very efficient, Newton-Raphson gives a better
> algorithm. For the square root of A this amounts to
> next approx = (prev approx + A/(prev approx))/2
> For nth roots, it amounts to new = (n-1)old/n + A/n / old^(n-1)
Newtons method might be fine *if* you have very fast multiplication and
division. Otherwise, using a long-division-like method will get the
cube root in a time that is similar to two divisions (implemented by
long division). Square root by a long-division-like method will use
about 1.5 times the time of a long division.
For those who do not know the long-division-like cube root method, I
will now show it and describe how (why) it works. We start out with a
simple and not very efficient algorithm, and then stepwise modify it
to a less obvious but much more efficient algorithm. The same method
can be used to develop the square root algorithm.
We assume that the number (a) which we want the cube root of is
between 0.125 and 1.0, as we can compute the cube root of floating
point numbers by this rule:
croot(x*2^3n) = croot(x)*2^n
croot(x*2^(3n+1)) = croot(x/4)*2^(n+1)
croot(x*2^(3n+2)) = croot(x/2)*2^(n+1)
If the exponent is small (8 bits) it is worthwhile to make the (mod 3)
and (div 3) into table lookups.
Assuming 0.5 <= x < 1.0 (as in most floating point formats), we get
0.125 <= a < 1.0. Now for the first version:
{ 0.125 <= a < 1.0 }
r := 0.0;
for n := 1 to number_of_bits do begin
if (r+2^-n)^3 <= a then
r := r+2^-n
end;
{ croot(a) - r < 2^-number_of_bits }
We repeatedly try to add a bit to the root and test if this doesn't
get us past the root. Note that invariants are enclosed in braces. We
then expand (r+2^-n)^3 by the cube-of-sums formula:
{ 0.125 <= a < 1.0 }
r := 0.0;
for n := 1 to number_of_bits do begin
if r^3+3*2^-n*r^2+3*2^-2n*r+2^-3n <= a then
r := r+2^-n
end;
{ croot(a) - r < 2^-number_of_bits }
We then do incremental computing of r^3:
{ 0.125 <= a < 1.0 }
r := 0.0;
r3 := 0.0;
for n := 1 to number_of_bits do begin
t := r3+3*2^-n*r^2+3*2^-2n*r+2^-3n;
if t <= a then begin
r3 = t;
r := r+2^-n
end
end;
{ croot(a) - r < 2^-number_of_bits }
We can actually avoid keeping r3, if we modify a instead:
{ 0.125 <= a < 1.0 }
r := 0.0;
for n := 1 to number_of_bits do begin
t := 3*2^-n*r^2+3*2^-2n*r+2^-3n;
if t <= a then begin
a := a-t;
r := r+2^-n
end
end;
{ croot(a0) - r < 2^-number_of_bits }
We then do incremental computing of r^2:
{ 0.125 <= a < 1.0 }
r := 0.0;
r2 := 0.0;
for n := 1 to number_of_bits do begin
t := 3*2^-n*r2+3*2^-2n*r+2^-3n;
if t <= a then begin
a := a-t;
r2 := r2+2^-(n-1)*r+2^-2n;
r := r+2^-n
end
end;
{ croot(a0) - r < 2^-number_of_bits }
Now there are only additions, multiplication by 3 or powers of 2,
which can be implemented by shifting. Note that a will be less than
2^-n and t will be less than 4*2^-n. To avoid loss of precision, it
might be beneficial to compute 2^n*a and 2^(n-2)*t instead, as these
will use more bits of precision and still be less than 1. This gives
the following version:
{ 0.125 <= a < 1.0 }
r := 0.0;
r2 := 0.0;
for n := 1 to number_of_bits do begin
t := .75*r2+3*2^-(n+2)*r+2^-(2n+2);
if 4*t <= 2*a then begin
a := 2*a-4*t;
r2 := r2+2^-(n-1)*r+2^-2n;
r := r+2^-n
end else
a := 2*a
end;
{ croot(a) - r < 2^-number_of_bits }
Note that we double a and quadruple t in the comparison (with the risk
of making them larger than one). We can halve a instead.
{ 0.125 <= a < 1.0 }
r := 0.0;
r2 := 0.0;
a := a/2;
for n := 1 to number_of_bits do begin
t := .75*r2+3*2^-(n+2)*r+2^-(2n+2);
if t <= a then begin
a := a-t;
r2 := r2+2^-(n-1)*r+2^-2n;
r := r+2^-n
end;
a := 2*a
end;
{ croot(a) - r < 2^-number_of_bits }
If you have a processor (like the ARM) which can do a free (constant)
barrel shift with every add, you would benefit from unfolding the loop
and use shifts of a constant number of bits to compute 2^-n, 2^-2n
etc. in each iteration.
I hope this will be of some help.
Yours,
Torben Mogensen (tor...@diku.dk)
loa...@vax.oxford.ac.uk
> In article <9...@se.alcbel.be>, mv...@se.alcbel.be (Marc Verbruggen) writes:
>>
>> Recently, someone explained me a method to calculate the square root of any
>> integer up to any required precision. In the following explanation I use
>> the following notation:
>>
>> ...method deleted...
>
> and b is the rest.
>
> 2. It can be done for cube roots but is rather messy.
>
> 3. Neither is very efficient, Newton-Raphson gives a better
> algorithm. For the square root of A this amounts to
> next approx = (prev approx + A/(prev approx))/2
you look carefully at the operations involved, it can be carried with only a
little more work than a _single_ division.
--
R Loader
1 . 4 1 4
------------------
2 . 00 00 00 00 00
1 1
------------------
1 . 00 2.4 * .4 = .96
. 96
------------------
. 04 00 2.81 * .01 = .0281
. 02 81
------------------
. 01 19 00 2.824 * .004 = .011296
. 01 12 96
------------------
. 00 07 04
etc.
Try finding a few digits of sqrt(2) by hand using N-R.
Noam Elkies
loa...@vax.oxford.ac.uk writes:
>Try finding a few digits of sqrt(2) by hand using Newton-Raphson.
Remember that you don't actually need to do the decimal divisions
until the last step. Starting with the approximation sqrt(2)~1, we get
3/2, 17/12, 577/408=((17/12)+(24/17))/2. Now 17/12=1.416666... is easy,
24/17=1.411765... only a bit worse, and the average 1.414216 is within
.000003 of the right answer.
--Noam D. Elkies (elk...@zariski.harvard.edu)
Dept. of Mathematics, Harvard University
Richard Harter
>In article <1992May29....@vax.oxford.ac.uk>
>loa...@vax.oxford.ac.uk writes:
>>Try finding a few digits of sqrt(2) by hand using Newton-Raphson.
>Remember that you don't actually need to do the decimal divisions
>until the last step. Starting with the approximation sqrt(2)~1, we get
>3/2, 17/12, 577/408=((17/12)+(24/17))/2. Now 17/12=1.416666... is easy,
>24/17=1.411765... only a bit worse, and the average 1.414216 is within
>.000003 of the right answer.
This is really only true for small integers.
Quite some time ago I posted what appears to be the most efficient algorithm
for the hand calculation of square roots. I had posted the formal details
but can't lay my hands on them. The general notion goes like this:
Let X.ABCDEF... be the number that we wish to take the square root of
and x.abcdef... be its root, where the letters indicate digits and the
leading letter is 0<X<100 and 0<x<10.
Now if we expand the square of the root in powers of 1/10 we get
x**2 + 2ax/10 + (2bx +a**2)/100 + (2cx + 2ab)/1000 + ...
If we were doing a formal series reversion this would give us
a = A/2x, b = (B - a**2)/2x, c = (C -2ab)/2x, ...
When we carry this out in powers of 1/10 we have to allow for the carry's.
This can be done formally -- the material I can't find shows how to do
it -- but it practice it can be easily estimated. To illustrate
N= 2.00000
x=1, sub x^2 1
__
Next digit 10
a=4, sub 2ax 8
--
Next digit 20
sub a^2 16
--
4
b=1, sub 2bx 2
__
Next digit 20
sub 2ab 8
__
12
c=4, sub 2cx 8 <--------- Setting c=5 will make the next term go neg
--
Next digit 40
sub 2ac+b^2 33
--
7
d=2, sub 2dx 4
Etc. The most striking disadvantage is that one tends to run out of
paper! The computational cost is approximately 1/4 that of the book
algorithm (counting multiplys and adds). The Newton Raphson iteration
is more painful if applied naively -- a judicious application of partial
truncation division brings it down to the book algorithm. Another problem
with the N-R scheme is that the number of iterations required to get a
specified accuracy. Besides being more efficient the above algorithm is
attractive in that it does not require any more digits in the original
number than there are in the answer.
Once upon a time I was a (moderately) fast calculator -- even today I mostly
never bother with mechanical aids. One of parlor tricks was to write a
number on the board and then write the square root immediately afterwards
to the same number of places as the original number. Five or six digits
was fairly easy. For multiplication I used cross products rather than
Trachtenberg, so the above was congenial since in either case one is taking
the running sum of products of two digits.
USENET News System
>
> In article <1992May30.1...@husc3.harvard.edu> elk...@ramanujan.harvard.edu (Noam Elkies) writes:
> >In article <1992May29....@vax.oxford.ac.uk>
> >loa...@vax.oxford.ac.uk writes:
> >>Try finding a few digits of sqrt(2) by hand using Newton-Raphson.
>
> for the hand calculation of square roots. I had posted the formal details
> but can't lay my hands on them. The general notion goes like this:
>
If you don't have a calculator, there are more efficient methods than
Newton Raphson. If you are writing multi-precision algorithms, the same
applies. But if you are programming or using a machine with built-in
multiplication and division the argument no longer applies.
Actually I use exp(ln(a)/n) the same as (almost) everyone else.
Terry Moore