Newton, Halley, Pell and the Optimal Iterative High-Order Rational Approximation of √N

Isaac Fried

Department of Mathematics, Boston University, Boston, MA, USA.

DOI: 10.4236/am.2018.97059                

Abstract

In this paper we examine single-step iterative methods for the solution of the nonlinear algebraic equation f (x) = x² - N = 0 , for some integer N, generating rational approximations p/q that are optimal in the sense of Pell’s equation p² - Nq² = k for some integer k, converging either alternatingly or oppositely.

Keywords

Iterative Methods, Super-Linear and Super-Quadratic Methods, Square Roots, Pell’s Equation, Optimal Rational Iterants, Root Bounds

1. Introduction

We present ever higher order single-point iterative methods for the numerical solution of the nonlinear equation $f\left(x\right)=0$. Then we show that for $f\left(x\right)=x^2-N$ these methods are optimal in the sense of Pell’s equation (see [1] [2] [3] [4] ), namely, that if the initial guess $x_0=p_0/q_0$ satisfies the diophantine Pell’s equation $p_0^2-N{q}_0^2=k$ , for some integer k, then the iterated value $x_1=p_1/q_1$ , obtained by a method of order n, satisfies the Pell equation $p_1^2-N{q}_1^2=k^n$.

Using a generalization of the recursive solution to Pell’s equation we generate super-linear and super-quadratic methods that converge alternatingly and oppositely to provide upper and lower bounds on the targeted root (see [5] [6] [7] ).

2. Pell’s Equation

Let N be a positive integer which is not a square. The pair of natural numbers p, q satisfying the general Pell’s equation (see [1] [2] [3] [4] ):

$p^2-N{q}^2=k$ (1)

are such that

${\left(\frac{p}{q}\right)}^2=N+\frac{k}{q^2}\text{or}\frac{p}{q}=\sqrt{N}\left(1+\frac{k}{2N{q}^2}\right)$ (2)

nearly, if $k/\left(N{q}^2\right)\ll 1$. If k > 0, then p/q is an overestimate of $\sqrt{N}$ , and if k < 0, then p/q is an underestimate of $\sqrt{N}$.

We verify (see [1] , Chapter 32) that a new solution pair $\left(p_1\mathrm{,}{q}_1\right)$ to the optimal, or minimal ( $k=\pm 1$ ), Pell’s equation is obtained from a known solution pair $\left(p_0\mathrm{,}{q}_0\right)$ by the expansion of

$p_1+\sqrt{N}{q}_1={\left(p_0+\sqrt{N}{q}_0\right)}^n\mathrm{,}$ (3)

where variable n is taken odd for $k=-1$.

For example, if we take in Equation (3) n = 2, then (see [8] [9] , and also [10] [11] )

$p_1+\sqrt{N}{q}_1=\left(p_0^2+N{q}_0^2\right)+\sqrt{N}\left(2p_0{q}_0\right)$ (4)

and

$p_1=p_0^2+N{q}_0^2,q_1=2p_0{q}_0,\text{or}{x}_1=\frac{x_0^2+N}{2x_0},x_1=\frac{p_1}{q_1},x_0=\frac{p_0}{q_0},$ (5)

which is Newton’s method, preferably written as

$x_1=x_0-\frac{x_0^2-N}{2x_0},x_1=x_0-\frac{f\left(x_0\right)}{f^{\prime }\left(x_0\right)},f\left(x\right)=x^2-N.$ (6)

Here

$p_1^2-N{q}_1^2={\left(p_0^2-N{q}_0^2\right)}^2\mathrm{.}$ (7)

3. Super-Linear Iterative Method

We start with the general power series expansion around $\sqrt{N}$

$\begin{array}{l}{x}_1=F\left(x_0\right)=F\left(\sqrt{N}\right)+F^{\prime }\left(\sqrt{N}\right)\left(x-\sqrt{N}\right)+\frac{1}{2!}{F}^{″}\left(\sqrt{N}\right){\left(x-\sqrt{N}\right)}^2\\ +O\left({\left(x-\sqrt{N}\right)}^3\right)\end{array}$ (8)

and ask that $\sqrt{N}=F\left(\sqrt{N}\right)$ , or that $\sqrt{N}$ is a fixed-point of iteration function $F\left(x\right)$. Now we pass from the general to the specific

$x_1=\frac{A{x}_0+B}{C{x}_0+D}=\frac{A\sqrt{N}+B}{C\sqrt{N}+D}+\frac{AD-BC}{{\left(C\sqrt{N}+D\right)}^2}\left(x_0-\sqrt{N}\right)+O\left({\left(x_0-\sqrt{N}\right)}^2\right)$ (9)

for parameters $A\mathrm{,}B\mathrm{,}C\mathrm{,}D$ , and we ask, here specifically, that

