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.