Next Article in Journal
Computational Modeling of Lymph Filtration and Absorption in the Lymph Node by Boundary Integral Equations
Next Article in Special Issue
Overlapping Grid-Based Optimized Single-Step Hybrid Block Method for Solving First-Order Initial Value Problems
Previous Article in Journal
Biomedical Image Classification via Dynamically Early Stopped Artificial Neural Network
Previous Article in Special Issue
Numerical Simulation of Micro-Bubbles Dispersion by Surface Waves
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Dynamics and Stability on a Family of Optimal Fourth-Order Iterative Methods

by
Alicia Cordero
1,*,
Miguel A. Leonardo Sepúlveda
2,3 and
Juan R. Torregrosa
1
1
Instituto de Matemática Multidisciplinar, Universitat Politècnica de València, Camino de Vera, s/n, 46022 Valencia, Spain
2
Ciencias Básicas y Ambientales (CBA), Instituto Tecnológico de Santo Domingo (INTEC), Area de Ciencia Básica y Ambiental, Av. Los Próceres, Gala, Apartado Postal 342-9 and 249-2, Santo Domingo 10602, Dominican Republic
3
Departamento de Matemática, Instituto Superior de Formación Docente Salomé Ureña, Recinto Félix Evaristo Mejía (ISFODOSU), Av. Caonabo con esq. Leonardo da Vinci, Urb. Renacimiento, Sector Mirador Sur, Santo Domingo 10114, Dominican Republic
*
Author to whom correspondence should be addressed.
Algorithms 2022, 15(10), 387; https://doi.org/10.3390/a15100387
Submission received: 30 September 2022 / Revised: 17 October 2022 / Accepted: 18 October 2022 / Published: 21 October 2022

Abstract

:
In this manuscript, we propose a parametric family of iterative methods of fourth-order convergence, and the stability of the class is studied through the use of tools of complex dynamics. We obtain the fixed and critical points of the rational operator associated with the family. A stability analysis of the fixed points allows us to find sets of values of the parameter for which the behavior of the corresponding method is stable or unstable; therefore, we can select the regions of the parameter in which the methods behave more efficiently when they are applied for solving nonlinear equations or the regions in which the schemes have chaotic behavior.

1. Introduction

One of the most important legacies in numerical analysis in the 20th century is the resolution of nonlinear equations by using high-order iterative methods, since we often encounter nonlinear mathematical models for which there are no analytical methods that allow us to find a solution of the nonlinear equation f ( x ) = 0 . This topic has attracted researchers from different areas, such as Traub in [1], who initiated the analysis of these methods.
Recently, Petković et al. (in [2]) and Amat et al. (in [3]) collected and updated the state of the art of multipoint methods. The advantage of this type of method is that it does not use higher-order derivatives, and it is of great practical importance because it overcomes the theoretical limitations of one-point methods. The latter are not good alternatives for increasing the order of convergence and the computational efficiency rate.
When the order of an iterative method increases, so does the number of functional evaluations per step. Ostrowski’s efficiency index (see [4]) gives a measure of the balance between these quantities, according to the formula p 1 d , where p is the order of convergence of the method and d is the number of functional evaluations per iteration. According to the Kung–Traub conjecture [5], the order of convergence of any multipoint method without memory cannot exceed the limit 2 n 1 , which is called the optimal order. Therefore, the optimal order of a method with three functional evaluations per step is 4. King’s family of methods [6], of which Ostrowski’s method [4] is a particular case, as is Jarratt’s method [7] and some of the families of multistep methods studied extensively in Traub’s book [8], are optimal fourth-order methods, since they perform only three functional evaluations per step. In recent years, as the extensive literature shows, there has been a renewed interest in the search for multistep methods in order to achieve optimal order convergence and, thus, better efficiency.
The rest of this paper is organized as follows. Section 2 is devoted to the presentation of a family of iterative methods that will be subjected to a dynamical analysis in Section 3. The stability of fixed points is studied, and the stable and unstable elements of the family are selected by means of the generation of the parameter plane. Some numerical tests are performed in Section 4, showing the performance of selected stable and unstable methods in comparison with other known ones. Finally, some conclusions are shown in Section 5.

2. The Proposed Family