$\sqrt{N}=\frac{A\sqrt{N}+B}{C\sqrt{N}+D},$ (10)

or, again, that $\sqrt{N}$ is a fixed-point of the rational iteration function in Equation (9). To satisfy Equation (10) we take $B=NC,D=A$ , and are left with

$x_1=\frac{A{x}_0+NC}{C{x}_0+A},x_1-\sqrt{N}=\frac{A/C-\sqrt{N}}{A/C+\sqrt{N}}\left(x_0-\sqrt{N}\right),A/C>0$ (11)

in which $A/C+\sqrt{N}$ in the second denominator is written for $A/C+x_0$.

Writing $x_0=p_0/q_0$ and $x_1=p_1/q_1$ , the iterative method assumes the form

$x_1=\frac{p_1}{q_1}=\frac{A{p}_0+CN{q}_0}{C{p}_0+A{q}_0}\text{or}{p}_1=A{p}_0+CN{q}_0,q_1=C{p}_0+A{q}_0,x_0=\frac{p_0}{q_0},C\ne 0$ (12)

for parameters A and C. Referring to Equation (12) we have

Lemma 1. If $p_0\mathrm{,}{q}_0$ are such that $p_0^2-N{q}_0^2=k$ , and $A^2-N{C}^2=m$ , then $p_1\mathrm{,}{q}_1$ are such that $p_1^2-N{q}_1^2=km$.

Proof. We verify that

$p_1^2-N{q}_1^2={\left(A{p}_0+NC{q}_0\right)}^2-N{\left(C{p}_0+A{q}_0\right)}^2=\left(A^2-N{C}^2\right)\left(p_0^2-N{q}_0^2\right)=km,$ (13)

and the result follows.

For instance, if $N=2,C=A=1,p_0=3$ and $q_0=2$ , then

$p_1=p_0+2q_0=7,q_1=p_0+q_0=5\text{and}{p}_1^2-2q_1^2=-1.$ (14)

Observe that the iterative method (11) converges linearly for any $A/C>0$ , since then

$-1<\frac{A/C-\sqrt{N}}{A/C+\sqrt{N}}<1$ (15)

and $\left|x_1-\sqrt{N}\right|<\left|x_0-\sqrt{N}\right|$.

4. Alternating Convergence

If in Equation (11), $A/C>\sqrt{N}$ , then $x_0-\sqrt{N}$ and $x_1-\sqrt{N}$ are of the same sign, but if $A/C<\sqrt{N}$ , they are of opposite signs. Also, the smaller $\left|A/C-\sqrt{N}\right|$ , the faster the convergence.

The method of Equation (11), as well as higher order methods, can be derived directly, in reverse, from the generalized Equation (3)

$p_1+\sqrt{N}{q}_1={\left(p+\sqrt{N}q\right)}^n{\left(p_0+\sqrt{N}{q}_0\right)}^m$ (16)

with n = 1 and m = 1. Indeed, expansion of Equation (16) brings it to the form

$p_1+\sqrt{N}{q}_1=\left(p{p}_0+Nq{q}_0\right)+\sqrt{N}\left(q{p}_0+p{q}_0\right),$ (17)

which elicits the pair of equations

$p_1=p{p}_0+Nq{q}_0,q_1=q{p}_0+p{q}_0$ (18)

with A = p and C = q in Equation (11).

For example, taking in Equation (11) $N=7$ , $p_0/q_0=8/3$ , $p_0^2-7q_0^2=1$ , $p/q=A/C=5/2$ , $p^2-7q^2=-3$ , we obtain from it the alternating sequence of convergents:

$x_1=\frac{p_1}{q_1}=\left\{\frac{8}{3},\frac{82}{31},\frac{844}{319},\frac{8686}{3283},\frac{89392}{33787},\frac{919978}{347719},\cdots \right\}$ (19)

with

${\left(\frac{p_1}{q_1}\right)}^2=\left\{7.1,6.997,7.00009,6.9999975,7.00000007,6.9999999979,\cdots \right\}$ (20)

5. The Method of Newton and Its Opposites

Taking in Equation (9) $A/C=x_0$ , or $A=p_0$ and $C=q_0$ , the linear method rises to become the quadratic method of Newton, otherwise directly obtainable from Equation (3) with n = 2 (see Equations (5)-(7)).

Here, for Newton’s method

