A Novel Family of Efficient Weighted-Newton Multiple Root Iterations

Deepak Kumar 1,∗, Janak Raj Sharma 2 and Lorentz Jăntschi 3,4,* 1 Department of Mathematics, Chandigarh University, Mohali 104413, India 2 Department of Mathematics, Sant Longowal Institute of Engineering and Technology, Longowal 148106, India; jrsharma@sliet.ac.in 3 Department of Physics and Chemistry, Technical University of Cluj-Napoca, 400114 Cluj-Napoca, Romania 4 Institute of Doctoral Studies, Babeş-Bolyai University, 400084 Cluj-Napoca, Romania * Correspondence: deepak.e9086@cumail.in (D.K.); lorentz.jantschi@chem.utcluj.ro (L.J.)


Introduction
Approximating the solution of nonlinear equations by numerical methods is an important problem in many branches of science and engineering. For example, problems from many areas such as Physics, Chemistry, Mathematical Biology, and Engineering science are reduced to finding solution of nonlinear equations [1][2][3][4]. In general, closed form solutions can not be obtained so researchers use iterative methods for approximating the solution. In this paper, our aim is to construct higher order multi-point iterative methods for the multiple zeros of the univariate function f (x), where f : C → C is analytic about the required zero. The advantages of multi-point methods over one-point methods are discussed in Traub's well-known book [4].
They have shown that this iterative scheme attains fourth-order convergence provided that the function G(u) satisfies the conditions In this article, we aim to develop multiple root solvers of high efficiency, meaning the methods with rapid convergence that require less computations. Proceeding in this way, we develop a family of seventh-order three-point Newton-type methods for computing multiple zeros. The proposed iterative scheme is the composition of three steps that uses the Liu-Zhou iteration (2) as the first two steps and a Newton-type iteration in the third step. The algorithm requires four function evaluations per iteration and, therefore, possesses the efficiency index 7 1/4 ≈ 1.627, which is better than the efficiency index 4 1/3 ≈ 1.587 of the basic method (2). In this sense, the proposed iteration is the modification over the iteration (2). Some special methods of the new family are established. The usefulness of the methods is demonstrated by performing numerical tests on several applied science problems. Thereby, we observe in each example that the new methods have far better numerical results than the existing methods. The convergence domains are also assessed using the graphical tool, namely, basins of attraction, which is an useful technique to check the convergence regions visually.
We summarize the contents of this article. In Section 2, the family of seventh-order iterative solvers is derived and its local convergence is studied. In order to check the convergence regions of the methods graphically, the basins of attractors are assessed in Section 3. Some numerical experiments are performed in Section 4 to verify the theoretical results and to compare the performance with the existing methods. In Section 5, concluding remarks are reported.

Development of Scheme
Our aim is to develop an iterative method for computing a multiple root with multiplicity m > 1, which accelerates the convergence rate of the two-step Liu-Zhou iteration (2) using a minimum number of functions and derivative evaluations. Thus, it will turn out to be judicious if we consider a somewhat complicated three-step iterative scheme of the following form, where u = f (y) m , a is scalar, and G, H : C → C are analytic functions about a neighborhood of 0. Note that the second and third steps are weighted by the factors G(u) and H(u), and, as a result, these are called weight factors or weight functions. Note that u and v are one-to-m − 1 and one-to-m multi-valued functions, respectively. Therefore, it is convenient to treat them as the principal root. As an example, we consider the case of v.

The principal root is given by
f (x) ≤ π; this convention of Arg(Z) for Z ∈ C agrees with that of Log[Z] command of Mathematica [21]. We employ this command in the later section of numerical simulation. Similarly, we treat for u.
In the following theorem, we shall prove the seventh-order convergence of the proposed scheme (3).