In this section, we present particular cases of optimal fourth-order iterative methods for solving the nonlinear equation f ( x ) = 0 , f : D R R , which was proposed in [9] by using the weight function procedure, the iterative expression of which is:
y n = x n f x n f x n x n + 1 = x n f x n f x n H u x n ,
where u ( x n ) = f ( y n ) f ( x n ) and H ( u ) is a real-valued weight function that achieves an optimal fourth-order convergence by using three functional evaluations. The following result indicates the necessary conditions that the weight function in question (1) has to fulfil so that the order of convergence reached is four [9].
Theorem 1.
Let f be a sufficiently differentiable function f : D R R with a single-zero ξ on the open interval D. Let H ( u ) be a real-valued differentiable function. If an initial approximation x 0 is sufficiently close to the required root ξ of equation f ( x ) = 0 , then the order of convergence of the scheme (1) is four when H ( u ) satisfies the following conditions:
H ( 0 ) = 1 , H ( 0 ) = 1 , H ( 0 ) = 4 a n d | H ( 0 ) | < + .
In this case, the error equation is
e n + 1 = 5 1 6 H ( 0 ) c 2 3 c 2 c 3 e n 4 + O e n 5 ,
where
e n = x n ξ a n d c k = 1 k ! f ( k ) ( ξ ) f ( ξ ) , k = 2 , 3 , 4 ,

Special Cases

In this section, we present some particular classes of iterative schemes obtained from expression (1) for different weight functions H ( u ) that satisfy the conditions of Theorem 1.
  • Firstly, we consider the well-known Ostrowski method, which is a particular case given by the following weight function:
    H ( u ) = 1 u 1 2 u .
    This weight function guarantees the fourth-order convergence of Ostrowski’s scheme, whose iterative expression is:
    y n = x n f x n f x n , x n + 1 = x n f x n f x n f x n f y n f x n 2 f y n .
  • The following weight function is now considered:
    H ( u ) = 1 + u 2 u 3 1 u ,
    which satisfies all of the conditions of Theorem 1. This scheme was proposed by Maheshwari in [10], and its iterative expression is:
    y n = x n f x n f x n , x n + 1 = x n f x n f x n f ( x n ) f ( x n ) f ( y n ) + f ( y n ) f ( x n ) 2 .
  • The following uniparametric family satisfies the conditions of Theorem 1:
    H ( u ) = α + u ( 1 + α + u + 2 α u ) α + u ,
    where α is a real parameter, and it provides a new iterative class of iterative schemes:
    y n = x n f x n f x n , x n + 1 = x n f x n f x n α f ( x n ) 2 + f ( x n ) f ( y n ) ( 1 + α ) + ( 1 + 2 α ) f ( y n ) 2 f ( x n ) [ f ( y n ) + α f ( x n ) ] .
  • Finally, using the convergence conditions to define the Taylor polynomial associated with the generic function H ( u ) , we obtain the last weight function presented here, thus satisfying Theorem 1.
    H ( u ) = H ( 0 ) 6 u 3 + 2 u 2 + u + 1 .
    We know that | H ( 0 ) | < , so if we substitute H ( 0 ) with β , we have a parametric family of iterative schemes of fourth-order convergence:
    y n = x n f x n f x n , x n + 1 = x n f x n f x n β 6 f ( y n ) f ( x n ) 3 + 2 f ( y n ) f ( x n ) 2 + f ( y n ) f ( x n ) + 1 .
Now, we will study the stability of the family given by Equation (11).

3. The Scaling Theorem

