Square Root Process Similar to Long Division

09172021, 11:00 PM
(This post was last modified: 09172021 11:20 PM by jeejohn.)
Post: #1




Square Root Process Similar to Long Division
I saw the thread "Third Order Convergence for Square Roots..."
I thought that I would offer up a process that allows calculations of square roots, similar to how long division is done: Using Albert Chan's 12345 to calculate the square root: 100 + 10 + 1 + 0.1 + 0.00 + 0.008 + .0000 + .00005 = 111.10805 100 _________ 12345 (200) (*100) = 10000 200 __________2345 210 _________2100 220 ___________245 221 __________221 222 ____________24 222.1 _________22.21 222.2 ___________1.79 222.20 __________0.0000 (222.20*.01=2.222 > 1.779) 222.20 __________1.79 222.208 ________1.777664 222.216 _________0.012336 222.2160 ________0.000000 (222.216*.0001=.00222216 > .001336 222.2160 ________0.012336 222.21605 _______0.001111080525 222.21610 _______0.0012251975 etc. Basically, we're keeping track of 2*X on the left, then we add a digit δ to get 2X+δ. Finally we multiple by δ to get 2Xδ+δ^2. We subtract this from the previous quotient to get a new quotient. Note that for the final quotient = 12345111.10805^2 = .0012251975 is the remainder left over from squaring 111.10805. I used the process to find square roots of numbers manually to 25+ digits fairly quickly. I have thought of somehow programming this on the HP48s, but I have not been able to get past the doubling of the one digit at each step in some quick ML way. I'm sorry about how this looks, but this forum appears to take out all of the extra spaces I had placed to make it more readable. I used underlines to space things out some, but it doesn't quite look right. Except for the first digit, the rest of the digits added are in bold (as δ), then multiplied by δ. In the next line the digit added is doubled (or add δ again?). John Edry 

09172021, 11:58 PM
Post: #2




RE: Square Root Process Similar to Long Division
(09172021 11:00 PM)jeejohn Wrote: I saw the thread "Third Order Convergence for Square Roots..." Nice. By the way, you can cut/paste table inside code block (like mine) I just deleted a post about oldschool square root from Third Order Convergence thread. It really should belong here. Here goes. For sqrt of value very close to 1, it is is even closer to 1. √(1+2ε) ≈ 1+ε We can even *assume* result is 1, then adjust for it. Example, √(0.9996) Here is old school method for square root. With nice divisor = 1, we do 3digits at a time. Numbers are scaled so that we work mostly with integer. 0.9996E6 = 1E3^2  400 // X = 1E3 Newtons correction (SX*X)/2/X, for 2nd 3digits. 400E3/2 = 200E3 = 200X + 0 Negative digits OK, we will fix later. It just meant √(0.9996) ≈ 0.999800 Other digits followed below patterns: 1st column = quotient 2nd column = remainder, then corrected (shelllike cross multiply, except center) Code: 1000 Normalize base1000 digits: √(.9996) = 0.999 799 979 995 998 999 719 915 973 591 ... Ref: book "Dead Reckoning: Calculating without instruments", chapter on roots 

09182021, 01:51 AM
Post: #3




RE: Square Root Process Similar to Long Division
If divisor is not 1, it work the same way
12345 = 111^2 + 24 // X = 111 111 24E3/2 = 12E3 = 108X + 12 108 12E3  108^2/2 = 6168 = 55X + 63 055 63E3  55*108 = 57060 = 514X + 6 √(12345) ≈ 111. 108 055 514 

09182021, 12:46 PM
Post: #4




RE: Square Root Process Similar to Long Division
(09172021 11:00 PM)jeejohn Wrote: Basically, we're keeping track of 2*X on the left, then we add a digit δ to get 2X+δ. Finally we multiple by δ to get 2Xδ+δ^2. We subtract this from the previous quotient to get a new quotient. Above is based on identity: (X+δ)^2 = X^2 + δ*(2X+δ) N = 12345 = (X+δ)^2 X = 100 R = NX^2 = 2345 R = δ*(2X+δ), but we cannot get δ without square root. Newton's method assumed (2X+δ) ≈ 2X, then solve for δ With this assumption, solved δ always overestmate. To compensate, we floored estimated result. R/2/X = 11.725 δ = 11 R = δ*(2X+δ) = 24 X = X+δ = 111 R/2/X ≈ 0.108108108 δ = 0.1081 R = δ*(2X+δ) = 0.00988561 X = X+δ = 111.1081 R/2/X ≈ 4.448645058E5 δ = 4.448646E5 // floor of negative number, same sig. digits as X √12345 ≈ 111.1081 + (4.448646E5) = 111.10805551354 

09182021, 01:07 PM
Post: #5