Theorem 1.
Let the function f : C → C be analytic in a domain containing multiple zero α of multiplicity m. Suppose that the starter x 0 is close enough to α. Then, the iterative technique expressed by (3) possesses seventh-order convergence, provided that the functions G(u) and H(u) satisfy the conditions Proof. Let the error at n-th stage be e n = x n − α. Using the Taylor's expansion of f (x n ) about α, we have that m! e m n 1 + C 1 e n + C 2 e 2 n + C 3 e 3 n + C 4 e 4 n + C 5 e 5 n + C 6 e 6 n + C 7 e 7 n + O(e 8 n ) (5) and where , k ∈ N. From (5), (6), and the first step of (3), we have that where ω i = ω i (m, C 1 , C 2 , . . . , C 7 ) are given in terms of m, C 1 , C 2 , . . . , C 7 with two explicitly written and ω 2 = 1 m 3 3m 2 C 3 + (m + 1) 2 C 3 1 − m(4 + 3m)C 1 C 2 . The remaining expressions of ω i (i = 3, 4, 5) are not being produced explicitly as these are very lengthy.
Expansion of f (y n ) and f (y n ) about α leads us to the expression and By using (6) and (9), we get expression of u n as where η i = η i (m, C 1 , C 2 , . . . , C 7 ) are given in terms of m, C 1 , C 2 , . . . , C 7 with one explicitly written coefficient, Expanding the function G(u) about origin using Taylor series By inserting the expressions (5), (6), and (11) in the second step of scheme (3) and simplifying, where In order to obtain higher order convergence, the coefficients of e n , e 2 n , and e 3 n should vanish. The resulting equations yield Combining Equations (12) and (13), we obtain that Developing Inserting (5) and (15) in the expression of v n , then where τ i = τ i (m, C 1 , C 2 , . . . , C 7 ). The Taylor expansion of H(u n ) about origin 0 is given by Therefore, by substituting (5), (6), (10), (16), and (17) into the last step of scheme (3), we obtain the error equation Error Equation (18) shows that at least fifth-order convergence is attained if H(0) = 1. Using this value in ξ 1 = 0, we will obtain H (0) = 2.

Remark 2.
Kung and Traub [23] have conjectured that the multi-point methods without memory requiring n functional evaluations can attain the maximum convergence order 2 n−1 . That means that with four functional evaluations one can develop a method of optimal order eight. The methods qualifying Kung-Traub hypothesis are also called optimal methods. Such methods are rare for multiple roots due to the complexity in finding the convergence order. Nevertheless, the proposed seventh methods are better than the existing sixth-order methods by Geum et al. [8] (see Formula (30) in the numerical section) in the sense that latter also require same number of evaluations, i.e., two functions and two derivatives.  [24]).

Some Special Cases
We can generate numerous special cases of the family (3) based on the forms of weight functions G(u) and H(u) that satisfy the conditions of Theorem 1. Some simple forms are given as follows, Case I: Let us describe the following polynomial forms of G and H directly from the proposed Theorem 1: Thus, the corresponding new seventh-order method (denoted by NM-I) is given by Case II: Let us describe the following functions that satisfy the conditions of Theorem 1: Thus, the corresponding seventh-order method (now denoted by NM-II) is expressed as Case III: Let us consider the following forms satisfying the conditions of Theorem 1: The corresponding new method (denoted by NM-III) is given by Case IV: Next, let us consider the following forms satisfying the conditions of Theorem 1: Then, the corresponding seventh-order iterative scheme (denoted by NM-IV) is given by

Complex Geometry of Methods
Here, we aim to assess the complex dynamics of new methods based on the geometrical technique, namely, basins of attraction, of the zeros of a polynomial f (z) in complex domain. Using this tool, one can get an important information about the stability and convergence of a method. The idea was introduced initially by Vrscay and Gilbert [25]. Recently, many authors have used this tool in their work, see, for example, in [2,26,27] and the references cited therein.
The initial point z 0 is chosen in a rectangular region R ∈ C containing all the zeros of f (z). Starting from the point z 0 , the method either converges to the zero of f (z) or eventually diverges. We choose the tolerance value 10 −3 up to maximum 25 iterations to stop the iteration process. If this tolerance is not attained in required iterations, then the method does not converge to any root. To plot the basins, we adopt the following strategy. A color is assigned to each point z 0 lying in the basin of corresponding root. Then, the point represents the attraction basin with that particular color provided that the method converges. Contrary to this, if the method fails to converge in specified iterations, then the point paints the black color.
In what follows we assess the basins of attraction by employing the methods NM-I-NM-IV on the following three polynomials. Problem 1. As the first example, consider the polynomial f 1 (z) = (z 3 + 4z) 3 which has zeros {0, ±2i} each with multiplicity 3. For drawing basins, we use a rectangle R of size [−3, 3] × [−3, 3] and fix the color green to each initial point in the basin of zero "−0", the color red to each point in the basin of zero "2i", and the color blue to every point in the basin of zero "−2i". Basins so drawn for the methods NM-I-NM-IV are shown in Figure 1. It is clear that the methods NM-IV and NM-II posses fewer divergent points (painted with black color), followed by NM-I and NM-III.  These graphics easily depict the convergence behavior of any method. If we select a value of z 0 in a place where different basins meet each other, it is difficult to guess which zero is going to be obtained by the method. Therefore, the selection of z 0 in such a region is not preferable. The zones with black color and with amalgam of different colors can not be suitable to choose the initial guess z 0 to acquire a unique root. The most attractive pictures are those with intricate boundaries between the basins of attraction. These boundaries have fractal-like shapes and correspond to the cases where the method is more demanding with choice of initial point. At such regions the dynamic behavior of initial guess is more unpredictable.

