Square Roots
June 1, 2012
In today’s exercise we look at four different algorithms for computing square roots of positive floating-point numbers.
The first method is bisection. For any x, the √x is between 1 and x. The bisection method repeatedly computes the square of the midpoint between the two and discards the half of the range in which the square root doesn’t appear; when the two endpoints are sufficiently close, their midpoint is returned as the result of the function.
In the first century, the Greek mathematician Hero of Alexandria (his name is often spelled Heron) devised an iterative method for calculating √x. If xk is a good approximation to √x, then the average of xk and x/xk is a better approximation; the accuracy doubles at each iterative step. The idea behind the algorithm is that √x is the geometric mean of 1 and x, but not knowing the geometric mean we use the arithmetic mean instead in a series of approximations.
The third method was invented by Sir Isaac Newton, and is again based on an iterative method. If xk is a good approximation to a function, then xk+1 = xk – f(xk) / f'(xk) is a better approximation. Knowing that 2x is the derivative of x2 makes it easy to compute the sequence of approximations.
The fourth method is an optimization that can be applied to either Heron’s or Newton’s methods. The input number x is reduced to the range [1,2) by successively multiplying or dividing by 2. Then, with the range so small, the loop is unrolled to a fixed number of iterations and the result divided or multiplied by the square root of 2 the same number of times as the original multiplication or division. This is fast because multiplying and dividing by 2 is easy for binary computers, and can be made even faster by table lookup of the starting point for the iteration, saving one or two steps.
Your task is to implement the four square root algorithms described above. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
Not sure if this is something I can fix client-side (please let me know if so), but I recommend for future posts using a more appropriate font (or some TeX variant like you are using on the solutions page for this exercise) when presenting numeric symbology.
For instance, the differentiation symbol used on the denominator of Newton’s method is barely visible. So at first glance, it looked like you made a typo and we were to perform “x_{k+1} = x_{k} – 1”.
Maybe its just the italics that makes it hard to read, but if future exercises use any equations not so familiar to us (unliked Newton’s method), it would be disappointing to work with incorrect equations simply because the typesetting was too difficult to read.
Just my two cents, and keep up the great exercises!
P.S. I would like to see some of these solutions come up with a clever approximation seed used for the Alexandrian/Babylonian method. Wikipedia has a good example, but I would be interested if there are any better solutions geared towards computer algebra.
Ardnew: Regarding math: WordPress, where I host the blog, does a bad job with math. They provide LaTeX with the dollarsign-latex command, but then it doesn’t resize with the text, so it looks more-or-less ugly depending on the text size. And it uses LaTeX inline mode even in displays, which is just wrong, as you will notice on the solution page of this exercise. Writing out math using HTML characters can work either well or poorly depending on the browser — it looks much better on my browser at home, where I generally write these exercises, than at work (both are Chrome, but they render quite differently); when I have trouble, I increase the text size in the browser. The proper solution is jsmath, but that requires me to host my own web site, which I prefer not to do, as I prefer to spend my limited time writing exercises than administering a web site (and since my site generates no revenue, I prefer to use the free hosting at WordPress). Thanks for the suggestion, but I think my best advice is just to increase the text size in your browser (or squint). If that doesn’t work, I have noticed that printing the exercise often provides a much better visual experience than reading it on-screen. And if you are still unclear, you can always view source (but be prepared to gag at all the junk that WordPress adds).
The clever approximation method that you mention usually takes the number of digits in x, divides by 2 (integer division) and subtracts 1, then sticks a 2 on front of that many zeroes if the original number of digits was even and a 6 on front of that many zeros if the original number of digits was odd. For instance, our example 125348 has six digits, then 6/2-1 = 2 so we have two zeroes, and the approximation is 200. A similar procedure is followed when x is less than 1. Fancier approximation methods have more starting points in the “lookup table” than 2 or 6. But the effect of a clever approximation is seldom more than a reduction of a single step in the iteration, and since the cost of computing the clever approximation is about the same as the cost of a single step, there is generally little saving. Where good approximations can work is where the range of x is sufficiently restricted (say, 32-bit integers) that a good approximation can save two or three steps in the iteration. If you use the “optimized” method mentioned in the exercise that restricts x to the range 1 ≤ x < 2 (or 4), it might be possible with some experimenting to find a clever approximation that actually does pay for itself; perhaps you would like to study this topic and report back to the rest of us.
Note that heron and newton are essentially the same (they generate the same iterates in exact arithmetic).
my version, i didn’t implement the last one because i was unclear of the algorithm.
http://pastebin.com/3YVeMVCY
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] 牛顿的方法在整数上工作得很好： class=”lang-python prettyprint-override”> def isqrt(n): x = n y = (x + 1) // 2 while y < x: x = y y = (x + n // x) // 2 return x [/code] 这将返回x * x不超过n的最大整数x。如果你想检查结果是否恰好是平方根，只需执行乘法来检查n是否是一个完美的正方形。 我讨论这个算法，以及另外三个计算平方根的算法my blog […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]
[…] Discuto este algoritmo, y otros tres algoritmos para calcular raíces cuadradas, en mi blog. […]
[…] I discuss this algorithm, and three other algorithms for calculating square roots, at my blog. […]