Chebyshev’s Method for Multiple Zeros of Analytic Functions: Convergence, Dynamics and Real-World Applications

Abstract

This paper deals with the convergence and dynamics of Chebyshev’s method for simple and multiple zeros of analytic functions. We establish a local convergence theorem that provides error estimates and exact domains of initial approximations to guarantee the Q-cubic convergence of the method right from the first iteration. Applications to some real-world problems as well as theoretical and numerical comparison with the famous Halley’s method are also provided.

1. Introduction

Suppose $f:\mathbb{C}\to \mathbb{C}$ is an arbitrary function. Finding the zeros of f is one of the main tasks that arise from real-world problems. Iteration methods are among the most efficient tools for approximating the zeros of f. Undoubtedly, the most famous iteration methods in the literature are Newton’s method, Halley’s method [1] and Chebyshev’s method [2]. It is well known that Newton’s method converges quadratically, while Halley and Chebyshev’s methods converge cubically to simple zeros, but all of them converge only linearly for multiple zeros. In 1870, Schröder [3] presented the following modification of Newton’s method:

$${\displaystyle z_{k+1}=z_k-m\frac{f\left(z_k\right)}{f^{\prime }\left(z_k\right)}}$$

(1)

which restores the quadratic convergence for zero with known multiplicity $m\ge 1$. In 1963, Obreshkov [4] developed the following modifications of Halley and Chebyshev’s methods:

$$z_{k+1}=z_k-{\left({\displaystyle \frac{m+1}{2m}}{\displaystyle \frac{f^{\prime }\left(z_k\right)}{f\left(z_k\right)}}-{\displaystyle \frac{1}{2}}{\displaystyle \frac{f^{″}\left(z_k\right)}{f^{\prime }\left(z_k\right)}}\right)}^{-1}(\mathrm{Halley}’\mathrm{s}\mathrm{method})$$

(2)

$$z_{k+1}=z_k-{\displaystyle \frac{m^2}{2}}{\displaystyle \frac{f\left(z_k\right)}{f^{\prime }\left(z_k\right)}}\left({\displaystyle \frac{3-m}{m}}+{\displaystyle \frac{f\left(z_k\right)}{f^{\prime }\left(z_k\right)}}{\displaystyle \frac{f^{″}\left(z_k\right)}{f^{\prime }\left(z_k\right)}}\right)(\mathrm{Chebyshev}’\mathrm{s}\mathrm{method}),$$

(3)

which converge cubically to zeros with multiplicity m. It is worth noting that, Traub [5] has suggested using any method for $f^{(m-1)}$ or $f^{1/m}$ instead of f to preserve its convergence order. For a more detailed historical survey about these famous methods, we refer the reader to the elegant papers of Ypma [6], Scavo and Thoo [7] and Ezquerro et al. [8].

In 2009, Proinov [9] established a local convergence theorem about Newton’s method (1) for multiple zeros of analytic functions. Very recently, based on Proinov’s general convergence theory developed in [9,10,11], Ivanov [12] has proven a general convergence theorem about Picard iteration

$$z_{k+1}=T{z}_k(k=0,1,2,\dots ),$$

(4)

where $T:D\subset \mathbb{C}\to \mathbb{C}$ is some iteration function and, as an application, has proven a local convergence result about Halley’s method (2) for simple and multiple zeros of analytic functions.

In this paper, using the approach of [12], we study the local convergence of Chebyshev’s method (3) for simple and multiple zeros of analytic functions. Thus, we obtain a local convergence theorem (Theorem 3) that provides initial conditions and error estimates to guarantee the Q-cubic convergence of Chebyshev’s method (3). More precisely, our result provides the following:

• Exact domains of initial approximations, that is, exact Q-convergence balls (Theorem 4), that ensure the Q-cubic convergence of Chebyshev’s method right from the first iteration;
• A priory and a posteriori error estimates;
• An assessment of the asymptotic error constant of the method.

To the best of our knowledge, there is no such convergence theorem about Chebyshev’s method (3), except Theorem 2.1 of [13], which, in fact, only detects the convergence order of the Chebyshev–Halley family (see [14]), which, in turn, includes Chebyshev’s method (3). So, Theorem 3 improves and complements the mentioned Osada’s theorem as regards Chebyshev’s method. Moreover, our results allow us to compare Chebyshev’s method with the famous Halley one, theoretically (see Remark 3). Note that, the present study is a kind of continuation of [15], where a detailed local convergence analysis of Chebyshev’s method (3) for simple and multiple polynomial zeros has been presented.

On the other hand, the dynamics of Chebyshev’s method has been widely studied for simple polynomial zeros (see [16,17,18,19] and references therein). Here, we present the basins of attraction of Chebyshev’s method (3) applied to polynomial and transcendental equations arising from some real-world problems. Our examples (see Examples 4 and 5) show that neither Chebyshev’s method nor Halley’s one satisfies the ‘nearest root principal’ when applied to irrational or transcendental equations instead to polynomials with only simple zeros. A comparison with Halley’s method (2) based on the convergence behavior and dynamics of the two methods is also provided.