$x_1-\sqrt{N}=\frac{1}{2\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^2$ (21)

nearly, if $x_0$ is close to $\sqrt{N}$.

The method

$\frac{p_1}{q_1}=2N\frac{p_0{q}_0}{p_0^2+N{q}_0^2},\text{or}{x}_1=\frac{2N{x}_0}{x_0^2+N}$ (22)

is such that

$p_1^2-N{q}_1^2=-N{\left(p_0^2-N{q}_0^2\right)}^2,$ (23)

or

$x_1-\sqrt{N}=-\frac{1}{2\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^2$ (24)

if $x_0$ is close to $\sqrt{N}$. Here, convergence is quadratic and from below. Compare Equations ((21) and (24)).

The average of methods (5) and (22)

$x_1=\frac{1}{2}\left(\frac{x_0^2+N}{2x_0}+\frac{2N{x}_0}{x_0^2+N}\right)$ (25)

is quartic

$x_1-\sqrt{N}=\frac{1}{8N\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^4.$ (26)

Or

$\frac{p_1}{q_1}=\frac{p_0^4+6N{p}_0^2{q}_0^2+N^2{q}_0^4}{4p_0{q}_0\left(p_0^2+N{q}_0^2\right)},\text{and}{p}_1^2-N{q}_1^2={\left(p_0^2-N{q}_0^2\right)}^4.$ (27)

For example, for N = 2 and $x_0=3/2$ we obtain from method (6) $x_1=17/12$ , from method (22) $x_1=24/17$ , and for their average $x_1=577/408$ , and

${17}^2-2\times {12}^2=1,{24}^2-2\times {17}^2=-2,{577}^2-2\times {408}^2=1.$ (28)

Here, ${\left(17/12\right)}^2=2.007,{\left(24/17\right)}^2=1.993,{\left(577/408\right)}^2=2.000006$.

The biased average method

$x_1=\frac{x_0^2+N}{2x_0}\left(\frac{1}{2}-ϵ\right)+\frac{2N{x}_0}{x_0^2+N}\left(\frac{1}{2}+ϵ\right),ϵ=\frac{1}{16N^2}{\left(x_0^2-N\right)}^2$ (29)

produces an oppositely converging quartic method such that, asymptotically

$x_1-\sqrt{N}=-\frac{1}{8N\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^4.$ (30)

Compare Equations ((26) and (30)).

The biased average method

$x_1=\frac{x_0^2+N}{2x_0}\left(\frac{1}{2}-ϵ\right)+\frac{2N{x}_0}{x_0^2+N}\left(\frac{1}{2}+ϵ\right),ϵ=\frac{1}{32N^2}{\left(x_0^2-N\right)}^2$ (31)

is a quintic method and such that

$x_1-\sqrt{N}=-\frac{1}{4N^2}{\left(x_0-\sqrt{N}\right)}^5+O\left({\left(x_0-\sqrt{N}\right)}^6\right),$ (32)

implying that the convergence of method (31) is alternating. Indeed, starting with $x_0=3/2$ we obtain from method (31)

$x_0^2=2+0.25,x_1^2=2-7.6\times {10}^{-7},x_1^2=2+2.5\times {10}^{-34}.$ (33)

6. More Convergence from Below

The noteworthy method

$x_1=\frac{3N-x_0^2}{2N}{x}_0$ (34)

converges to $\sqrt{N}$ quadratically and from below,

$x_1-\sqrt{N}=-\frac{3}{2\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^2.$ (35)

We write $x_0=p_0/q_0$ and $x_1=p_1/q_1$ and have for Equation (34) that

$p_1^2-N{q}_1^2=\left(p_0^2-4N{q}_0^2\right){\left(p_0^2-N{q}_0^2\right)}^2.$ (36)

7. Super-Linear Alternating Methods

We put in Equation (11)

$A=x_0\left(1+2ϵ\right),C=1,$ (37)

and obtain

$x_1-\sqrt{N}=ϵ\left(x_0-\sqrt{N}\right)+\frac{1}{2\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^2$ (38)

nearly, if $x_0$ is close to $\sqrt{N}$ and $ϵ\ll 1$ , the super-linear method

$x_1=x_0-\frac{x_0^2-N}{2x_0\left(1+ϵ\right)}\text{or}{x}_1=x_0-\frac{x_0^2-N}{2x_0}\left(1-ϵ\right).$ (39)

A small negative $ϵ$ causes method (39) to ultimately oscillate, or alternate.

With

$ϵ=-\frac{1}{4N}\left(x_0^2-N\right)$ (40)

method (39) becomes cubic and of alternating convergence

$x_1-\sqrt{N}=-\frac{3}{4N}{\left(x_0-\sqrt{N}\right)}^3+O\left({\left(x_0-\sqrt{N}\right)}^4\right).$ (41)

8. Stacked Methods

From

$\begin{array}{c}{p}_2+\sqrt{N}{q}_2=\left(p_0+\sqrt{N}{q}_0\right)\left(p_1+\sqrt{N}{q}_1\right)\\ =\left(p_0{p}_1+N{q}_0{q}_1\right)+\sqrt{N}\left(p_0{q}_1+p_1{q}_0\right)\end{array}$ (42)

we have the stacked method

$\frac{p_2}{q_2}=\frac{p_0{p}_1+N{q}_0{q}_1}{p_0{q}_1+p_1{q}_0},$ (43)

or

$x_2=\frac{x_0{x}_1+N}{x_0+x_1}.$ (44)

It is such that if

$x_0=\sqrt{N}+{ϵ}_0,\text{and}{x}_1=\sqrt{N}+{ϵ}_1,$ (45)

then

$x_2=\sqrt{N}+\frac{{ϵ}_0{ϵ}_1}{2\sqrt{N}}=\sqrt{N}+\frac{1}{2\sqrt{N}}\left(x_0-\sqrt{N}\right)\left(x_1-\sqrt{N}\right)$ (46)

nearly, if both epsilons are small compared with $\sqrt{N}$.

If ${ϵ}_0{ϵ}_1<0$ , then $x_2<\sqrt{N}$ , and if ${ϵ}_0{ϵ}_1>0$ , then $x_2>\sqrt{N}$. For example, for N = 2 we obtain from the stacked method of Equation (44) the alternatingly converging sequence

$x_2=\left\{\frac{1}{1},\frac{3}{2},\frac{7}{5},\frac{41}{29},\frac{577}{408},\frac{47321}{33461},\cdots \right\}$ (47)

with

$x_2^2=\left\{1,2.25,1.96,1.9988,2.000006,1.9999999991,\cdots \right\}$ (48)

9. Halley’s Third-Order Method

Halley’s cubic iterative method

$x_1=x_0-\frac{\det \left[\begin{array}{cc}1& f_0\\ 0& 2{f^{\prime }}_0\end{array}\right]}{\det \left[\begin{array}{cc}{f^{\prime }}_0& f_0\\ {f^{″}}_0& 2{f^{\prime }}_0\end{array}\right]}\cdot f_0$ (49)

becomes for $f\left(x\right)=x^2-N$ and $x_0=p_0/q_0$

$x_1=\frac{p_1}{q_1},p_1=p_0\left(p_0^2+3N{q}_0^2\right),q_1=q_0\left(3p_0^2+N{q}_0^2\right),$ (50)

and is verified to be such that

$p_1^2-N{q}_1^2={\left(p_0^2-N{q}_0^2\right)}^3=k^3\text{if}{p}_0^2-N{q}_0^2=k,$ (51)

implying that if $p_0/q_0$ is an underestimate (k < 0), then so is $p_1/q_1$ , and if $p_0/q_0$ is an overestimate (k > 0), then so is $p_1/q_1$.

Otherwise, here

$x_1-\sqrt{N}=\frac{1}{4N}{\left(x_0-\sqrt{N}\right)}^3$ (52)

nearly, if $p_0/q_0$ is close to $\sqrt{N}$.

10. Fourth-Order Method

The quartic method (see [12] [13] for higher order methods):

$x_1=x_0-\frac{\det \left[\begin{array}{ccc}1& f_0& 0\\ 0& 2{f^{\prime }}_0& f_0\\ 0& 3{f^{″}}_0& 3{f^{\prime }}_0\end{array}\right]}{\det \left[\begin{array}{ccc}{f^{\prime }}_0& f_0& 0\\ {f^{″}}_0& 2{f^{\prime }}_0& f_0\\ {f^{‴}}_0& 3{f^{″}}_0& 3{f^{\prime }}_0\end{array}\right]}\cdot f_0$ (53)

becomes for $f\left(x\right)=x^2-N$ and $x_0=p_0/q_0$

$x_1=p_1{q}_1,p_1=p_0^4+6N{p}_0^2{q}_0^2+N^2{q}_0^4,q_1=4p_0{q}_0\left(p_0^2+N{q}_0^2\right),$ (54)

observed to be a repeated second order method and such that

$p_1^2-N{q}_1^2={\left(p_0^2-N{q}_0^2\right)}^4\mathrm{.}$ (55)

Otherwise, here

$x_1-\sqrt{N}=\frac{1}{8N\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^4$ (56)

if $x_0$ is close to $\sqrt{N}$. Convergence here is from above.

11. Fifth-Order Method

The quintic method

$x_1=x_0-\frac{\det \left[\begin{array}{cccc}1& f_0& & \\ 0& 2{f^{\prime }}_0& f_0& \\ 0& 3{{f^{\prime }}^{\prime }}_0& 3{f^{\prime }}_0& f_0\\ 0& 4{{{f^{\prime }}^{\prime }}^{\prime }}_0& 6{{f^{\prime }}^{\prime }}_0& 4{f^{\prime }}_0\end{array}\right]}{\det \left[\begin{array}{cccc}{f^{\prime }}_0& f_0& & \\ {{f^{\prime }}^{\prime }}_0& 2{f^{\prime }}_0& f_0& \\ {{{f^{\prime }}^{\prime }}^{\prime }}_0& 3{{f^{\prime }}^{\prime }}_0& 3{f^{\prime }}_0& f_0\\ {{{{f^{\prime }}^{\prime }}^{\prime }}^{\prime }}_0& 4{{{f^{\prime }}^{\prime }}^{\prime }}_0& 6{{f^{\prime }}^{\prime }}_0& 4{f^{\prime }}_0\end{array}\right]}\cdot f_0$ (57)

becomes for $f\left(x\right)=x^2-N$ and $x_0=p_0/q_0$

$x_1=\frac{p_1}{q_1},p_1=p_0\left(p_0^4+10N{p}_0^2{q}_0^2+5N^2{q}_0^4\right),q_1=q_0\left(5p_0^4+10N{p}_0^2{q}_0^2+N^2{q}_0^4\right),$ (58)

and happens to be such that

$p_1^2-N{q}_1^2={\left(p_0^2-N{q}_0^2\right)}^5\mathrm{.}$ (59)

Otherwise, here

$x_1-\sqrt{N}=\frac{1}{16N^2}{\left(x_0-\sqrt{N}\right)}^5$ (60)

if $x_0$ is close to $\sqrt{N}$.

12. A Rational Quadratic Method

The method

$x_1=\frac{2\sqrt{2}+4x_0+\sqrt{2}{x}_0^2}{2+2\sqrt{2}{x}_0+x_0^2}$ (61)

is merely $x_1=\sqrt{2}$ in disguise. Replacement of $\sqrt{2}$ by the good rational approximation p/q turns the scheme into

$x_1=\frac{2p+4q{x}_0+p{x}_0^2}{2q+2p{x}_0+q{x}_0^2},$ (62)

and for the specific $p/q=7/5,7^2-2\times 5^2=-1$ , it becomes

$x_1=\frac{p_1}{q_1}=\frac{14+20x_0+7x_0^2}{10+14x_0+5x_0^2},x_0=\frac{p_0}{q_0},p_1^2-2q_1^2=-{\left(p_0^2-2q_0^2\right)}^2.$ (63)

Starting with $x_0=7/5$ we obtain $x_1=239/169$ , $x_1^2=1.999965$. Starting with $x_0=17/12$ we obtain $x_1=8119/5741$ , $x_1^2=1.99999997$. Then

$\begin{array}{l}x=1855077841/1311738121,{1855077841}^2-2\times {1311738121}^2=-1\\ {\left(1855077841/131173812\right)}^2=\mathrm{1.99999999999999999942.}\end{array}$ (64)

From

$x_1=\frac{p_1}{q_1}=\frac{6+8x_0+3x_0^2}{4+6x_0+2x_0^2},x_0=\frac{p_0}{q_0},p_1^2-2q_1^2={\left(p_0^2-2q_0^2\right)}^2,$ (65)

obtained from Equation (62) with $p/q=3/2,3^2-2\times 2^2=1$ , we compute

$x_1=\left\{3/2,99/70,114243/80782,152139002499/107578520350\right\}$ (66)

with

$\begin{array}{l}{152139002499}^2-2\times {107578520350}^2=1\\ {\left(152139002499/107578520350\right)}^2=\mathrm{2.000000000000000000000086.}\end{array}$ (67)

13. The General Rational Super-Quadratic Method

We start by writing

$x_1=\frac{p_1}{q_1}=\frac{A{x}_0^2+B{x}_0+C}{P{x}_0^2+Q{x}_0+R}$ (68)

to have

$x_1-\sqrt{N}=\frac{p_1\left(x_0\right)-q_1\left(x_0\right)\sqrt{N}}{q_1\left(x_0\right)}.$ (69)

To have a factor ${\left(x_0-\sqrt{N}\right)}^2$ in the numerator of the right-hand side of Equation (69), we ask that

$p_1\left(x\right)-q_1\left(x\right)\sqrt{N}=0,\text{and}\text{that}{\left(p_1\left(x\right)-q_1\left(x\right)\sqrt{N}\right)}^{\prime }=0,\text{at}x=\sqrt{N}$ (70)

resulting in

$P=1,B=2N,R=N,C=AN,Q=2A,$ (71)

and the method

$x_1-\sqrt{N}=\frac{A-\sqrt{N}}{x_0^2+2A{x}_0+N}{\left(x_0-\sqrt{N}\right)}^2$ (72)

that can be raised to cubic with the choice $A=x_0$.

Instead, we leave $A=p/q,x_0=p_0/q_0$ to have the method

$x_1=\frac{p_1}{q_1}=\frac{A{x}_0^2+2N{x}_0+AN}{x_0^2+2A{x}_0+N},A=\frac{p}{q},x_0=\frac{p_0}{q_0}$ (73)

such that

$p_1^2-N{q}_1^2=\left(p^2-N{q}^2\right){\left(p_0^2-N{q}_0^2\right)}^2.$ (74)

For example, for $N=7,A=8/3,x_0=8/3,8^2-7\times 3^2=1$ , we obtain from Equation (73)

$\begin{array}{l}{x}_1=\frac{2024}{765},x_1^2=7.0000017,{2024}^2-7\times {765}^2=1,x_1=\frac{130576328}{49353213},\\ x_1^2=7.0000000000000004,{130576328}^2-7\times {49353213}^2=1.\end{array}$ (75)

For $N=7,A=5/2,x_0=5/2,5^2-7\times 2^2=-3$ , we obtain from Equation (73)

$\begin{array}{l}{x}_1=\frac{545}{206},x_1^2=6.99936,{545}^2-7\times {206}^2=-27,x_1=\frac{6113945}{2310854},\\ x_1^2=6.9999999996,{6113945}^2-7\times {2310854}^2=-2187.\end{array}$ (76)

Equation (73), as well as higher order methods, could have been derived directly, in reverse, from

$p_1+\sqrt{N}{q}_1={\left(p+\sqrt{N}q\right)}^n{\left(p_0+\sqrt{N}{q}_0\right)}^m$ (77)

with $n=1,m=2$.

14. The Direct Construction of a Super-Quadratic Method

To locate root a of $f\left(x\right)$ , $f\left(a\right)=0$ , we start by writing the fixed-point iterative method

$x_1=F\left(x_0\right),F\left(x\right)=x+Af\left(x\right)+B{f}^2\left(x\right)$ (78)

for constants A and B. Then we require that

$F^{\prime }\left(a\right)=0,F^{″}\left(a\right)=ϵ,$ (79)

where $ϵ$ is any parameter.

Differentiating $F\left(x\right)$ once and twice, the previous system of two equations in the two unknowns A and B becomes

$\left[\begin{array}{cc}{f}^{\prime }& 2f{f}^{\prime }\\ f^{″}& 2\left({f^{\prime }}^2+2f{f}^{″}\right)\end{array}\right]\left[\begin{array}{c}A\\ B\end{array}\right]=\left[\begin{array}{c}-1\\ ϵ\end{array}\right],$ (80)

which we solve to have

$A=-\frac{1}{{f^{\prime }}^3}\left(ϵf{f}^{\prime }+{f^{\prime }}^2+f{f}^{″}\right),B=\frac{1}{2{f^{\prime }}^3}\left(ϵf^{\prime }+f^{″}\right).$ (81)

Since root a of $f\left(x\right)$ is unknown we replace a by $x_0$ to have the method

$x_1=x_0-\frac{f_0}{{f^{\prime }}_0}-\frac{1}{2}\frac{f_0^2}{{f^{\prime }}_0^2}\left(ϵ+\frac{{f^{″}}_0}{{f^{\prime }}_0}\right)\mathrm{,}$ (82)

where $f_0=f\left(x_0\right)$ etc. Here

$x_1-a=-\frac{1}{2}ϵ{\left(x_0-a\right)}^2+O\left({\left(x_0-a\right)}^3\right),$ (83)

and convergence is from below if $ϵ>0$ , while convergence is from above if $ϵ<0$.

For

$f\left(x\right)=x^2-N$ (84)

the method becomes

$x_1=\frac{1}{8x_0^3}\left(-N^2+6N{x}_0^2+3x_0^4-ϵ{\left(x_0^2-N\right)}^2{x}_0\right).$ (85)

For example, for $N=2,ϵ=1/25$ , and $x_0=1.5$ we have $x_1=1.414213$ and $x_1^2=1.9999983$. For $ϵ=0$ we have $x_1=1.414352$ and $x_1^2=2.00039$.

The choice

$ϵ=\frac{1}{2N\sqrt{N}}\left(x_0^2-N\right)$ (86)

makes method (82) the quartic

$x_1-\sqrt{N}=-\frac{7}{8N\sqrt{N}}{\left(x_0-\sqrt{N}\right)}^4+O\left({\left(x_0-\sqrt{N}\right)}^5\right).$ (87)

15. The Simplest of All Methods

A simple routine for constructing a rational approximation to an irrational number consists of starting with any good rational approximation p/q to, say, $\sqrt{2}$ , then adding one to p if ${\left(p/q\right)}^2<2$ , or adding one to q if ${\left(p/q\right)}^2>2$. Starting with 3/2 we obtain this way the alternating sequence

$\frac{3}{2},\frac{3}{3},\frac{4}{3},\frac{5}{3},\frac{5}{4},\frac{6}{4},\frac{6}{5},\frac{7}{5},\frac{8}{5},\frac{8}{6},\frac{9}{6},\frac{9}{7},\frac{10}{7},\frac{10}{8},\frac{11}{8},\frac{12}{8},\frac{12}{9},\frac{13}{9},\cdots ,$ (88)

where ${\left(12/9\right)}^2=1.78,{\left(13/9\right)}^2=2.086$.

The method is sluggish, yet we can glean from this long sequence some very good Pell approximations to $\sqrt{2}$ , such as $1/1,k=-1$ ; $3/2,k=1$ ; $7/2,k=-1$ ; $17/12,k=1$ ; $41/29,k=-1$ ; $99/70,k=1$ ; $239/169,k=-1$ ; $577/408,k=1$. Number $k=p^2-N{q}^2$.

Going up to 4-digit approximations we find ${\left(3363/2378\right)}^2=2.000000180$ , ${3363}^2-2\times {2378}^2=1$ , and then ${\left(8119/5741\right)}^2=1.999999970$ , ${8119}^2-2\times {5741}^2=-1$. Among the 5-digit approximations we find ${\left(19601/13860\right)}^2=2.000000005$ , ${19601}^2-2\times {13860}^2=1$ and ${\left(47321/33461\right)}^2=1.999999999$ , ${47321}^2-2\times {33461}^2=-1$.

Thus, the alternating sequence of rational approximations to $\sqrt{2}$

$\frac{1}{1}\mathrm{,}\frac{3}{2},\frac{7}{5},\frac{17}{12},\frac{41}{29}\mathrm{,}\frac{99}{70}\mathrm{,}\frac{239}{169}\mathrm{,}\frac{577}{408}\mathrm{,}\frac{1393}{985}\mathrm{,}\frac{3363}{2378}\mathrm{,}\frac{8119}{5741}\mathrm{,}\frac{19601}{13860}\mathrm{,}\frac{47321}{33461}\mathrm{,}\frac{114243}{80782}$ (89)

is of excellent p/q rational approximations to $\sqrt{2}$ such that $p^2-2q^2=-1$ if $p/q<\sqrt{2}$ , and $p^2-2q^2=1$ if $p/q>\sqrt{2}$.

For N = 7 we find this way $8/3,k=1$ ; $127/48,k=1$ ; $2024/765,k=1$ for the upper bounds, and $2/1,k=-3$ ; $5/2,k=-3$ ; $37/14,k=-3$ ; $82/31,k=-3$ for the lower bounds.

To understand the convergence mechanism of this algorithm, let p/q be the last fraction less then $\sqrt{2}$ , namely, such that $p/q<\sqrt{2}$ , but $\left(p+1\right)/q>\sqrt{2}$. Then

$\frac{p}{q}<\sqrt{2}<\frac{p}{q}+\frac{1}{q},$ (90)

and the bounds on $\sqrt{2}$ become tighter as q increases by the repeated addition of 1 to it.

16. Bisection by Mediants

Mediant m of the two nonzero rationals $a/b<c/d$ is

$m=\frac{a+c}{b+d}.$ (91)

Lemma 2. We have

$\frac{a}{b}<m<\frac{c}{d}.$ (92)

Proof. Since $a/b<c/d$ , $bc-ad>0$ , and the result follows.

Lemma 3. We have

If $bc-ad=k$ , then $m-\frac{a}{b}=\frac{k}{b\left(b+d\right)}$ , and $\frac{c}{d}-m=\frac{k}{d\left(b+d\right)}$. (93)

Proof. The result follows by some simple algebra.

For example, from Equation (89) we have that $7/5<\sqrt{2}<3/2$ with $3/2-7/5=1/10$. Here the mediant $m=10/7$ , and $7/5\mathrm{<mo>\; <\; </mo>}\sqrt{2}\mathrm{<mo>\; <\; </mo>}10/7$ with $10/7-7/5=1/35$. The next $m=17/12$ , and $7/5\mathrm{<mo>\; <\; </mo>}\sqrt{2}\mathrm{<mo>\; <\; </mo>}17/12$ with $17/12-7/5=1/60$ ; all spreads between the upper and lower bounds having a numerator equal to one.

Unlike ordinary bisections, bisection by mediants converges to a rational number in a finite number of steps. For example, by mediants

$\frac{1}{1}<2<\frac{4}{1},\frac{1}{1}<2<\frac{5}{2},\frac{6}{3}\le 2\le \frac{6}{3},$ (94)

while by ordinary bisection

$\frac{1}{1}<2<\frac{4}{1},\frac{2}{2}<2<\frac{5}{2},\frac{7}{4}<2<\frac{10}{4},\frac{14}{8}<2<\frac{17}{8},\frac{28}{16}<2<\frac{31}{16},\cdots .$ (95)

17. Root Bracketing

We start with the following result.

Lemma 4. Let the integer pair $\left(p_0\mathrm{,}{q}_0\right)$ satisfy Pell’s equation $p_0^2-N{q}_0^2=1$ , and let $p_1=A{p}_0+CN{q}_0$ , $q_1=C{p}_0+A{q}_0$. Then

$\frac{p_0}{q_0}-\frac{p_1}{q_1}=\frac{C}{q_0{q}_1}.$ (96)

Proof. The result follows by common denominator.

Numerical example. For $N=7,A=2,C=1$ we have that $A^2-7C^2=-3$. Hence, in accordance with Lemma 1

$p_1=2p_0+7q_0,q_1=p_0+2q_0,\text{are}\text{such}\text{that}{p}_1^2-7q_1^2=-3\left(p_0^2-7q_0^2\right)=-3.$ (97)

Choosing the Pell (k = 1) pair $\left(p_0,q_0\right)=\left(127/48\right)$ , ${127}^2-7\times {48}^2=1$ , we obtain the Pell ( $k=-3$ ) pair $\left(p_1,q_1\right)=\left(590,223\right)$ , ${590}^2-7\times {223}^2=-3$ , and

$\frac{590}{223}<\sqrt{7}<\frac{127}{48}\text{of}\text{spread}\frac{127}{48}-\frac{590}{223}=\frac{1}{10704}$ (98)

of a numerator equal to one.

Similarly, choosing the Pell (k = 1) pair $\left(p_0,q_0\right)=\left(2024,765\right)$ we obtain the Pell ( $k=-3$ ) pair $\left(p_1,q_1\right)=\left(9403,3554\right)$ , and

$\frac{9403}{3554}<\sqrt{7}<\frac{2024}{765}\text{of}\text{spread}\frac{2024}{765}-\frac{9403}{3554}=\frac{1}{2718810}.$ (99)

The mediant in Equation (99) is $m=\left(9403+2024\right)/\left(3554+765\right)=11427/4319$ , and with it

$\frac{9403}{3554}<\sqrt{7}<\frac{11427}{4319}\text{of}\text{spread}\frac{11427}{4319}-\frac{9403}{3554}=\frac{1}{15349726}.$ (100)

18. Conclusion

In this paper we have examined single-step iterative methods for the solution of the nonlinear algebraic equation $f\left(x\right)=x^2-N=0$ , for some integer N, which produce rational approximations p/q that are optimal in the sense of Pell’s equation $p^2-N{q}^2=k$ for some integer k. We have also considered the most elementary bisection method for iteratively creating upper and lower bounds on the targeted root.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Silverman, J.H. (2006) A Friendly Introduction to Number Theory, Chapter 12. Prentice Hall, New Jersey.
[2]	Whitford, E.E. (1912) The Pell Equation. E. E. Whitford, New York.
[3]	Barbeau, E.J. (2000) Pell’s Equation. Springer, New York.
[4]	Jacobson Jr., M.J. and Williams, H.C. (2009) Solving the Pell Equation. Springer, New York.
[5]	Fried, I. (2009) Oppositely Converging Newton-Raphson Method for Nonlinear Equilibrium Problems. International Journal for Numerical Methods in Engineering, 79, 375-378. https://doi.org/10.1002/nme.2574
[6]	Fried, I. (2013) High-Order Iterative Bracketing Methods. International Journal for Numerical Methods in Engineering, 94, 708-714. https://doi.org/10.1002/nme.4467
[7]	Fried, I. (2014) Effective High-Order Iterative Methods via the Asymptotic Form of the Taylor-Lagrange Remainder. Journal of Applied Mathematics, 2014, Article ID: 108976. https://doi.org/10.1155/2014/108976
[8]	Brown, L.M. (1997) An Algorithm for Square Roots: An Episode in the Campaign Against Dotage. The Mathematical Gazette, 81, 428-429. https://doi.org/10.2307/3619622
[9]	McBride, A. (1999) Remarks on Pell’s Equation and Square Root Algorithms. The Mathematical Gazette, 83, 47-52. https://doi.org/10.2307/3618682
[10]	Dunkel, O. (1927) A Note on the Computation of Arithmetic Roots. The American Mathematical Monthly, 34, 366-368. https://doi.org/10.2307/2300041
[11]	Uspensky, J.V. (1927) Note on the Computation of Roots. The American Mathematical Monthly, 34, 130-134. https://doi.org/10.1080/00029890.1927.11986665
[12]	Neta, B. (1979) A Sixth-Order Family of Methods for Nonlinear Equations. International Journal of Computer Mathematics, Section B, 7, 157-161. https://doi.org/10.1080/00207167908803166
[13]	Fried, I. (2016) A Remarkable Chord Iterative Method for Roots of Uncertain Multiplicity. Applied Mathematics, 7, 1207-1214. https://doi.org/10.4236/am.2016.711106