An iterative algorithm for approximating roots of integers

Fixed-point iteration is a technique commonly used for approximating irrational constants, such as calculating the golden ratio $\phi\approx 1.618$ or finding square roots via the Babylonian method. More generally, the Newton-Raphson method is an iterative process used for estimating the roots of functions, which can be applied to find the roots of certain integers. For example, searching for the real root of $f(x)=x^{3}-5$ will yield $x=\sqrt[3]{5}$.

In this paper, we will introduce an algorithm for estimating rational powers of integers that employs fixed-point iteration, provide examples of its use, and compare it to the existing methods mentioned above.

It is important to note the definition of a “fixed point” in the context of functions and iterative processes, as well as what it means for a fixed point to be “stable.”

In other words, a fixed point of a function occurs wherever the input and output values are equal.

Not all fixed points are of equal importance, however. One important distinction to be made for a given fixed point is whether it can be classified as “stable” or “unstable,” which describes how the function acts in the neighborhood of this point.

The stability of fixed points of a function $f(x)$ can be identified by the magnitude of the derivative at those points [1].

When a fixed point is “superstable,” it is said that the algorithm “converges (at least) quadratically” because the error in the approximation $x_{n+1}$ must be some quadratic function of the error in $x_{n}$.

In the Newton-Raphson method, we see an example of how fixed-point iteration can be used to narrow in on an estimated value for the root of a given integer. Suppose that we begin the process in a slightly different way.

For clarity of this example, suppose we are trying to find an estimated value of $\sqrt{7}$. Rather than defining the function $f(x)=x^{2}-7$ and trying to find the positive root, we simply consider the equation $r^{2}=7$ for some $r$. Ideally, we could divide each side of the equation by $r$ to get $r=7/r$ and iterate as $r_{n+1}=7/r_{n}$ with arbitrary $r_{1}$ to get closer and closer to $r_{n}=\sqrt{7}$ as $n\to\infty$. There is, however, a fatal flaw with this approach: the fixed point of $r=\sqrt{7}$ has only neutral stability, so the value does not converge as we might hope. Instead, it can be easily seen that any nonzero starting value of $r_{1}$ will simply yield values alternating between $r_{1}$ and $7/r_{1}$.

This flaw can be remedied if we offset $r$ by some constant value. Suppose that we rewrite $r=2-c$. (As for where the 2 and the negative sign come from, that will be explained shortly.) Then, we have the initial equation $(2-c)^{2}=4-4c+c^{2}=4+c(c-4)=7\iff c=\frac{7-4}{c-4}$ and can turn this into the fixed-point iteration with

Taking the derivative of the right-hand side with respect to $c$ at the fixed point of $c=2-\sqrt{7}$ gives

The absolute value of the derivative at the fixed point is less than 1, so it is stable! This means that repeatedly iterating on the value of $c$ will cause the value to tend towards $2-\sqrt{7}$, at which point we may subtract $2-c$ to recover our approximation for $\sqrt{7}$ as desired.

This is where the earlier rewriting of $r=2-c$ is important. Had we chosen a different constant or sign, it is entirely possible that the fixed point would not have been stable, ruining the results. In the following section, we generalize our results and make clear exactly how to choose this constant to make the algorithm work (and optimize its rate of convergence).

First, we want to generalize to arbitrary input $x\in\mathbb{R}^{+}$, not just the number 7 used in the previous example. (Note that each of the examples given in this paper use integer $x$, but there is no limitation that prevents calculating $\sqrt{\pi}$, for example. The main drawback is simply that it is more difficult to calculate by hand and of less interest to the general reader.) Start with the quadratic equation

where $b/2$ is our “offset” value and we can solve to get $\frac{b}{2}-c=\pm\sqrt{x}$. Next, rewrite equation (1) in expanded form as

and perform algebraic manipulation to obtain

which can be easily transformed into the fixed-point iteration of

From here, it is easily verified that the two fixed points for $c$ are given by $c_{-}=\frac{b}{2}-\sqrt{x}$ and $c_{+}=\frac{b}{2}+\sqrt{x}$. Determining which of these is stable is fundamental in knowing how to proceed.

We define $f(c)=\frac{x-(b^{2}/4)}{c-b}$ and take the absolute value of the derivative to get