If the conjugation classes and the scaling theorem are verified for the method under study, it is possible to generalize a dynamical study on the polynomial of the second degree on which this study is performed. In this case, by applying this to a generic polynomial, such as p ( z ) = ( z a ) ( z b ) , we can assume that the result obtained will be equivalent to the result that we will obtain with any quadratic polynomial, since there will always be an affine transformation that takes us from this generic polynomial to any quadratic polynomial (see, for example, [11]). We emphasize that this result is only fulfilled in methods that contain derivatives, which is the justification of this study.
Let C ^ be the Riemann sphere. Let R f : C ^ C ^ be the fixed-point operator of iterative scheme (11) on a function f, that is,
R f ( z ) = z β 6 u 3 + 2 u 2 + u + 1 f ( z ) f ( z ) ,
where
u = f ( y ) f ( z ) and y = z f ( z ) f ( z ) .
In order to generalize the dynamical analysis of family (11), we prove the following result.
Theorem 2
(Scaling Theorem). Let A ( z ) = α z + μ , α 0 be an affine mapping. In addition, let f ( z ) be a holomorphic function and h ( z ) = λ ( f A ) ( z ) with λ 0 ; then, the fixed-point operator R ¯ f is an affine conjugate to R ¯ h by A, i.e.,
A R h ¯ A 1 ( z ) = R f ¯ ( z ) .
Proof. 
We will prove the following identity:
R f ¯ A ( z ) = A R h ¯ ( z ) , for all z C ^ .
R f ( z ) = z β 6 u 3 + 2 u 2 + u + 1 f ( z ) f ( z ) = z β 6 f ( y ) f ( z ) 3 + 2 f ( y ) f ( z ) 2 + f ( y ) f ( z ) + 1 f ( z ) f ( z ) .
By developing the left side,
R f ¯ A ( z ) = A ( z ) β 6 f ( A ( y ) ) f ( A ( z ) ) 3 + 2 f ( A ( y ) ) f ( A ( z ) ) 2 + f ( A ( y ) ) f ( A ( z ) ) + 1 f ( A ( z ) ) f ( A ( z ) ) .
Furthermore, we know that the fixed-point operator related to function h, for α 0 , is:
R h ( z ) = z β 6 u 3 + 2 u 2 + u + 1 h ( z ) h ( z ) = z β 6 λ f ( A ( y ) ) λ f ( A ( z ) ) 3 + 2 λ f ( A ( y ) ) λ f ( A ( z ) ) 2 + λ f ( A ( y ) ) λ f ( A ( z ) ) + 1 λ f ( A ( z ) ) λ α f ( A ( z ) ) = z 1 α β 6 f ( A ( y ) ) f ( A ( z ) ) 3 + 2 f ( A ( y ) ) f ( A ( z ) ) 2 + f ( A ( y ) ) f ( A ( z ) ) + 1 f ( A ( z ) ) f ( A ( z ) ) .
Considering that A ( ξ 1 ξ 2 ) = A ( ξ 1 ) A ( ξ 2 ) + μ , for all ξ 1 , ξ 2 C ,
A R h ¯ ( z ) = A z 1 α β 6 f ( A ( y ) ) f ( A ( z ) ) 3 + 2 f ( A ( y ) ) f ( A ( z ) ) 2 + f ( A ( y ) ) f ( A ( z ) ) + 1 f ( A ( z ) ) f ( A ( z ) ) = A ( z ) A 1 α β 6 f ( A ( y ) ) f ( A ( z ) ) 3 + 2 f ( A ( y ) ) f ( A ( z ) ) 2 + f ( A ( y ) ) f ( A ( z ) ) + 1 f ( A ( z ) ) f ( A ( z ) ) + μ = A ( z ) β 6 f ( A ( y ) ) f ( A ( z ) ) 3 + 2 f ( A ( y ) ) f ( A ( z ) ) 2 + f ( A ( y ) ) f ( A ( z ) ) + 1 f ( A ( z ) ) f ( A ( z ) ) .
From the results obtained in (14) and (15),
R f ¯ A ( z ) = A R h ¯ ( z ) .
Therefore, the method in (11) satisfies the scaling theorem. □
Now, we analyze the dynamical behavior of this class in terms of parameter β in order to determine the methods with good stability properties and avoid the members of the family whose behavior is unstable.

4. Dynamical Analysis

In this part, we use several tools of complex dynamics. Thus, we need the function f to be defined on the Riemann sphere C ^ .
In the following, we recall some concepts of complex dynamics. Let us suppose a fixed-point function on an arbitrary polynomial p ( z ) = ( z a ) ( z b ) ; this gives a rational function, R .
Indeed, if given a rational function R : C ^ C ^ , where C ^ represents the Riemann sphere, the orbit at a point z 0 can be defined as
{ z 0 , R ( z 0 ) , R 2 z 0 , R 3 z 0 , , R n z 0 , } .
It is necessary to study the phase plane of R by classifying the fixed points according to the asymptotic behavior related to their orbits. A z 0 C ^ is said to be a fixed point if R z 0 = z 0 . A periodic point z 0 , of period p > 1 , is a point where R p z 0 = z 0 and R k z 0 z 0 , where k < p . A preperiodic point is a point z 0 that is not periodic, but there exists a k > 0 such that R k z 0 is periodic. A point z 0 is a critical point of the rational function R if R ( z 0 ) = 0 . If a critical point is different from any of those related to the roots of polynomial p ( z ) , it is called a free critical point. In a similar way, the fixed points that are different from the roots of p ( z ) are called strange fixed points.
Moreover, a fixed point z 0 is called an attractor if R z 0 < 1 , a superattractor if R z 0 = 0 , repulsive if R z 0 > 1 , and parabolic if R z 0 = 1 . Fixed points other than those associated with the roots of the polynomial p ( z ) are called strange fixed points. The basin of attraction of an attractor α is defined as:
A ( α ) = z 0 C ^ : R n z 0 α , n .
For the rational function R , the Fatou set of R , which is denoted by F ( R ) , is the set of points z C ^ whose orbits tend to one attractor (a fixed point or periodic orbit). Its complement in C ^ is the Julia set, J ( R ) . That is, the basin of attraction from every fixed or periodic point belongs to the Fatou set, while the boundaries of these basins of attraction are in the Julia set (for more explanation, see [12,13,14]).
The following result justifies the interest in the parameter planes that we will introduce later. The connected component of the basin of attraction, which contains the point or periodic orbit of attraction, is called the immediate basin of attraction.
Theorem 3
(Fatou–Julia). Let R be a rational function. An immediate basin of attraction of any attractor has at least one critical point.
To any rational operator that satisfies the scaling theorem, we can apply a Möbius transformation, which is stated as follows. If p ( z ) = ( z a ) ( z b ) , then R p depends on parameters a and b. To eliminate this dependency, we apply a Möbius transformation h as follows:
O p ( z ) = h R p h 1 ( z ) ,
where
h ( z ) = z a z b ,
which satisfies
  • h ( ) = 1 ,
  • h ( a ) = 0
  • h ( b ) = .
