How to calculate a square root

I was thinking about on how to calculate a square root.

Usually you are told "use this formula" and although most people can follow simple formulas, I always get that feeling "why are we doing this?". You can always use a calculator and your calculations will often be wrong, they will be just approximations...

Ok approximations is the name of the lesson today, next week we will calculate square roots up to any number of decimals we want, or infinite precision.

How can we do those approximations?

If we are doing approximations I can always give you the next integer or the lower integer of the real answer and you could go happy with it.

Now let us suppose that those numbers are a and b respectively. Given x, we can easily find numbers a and b such that: a <= sqrt(x) <= b. That is almost trivial.

Now we know for sure that sqrt(x) * sqrt(x ) = x.

Let us suppose we take a and do the following:

c = x / a => a * c = x

We don't have sqrt(x), but we have a and c which are close approximations.

And we also know:

a * a < x
b * b > x
c * c > x
a < c < b

What we need now is intuition.

We know sqrt(x) is between a and c.

What would happen if we take a1= (a+c)/2?

If a1 * a1 == x, then sqrt(x) = a1
If a1 * a1 < x, then a = a1, c = x/a, repeat the algorithm
If a1 * a1 > x, then c = a1, a = x/c, repeat the algorithm

So you can approximate easily sqrt(x) just by doing divisions. Actually if you keep finding the median it should converge the double of digits in each step (since you are dividing by 2).

The problem is how close are we going?

When you are calculating sqrt(1), a = c, so this question is irrelevant.

When you are calculating sqtr(10), a = 3, b = 4, c = 10 / 3 = 3.33333..., so a1 = 3.166665, the actual sqrt(10) = 3.16227, pretty close to be the first approximation, but a1 > sqrt(x) even for the first iteration. The problem here is that 10 is a number to small.

When you are calculating sqrt(10000000), a = 3162, c = 3162.55534471853257432005..., so a1 = 3162,277672359266287160025, the actual sqrt = 3162,2776601683793319988935, pretty close to be the first approximation, and still the approximation is higher... (red numbers are the ones we didn't get right).

So calculating the integer part now seems to be the hardest part, but it is actually the simplest.

First remove the border conditions like negative numbers, zero, one, all that can be solved without doing any math.

Then calculate a number a so that a * a <= x and (a+1) * (a+1) > x. To do this first start with a = x / 2 and iterate using binary search.

For example if x = 2, start with a = 1. So a = 1, x = 2, a * a = 1, (a+1) * (a+1) = 4 => Found.

If x = 19, start with a = 9 (the integer part of x/2).

So a = 9, x = 19, a * a = 81, too big => new a = a /2 = 4 (the integer part of a /2)

So a = 4, x = 19, a * a = 16, (a+1) * (a+1) = 25 => Found.

Simple, isn't it?