We first check the fixed point of $c_{+}=\frac{b}{2}+\sqrt{x}$. Substituting yields

Consider the case in which we have picked $b$ such that $0<b<2\sqrt{x}$. We can rewrite $\sqrt{x}=b/2+\varepsilon$ for some $\varepsilon>0$. Then, we can simplify to

In this case, then, the fixed point $c_{+}$ is unstable. Consider instead the case of $0<2\sqrt{x}<b$. We rewrite $\sqrt{x}=b/2-\varepsilon$ for some $0<\varepsilon<b/2$. Then, we can simplify to

where the final inequality holds true because of the bounds on $\varepsilon$. This is unstable once more.

Next, let us consider the case of $-2\sqrt{x}<b<0$. Rewriting $\sqrt{x}=-b/2+\varepsilon$ for some value $\varepsilon>0$ gives

Unlike the first two cases, this shows that the fixed point $c_{+}$ is stable for such a selection of $b$.

Fourth, we consider the case of $b<-2\sqrt{x}<0$. Rewriting $\sqrt{x}=-b/2-\varepsilon$ for some value $0<\varepsilon<-b/2$ gives

Once again, the fixed point is stable.

The last inequality here may not be immediately clear, but it can be observed that the value of the denominator lies within the open interval of $b$ (which is greater in magnitude than any possible value of $\varepsilon$) to $b/2$, at which point the numerator and denominator are equal in magnitude.

Now, we address the few remaining special cases. When we have the exact equality $b=2\sqrt{x}$, we see $|f^{\prime}(c_{+})|\to\infty$, making it unstable. When $b=0$, we have that $|f^{\prime}(c_{+})|=1$, indicating a neutral stability (the same situation described as a “flaw” in the beginning of section 2). Finally, when we have $b=-2\sqrt{x}$, we see $|f^{\prime}(c_{+})|=0$, making the point superstable.

Nearly identical calculations reveal that whenever the point $c_{+}$ is stable, the point $c_{-}$ is unstable, and vice versa. A summary table of all the stability results above for the different possible values of $b$ is provided below:

Essentially, the key takeaway is that there is a kind of symmetry about $b=0$ with respect to the stability of these fixed points. To make things simpler, we can simply choose to consider strictly positive $b$ and deal only with $c_{-}$, or we may choose to consider only negative $b$ and the corresponding $c_{+}$; either is equally valid. In this analysis, we opt with the former.

There are two connections we will draw. Namely, we will discuss connections to continued fractions and the Newton-Raphson method mentioned before.

We recall that the special case $d=2$ in equation (4) results in (2):

In this case, we see that the recurrence relation allows for an infinite generalized continued fraction to be constructed [2].

In 3.2, the numerators were each 1, making a “simple” continued fraction.

Next, the Newton-Raphson method is another technique used for approximating roots of positive real numbers, and it also relies on fixed-point iteration in its calculations.

Suppose we want to find an approximation for $\sqrt{7}$. We know that we can achieve this by searching for the positive root of $f(x)=x^{2}-7$. By setting $x_{1}=1$ and performing successive linear approximations with

we will have $x_{n}$ converge to the desired value. (This is equivalent to the Babylonian method, after algebraic simplification.)

In fact, we see that at the fixed point $x_{n}=\sqrt{7}$, the derivative of the right-hand side with respect to $x_{n}$ is zero. This indicates that the Newton-Raphson method converges at least quadratically (the number of accurate decimal places roughly doubles in each iteration), which is better than the method explored in this paper, which only converges linearly (since the derivative is generally nonzero, as discussed in the section on stability). This can be explained by the theory of Taylor expansion [3].

Thus, the Newton-Raphson method and higher-order analogues like Halley’s method experience a much faster rate of convergence, with the trade-off being that elementary differentiation of the polynomial $f(x)$ is necessary for their computation. Since this polynomial is always of the form $f(x_{n})=x_{n}^{d}-k$ with degree $d$ and constant $k$, we know that $f^{\prime}(x_{n})$ in right-hand side of (5) will always be equal to $dx_{n}^{d-1}$ by the power rule. For Halley’s method and higher-order root-finding algorithms, the higher-order derivatives can be similarly computed without much effort.