2. Preliminaries

Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$ be analytic in a neighborhood of $\xi \in \mathcal{D}$, which is a zero of f with multiplicity $m\ge 1$, then we define the function $E:\mathbb{C}\to {\mathbb{R}}_{+}$ with

$$E\left(z\right)=\gamma |z-\xi |,$$

(5)

where $\gamma$ is defined by (see Yakoubsohn [20])

$$\gamma =\gamma (f,\xi )=\underset{k>m}{sup}{\left|{\displaystyle \frac{m!f^{\left(k\right)}\left(\xi \right)}{k!f^{\left(m\right)}\left(\xi \right)}}\right|}^{1/(k-m)}.$$

(6)

The following definition introduces a wide class of iteration functions associated with the function E defined by (5).

Definition 1

([12] Definition 3.1). The function $T:D\subset \mathbb{C}\to \mathbb{C}$ is called gamma iteration function at the point $\xi \in \mathbb{C}$ if there exists a quasi-homogeneous function $\varphi :J\to {\mathbb{R}}_{+}$ of degree $p\ge 0$, such that for each $z\in \mathbb{C}$ with $E\left(z\right)\in J$, it follows that

$$z\in D\mathit{and}\left|T\right(z)-\xi |\le \varphi \left(E\right(z\left)\right)|z-\xi |,$$

(7)

where the function E is defined by (5). The function ϕ shall be called a control function of T.

The following general convergence theorem shall be used for the proof of our main result:

Theorem 1

([12] Theorem 3.4). Let $T:D\subset \mathbb{C}\to \mathbb{C}$ be a gamma iteration function at the point $\xi \in D$ with control function ϕ of the exact degree $p\ge 0$, and $z_0\in \mathbb{C}$ be an initial guess satisfying

$$E\left(z_0\right)\in \{t\in J:\varphi \left(t\right)<1\}.$$

(8)

Then, ξ is a unique fixed point of T, Picard iteration (4) is well defined and converges to ξ with Q-order $r=p+1$, and for all $k\ge 0$, we have the following error estimates:

$$|z_{k+1}-\xi |\le {\lambda }^{r^k}|z_k-\xi |\mathit{and}|z_k-\xi |\le {\lambda }^{S_k\left(r\right)}|z_0-\xi |,$$

(9)

where $\lambda =\varphi \left(E\left(z_0\right)\right)$. In addition, if $z_k\ne \xi$ for sufficiently large k, then we have the following estimate of the asymptotic error constant:

$$\underset{k\to \infty }{lim\; sup}{\displaystyle \frac{|z_{k+1}-\xi |}{|z_k{-\xi |}^r}}\le {\gamma }^p\underset{t\to 0^{+}}{lim}{\displaystyle \frac{\varphi \left(t\right)}{t^p}}.$$

Definition 2

([12] Definition 3.5). Let $\xi \in \mathbb{C}$ and $\rho >0$. The set $B(\xi ,\rho )=\{z\in \mathbb{C}:|z-\xi |<\rho \}$ is said to be a Q-convergence ball (with center ξ and radius ρ) of some iterative method if for any initial point $z_0\in B(\xi ,\rho )$, the iterative sequence ${\left(z_k\right)}_{k=0}^{\infty }$ generated by the method is well defined and converges to ξ with a Q-order of at least $r\ge 1$.

The next theorem provides a Q-convergence ball of Picard iteration (4).

Theorem 2

([12] Theorem 3.6). Let $T:D\subset \mathbb{C}\to \mathbb{C}$ be a gamma iteration function at the point $\xi \in D$ with a control function ϕ of the exact degree $p\ge 0$, and R be the minimal positive solution of the equation $\varphi \left(t\right)=1$. Then, $B(\xi ,R/\gamma )$ is a Q-convergence ball of Picard iteration (4) with Q-order $r=p+1$ and with error estimates (9). Also, we have the following estimate of the asymptotic error constant:

$$\underset{k\to \infty }{lim\; sup}{\displaystyle \frac{|z_{k+1}-\xi |}{|z_k{-\xi |}^r}}<{\left({\displaystyle \frac{\gamma }{R}}\right)}^{r-1}.$$

We end this section with two technical lemmas that play a crucial role in the proofs of our results.

Lemma 1

([12] Lemma 4.1). Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, $\xi \in \mathcal{D}$ be a zero of f with multiplicity $m\ge 1$, and $E\left(z\right)<1$, where the function $E:\mathbb{C}\to {\mathbb{R}}_{+}$ is defined by (5). Then, the following estimates hold true:

$$\begin{array}{ccc}\hfill \left|a\right|& \le & {\displaystyle \frac{E\left(z\right)}{1-E\left(z\right)}},\left|b\right|\le {\displaystyle \frac{(m+1)E\left(z\right)-mE{\left(z\right)}^2}{{(1-E\left(z\right))}^2}},\left|l\right|\le {\displaystyle \frac{2E{\left(z\right)}^2}{{(1-E\left(z\right))}^3}}\mathit{and}\hfill \\ \hfill \left|c\right|& \le & {\displaystyle \frac{m(m-1){(1-E\left(z\right))}^2-2(m-1)E\left(z\right)+2m}{{(1-E\left(z\right))}^3}}E\left(z\right),\hfill \end{array}$$

