Calculating Square Roots - Matrix Algorithm

[danago]

Hi. I only just recently found out about an algorithm for calculating the square roots of a number.

Lets say i want to evaluate $$\sqrt {n}$$. I can make an approximation by inspection, and say $$\sqrt n \approx \frac{a}{b}$$. Now, using this approximation, i can write:

$$\left[ {\begin{array}{*{20}c} 1 & n \\ 1 & 1 \\ \end{array}} \right]\left[ {\begin{array}{*{20}c} a \\ b \\ \end{array}} \right] = \left[ {\begin{array}{*{20}c} {a + bn} \\ {a + b} \\ \end{array}} \right]$$

Treating the resultant matrix as a fraction ($$\frac{{a + bn}}{{a + b}}$$
), i have a better approximation of $$\sqrt {n}$$. If i keep repeating this method with the new approximation, over and over again, i get a more accurate answer. So the next step would be:

$$\left[ {\begin{array}{*{20}c} 1 & n \\ 1 & 1 \\ \end{array}} \right]\left[ {\begin{array}{*{20}c} {a + bn} \\ {a + b} \\ \end{array}} \right] = \left[ {\begin{array}{*{20}c} {a + an + 2bn} \\ {2a + b + bn} \\ \end{array}} \right]$$

And $$\frac{{a + an + 2bn}}{{2a + b + bn}}$$ would be an even better approximation to $$\sqrt {n}$$.

Even if the starting approximation is way off, if a lot of iterations are complete, the answer will still be accurate.


Im just wondering, why does this method work? Is there a name for this method, or anywhere i can look to find more information?

Thanks,
Dan.

 

[Integral]

That sure seems like a complex method.

A Newton's method root finder yields a very simple iterative method:

$$x_{n+1} = \frac 1 2 (x_n + \frac A {x_n})$$

So A is the number you are computing the root of and $x_0$ is your inital guess, it doesn't have to be very good. Compute $x_1$ with the formula, continue.

This procedure converges quadratically.

 

[danago]

Hmm that method does seem good. I am still interested in the matrix method though. I am not interested in it as a way to actually calculate square roots, but as a mathematical phenomenon. Just curious as to why it works.

Ofcouse I am still always going to resort to my calculator to evaluate them :)

 

[tim_lou]

well, you are essentially taking an initial value of x, and then iterate through
$$x_{i+1}=\frac{x_{i}+n}{x_{i}+1}$$
so, let
$$f(x)=\frac{x+n}{x+1}$$
the fix point of f
$$f(x)=x$$
is sqrt(n)

you can check derivatives and find out when this fix point is asymptotically stable.

 

[danago]

well, you are essentially taking an initial value of x, and then iterate through
$$x_{i+1}=\frac{x_{i}+n}{x_{i}+1}$$
so, let
$$f(x)=\frac{x+n}{x+1}$$
the fix point of f
$$f(x)=x$$
is sqrt(n)

you can check derivatives and find out when this fix point is asymptotically stable.

Yep, i understand that i am essentially iterating through $$x_{i+1}=\frac{x_{i}+n}{x_{i}+1}$$, but what do you mean when you say 'fix point'?

 

[danago]

Just researched it on wikipedia, and now understand what a fixed point is.

$$\begin{array}{l} f(x) = \frac{{x + n}}{{x + 1}} \\ \frac{{x + n}}{{x + 1}} = x \\ x^2 + x = x + n \\ x^2 = n \\ x = \sqrt n \\ \end{array}$$

That pretty much exactly answers my question. Thanks very much for that

 

[tim_lou]

hehe, if you want to dig deep. Do you know why this works for All n?

well, the reason, is, the absolute value of f'(x) at the fixed point is less than 1 for all n except n=0. there is a theorem that says (from my difference equation book)

if x* is a fix point for f and f is continuously differentiable at x* then
(1) if |f'(x*)| < 1, then x* is asymptotically stable.
(2) if |f'(x*)| > 1, then x* is unstable.

so that in your case, x* is asymptotically stable.

what if |f'(x*)|=0? then you need to do a lot more tests... if you want to know more. go to the library and pick up a difference equation book...

 

[Hurkyl]

The analysis is a lot easier in the case of multiplying by a matrix.

The eigenvalues of the matrix

$$\left[ {\begin{array}{*{20}c} 1 & n \\ 1 & 1 \\ \end{array}} \right]$$

are

$$1 \pm \frac{\sqrt{n}}{2}$$

If you take (almost) any vector and repeatedly multiply by the matrix A, the result will approach the eigenspace associated to the largest eigenvalue in magnitude. Here, that's $1 + \sqrt{n} / 2$.


Given the rest of the discussion, I expect that the eigenspace associated to the eigenvalue $1 + \sqrt{n} / 2$ is the set of all:

$$a \left[ {\begin{array}{c} \sqrt{n} \\ 1 \end{array}} \right]$$

In fact, I would go so far as to expect A to factor as:

$$\left[ {\begin{array}{cc} 1 & n \\ 1 & 1 \\ \end{array}} \right] = \left[ {\begin{array}{cc} \sqrt{n} & -\sqrt{n} \\ 1 & 1 \\ \end{array}} \right] \left[ {\begin{array}{cc} 1 + \frac{\sqrt{n}}{2} & 0 \\ 0 & 1 - \frac{\sqrt{n}}{2} \\ \end{array}} \right] \left[ {\begin{array}{cc} \sqrt{n} & -\sqrt{n} \\ 1 & 1 \\ \end{array}} \right]^{-1}$$