Operators O p and R p have the same qualitative properties of stability.
With Theorem 2, we may find all of the stable behaviors of a rational function R by analyzing its performance on the set of critical points, as we see below.
Using polynomial p ( z ) = ( z a ) ( z b ) , the iterative scheme (11) on this polynomial provides the following rational operator:
M p ( z , β , a , b ) = z + ( z a ) ( z b ) 1 + ( z a ) ( z b ) 6 ( 2 z a b ) 2 q ( z ) + ( z a ) 2 ( z b ) 2 β 6 ( 2 z a b ) 6 2 z a b ,
where q ( z ) = a 2 + 4 a b + b 2 6 ( a + b ) z + 6 z 2 . By applying the Möbius transformation on M p ( z , β , a , b ) , we obtain the rational operator to work, which is an operator free of parameters a and b.
h ( u ) = u a u b ,
and its inverse function is
h 1 ( u ) = u b a u 1 .
To get the fixed-point operator,
R p ( z , β ) = h M p h 1 ( z ) = z 4 6 ( z + 1 ) 2 ( 5 + z ( z + 4 ) ) + β 6 ( z + 1 ) 2 ( 1 + z ( 5 z + 4 ) ) + z 4 β .
So, M p becomes a conjugate using the Möbius transformation for R p ( z , β ) in order to simplify the dynamical analysis. Note that the rational operator R p ( z , β ) is not reliant upon the roots of p ( z ) , a, and b.

4.1. Stability of the Fixed Points

In this section, we analyze the stability of fixed points. We give some important results on the stability of strange fixed points, and we analyze the behavior of the rational operator on free critical points. The strange fixed points are z = 1 , and the roots of the polynomial are
q ( z , β ) = 6 z 6 + 42 z 5 + 126 z 4 + ( 180 + β ) z 3 + 126 z 2 + 42 z + 6 ,
which is denoted by s j ( β ) , j = 1 , 2 , , 6 .
The stability of the strange fixed points is analyzed in the following result.
Theorem 4.
The fixed points z = 0 and z = , which are associated with the roots a and b, are superattractors regardless of the value of β. z = 1 is a strange fixed point when β 240 , and its stability is determined as follows:
(a) 
If | β 240 | > 768 , then the fixed point z = 1 is an attractor,
(b) 
If | β 240 | < 768 , then the fixed point z = 1 is repulsive,
(c) 
If | β 240 | = 768 , then the fixed point z = 1 is neutral.
Proof. 
By evaluating the rational operator R p ( z , β ) , as well as its derivative at z = 0 , we conclude that z = 0 is a superattractor fixed point. For the punto z = , we derive the inverse of the operator R p ( z , β ) , i.e., 1 R p 1 z , and we evaluate it in z = 0 . We can affirm that z = is a superattractor fixed point.
Now, we evaluate R p ( z , β ) on the fixed point z = 1 in order to determine its stability.
(a)
R p ( 1 , β ) = 768 β 240 < 1 , which implies that 768 < | β 240 | , and therefore R p ( 1 , β ) < 1 if and only if | β 240 | > 768 and β 240 .
(b)
R p ( 1 , β ) = 768 β 240 > 1 , which implies that 768 < | β 240 | , and, therefore, R p ( 1 , β ) > 1 if and only if | β 240 | < 768 and β 240 .
(c)
Finally, R p ( 1 , β ) = 768 β 240 = 1 , which implies that 768 = | β 240 | , and, therefore, R p ( 1 , β ) = 1 if and only if | β 240 | = 768 and β 240 .
This completes the proof. □
In Figure 1, the stability of the strange fixed point z = 1 is presented.
Theorem 4 shows that we have two fixed points that are superattractors, z = 0 and z = , and the strange fixed point z = 1 , for which we find different regions where this behaves as an attractor, repulsor, and parabolic for β 240 . However, when we solve the equation R p ( z , β ) = z , six other solutions appear; these are strange fixed points that are solutions of the polynomial q ( z , β ) (18), which is denoted by s j ( β ) , j = 1 , 2 , , 6 . In the following result, we analyze the stability of these strange fixed points.
Theorem 5.
The strange fixed points of the rational operator R p ( z , β ) , s j ( β ) , j = 1 , 2 , are all repulsors in a complex plane.
On the other hand, the strange fixed points s j ( β ) , j = 3 , 4 are:
1. 
attractors in the open interval 32 9 < β < 3.55728 . In particular, they are superattractors for β 3.55728 4.9603 i .
2. 
parabolics for β = 4 9 ( 8 5 i 5 ) .
3. 
repulsors in the rest of the complex plane.
Finally, the strange fixed points s j ( β ) , j = 5 , 6 are:
1. 
attractors in 524.182 < β < 289.865 and 289.865 < β < 189.865 ) . In particular, they are superattractors for β 289.865 .
2. 
parabolic at all points on the boundary of the circumference of the cone—in particular, for the complex numbers β = 4 9 ( 8 5 i 5 ) .
3. 
repulsors in the rest of the complex plane.
Theorem 5 cannot be analytically proven, since the substitution of these strange fixed points into the stability function results in roots of higher-order polynomials. Still, the stability regions in the complex plane of the strange fixed points s 1 ( β ) and s 2 ( β ) can be visualized in Figure 2a, the stability regions of the fixed points s 3 ( β ) and s 4 ( β ) can be visualized in the complex plane in Figure 2b, and the stability regions of the fixed points s 5 ( β ) and s 6 ( β ) can be visualized in the complex plane in Figure 2c. Note that in the foreground, no stability region is visualized, which justifies the behavior of such strange fixed points; they are repulsors in the whole complex plane, and the strange fixed points s 3 ( β ) and s 4 ( β ) are attractors in the small region described in Theorem 5, while the strange fixed points s 5 ( β ) and s 6 ( β ) , are parabolic in the boundary of the circumference of the cone and are attracting inside the cone, being superattracting at the vertex of the cone. They are repulsive in the rest of the complex plane.

