Abstract
For solving nonlinear systems of big size, such as those obtained by applying finite differences for approximating the solution of diffusion problem and heat conduction equations, three-step iterative methods with eighth-order local convergence are presented. The computational efficiency of the new methods is compared with those of some known ones, obtaining good conclusions, due to the particular structure of the iterative expression of the proposed methods. Numerical comparisons are made with the same existing methods, on standard nonlinear systems and a nonlinear one-dimensional heat conduction equation by transforming it in a nonlinear system by using finite differences. From these numerical examples, we confirm the theoretical results and show the performance of the presented schemes.
1. Introduction
The design of fixed point iterative methods for solving equations and systems of nonlinear equations is an important and challenging task in the field of numerical analysis. Nonlinearity is ubiquitous in physical phenomena. Fluid and plasma mechanics, gas dynamics, elasticity, relativity, chemical reactions, combustion, ecology, biomechanics, economics modeling problems, transport theory, and many other phenomena are all governed inherently by nonlinear equations. For this reason, an ever increasing proportion of modern mathematical research is devoted to the analysis of nonlinear systems and nonlinear phenomena. These and other more examples allow us to affirm that finding the solution of a nonlinear system is a classical and difficult problem with many applications in sciences and engineering, wherein is a sufficiently Fréchet differentiable function in an open convex set . In the last years, many iterative methods have been constructed for solving nonlinear systems; see, for example, [1–6] and the references therein. The best known method for finding a solution is Newton’s scheme:where is the Jacobian matrix of function evaluated in the th iteration.
Based on Newton’s or Newton-like iterations, some high-order methods for computing a solution of nonlinear system have been proposed in the literature. The aim of these new schemes is to accelerate the convergence or to improve the computational efficiency. For example, among other authors, Soleymani et al. in [7] constructed a multistep class of sixth-order iterative method for solving nonlinear systems; Hueso et al. in [4] developed sixth-order iterative methods requiring two evaluations of function and two of Jacobian per iteration; Sharma and Arora in [6] designed a sixth-order method which requires three functional and two Jacobian evaluations per iteration, and Xiao and Yin [8] developed a three-step iterative scheme of order five. On the other hand, Wang et al. in [9] constructed a seventh-order derivative free iterative method by evaluating the first order divided difference operator three times per iteration.
Specifically, the sixth-order scheme designed by Soleymani et al. [7], that we denote by SLB, is a three-step Jarratt-type method whose iterative expression is
Hueso et al. in [4] developed several iterative schemes of order six; we use one of them in the numerical section for comparing it with our proposed schemes on different test problems. In particular, in [4] the authors present the following method that we denote by HMT: where denotes the identity matrix of size .
On the other hand, Xiao and Yin in [8] described the following three-step fifth-order scheme that we denote by XY:
In order to compare the different methods under the point of view of the computational cost, Ostrowski in [10] defined the efficiency index as , where is the order of convergence and is the number of functional evaluations per iteration. Let us remark that for evaluating function we need scalar functional evaluations (the coordinate functions of ), whilst for evaluating Jacobian it is necessary to evaluate functions (all the entries of matrix ). On the other hand, all the iterative methods for solving nonlinear systems require one or more matrix inversion; that is, one or more linear systems must be solved. So, the number of operations needed for solving a linear system plays in this context an important role. For this reason, the authors introduced in [2] the computational efficiency index, , which combines the efficiency index defined by Ostrowski and the number of products-quotients required per iteration. We define this index as , where is the number of products-quotients per iteration.
We recall that the number of products and quotients required for solving a linear system by Gaussian elimination is , where is the size of the system. In addition, for solving linear systems with the same matrix of coefficients, by using LU decomposition, we need products-quotients. By using this information, in Section 3, we compare the computational efficiency indices of the different methods used in this manuscript.
The main objective of this paper is to develop high-order iterative methods in such a way that they involve as lower computational cost as possible and they have good stability properties on big-sized systems. From Newton’s method and by using variation parameters procedure we present a two-step iterative method of order five. This order can be increased in three units by adding a new step with the same structure as the second one. This idea is generalized for obtaining an iterative method of arbitrary order increasing in three units the order each time that we add a new step with the same structure as the previous one. In each new step we only need a new functional evaluation.
In order to analyze the convergence properties of the different new schemes that will be introduced in this paper, we need to recall several concepts and tools, some of them introduced by the authors in [2].
1.1. Basic Definitions
Let be a sequence in which converges to . Then, convergence is called of order , , if there exists ( if ) and such thatorwhere .
The following notation was introduced in [2], but we present it for completeness. Let be sufficiently Fréchet differentiable in . The th derivative of at , , is the -linear function such that . It is easy to observe the following:(1).(2), for all permutation of .
From the above properties we can use the following notation:(a).(b).
On the other hand, for lying in a neighborhood of a solution of , we can apply Taylor’s expansion and, assuming that the Jacobian matrix is nonsingular, we havewhere . We observe that since and .
In addition, we can express aswhere is the identity matrix. Therefore, . From (8), we obtainwhereThe equationwhere is a -linear function , is called the error equation and is the order of convergence. Observe that is .
We summarize the contents of this paper. In Section 2 we describe new iterative methods for solving nonlinear systems and show a procedure for constructing schemes with arbitrary order of convergence. The efficiency index and the computational efficiency index of our methods and of other known ones are analyzed in Section 3. Section 4 is devoted to numerical tests for confirming the theoretical results and to show the performance of the presented schemes. The examples used are a nonlinear one-dimensional heat conduction equation by transforming it in a big-sized nonlinear system by using finite differences and standard nonlinear systems. The paper finishes with some conclusions and the references used.
2. Development and Convergence of the Method
By adding a new step to Newton’s method, we construct the following two-step scheme with fifth-order of convergence:where , , and are arbitrary parameters and is the identity matrix of size .
The following result establishes the convergence of iterative method (12).
Theorem 1. Let be sufficiently differentiable at each point of an open neighborhood of , that is, a solution of the system , and the initial estimation is close enough to . Let us suppose that is continuous and nonsingular in . Then, sequence obtained from expression (12) converges to with order 5 if the parameters take values , with the error equation in this case beingwhere .
Proof. By using Taylor expansion of and around , From the above expression, we havewhereThen, So,On the other hand,In a similar way as before, we obtain where soTherefore, we obtain If we replace, in the previous expression, the values of parameters , , and , we obtain the error equation: and the proof is finished.
If we add in expression (12) a third step with the same structure as the second one, we haveIt can be proved that method (25) has eighth-order of convergence if the parameters take the values , and and, for these values of parameters, the method is denoted as CCGT1. This idea can be generalized for obtaining an iterative method of arbitrary order of convergence. With each new step the order of convergence increases in three units, needing only one new functional evaluation.
Theorem 2. Let be the iterative expression of a method of order , with asymptotic error constant , where the two first steps are those of (12). By adding a new step in the form then the order of the resulting method is , if , with its error equation beingwhere .
Proof. Let us suppose that where , .
Then,Therefore,and, if we take , we have order . In that case, the error equation is and the proof is finished.
A small variation in the iterative expression of (12) allows us to construct another new scheme:whose order of convergence is described in the following result.
Theorem 3. Let be sufficiently Fréchet differentiable at each point of an open neighborhood of , that is, a solution of the system , and the initial estimation is close enough to . Let us suppose that is continuous and nonsingular in . Then, sequence obtained from expression (32) converges to with order 5 if the parameters take the values , and . In this case, the error equation is
In a similar way as before, this structure can be extended in order to construct an iterative scheme of arbitrary order.
Theorem 4. Let be the iterative expression of a method of order , with asymptotic error constant , where the two first steps are those of (32). By adding a new step in the form then the order of the resulting method is , if , with the error equation being
We will denote by CCGT2 the following three-step eighth-order iterative method:
3. Computational Efficiency
We are going to use two indices for comparing the different iterative schemes for solving nonlinear systems: the multidimensional extension of the efficiency index defined by Ostrowski as and the computational efficiency index defined in the Introduction as , where is the order of convergence, is the number of functional evaluations per iteration, and is the number of products-quotients per iteration. To compute in any iterative method we need to calculate scalar functions. The number of scalar functional evaluations is for any new evaluation of derivative .
In Table 1, the efficiency indices of methods CCGT1, CCGT2, SLB, XY, and HMT are presented. The number of Jacobian evaluations is the same in all these schemes, but the number of functional evaluations and the order of convergence are different.
In Figure 1 we show the classical index defined by Ostrowski for the mentioned methods and systems of sizes from 2 to 50. We can observe that the index of CCGT1 and CCGT2 is the same as well as the index of SLB and HMT. For every size, except , the best index corresponds to method CCGT1.
(a) Sizes from 2 to 20
(b) Sizes from 21 to 50
On the other hand, in order to compute an inverse linear operator we solve a linear system where we have to do products-quotients for obtaining decomposition and solving two triangular linear systems. In addition, we need products for matrix-vector multiplication.
Taking into account the previous considerations, we calculate of method CCGT1. For each iteration, we need to evaluate function three times and twice Jacobian , so functional evaluations are needed. In addition, we must solve one linear system with as coefficients matrix (i.e., products-quotients), six linear systems with as coefficients matrix (i.e., products-quotients), and four matrix-vector products ( products-quotients). Therefore, the value of index for method CCGT1 on a nonlinear system of size is
In Table 2, we show index of schemes CCGT1, CCGT2, SLB, XY, and HMT. In it, NFE is the number of functional evaluations, NLS1 denotes the number of linear systems with the matrix of coefficients to be solved, NLS2 is the number of linear systems with another matrix of coefficients that are solved, and denotes the number of products matrix-vector.
Let us observe that, although the classical index is similar in all these cases, it is not the case of the computational efficiency index since the number of inverse linear operators is different for each scheme. In Figure 2 the computational efficiency index for several methods and systems of sizes from 2 to 50 is shown. We can observe that until the best index corresponds to method XY; meanwhile, for the best index is the one corresponding to CCGT2.
(a) Sizes from 2 to 20
(b) Sizes from 21 to 50
4. Numerical Results
We check the numerical behavior of our method on a nonlinear one-dimensional heat conduction equation. A heat transfer problem is said to be one-dimensional if the temperature in the medium varies in one direction only and thus heat is transferred in one direction, and the variation of temperature and thus heat transfer in other directions are negligible or zero. For example, heat transfer through the glass of a window can be considered to be one-dimensional since heat transfer through the glass occurs predominantly in one direction (the direction normal to the surface of the glass) and heat transfer in other directions (from one side edge to the other and from the top edge to the bottom) is negligible.
To describe a heat transfer problem completely, an initial condition () and two boundary conditions must be given for each direction of the coordinate system along which heat transfer is significant. Therefore, we need to specify two boundary conditions for one-dimensional problems, four boundary conditions for two-dimensional problems, and six boundary conditions for three-dimensional problems. Different authors have approximated the solution of these problems by means of numerical techniques; see, for example, [11, 12] and the references therein.
In our study a particular case is used, corresponding to the following heat conduction equation:where . The initial condition is and the boundary conditions areBy applying an implicit method of finite differences we can transform problem (38) in a family of nonlinear systems, which provides the approximated solution in a time from the approximated solution in . We choose the spacial step and the temporal step , where and are the number of -subintervals and -subintervals, respectively, and is the final instant of our study, so we have selected a grid of domain with points ,We want to estimate the solution of (38) at these point, by transforming it in many nonlinear systems, as much as the number of . To do that, we use the following approximations:Denoting by the estimation of the solution at and by replacing them in (38), we construct the following nonlinear system:for and . Some algebraic manipulations allow us to transform this system in for and . For a fixed , we have the following nonlinear system of size :The unknowns of this system are , that is, the approximations of the solution in . We observe that for solving this system we need the solution in .
We are going to solve this system for different values of , using in each case or higher, , and methods CCGT1 and CCGT2. As initial guess, we use the solution at . We compare the obtained values with the exact solution , in order to analyze the stability and consistence of the new method.
All computations are performed in the programming package MATLAB using variable precision arithmetic with 50 digits of mantissa. For every value of , we analyze the mean number of iterations (iter) needed to converge to the solution such that is satisfied, where denotes the Euclidean norm. Let us recall that the iterative methods are applied for solving the nonlinear system that involves each column of the solution matrix; the mean of the number of iterations needed when all the columns have been calculated is shown in the tables. The processor of the machine used is Intel(R) Xeon(R) CPU E5-2420 v2 @ 2.20 GHz, with 64 GB of RAM. The results are presented in Tables 3 and 4, where the good performance of both methods in terms of exact error and CPU time (in seconds) is observed, even for .
In Figure 3(a) we show the approximated solution of the problem when , by taking and and using method CCGT1. It is a good approximation if we compare it with the exact solution (Figure 3(b)). In Figure 4 the absolute value of the exact error is presented.
(a) Approximated
(b) Exact
In the rest of the section, the new method CCGT1 is compared with methods SLB, HMT, and XY for solving some academical nonlinear systems. The numerical results are shown in Tables 5, 6, and 7. All experiments have been carried out on MATLAB using variable precision arithmetic with 2000 digits of mantissa. To verify the theoretical order of convergence , we calculate the approximated computational order of convergence (ACOC) introduced in [13] as
Example 1. The first nonlinear system is defined by (see [9])
In this test example we use and the initial estimation , with the solution being . In Table 5 we show the values of and for as well as the value of ACOC.
Example 2. The second nonlinear system is also defined in [9]:
The numerical results are displayed in Table 6. The initial estimation is and the size of the system is , with the solution being . We show the same information as in the previous example.
Example 3. Finally, the third example is
The solution of this system is . The numerical results of this example by using the initial estimation and are displayed in Tables 7 and 8, respectively. For the first initial guess, the results are very bad, and only method CCGT1 is convergent. For the second one, the numerical solutions are standard and confirm the theoretical results.
5. Conclusions
We have designed two three-step iterative methods of order eight for solving nonlinear systems with a competitive efficiency index. Moreover, higher order methods can be constructed by adding new steps with the same structure involving one new functional evaluation (per step) and increasing the order of convergence in three units per step. These methods are specially useful in nonlinear systems with a big size. Their performance has been checked by means of a heat transfer problem showing excellent results and also has been compared with recent high-order methods on academical examples.
Competing Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This research was partially supported by Ministerio de Economía y Competitividad MTM2014-52016-C2-2-P and by Generalitat Valenciana PROMETEO/2016/089.