(10)

where the numbers $a,b,c$ and l are defined by

$$a=\sum _{k=m+1}^{\infty }{A}_k,b=\sum _{k=m+1}^{\infty }k{A}_k,c=\sum _{k=m+1}^{\infty }k(k-1)A_k\mathit{and}l=\sum _{k=m+2}^{\infty }(k-m)(k-m-1)A_k,$$

(11)

with $A_k$ defined by

$$A_k=\left({\displaystyle \frac{m!f^{\left(k\right)}\left(\xi \right)}{k!f^{\left(m\right)}\left(\xi \right)}}\right){(z-\xi )}^{k-m}.$$

Remark 1.

As for the numbers defined by (11), it is not hard to verify that $l=c+m(m+1)a-2mb$.

Lemma 2

([12] Lemma 4.2). Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$ and $\xi \in \mathcal{D}$ be a zero of f with multiplicity $m\ge 1$. If f is analytic around ξ, then we can define $a,b$ and c with (11) and the following claims hold true:

(i)

If $z\in \mathcal{D}$ is not a zero of f and $a\ne -1$, then we have

$${\displaystyle \frac{f^{\prime }\left(z\right)}{f\left(z\right)}}={\displaystyle \frac{m+b}{(1+a)(z-\xi )}}.$$

(ii)

If $z\in \mathcal{D}$ is not a zero of $f^{\prime }$ and $b\ne -m$, then we have

$${\displaystyle \frac{f^{″}\left(z\right)}{f^{\prime }\left(z\right)}}={\displaystyle \frac{m^2-m+c}{(m+b)(z-\xi )}}.$$

3. Local Convergence Analysis

In this section, we present a comprehensive local convergence analysis of Chebyshev’s method for simple and multiple zeros of analytic functions,

$$z_{k+1}=C_f\left(z_k\right),$$

(12)

where Chebyshev’s iteration function $C_f:D\subset \mathbb{C}\to \mathbb{C}$ is defined by