4.2. Critical Points

The critical points are obtained by solving equation R p ( z , β ) = 0 . We obtain the points z = 0 , z = , which correspond to the roots of polynomial p ( x ) , and the free critical points described in the following results.
Theorem 6.
The free critical points of operator R p ( z , β ) are:
1. 
c r j ( β ) = 120 + 3 β 7 240 β β 2 4 ( β 30 ) , j = 1 , 2 .
2. 
z = 1 .
Let us also remark that the critical points satisfy c r 1 ( β ) = 1 c r 2 ( β ) . So, they are not independent, and they have the same parameter plane.
There are also values of parameter β for which the critical points coincide:
(a)
The critical points c r 1 ( β ) = c r 2 ( β ) = c r 3 ( β ) match for β = 0 .
(b)
The points c r 1 ( β ) = c r 2 ( β ) match if β = 240 .
The free critical point z = 1 is a pre-image of the strange fixed point z = 1 because R ( 1 , β ) = 1 ; therefore, there is no need to obtain its related parameter plane because the information that would be obtained is contained in the stability function of z = 1 .

5. Parameter and Dynamical Planes

In this section, our goal is to determine the regions of the parameter plane for which we can obtain the most stable members of the family, since the values of β in those regions will provide us with the best methods in terms of numerical stability.

5.1. Parameter Planes

As c r 1 ( β ) and c r 2 ( β ) are conjugates, there exists just one free independent critical point, the related parameter plane of which can be seen in Figure 3. This is focused in a black region that mostly corresponds to regions shown in the stability planes of the strange fixed points s 3 ( β ) and s 4 ( β ) and the bases of the cones defining the stability functions of the strange fixed points s 5 ( β ) and s 6 ( β ) (see Figure 2b,c). Let us assume that there are free independent critical points acting in the capacity of the starting point of the iterative schemes of the family associated with each complex value of β ; such points are shown in red on the complex plane if the method converges on any of the roots (zero and infinity) and in black in other cases. Following this procedure, we obtain the parameter plane presented in Figure 3 by using the processes shown in [10]. For this representation, we used a mesh of 1500 × 1500 points, a maximum of 500 iterations, and a tolerance of 10 3 in the stopping criterion. We also show a close-up of this parameter plane in Figure 3b, focusing on the larger red area; there, the family members are, in general, very stable. These will be the best values to use for iterative methods in terms of their stability (see [15]).
We emphasize how important the parameter planes are, since they allow us to precisely identify the values of the β parameter for which the family of iterative methods have stable behavior or the members have chaotic behavior. In these, we can observe two black- and red-colored regions; the red-colored region allows us to obtain dynamic planes whose behavior is stable, and the black-colored regions imply methods that mostly correspond to values of the parameter in which there are strange fixed points that are attractors or have periodic attractor orbits; this is guaranteed by the (Fatou–Julia), Theorem 3.

5.2. Dynamical Planes

The dynamical planes appearing in Figure 4 and Figure 5 were obtained using a mesh of 2000 × 2000 points of the complex plane, as well as a maximum of 1000 iterations. These show the basins of attraction that correspond to different values that we selected from the parameter plane—both from the red-colored zone, where we had more stable members of the family, and the black-colored zones of the parameter plane, where we had unstable members of the family. In both cases, different stable behaviors, unstable behaviors, and periodic basins of attraction, among others, are shown.
In each of the dynamical planes, the circles represent the fixed points of R ( z , β ) (with the independence of their stability), the squares represent the critical points, and the starred point represents the superattracting fixed points.