Numerical Examples
In this section, we employ the special cases NM-i, i = I, II, III, IV of family (3) on some nonlinear equations to test the validity of theoretical results derived in previous sections. Performance is compared with some well-known sixth-order methods such as the two-and three-point methods by Geum et al. [8,9]. The two-point method [8], applicable for m > 1, is given as where u = f (y n ) f (x n ) 1 m and s = f (y n ) f (x n ) 1 m−1 , and Q f : C 2 → C is a holomorphic function in neighborhood of origin (0, 0). Considering the following cases of function Q f (u, s) in the Formula (30) and denoting the corresponding iterative method by GKN-1(j), j = a, b, c, d: 4m 2 −8m+7 and d 1 = 2(m − 1). Next, the three-point method [9] for m ≥ 1 is expressed as wherein m . Functions Q f and K f are analytic in a neighborhood of 0 and (0, 0), respectively. We consider the following combinations of Q f (u) and K f (u, v), and denote the corresponding iterative schemes by GKN-2(j), j = a, b, c, d: Computations are performed in the programming package of Mathematica software using multiple-precision arithmetic. The results displayed in Tables 1-4 contain (i) the iteration number (n) in which the solution is obtained with required accuracy, (ii) the last three errors e n = |x n+1 − x n |, (iii) the computational order of convergence (COC), and (iv) the CPU time (CPU-time) elapsed during the execution of program. The required iteration number (n) and elapsed time are recorded when the criterion |x n+1 − x n | + | f (x n )| < 10 −350 is satisfied. Computational order of convergence (COC) is calculated by the formula (see [28])   For numerical testing we choose the following problems.

Example 1.
Finding the eigenvalues of the characteristic equation of a square matrix of order greater than 4 is a big problem [29]). We consider the following 9 × 9 matrix.
The characteristic polynomial of the matrix (M) is given by One of the eigen values is 3 with multiplicity 4. Numerical results produced for this example, taking initial guess x 0 = 2.25, are shown in Table 1.

Example 2.
We consider the manning problem arising in isentropic supersonic flow around a sharp expansion corner. Let b = γ+1 γ−1 , where γ is the specific heat ratio of the gas. The following relation is developed between the Mach numbers before the corner and after the corner denoted by M 1 and M 2 , respectively (see [3]), Let us consider a particular case: To solve the equation for M 2 given that M 1 = 1.5, γ = 1.4, and δ = 10 0 . Then, we have that Now considering this equation for three times, the required nonlinear function is given as This function has one zero 1.8411027704 . . . with multiplicity 3. To find the zero, we choose initial approximation x 0 = 1.50. The obtained results are shown in Table 2.
Example 3. Next we assume nonlinear test function (see [8]) The multiple zero of function f 3 is 2 with multiplicity 5, which is calculated using initial approximation x 0 = 1.5. The obtained results are shown in Table 3.

Example 4.
Consider the isothermal continuous stirred tank reactor (CSTR) problem [30]. Components A and R are fed to the reactor at the corresponding rates of Q and q − Q. Then, the following scheme develops in the reactor (see [30]), This reaction scheme was analyzed in [30] to design simple feedback control system. In the analysis, the following equation was given for the transfer function of the reactor, The results so produced by the methods using x 0 = −2.80 are shown in Table 4. We observe from the numerical results that the errors become smaller as the iterations proceed, which points to the increasing accuracy in the values of successive approximations. Per iteration, the number of significant figures gained by the proposed methods is larger than the existing methods because of the higher order. The reading '0 of error |e n | indicates that at this stage the stopping condition |x n+1 − x n | + | f (x n )| < 10 −350 has been reached. The results of penultimate column of each table support the theoretical seventh-order of convergence. This shows that the convergence order is preserved. However, this is not true for the existing sixth-order methods GKN-1(j), j= a, b, c, d, as the sixth-order convergence is not preserved in last problem. The computational efficiency can be observed by the readings of elapsed CPU-time displayed in the last column of each table. Indeed, the new methods are more efficient since they consume less execution time than the existing ones. We have also applied the methods on other different problems to confirm the accuracy and efficiency and results are found on a par with the above conclusions.

Conclusions
A class of seventh-order numerical methods has been designed for computing multiple zeros of nonlinear functions. Local convergence analysis has been shown under standard assumptions which proves the convergence order seven. Some particular cases have been explored and their performance has been checked by using two different ways viz. by numerical testing and by graphical tool of attraction basins. Comparison of performance of the methods with existing methods has also been shown. In addition, a comparison of estimated CPU-time has been performed in order to rank the algorithms. We emphasize that the ranking obtained in this way matches well with the ranking obtained from the computational efficiency. As remarked earlier that according to Kung-Traub conjecture one can develop a method with optimal eighth convergence using four function evaluations. Therefore, this will be a motivational factor for us in future endeavor to develop such methods