An iteration for the k -th root with cubic convergence

An iteration for the $k$ -th root with cubic convergence

The other day I encountered a formula for $f$ such that the iteration $y = f.k.a.y$ converged cubically to $y = k$ . Cubic convergence means that in the transition from one estimate to the next, the number of correct leading digits is roughly tripled: a relative error $\delta$ leads to a relative error of the order of $\delta^3$ for the next approximation. This note is not about the merits of this iteration or the lack thereof, but on how it could have been derived. When I saw the formula, it was patently obvious that $V_{a}$ was a fixpoint of $f.k.a$ , the cubic convergence was not obvious at all, and the formula was as far as I was concerned, a Big Rabbit. The purpose of this note is to clo something about these latter issues.

Let $\alpha$ denote our target value for $y$ , i.e. $\alpha = {V}_{a}^{k}$ . Then,as remarked, $\alpha$ has to be a fixpoint of f.k.a. To eliminate the unknown $\alpha$ from that requirement,we observe for any $f,k,a,y$ and $x$ :

$\alpha$ is fix point of $f.k.a$

$$ \begin{array}{l} \begin{array}{r l} {\alpha} & {= \mathrm {f . k . a . a}} \ {\Leftrightarrow} & {{\mathrm {L e i b n i z P r i n c i p l e} }} \end{array} \ y = a \wedge y = f. k. a. y \ \begin{array}{r l} {\equiv} & {{\mathrm {b e c a u s e w e d o n t c o n s i d e r n e g a t i v e o r}} \ {\mathrm {c o m p l e x v a l u e s ,}} & {y = \alpha \equiv y ^ {k} = a }} \end{array} \ y ^ {k} = a \quad \wedge \quad y = f. k. a. y \ \end{array} $$

$$ \begin{array}{l} \mathrm {W h e n} y ^ {k} = a, \mathrm {t h e a b o v e s a y s t h a t} \ y := f. k. a. y \ \end{array} $$

should leave y unchanged , or -in terms of familiar operators plus and times- should increase $y$ by 0 or should multiply $y$ by 1 . For positive y these alternatives are equivalent , we choose the multiplicative version, more precisely we rewrite

$$ f. k. a. y = y \cdot g. k. a. y ^ {k} $$

where $g$ has the property

$$ a = b \Rightarrow g. k. a. b = 1 $$

In this last requirement the constant 1 was introduced because it is the neutral element of the multiplication. For the same reason, the requirement is trivially met by choosing

$$ g. k. a. b = \frac {a}{b} \quad \mathrm {o r} \quad g. k. a. b = \frac {b}{a}, $$

but neither gives a converging iteration: for an approximation $y = \infty \left( {1 + \delta }\right)$ they lead to next approximations

$$ y = \alpha \cdot (1 - (k + 1) \cdot \delta + \dots) \quad a n d \quad y = \alpha \cdot (1 + (k + 1) \cdot \delta + \dots) $$

which are not actually improvements. (With either choice, $y = \alpha$ would be a fixpoint, but the iteration would be unstable.)

of $\delta$ ", which itself is short for something else. (End of Note.)

The question is whether we can reach our goal by suitably "averaging" the last two choices for $g$ considered. To get enough freedom I propose to "average" denominators and numerators separately, and to consider

$$ g \cdot k \cdot a \cdot b = \frac {p \cdot a + q \cdot b}{r \cdot a + s \cdot b} $$

for suitable $p,q,r,s\ldots$ Since these are 4 homogenous parameters, they give us only 3 degrees of freedom, but that is enough to aim for cubic convergence.

More precisely, we are now considering

$$ y := y \cdot \frac {p \cdot a + q \cdot y ^ {k}}{r \cdot a + s \cdot y ^ {k}} $$

To study its convergence we substitute $y = \alpha \cdot (1 + \delta)$ in the right-hand side and remove all canceling factors $\alpha$ . One ends up studying

$$ C _ {0} + C _ {1} \cdot S + C _ {2} \cdot S ^ {2} + \dots $$

when these are the first terms of the Taylor expansion in $s$ of

$$ \frac {(1 + \delta) \cdot (p + q \cdot (1 + \delta) ^ {k})}{r + s \cdot (1 + \delta) ^ {k}} $$

The requirements ${C}{0} = 1,{c}{1} = 0,{c}_{2} = 0$

$\left( 0\right) \;p + q = r + s \neq 0$ (1) $p + q + {qk} - {sk} = 0$ (using (0)) $\left( 2\right) \;{9k}\left( {k + 1}\right) - {sk}\left( {k - 1}\right) = 0\;\left( \text{ using }\left( 0\right) \right)$

With $k \geq 1$ ,(2) gives $q = k - 1,s = k + 1$ , (1) then gives $p = k + 1$ ,and(2) $r = k - 1$ . The iteration with cubic convergence is

$$ y := y \cdot \frac {(k + 1) \cdot a + (k - 1) \cdot y k}{(k - 1) \cdot a + (k + 1) \cdot y k} $$

Looking for other things, I encountered this iteration at the end of a letter to my parents of August 1951. I had added that I was delighted by this discovery, but that was all: no hint of a proof of the cubic convergence and no indication of how (or why) I had derived this formula. So I set myself the task of finding a way in which that could be done. I am sure that at the time I did it differently.

Austn, 29 January 1999

Prof. Dr Edsger W. Dgkstra Department of Computer Sciences The University of Texas at Austin Austin, TX 78712 - 1188 USA