Unstable and Stable Dynamic Planes

For the values β = 0 , β = 1 , β = 10 , and β = 100 (Figure 5a–d, respectively), we observe stable behavior corresponding to the red-colored regions in the parameter plane. The orange regions in the dynamical planes come from the superattracting fixed point z = 0 , and the blue area comes from the fixed point z = . Finally, Figure 4a–d correspond to regions of the parameter plane where we obtain families whose behavior is unstable. Moreover, for β = 350 , we obtain a dynamical plane with two periodic orbits at points z = 0.254854 0.96698 i and z = 0.897069 + 0.44189 i . It is worth noting that the dynamical plane corresponding to z = 200 650 i has three basins of attraction, which are shown in green, orange, and blue colors; for z = 300 , there are four basins of attraction, which are shown in green, orange, blue, and red; for z = 700 , there are three basins of attraction, unlike the stable dynamical planes, which only have two.

6. Numerical Results

We use the functions proposed in Table 1 to perform some numerical tests in order to analyze the results obtained with different known methods in comparison with the one proposed in the iterative family (11).
We call the methods used for the numerical tests M1, M2, and M3. M3 is the method that we worked with in this manuscript, and M1 and M2 are two known fourth-order methods that we use to compare the numerical results; these are the scheme of C. Chun and the BCMT family, which can be found in [16,17].
M1 = y n = x n f x n f x n x n + 1 = y n f x n + 2 f y n f x n f y n f x n ,
where n = 0 , 1 , 2 .
M2 = y n = x n f x n f x n x n + 1 = x n f x n f x n α + 2 3 α 2 f x n f x n + α f y n 3 2 ( 1 + α ) α 2 f x n f x n + α f y n + 3 α 2 + 5 α + 4 3 α 2 .
where α = 1 and n = 0 , 1 , 2 , .
M3 = y n = x n f x n f x n , x n + 1 = x n f x n f x n β 6 f ( y n ) f ( x n ) 3 + 2 f ( y n ) f ( x n ) 2 + f ( y n ) f ( x n ) + 1 .
where β = 10 and n = 0 , 1 , 2 , .
The following expression can be used to estimate the theoretical order of convergence:
p C O C = ln x n + 1 ξ / x n ξ ln x n ξ / x n 1 ξ ,
as presented by S. Weerakoon and T.G.I. Fernando in [18]. However, the principal inconvenience of this COC is that it requires knowledge of the precise root ξ . In order to solve such a problem, in [19], Cordero and Torregrosa overhauled the definition of COC so that it would not need the root to be exactly known at all, as shown in the following:
p A C O C = ln x n + 1 x n / x n x n 1 ln x n x n 1 / x n 1 x n 2 .
Calculations were performed by using the Matlab programming software with variable-precision arithmetic, with 200 digits of mantissa and ϵ = 10 8 as the error tolerance. For the computer programs, the following stopping criteria were used:
(i)
x n + 1 x n < ϵ or
(ii)
f x n + 1 < ϵ .
Table 2 displays the best outcomes for each experiment. We can see how the following method is highly efficient compared to schemes such as Chun’s and the BCMT, which we used for comparison.

7. Conclusions

In our manuscript, we found a large, optimal, and general kind of fourth-order scheme that is free of second derivatives. The dynamical analysis that was performed on the proposed class of methods allowed us to choose family members that were particularly stable, discarding the ones with unstable performance, where strange fixed points and some periodical orbits were present. Lastly, we refer to Table 2 to show that the suggested family is equal to or more competitive than the classical method, the recognized and efficient methods of Ostrowski and Maheshwari, and other recognized and efficient methods found in the literature. The dynamic study shows us that there are values of the parameter β for which we can obtain highly efficient iterative families of methods; some of the values that we can use are β = 0 , β = 1 , β = 10 , and β = 100 .
This behavior is justified in the dynamic planes given in Figure 5a–d, since there is global convergence to the roots of a second-degree polynomial.
In addition, if we use parameter values corresponding to regions of instability, one of the consequences of using these schemes is that these do not converge to the root, and if they do converge, they would do so with a large number of iterations, that is, these members are unreliable. This fact can be confirmed in Figure 4a–d.
The results shown in Table 2 and Table 3 are satisfactory from the point of view of numerical calculations. It can be observed that the number of iterations, computational time, and theoretical convergence order reached by our family, in comparison with the known methods of Chun and the BCMT, are satisfactory, as we expected, since we took one of the values of the b e t a parameter for which it was theoretically shown that we would have efficient iterative methods. In addition, another piece of evidence that shows the degree of accuracy with which our method converges to the root is the computation of the absolute errors x n + 1 x n and f ( x n + 1 ) , which are shown in Table 3.