RE: Square Root Process Similar to Long Division
Newton's method converge fast, but required expensive mul/div.
HP has a digitbydigit method for squareroot, that is very efficient. First, we halved R R/2 = δ*(X+δ/2) if δ=1, δ*(X+δ/2) = (X+0.5) if δ=2, δ*(X+δ/2) = 2*(X+1) = (X+0.5) + (X+1.5) if δ=3, δ*(X+δ/2) = 3*(X+1.5) = (X+0.5) + (X+1.5) + (X+2.5) ... This way, we avoided division to get δ. When R/2  δ*(X+δ/2) < 0, we gone too far. √12345, digitbydigit: 1.2345  1^2 = .2345 R/2 = .11725 11.725  10.5 = 1.225 1.225  11.5 = 10.275 122.5  110.5 = 12 12  111.5 = 99.5 1200  1110.5 = 89.5 89.5  1111.5 = 1022 8950  11110.5 = 2160.5 895000  111100.5 = 783899.5 783899.5  111101.5 = 672798 672798  111102.5 = 561695.5 561695.5  111103.5 = 450592 450592  111104.5 = 339487.5 339487.5  111105.5 = 228382 228382  111106.5 = 117275.5 117275.5  111107.5 = 6168 6168  111108.5 = 104940.5 616800  1111080.5 = 494280.5 61680000  11110800.5 = 50569199.5 50569199.5  11110801.5 = 39458398 39458398  11110802.5 = 28347595.5 28347595.5  11110803.5 = 17236792 17236792  11110804.5 = 6125987.5 6125987.5  11110805.5 = 4984818 √12345 ≈ 111.10805 (underestimated due to positive remainder) Except for initial setup, algorithm only use subtraction Quote:Note that for the final quotient = 12345111.10805^2 = .0012251975 is the remainder We confirm by undo scaling and halving R first digit 1.2345/12345 = 1e4, other 7 digits 100^7 = 1e14 12345  111.10805^2 = 6125987.5 * 1e10 * 2 = 0.0012251975 

09182021, 08:55 PM
(This post was last modified: 09182021 09:00 PM by jeejohn.)
Post: #6




RE: Square Root Process Similar to Long Division
I had never taken the time to "reverse engineer" the HP square root algorithm. I knew it had to be something similar to what I posted as it calculated square roots only a little bit slower than doing divisions.
I take it that dividing R by 2 in DEC is still fast? In HEX it would just be a bit shift. I see that this allows each subtraction step to be incremented by 1. Nice. Back in the days of the HP48 newsgroup, I had posted a challenge over redoing the HP48 ML multiplication in a different possibly faster way. It used one of the multipliers  say X and calculated 4X. Then using a bunch of GOTOs, if digit=1 > +X digit=2 > +X+X digit=3 > +4XX (instead of +X+X+X) digit=4 > +4X (instead of +X+X+X+X) etc. It was slightly faster, but not really worth the effort. All the tests and GOTOs steps added a lot of clock cycles. For division, something similar can be done and also an idea was to add a second loop to follow a negative remainder to "Add Up" to a positive remainder. Something similar could be done for the square root algorithm. When I did my manual square root calculations, I would keep a table of 2Xδ updated to grab quickly to add at the bottom and subtract to get the next remainder. Interesting how the processor capabilities affects which algorithm is used to calculate irrational functions. HP finding out about the Cordic algorithm galvanized HP into creating the HP9100 (the first desktop scientific calculator), which in turn, encourage Bill Hewlett to challenged his workers to create the HP35. John Edry 

11062022, 04:33 PM
Post: #7




RE: Square Root Process Similar to Long Division
(09182021 01:07 PM)Albert Chan Wrote: Newton's method converge fast, but required expensive mul/div. We can do as Newton did and use this equation instead: \( \begin{align} \frac{1}{x^2} = d \end{align} \) Using Newton's method we get the following recurrence: \( \begin{align} x \leftarrow x + \frac{x(1  d \cdot x^2)}{2} \end{align} \) So we can get rid of the expensive division and only have to divide by 2. To get \(\sqrt{d}\) we have to multiply the final result \(x\) by \(d\). Program This is the corresponding Python program: Code: def newton(d, x, n): Examples d = Decimal('.9996') x = Decimal('1') newton(d, x, 7) 0.9996 0.99979992 0.99979997999599359936 0.99979997999599899971991597354766255640239090891735007232 0.9997999799959989997199159735914171390272639623863827366463080538247902204151860911590260475681681352 0.9997999799959989997199159735914171390272639623863827366491803235872148569016030640762583498945948583 0.9997999799959989997199159735914171390272639623863827366491803235872148569016030640762583498945948583 d = Decimal('12345') x = Decimal('0.009') newton(d, x, 7) 111.105 111.1080553875 111.1080555135405110305346994140625 111.1080555135405112450044387430752408689311675830734154751136313016974218865863978862762451171875 111.1080555135405112450044387430752414899113774596977299764856731617825851115519126275283305922167024 111.1080555135405112450044387430752414899113774596977299764856731617825903175167629180920706522637384 111.1080555135405112450044387430752414899113774596977299764856731617825903175167629180920706522637384 Depending on the use case only a few iterations might be needed. Still the multiplications might get tedious after a while. 

11062022, 07:25 PM
Post: #8




RE: Square Root Process Similar to Long Division
I learned the longdivisionlike square root algorithm from a book called Mathematics for Statistics, by W. L. Bashaw. The book was in English, a language I had only just started learning at the time, but the algorithm was explained with a workedout example that made it easy to follow. I was a bit surprised that no one else seemed to know this algorithm and it wasn't taught in school, even though calculators were not standard equipment yet at the time.
Regarding the numerical approach, I knew the r > (r + x/r)/2 iteration from applying Newton's method to r^2x=0. It does require one division per step, but since it converges pretty quickly, that didn't seem like a big issue to me. I always assumed that the fourbanger calculators back then used this algorithm, which would make sense because it's easy to implement, and it would be consistent with √ typically being slower than the arithmetic functions, but not very slow. 

« Next Oldest  Next Newest »

User(s) browsing this thread: 1 Guest(s)