$$C_f\left(z\right)=\left\{\begin{array}{cc}{\displaystyle z-\frac{m^2}{2}\frac{f\left(z\right)}{f^{\prime }\left(z\right)}\left(\frac{3-m}{m}+\frac{f\left(z\right)}{f^{\prime }\left(z\right)}\frac{f^{″}\left(z\right)}{f^{\prime }\left(z\right)}\right)}\hfill & \mathrm{if}{f}^{\prime }\left(z\right)\ne 0,\hfill \\ z\hfill & \mathrm{if}{f}^{\prime }\left(z\right)=0.\hfill \end{array}\right.$$

(13)

In what follows, for $m\ge 1$, we define the number $\tau$ with

$$\tau ={\displaystyle \frac{2m}{3m+1+\sqrt{m^2+6m+1}}}$$

(14)

and the real functions $G,H,F$ and $\varphi$ with

$$G\left(t\right)=2{(m+1-mt)}^3{t}^3+m(1-t)t^{2}F\left(t\right),H\left(t\right)=2{(2m{t}^2-(3m+1)t+m)}^3,$$

(15)

$$\begin{array}{ccc}\hfill F\left(t\right)& =& (6(1-t)+|m-3\left|\right)(m+1-mt)+m(2-t)(m(m-1){(1-t)}^2-2(m-1)t+2m)\hfill \\ & +& 2m\left(\right|m-3|(m+1-mt)+1){(1-t)}^2+m^2(m-1){(1-t)}^3.\hfill \end{array}$$

(16)

and

$$\varphi \left(t\right)={\displaystyle \frac{G\left(t\right)}{H\left(t\right)}}.$$

(17)

Note that $\tau$ is the unique solution of the equation $2m{t}^2-(3m+1)t+m=0$ in the interval $[0,1/2)$, and the function $\varphi$ is quasi-homogeneous of exact second degree ([21] [Definition 8]) in $[0,\tau )$.

The following is our main lemma in this paper.

Lemma 3.

Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, $\xi \in \mathcal{D}$ be a zero of f with multiplicity $m\ge 1$, and f be analytic around ξ. Then Chebyshev’s iteration function (13) is a gamma iteration function at the point ξ with a control function ϕ defined by (17).

Proof. 

Define the function E with (5) and the number $\tau$ with (14). Let $z\in \mathcal{D}$ be such that $E\left(z\right)<\tau$; then, we have

$$2mE{\left(z\right)}^2-(3m+1)E\left(z\right)+m>0\Rightarrow H\left(E\left(z\right)\right)>0,$$

(18)

where the function H is defined by (15). Now, we need to prove that

$$|C_f\left(z\right)-\xi |\le \varphi \left(E\left(z\right)\right)|z-\xi |,$$

(19)

where the function $\varphi$ is defined by (17). Let $f\left(z\right)\ne 0$. Define the numbers a, b, c and l with (11). Then, from the triangle inequality, the second inequality of (10) and (18), we obtain

$$|m+b|\ge m-\left|b\right|\ge {\displaystyle \frac{2mE{\left(z\right)}^2-(3m+1)E\left(z\right)+m}{{(1-E\left(z\right))}^2}}>0$$

(20)

which means that $b\ne -m$, and so, according to Lemma 2 (i), we deduce that $f^{\prime }\left(z\right)\ne 0$.

Now, it remains to prove (19). From the definition of Chebyshev’s iteration function (13) and Lemma 2, we obtain

$$C_f\left(z\right)-\xi =z-\xi -{\displaystyle \frac{m^2(1+a)(z-\xi )}{2(m+b)}}\left({\displaystyle \frac{3-m}{m}}+{\displaystyle \frac{(1+a)(m^2-m+c)}{{(m+b)}^2}}\right)=\sigma (z-\xi ),$$

(21)

where (see Remark 1)

$$\begin{array}{ccc}\hfill \sigma & =& {\displaystyle \frac{2{(m+b)}^3-m(3-m)(1+a){(m+b)}^2-m^2{(1+a)}^2(m^2-m+c)}{2{(m+b)}^3}}\hfill \\ & =& {\displaystyle \frac{2b^3+m(6+(m-3)(1+a))b^2+m^2[2(m-3)ab-(2+a)ac-m(m-1)a^2-l]}{2{(m+b)}^3}}.\hfill \end{array}$$

(22)

From this, the triangle inequality and the estimates (10) and (20), we obtain

$$\begin{array}{ccc}\hfill \left|\sigma \right|& \le & {\displaystyle \frac{{2\left|b\right|}^3{+m(6+|m-3\left|\right|1+a\left|\right)\left|b\right|}^2+m^2{\left[2\right|m-3\left|\right|a\left|\right|b|+|2+a\left|\right|a\left|\right|c|+m(m-1)|a|}^2+\left|l\right|]}{2{(m-|b\left|\right)}^3}}\hfill \\ & \le & {\displaystyle \frac{2{(m+1-mE\left(z\right))}^{3}E{\left(z\right)}^3+m(1-E\left(z\right))E{\left(z\right)}^{2}F\left(E\left(z\right)\right)}{H\left(E\right(z\left)\right)}}={\displaystyle \frac{G\left(E\right(z\left)\right)}{H\left(E\right(z\left)\right)}}=\varphi \left(E\left(z\right)\right).\hfill \end{array}$$

(23)

which, together with (21), proves (19) and completes the proof. □

The following is the main theorem of this paper:

Theorem 3.

Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, $\xi \in \mathcal{D}$ be a zero of f with multiplicity $m\ge 1$, and f be analytic in a neighborhood of ξ. Suppose $z_0\in \mathcal{D}$ is an initial approximation satisfying

$$\Psi \left(E\left(z_0\right)\right)>0,$$

(24)

where the function E is defined by (5) and $\Psi =H-G$, with H and G defined by (15). Then, Chebyshev’s iteration (12) is well defined and converges Q-cubically to ξ, with the following error estimates for all $k\ge 0$:

$$|z_{k+1}-\xi |\le {\lambda }^{3^k}|z_k-\xi |\mathit{and}|z_k-\xi |\le {\lambda }^{(3^k-1)/2}|z_0-\xi |,$$

(25)

where $\lambda =\varphi \left(E\left(z_0\right)\right)$ and the function ϕ is defined by (17). Additionally, if $z_k\ne \xi$ for sufficiently large k, then the following estimate of the asymptotic error constant holds:

$$\underset{k\to \infty }{lim\; sup}{\displaystyle \frac{|z_{k+1}-\xi |}{{|z_k-\xi |}^3}}\le {\displaystyle \frac{m(3m^2+1)+(m+1)(m+6+|m-3|(2m+1))}{2m^2}}{\gamma }^2.$$

(26)

Proof. 

Since $G\left(t\right)>0$ in the interval $(0,1/2)$, from $\Psi \left(t\right)>0$, we obtain $H\left(t\right)>G\left(t\right)>0$, which implies that $t<\tau$ and $\varphi \left(t\right)<1$, where $\tau$ is defined by (14). This means that the initial condition (24) implies that (8) with $J=[0,\tau )$, and so, the proof follows immediately from Lemma 3 and Theorem 1. □

Remark 2.

The initial condition (24) can be represented in the form $E\left(z_0\right)<R$, where R is the minimal positive solution of the equation $\Psi \left(t\right)=0$ in the interval $[0,\tau )$, where τ is defined by (14). Note that such R exists because Ψ is continuous, $\Psi \left(0\right)=2m^3>0$ and $\Psi \left(\tau \right)<0$.

Mindful of Theorem 2, we formulate Theorem 3 in the following equivalent form:

Theorem 4.

Let $\xi \in \mathcal{D}$ be a zero of $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$ with multiplicity $m\ge 1$, and f be analytic in a neighborhood of ξ. If R is the minimal positive solution of the equation $\Psi \left(t\right)=0$, where the function Ψ is the same as in Theorem 3, then $B(\xi ,R/\gamma )$ is a Q-convergence ball that guarantees Q-cubical convergence of Chebyshev’s method with estimates (25) and (26).

In the case of $m=1$, i.e., when $\xi$ is a simple zero of f, from Theorem 3, we obtain the following result:

Theorem 5.

Let $\xi \in \mathcal{D}$ be a simple zero of $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, which is supposed to be analytic in a neighborhood of ξ. Suppose $z_0\in \mathcal{D}$ is an initial approximation satisfying

$$E\left(z_0\right)<\mathcal{R}=0.119\dots ,$$

(27)

where the function E is defined by (5) with $\gamma ={sup}_{k>1}{|f^{\left(k\right)}\left(\xi \right)/(k!f^{\prime }\left(\xi \right))|}^{1/(k-1)}$ and $\mathcal{R}$ is the minimal positive solution of the equation

$$7t^6-40t^5+85t^4-82t^3+39t^2-12t+1=0.$$

(28)

Then, Chebyshev’s iteration (12) (for $m=1$) is well defined and converges Q-cubically to ξ with error estimate (25), where the function ϕ is defined by

$$\varphi \left(t\right)={\displaystyle \frac{{(2-t)}^3{t}^3-(1-t)(2t^3-12t^2+23t-15)t^2}{{(2t^2-4t+1)}^3}}.$$

(29)

Additionally, if $z_k\ne \xi$ for sufficiently large k, then the following estimate of the asymptotic error constant holds:

$$\underset{k\to \infty }{lim\; sup}{\displaystyle \frac{|z_{k+1}-\xi |}{|z_k{-\xi |}^3}}\le 15{\gamma }^2.$$

(30)

Proof. 

This proof follows directly from Theorem 3, since (27) and (24) are equivalent for $m=1$. □

Analogically to Theorem 4, we formulate Theorem 5 as a ball convergence one.

Theorem 6.

Let $\xi \in \mathcal{D}$ be a simple zero of $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, which is supposed to be analytic in a neighborhood of ξ. Then, $B(\xi ,\mathcal{R}/\gamma )$ with $\mathcal{R}$ defined by (27) is a Q-convergence ball that guarantees Q-cubical convergence of Chebyshev’s method with estimate (25), where ϕ is defined by (29), and the estimate of the asymptotic error constant (30).

Remark 3.

What is interesting is that in the cases of $m=2,\dots ,13$, from Theorem 3, we obtain theorems with larger convergence domains (balls) and better error estimates than Theorems 5 and 6. In fact, Theorem 3 gives the largest convergence domain when ξ is a zero of f with multiplicity $m=3$. Bearing in mind Remark 2, the values of R for $m=1,\dots ,13$ are given in

Table 1. The values of the convergence radii R and $R_H$ for different multiplicities m.

m	1	2	3	4	5	6	7	8	9	10	11	12	13
R	$0.119$	$0.154$	$0.169$	$0.163$	$0.156$	$0.149$	$0.144$	$0.139$	$0.134$	$0.130$	$0.126$	$0.123$	$0.120$
$R_H$	$0.143$	$0.121$	$0.100$	$0.084$	$0.073$	$0.064$	$0.056$	$0.051$	$0.046$	$0.042$	$0.039$	$0.036$	$0.033$

with three decimal places.

In light of Smale’s concept of approximate zero of the second kind ([22,23]), we give the following definition associated with Chebyshev’s method (12).

Definition 3.

Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$ and $\xi \in \mathcal{D}$ be a zero of f. The number $z_0\in \mathcal{D}$ is said to be an approximate zero of the second kind for Chebyshev’s method if the iterative sequence (12) is well defined and

$$|z_k-\xi |\le {\left({\displaystyle \frac{1}{2}}\right)}^{(3^k-1)/2}|z_0-\xi |\mathit{for}\mathit{all}k\ge 0.$$

(31)

According to Definition 3, the next consequence of Theorem 4 provides a ball of approximate zeros of the second kind for the Chebyshev’s method in the case of $m=1$.

Corollary 1.

Let $\xi \in \mathcal{D}$ be a simple zero of a function $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, which is analytic in a neighborhood of ξ. Then, any $z_0\in B(\xi ,\mathcal{R}/\gamma )$ is an approximate zero of the second kind for Chebyshev’s method, where $\mathcal{R}=0.099\dots$ is the minimal positive solution of the equation

$$6t^6-32t^5+62t^4-52t^3+24t^2-12t+1=0.$$

Hereinafter, we aim to compare the convergence properties of Chebyshev’s method (12) and Halley’s one (2) using Theorem 5 and Theorem 4.4 of [12], which we formulate in the following way:

Theorem 7

([12] Theorem 4.4). Let $f:\mathcal{D}\subset \mathbb{C}\to \mathbb{C}$, $\xi \in \mathcal{D}$ be a zero of f with multiplicity $m\ge 1$ and f be analytic in a neighborhood of ξ. Suppose $z_0\in \mathcal{D}$ is an initial approximation, such that

$$E\left(z_0\right)<R_H,$$

(32)

where the function E is defined by (5) and $R_H$ is the minimal positive solution of the equation $h\left(t\right)=g\left(t\right)$, with h and g defined by

$$\begin{array}{ccc}\hfill h\left(t\right)& =& [(3m+1)(m+1-mt)+m^2(m-1){(1-t)}^2]t^2,\hfill \\ \hfill g\left(t\right)& =& (m+1){(2m{t}^2-(3m+1)t+m)}^2-m^2(2t+(m-1){(1-t)}^2).\hfill \end{array}$$

Then the Halley’s iteration (2) is well defined and converges Q-cubically to ξ with error estimate (25), where $\varphi =h/g$, and with the following estimate of the asymptotic error constant:

$$\underset{k\to \infty }{lim\; sup}{\displaystyle \frac{|z_{k+1}-\xi |}{{|z_k-\xi |}^3}}\le {\displaystyle \frac{{(m+1)}^3-m(m-1)}{2m^3}}{\gamma }^2.$$

(33)

Remark 4.

One can see from Table 1 that, in the case of simple zeros, Theorem 7 gives a larger convergence domain and therefore better error estimates for Halley’s method than Theorem 5 for Chebyshev’s one. For all other multiplicities ($m\ge 2$), Chebyshev’s method behaves better than Halley’s one. On the other hand, it is easy to check that the estimate of the asymptotic error constant (33) is better than (26) for all $m\ge 1$.

4. Convergence and Dynamics in Some Real-World Applications

In this section, we apply Chebyshev’s method (3) and Halley’s method (2) for the approximation of the roots of six nonlinear equations arising from different real-world problems. A numerical comparison based on the convergence behavior and dynamics of the two methods is conducted within the examples. In addition, we wish to note that we are in the process of developing integrated physics, mathematics and art courses within STEAM (Science, Technology, Engineering, the Arts and Mathematics) education at the Bulgarian high-school level. More precisely, we propose implementing the Polynomiography approach (see [24]) to high-school physics education in Bulgaria to increase students’ interest, demonstrating the beauty of science through the ‘chain’ Physics–Mathematics–Art. So, we intend to use the following examples within these courses.

Suppose $f:\mathbb{C}\to \mathbb{C}$, and $z_0$ is an initial approximation. In the examples below, we apply Chebyshev and Halley’s methods for the approximation of a zero $\xi$ of f and for obtaining their basins of attraction. For any example, the dynamical planes are plotted using a mesh of $500\times 500$ points and 50 iterations. The yellow points in the figures depict the roots of the corresponding equation. For the convergence tests, we define the following error estimate

$${\epsilon }_k=|z_k-z_{k+1}|k\ge 0.$$

In the tables below, for any example, we present the approximated root $\xi$ and the estimates ${\epsilon }_k$ and ${\epsilon }_{k+1}$ for the two methods. All results in the tables are given with at least six decimal places. All calculations and figures are performed in CAS Wolfram Mathematica 8.0.

Example 1

(Adiabatic flame temperature [25]). It is known that the adiabatic energy balance can be presented by the following cubic equation:

$$9.433\times {10}^{-8}{z}^3+0.001z^2+7.256z-2324.618=0,$$

(34)

where the variable z represents the adiabatic flame temperature in Kelvin. Note that Equation (34) possesses the following roots: ${\xi }_{1,2}=-6242.729\pm 6462.059i$ and ${\xi }_3=305.247$. The basins of attraction of the roots of (34) are obtained by the Chebyshev and Halley’s methods in the square

$$\{z\in \mathbb{C}:|Re\left(z\right)|\le 7000\mathit{and}|Im\left(z\right)|\le 7000\}$$

and the initial guess $z_0=799.085274+2546.595701i$ is randomly taken from the same square.

Despite the ‘rough’ choice of the initial guess, one can see from

Table 2. Convergence results for Example 1.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	$305.247507$	3	$6.085157\times {10}^{-6}$	$1.139171\times {10}^{-13}$
Halley’s	$305.247507$	3	$1.251784\times {10}^{-10}$	$5.695856\times {10}^{-14}$

that the root ${\xi }_3$ is found by the two methods only in the third iteration, with an error estimation of less than ${10}^{-5}$ for Chebyshev’s method and less than ${10}^{-9}$ for Halley’s one, while at the next iteration, the root is found with accuracy of ${10}^{-12}$ and ${10}^{-13}$, respectively. The basins of attraction for this example are plotted in Figure 1.

Example 2

(Chemical equilibrium [25]). The equilibrium conversion of nitrogen and hydrogen to ammonia starting with $3mol$ of hydrogen and 1 mol of nitrogen in reaction conditions of 250 atm and $500$ °C is described by the following polynomial equation:

$$z^4-7.790z^3+14.744z^2+2.511z-1.674=0,$$

(35)

where z represents the moles of nitrogen reacted at equilibrium, which, for $1mol$ of nitrogen, also represents the fractional conversion of nitrogen. Equation (35) has four simple roots: ${\xi }_1=-0.384094$, ${\xi }_2=0.277759$ and ${\xi }_{3,4}=3.948542\pm 0.316123i$. The basins of attraction of the roots of (35) are obtained by Chebyshev and Halley’s methods in the square

$$\{z\in \mathbb{C}:|Re\left(z\right)|\le 5\mathit{and}|Im\left(z\right)|\le 5\}.$$

The initial guess $z_0=-2.736516-0.029370i$, again, is randomly taken from the same square.

The numerical results for this example are given in

Table 3. Convergence results for Example 2.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	$-0.384094$	5	$1.293772\times {10}^{-8}$	$1.240770\times {10}^{-23}$
Halley’s	$-0.384094$	5	$1.504136\times {10}^{-15}$	$9.860761\times {10}^{-32}$

and the basins of attraction are plotted in Figure 2. It is seen in the table that Halley’s method behaves slightly better that Chebyshev’s one.

In the following example, we shall consider a well-known engineering problem that leads to a polynomial equation with multiple roots.

Example 3

(Deflection of horizontal beam [26] [Section 5.1]). Consider a 10-meter beam with embedded ends, i.e., the ends are not only supported at both ends, but are also clamped horizontally to something (for example, walls). The deflection of the beam is described by the following Euler–Bernoulli equation:

$$EI{\displaystyle \frac{d^{4}f\left(z\right)}{d{z}^4}}=w\left(z\right),$$

where E is Young’s modulus of elasticity for the material, I is the moment of inertia of a cross-section of the beam and $w\left(z\right)$ is the load per unit length, which is constant when the beam has a uniform cross-section and is supporting its own weight. Taking $E=10$, $I=5$ and $w\left(z\right)=100$, as in [26] [Example 5.1], we obtain the following boundary value problem:

$${\displaystyle \frac{d^{4}f\left(z\right)}{d{z}^4}}=2,f\left(0\right)=f^{\prime }\left(0\right)=f\left(10\right)=f^{\prime }\left(10\right)=0,$$

with the solution

$$f\left(z\right)={\displaystyle \frac{1}{12}}{z}^4-{\displaystyle \frac{5}{3}}{z}^3-{\displaystyle \frac{20}{3}}{z}^2.$$

(36)

Note that f has two simple zeros of ${\xi }_1=-3.416407$ and ${\xi }_2=23.416407$, and a zero of ${\xi }_3={\xi }_4=0$ with multiplicity $m=2$.

For this example, the basins of attraction of the zeros of (36) are obtained by Chebyshev and Halley’s methods in the rectangle

$$\{z\in \mathbb{C}:-5\le Re(z)\le 25\mathit{and}|Im\left(z\right)|\le 25\}.$$

The initial guess $z_0=3.188976+0.321727i$ is randomly taken from the same rectangle.

The numerical results for this example are given in

Table 4. Convergence results for Example 3.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	0	3	$4.040452\times {10}^{-7}$	$7.386719\times {10}^{-21}$
Halley’s	0	3	$1.353083\times {10}^{-9}$	$1.105267\times {10}^{-25}$

and the basins of attraction are plotted in Figure 3.

Example 4

(Blood rheology [27,28]). It is well known that the blood is modeled as Casson Fluid. The Casson fluid model says that basic fluids such as blood and water flow through a tube such that the velocity gradient near the wall of the tube and the center core of the fluids travel as plugs with minimal deformation. So, to calculate the plug flow of Casson fluids, we use the following irrational equation

$${\displaystyle \frac{1}{21}}{z}^4-{\displaystyle \frac{4}{3}}z+{\displaystyle \frac{16}{7}}\sqrt{z}+0.6=0$$

(37)

with zeros of ${\xi }_{1,2}=-2.400167\pm 2.205360i$ and ${\xi }_{3,4}=1.919730\pm 1.470124i$. The basins of attraction for this example are obtained in the square

$$\{z\in \mathbb{C}:|Re\left(z\right)|\le 3\mathit{and}|Im\left(z\right)|\le 3\}$$

and the initial approximation $z_0=-1.664278-0.570523i$ is randomly taken from the same square.

The numerical results for this example are given in

Table 5. Convergence results for Example 4.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	$1.919730-1.470124i$	4	$8.807333\times {10}^{-7}$	$2.220446\times {10}^{-16}$
Halley’s	$-2.400167+2.205360i$	4	$4.09973\times {10}^{-3}$	$1.374906\times {10}^{-9}$

and the basins of attraction are plotted in Figure 4. What is interesting is that from one and the same initial approximation, the two methods converge to different roots. Taking a closer look at the initial guess and the approximated roots, one can see that neither Chebyshev’s method nor Halley’s satisfy the ‘nearest root principal’ (see [17] and references therein) when applied to irrational equations. In fact, this is easily observable from Figure 4, where the white points depict the initial guess $z_0$ while the black ones are the approximated zeros.

In the rest of the examples, we shall consider two problems that lead to transcendental equations.

Example 5

(Fractional conversion in a reactor [25]).  The fractional chemical conversion in a reactor is proposed to be determined by the equation

$$z-exp(0.8z-1.691954)=0,$$

(38)

with $z=(1-1.25X)/(1-X)$, where X denotes the fractional conversion of the limiting reactant. Equation (38) has two real roots, ${\xi }_1=0.219513$ and ${\xi }_2=3.775662$, which can be easily used for finding X.

For this example, the basins of attraction of the zeros of (38) are obtained in the rectangle

$$\{z\in \mathbb{C}:-1\le Re(z)\le 20\mathit{and}|Im\left(z\right)|\le 20\}.$$

The initial approximation $z_0=11.223362-13.190431i$ is randomly taken from the same rectangle.

The numerical results and the dynamical planes for this example are given in

Table 6. Convergence results for Example 5.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	$0.219513$	9	$4.392337\times {10}^{-13}$	$5.048709\times {10}^{-29}$
Halley’s	$0.219513$	6	$5.801288\times {10}^{-10}$	$2.775557\times {10}^{-26}$

and Figure 5. One can see from the table that Halley’s method behaves considerably better than Chebyshev’s, while the figure shows that both methods do not satisfy the ‘nearest root principal’ when applied to transcendental equations.

Example 6

(Torsional vibrations of a ball [29]).  The natural vibrations of a spherical cavity located in an infinitely elastic medium are considered in this example. The transcendental equation

$$z/tanz=1-z^2/3,$$

(39)

describes the complex spectrum of the natural frequencies z of the torsional vibrations of a spherical cavity in a continuous elastic medium. It is easy to see that Equation (39) has an infinite number of real roots and singularities. Since the main difficulty of applying the two methods is the need for the existence of a second derivative, we enlarge the domain of (39), and so we shall consider the equation

$$(1-z^2/3)sinz-zcosz=0,$$

(40)

which has the same roots as (39) but no singularities. We shall seek the roots of (40) in the interval $[-6,6]$ which are ${\xi }_{1,2}=\pm 5.763459$ and ${\xi }_3=0$. The basins of attraction of the zeros of (40) are obtained in the square

$$\{z\in \mathbb{C}:|Re\left(z\right)|\le 10\mathit{and}|Im\left(z\right)|\le 10\}$$

while the initial approximation $z_0=5.609641-2.798980i$ is randomly taken from the square

$$\{z\in \mathbb{C}:|Re\left(z\right)|\le 6\mathit{and}|Im\left(z\right)|\le 6\}$$

because we wish to approximate the roots inside it but not outside.

The numerical results and the dynamical planes for this example are given in

Table 7. Convergence results for Example 6.

Method	Root	k	${\mathit{\epsilon }}_{\mathit{k}}$	${\mathit{\epsilon }}_{\mathit{k}+1}$
Chebyshev’s	$5.763459$	8	$2.675586\times {10}^{-7}$	$6.652038\times {10}^{-21}$
Halley’s	$5.763459$	5	$2.790587\times {10}^{-10}$	$4.930320\times {10}^{-30}$

and Figure 6. It is easy to see that for this example, Halley’s method again behaves better than Chebyshev’s method.

5. Conclusions

A detailed local convergence analysis of Chebyshev’s method (12) for simple and multiple zeros of analytic functions has been provided. As a result, a local convergence theorem (Theorem 3) that provides the exact domains of the initial conditions, together with a priori and a posteriori error estimates to guarantee the Q-cubic convergence of Chebyshev’s method (12) right from the first iteration, has been obtained. To the best of the authors’ knowledge, there are no existing results pertaining to this method, except the theorem of Osada [13], which only detects the convergence order of the method. As consequences, two ball Q-convergence theorems (Theorems 4 and 6) as well as a Smale-type result (Corollary 1) that provides a ball of approximate zeros of the second kind for Chebyshev’s method in the case of $m=1$, have been obtained.

Furthermore, theoretical and numerical comparisons between Chebyshev’s method and the famous Halley’s one have been conducted on the bases of Theorem 3 and the convergence behavior and dynamics of the two methods when applied to six real-world problems. The comparison shows that Chebyshev’s method (12) has some theoretical advantages (see Remark 4) over Halley’s method (2), but this is not the case in the numerical comparison, where, in almost all examples, Halley’s method behaves better than Chebyshev’s one. In addition, our examples (Examples 4 and 5) have shown that neither Chebyshev’s method nor Halley’s one satisfies the ‘nearest root principal’ when applied to irrational or transcendental equations. As for the basins of attraction, one can see that the two methods produce their own specific shapes in the plain, which hints that all iteration methods have their own ‘art style’. For example, in all figures produced by Chebyshev’s method, clover-like shapes are observed, while in those produced by Halley’s one, balloon-like shapes prevail. Along this line, we intend to further implement the provided examples, and others, in some new integrated physics, mathematics and art courses within STEAM education in Bulgarian high-schools in order to increase students’ interest by revealing the art in the science.

Finally, we can further extend the convergence analysis of Chebyshev’s method for operator equations in Banach spaces.