Author Contributions

Conceptualization, A.C. and J.R.T.; software, M.A.L.S.; validation, A.C. and J.R.T.; formal analysis, J.R.T.; investigation, A.C. and M.A.L.S.; writing—original draft preparation, M.A.L.S.; writing—review and editing, A.C. and J.R.T.; visualization, M.A.L.S.; supervision, A.C. and J.R.T. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by PGC2018-095896-B-C22 and by MCIN/AEI/10.13039/5011000113033 as part of “ERDF: A way to making Europe”.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Acknowledgments

The authors would like to thank the anonymous reviewers for their comments and suggestions, as they have improved the final version of this manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Traub, I.F. Iterative Methods for the Solution of Equations; Prentice-Hall: Hoboken, NJ, USA, 1964. [Google Scholar]
  2. Petković, M.S.; Neta, B.; Petković, L.D.; Džunić, J. Multipoint Methods for Solving Nonlinear Equations; Academic Press: Amsterdam, The Netherlands, 2013. [Google Scholar]
  3. Amat, S.; Busquier, S. Advances in Iterative Methods for Nonlinear Equations; SEMA SIMAI Springer Series; Springer: Cham, Switzerland, 2016; Volume 10. [Google Scholar]
  4. Ostrowski, A.M. Solutions of Equations and System of Equations; Academic Press: New York, NY, USA, 1960. [Google Scholar]
  5. Kung, H.T.; Traub, J.F. Optimal order of one-point and multi-point iteration. ACM 1974, 21, 643–651. [Google Scholar]
  6. King, R. A family of fourth order methods for nonlinear equations. SIAM J. Numer. Anal. 1973, 10, 876–879. [Google Scholar] [CrossRef]
  7. Jarratt, P. Some fourth order multipoint methods for solving equations. Math. Comput. 1966, 20, 434–437. [Google Scholar] [CrossRef]
  8. Traub, J.F. Iterative Methods for the Solution of Equations; Chelsea Publishing Company: New York, NY, USA, 1982. [Google Scholar]
  9. Chicharro, F.I.; Cordero, A.; Garrido, N.; Torregrosa, J.R. Wide stability in a new family of optimal fourth-order iterative methods. Comput. Math. Methods 2019, 1, e1023. [Google Scholar] [CrossRef] [Green Version]
  10. Maheshwari, A.K. A fourth order iterative method for solving nonlinear equation. Appl. Math. Comput. 2009, 211, 383–391. [Google Scholar] [CrossRef]
  11. Wang, X.; Chen, X. The Dynamical Analysis of a Biparametric Family of Six-Order Ostrowski-Type Method under the Möbius Conjugacy Map. Fractal Fract. 2022, 6, 174. [Google Scholar] [CrossRef]
  12. Fatou, P. Sur les substitutions rationnelles. Comp. Rend. Heb. S. Acad. Sci. 1917, 164, 806–808, 1917, 165, 992–995. [Google Scholar]
  13. Julia, G. Mémoire sur l’iteration des fonctions rationnelles. Mat. Pur. Appl. 1918, 8, 47–245. [Google Scholar]
  14. Blanchard, P. Complex analytic dynamics on the riemann sphere. Bull. Am. Math. Soc. 1984, 11, 85–141. [Google Scholar] [CrossRef] [Green Version]
  15. Chicharro, F.I.; Cordero, A.; Torregrosa, J.R. Drawing Dynamical and Parameters Planes of Iterative Families and Methods. Sci. World J. 2013, 2013, 780153. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Behl, R.; Cordero, A.; Motsa, S.S.; Torregrosa, J.R. An optimal fourth-order family of methods for multiple roots and its dynamics. Numer. Algor. 2016, 71, 775–796. [Google Scholar] [CrossRef]
  17. Chun, C. Construction of Newton-like iterative methods for solving nonlinear equations. Numer. Math. 2006, 104, 297–315. [Google Scholar] [CrossRef]
  18. Weerakoon, S.; Fernando, T.G.I. A variant of Newton’s method with accelerated third-order convergence. Appl. Math. Letters. 2000, 13, 87–93. [Google Scholar] [CrossRef]
  19. Cordero, A.; Torregrosa, J.R. Variants of Newton’s method using fifth-order quadrature formulas. Appl. Math. Comput. 2007, 190, 686–698. [Google Scholar]
