9
$\begingroup$

I found a picture of Evan O'Dorney's winning project that gained him first place in the Intel Science talent search. He proposed a numerical method to find the square root, that gained him $100,000 USD.

Below are some links of pictures of the poster displaying the method.

How does this numerical method work and what is the proof?

His method makes use of Moebius Transformation.

  • 0
    @Lembik: Of course not. It is obvious (and thus long known) that the continued fraction convergents follow a recurrence relation, and clearly that can be computed by repeated squaring of the recurrence matrix. I found this by myself when I was about 18, so it's nothing great if Evan found it on his own. However, I'm sure that he knew it was well known but chose to portray otherwise to the judges, and I find it fishy that the professor who gave him the problem and obviously knew the solution let him do that.2015-07-24

2 Answers 2

8

Essentially if you are interesting in evaluating $\sqrt{a}$, the idea is to first find the greatest perfect square less than or equal to $a$. Say this is $b^2$ i.e. $b = \lfloor \sqrt{a} \rfloor \implies b^2 \leq a < (b+1)^2$. Then consider the function $f(x) = b + \dfrac{a-b^2}{x+b}$ $f(b) = b + \underbrace{\dfrac{a-b^2}{2b}}_{\in [0,1]} \in [b,b+1]$ $f(f(b)) = b + \underbrace{\dfrac{a-b^2}{f(b) + b}}_{\in [0,1]} \in [b,b+1]$ In general $f^{(n)}(b) = \underbrace{f \circ f \circ f \circ \cdots f}_{n \text{times}}(b) = b + \dfrac{a-b^2}{f^{(n-1)}(b)+b}$ Hence, $f^{(n)}(b) \in [b,b+1]$ always.

If $\lim\limits_{n \to \infty}f^{(n)}(b) = \tilde{f}$ exists, then $\tilde{f} = b + \dfrac{a-b^2}{\tilde{f}+b}$ Hence, $\tilde{f}^2 + b \tilde{f} = b \tilde{f} + b^2 + a - b^2 \implies \tilde{f}^2 = a$

To prove the existence of the limit look at $(f^{(n)}(b))^2 - a = \left(b + \dfrac{a-b^2}{f^{(n-1)}(b)+b} \right)^2 - a = \dfrac{(a-b^2)(a-(f^{(n-1)}(b))^2)}{(b+f^{(n-1)}(b))^2} = k_{n-1}(a,b)((f^{(n-1)}(b))^2-a) $ where $\vert k_{n-1}(a,b) \vert \lt1$. Hence, convergence is also guaranteed.

EDIT

Note that $k_{n-1}(a,b) = \dfrac{(a-b^2)}{(b+f^{(n-1)}(b))^2} \leq \dfrac{(b+1)^2 - 1 - b^2}{(b+b)^2} = \dfrac{2b}{(2b)^2} = \dfrac1{2b}$. This can be interpreted as larger the number, faster the convergence.

Comment: This method works only when you want to find the square of a number $\geq 1$.

EDIT

To complete the answer, I am adding @Hurkyl's comment. Functions of the form $g(z) = \dfrac{c_1z+c_2}{c_3z+c_4}$are termed Möbius transformations. With each of these Möbius transformations, we can associate a matrix $M = \begin{bmatrix} c_1 & c_2\\ c_3 & c_4\end{bmatrix}$ Note that the function, $f(x) = b + \dfrac{a-b^2}{x+b} = \dfrac{bx + a}{x+b}$ is a Möbius transformation.

Of the many advantages of the associated matrix, one major advantage is that the associate matrix for the Möbius transformation $g^{(n)}(z) = \underbrace{g \circ g \circ \cdots \circ g}_{n \text{ times}} = \dfrac{c_1^{(n)} z + c_2^{(n)}}{c_3^{(n)} z + c_4^{(n)}}$ is nothing but the matrix $M^n = \begin{bmatrix}c_1 & c_2\\ c_3 & c_4 \end{bmatrix}^n = \begin{bmatrix}c_1^{(n)} & c_2^{(n)}\\ c_3^{(n)} & c_4^{(n)} \end{bmatrix}$ (Note that $c_k^{(n)}$ is to denote the coefficient $c_k$ at the $n^{th}$ level and is not the $n^{th}$ power of $c_k$.)

Hence, the function composition is nothing but raising the matrix $M$ to the appropriate power. This can be done in a fast way since $M^n$ can be computed in $\mathcal{O}(\log_2(n))$ operations. Thereby we can compute $g^{(2^n)}(b)$ in $\mathcal{O}(n)$ operations.

  • 0
    Marvis - Appreciate for your edit, it will take me a extra few hour without your edit!2012-10-28
4

The iteration to find $\sqrt k$ is $f(x) = \frac{d x+k}{x+d}$ where $d = \lfloor \sqrt k \rfloor$. The iterations start with $x = d$.

If $x$ is a fixed point of this, $x = \frac{d x+k}{x+d}$, or $x(x+d) = dx + k$ or $x^2 = k$, so any fixed point must be the square root.

Now wee see if the iteration increases or decreases. If $y = \frac{d x+k}{x+d}$, $y - x = \frac{d x+k}{x+d} - x = \frac{d x+k - x(x+d)}{x+d} = \frac{k - x^2}{x+d} $ so if $x^2 < k$, $y > x$ and if $x^2 > k$, $y < x$.

Also, proceeding like analyses of Newton's method, $y^2-k = \frac{(d x+k)^2}{(x+d)^2} - k = \frac{d^2 x^2 +2 d x k + k^2 - k(x+d)^2}{(x+d)^2} = \frac{d^2 x^2 +2 d x k + k^2 - k(x^2 + 2dx +d^2)}{(x+d)^2} = \frac{d^2 x^2 +2 d x k + k^2 - kx^2 - 2dkx -kd^2)}{(x+d)^2} = \frac{d^2 x^2 + k^2 - kx^2 -kd^2)}{(x+d)^2} = \frac{d^2 (x^2-k) + k^2 - kx^2)}{(x+d)^2} = \frac{d^2 (x^2-k) - k(x^2-k))}{(x+d)^2} = \frac{(d^2-k) (x^2-k)}{(x+d)^2} = (x^2-k)\frac{d^2-k}{(x+d)^2} $.

Since $d = \lfloor \sqrt k \rfloor$, $d < \sqrt k < d+1$ or $d^2 < k < d^2 + 2d +1$ or $-2d - 1 < d^2 - k < 0$, so $|d^2-k| < 2d+1$. Using this, $|y^2-k| < |x^2-k|\frac{2d+1}{(x+d)^2}| = |x^2-k|\frac{2d+1}{x^2+2dx+d^2} $, so $|y^2-k|< |x^2-k|$, and each iteration gets closer to the square root.

Since the starting iterate is $d$, all following iterates exceed $d$ so $|y^2-k| < |x^2-k|\frac{2d+1}{(d+d)^2}| < |x^2-k|\frac{2d+1}{4d^2}| < |x^2-k|\frac{1+1/(2d)}{2d}| \le 3|x^2-k|/4$ since $d \ge 1$.

This show that the iteration converges. However, this does not show that it converges quadratically like Newton's, only that it converges linearly.