Figure 1. Strange fixed point stability region where z = 1 .
Figure 1. Strange fixed point stability region where z = 1 .
Algorithms 15 00387 g001
Figure 2. Strange fixed points s j ( β ) , j = 1 , , 6 . (a) s 1 ( β ) and s 2 ( β ) . (b) s 3 ( β ) and s 4 ( β ) . (c) s 5 ( β ) and s 6 ( β ) .
Figure 2. Strange fixed points s j ( β ) , j = 1 , , 6 . (a) s 1 ( β ) and s 2 ( β ) . (b) s 3 ( β ) and s 4 ( β ) . (c) s 5 ( β ) and s 6 ( β ) .
Algorithms 15 00387 g002
Figure 3. Parameter planes ( β ) for free critical points. (a) Parameter plane with z = c r 1 ( β ) . (b) Zoomed-in view of the parameter plane.
Figure 3. Parameter planes ( β ) for free critical points. (a) Parameter plane with z = c r 1 ( β ) . (b) Zoomed-in view of the parameter plane.
Algorithms 15 00387 g003
Figure 4. Dynamical planes whose behavior is unstable. (a) β = 200 650 i . (b) β = 300 . (c) β = 700 . (d) β = 350 .
Figure 4. Dynamical planes whose behavior is unstable. (a) β = 200 650 i . (b) β = 300 . (c) β = 700 . (d) β = 350 .
Algorithms 15 00387 g004
Figure 5. Dynamical planes whose behavior is stable. (a) β = 1 . (b) β = 10 . (c) β = 100 200 i . (d) β = 0 .
Figure 5. Dynamical planes whose behavior is stable. (a) β = 1 . (b) β = 10 . (c) β = 100 200 i . (d) β = 0 .
Algorithms 15 00387 g005
Table 1. Test functions.
Table 1. Test functions.
f ( x ) ξ
f 1 ( x ) = e x + cos x 1.7461395304080124176507030889537802
f 2 ( x ) = x 4 + sin π x 2 5 1.4142135623730950488016887242096980
f 3 ( x ) = e x 1.5 tan 1 x 0.7676532662012788981900298911398052
f 4 ( x ) = 10 x e x 2 1 1.6796306104284499406749203388379704
Table 2. Benchmark comparison of fourth-order methods in relation to the convergence order and computational time on the function f ( x ) .
Table 2. Benchmark comparison of fourth-order methods in relation to the convergence order and computational time on the function f ( x ) .
f ( x ) IterationsComputational Time in (s) ACOC
M1M2M3M1M2M3M1M2M3
f 1 ( x ) , 2 3330.0502110.0446470.0539684.14973.47474.1213
f 2 ( x ) , 2 5440.1856600.1033320.0962033.99283.95844.0173
f 3 ( x ) , 1 3330.0548890.0585690.0610123.76523.66403.8134
f 4 ( x ) , 1.9 4540.0748890.1237470.0983243.96303.98853.9945
Table 3. Benchmark comparison of fourth-order methods in relation to the absolute error x n + 1 x n and f ( x n + 1 ) .
Table 3. Benchmark comparison of fourth-order methods in relation to the absolute error x n + 1 x n and f ( x n + 1 ) .
f ( x ) x n + 1 x n f ( x n + 1 )
M1M2M3M1M2M3
f 1 ( x ) 7.3 × 10 16 2.59 × 10 16 1.67 × 10 16 1.16 × 10 62 2.87 × 10 64 2.58 × 10 65
f 2 ( x ) 1.16 × 10 27 1.35 × 10 13 2.31 × 10 19 2.19 × 10 207 2.17 × 10 50 3.84 × 10 74
f 3 ( x ) 5.69 × 10 10 9.24 × 10 9 1.33 × 10 10 5.64 × 10 27 8.89 × 10 32 1.11 × 10 39
f 4 ( x ) 3.7 × 10 19 1.3 × 10 25 2.18 × 10 28 2.22 × 10 73 7.57 × 10 99 1.79 × 10 110
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Cordero, A.; Leonardo Sepúlveda, M.A.; Torregrosa, J.R. Dynamics and Stability on a Family of Optimal Fourth-Order Iterative Methods. Algorithms 2022, 15, 387. https://doi.org/10.3390/a15100387

AMA Style

Cordero A, Leonardo Sepúlveda MA, Torregrosa JR. Dynamics and Stability on a Family of Optimal Fourth-Order Iterative Methods. Algorithms. 2022; 15(10):387. https://doi.org/10.3390/a15100387

Chicago/Turabian Style

Cordero, Alicia, Miguel A. Leonardo Sepúlveda, and Juan R. Torregrosa. 2022. "Dynamics and Stability on a Family of Optimal Fourth-Order Iterative Methods" Algorithms 15, no. 10: 387. https://doi.org/10.3390/a15100387

APA Style

Cordero, A., Leonardo Sepúlveda, M. A., & Torregrosa, J. R. (2022). Dynamics and Stability on a Family of Optimal Fourth-Order Iterative Methods. Algorithms, 15(10), 387. https://doi.org/10.3390/a15100387

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop