On a Family of High-Order Iterative Methods under Kantorovich Conditions and Some Applications

Abstract

This paper is devoted to the study of a class of high-order iterative methods for nonlinear equations on Banach spaces. An analysis of the convergence under Kantorovich-type conditions is proposed. Some numerical experiments, where the analyzed methods present better behavior than some classical schemes, are presented. These applications include the approximation of some quadratic and integral equations.

1. Introduction

This paper deals with the approximation of nonlinear equations as follows 𝐹(𝑥)=0, (1.1) where 𝐹∶Ω⊆𝑋→𝑌 is a nonlinear operator between Banach spaces, using the following family of high-order iterative methods: &𝑦_𝑛=𝑥_𝑛−𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛,&𝑥_𝑛+1=𝑦_𝑛−𝐼+𝐿_𝐹𝑥_𝑛+𝐿_𝐹𝑥_𝑛²𝐺_𝐹𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑦_𝑛, (1.2) where 𝐼 is the identity operator on 𝑋 and for each 𝑥∈𝑋, 𝐿_𝐹(𝑥) is the linear operator on Ω⊆𝑋 defined by the following 𝐿_𝐹(𝑥)=𝐹^′(𝑥)⁻¹𝐹^′′(𝑥)𝐹^′(𝑥)⁻¹𝐹(𝑥), (1.3) assuming that 𝐹′(𝑥)⁻¹ exists and 𝐺_𝐹∶Ω⊆𝑋→ℒ(𝑋,𝑋) is a given nonlinear operator (usually depending on the operator 𝐹 and its derivatives), here ℒ(𝑋,𝑋) denotes the space of bounded linear operators from 𝑋 to 𝑋.

The second step can be interpreted as an acceleration of the initial one (in our case Newton’s method). Indeed, this family was introduced for scalar equations 𝑓(𝑡)=0 in [1], for any initial scheme, Traub’s theorem reads:

Theorem 1.1. For all sufficiently smooth function 𝑔_𝑓(𝑥), the following iterative method 𝑦_𝑛=Φ𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−1+𝐿_𝑓𝑥_𝑛+𝐿_𝑓𝑥_𝑛²𝑔_𝑓𝑥_𝑛𝑓𝑦_𝑛𝑓^′𝑥_𝑛, (1.4) has order of convergence min{𝑝+2,2𝑝}, where 𝑝 is the order of Φ(𝑥).

In this paper, we consider as the function Φ(𝑥) the classical Newton method. We have mainly three reasons. First, because we can recover many well-known high-order iterative methods. Second, because the domain of convergence of Newton’s method is bigger than high order schemes [2]. Finally, since in practice it is a good strategy to start with a simple method when we are not sufficiently close to the solution [3].

On the other hand, conditions are imposed on 𝑥₀ and on 𝐹 in order to assure the convergence of {𝑥_𝑛}_𝑛 to a solution 𝑥^∗ of 𝐹(𝑥)=0. This analysis, usually known as Kantorovich type, are based on a relationship between the problem in a Banach space and a single nonlinear scalar equation which leads the behavior of the problem. A priori error estimates, depending only on the initial conditions, and, hence, the order of convergence can be obtained by using Kantorovich type theorems.

A review to the amount of literature on high-order iterative methods in the two last decades (see for instance [4] and its references, or this incomplete list of recent papers [5–16]) may reveal the importance of high-order schemes. The main practical difficulty related to the classical third-order iterative methods is the evaluation of the second-order derivative. For a nonlinear system of 𝑚 equations and 𝑚 unknowns, the first Fréchet derivative is a matrix with 𝑚² entries, while the second Fréchet derivative has 𝑚³ entries. This implies a huge amount of operations in order to evaluate every iteration. However, in some cases, the second derivative is easy to evaluate. Some clear examples of this case are the approximation of Hammerstein equations where the second Fréchet derivative is diagonal by blocks or quadratic equations where it is constant.

The structure of this paper is as follows: in Section 2 we present some particular examples of methods included in the family and in Section 3, we assert convergence and uniqueness theorems (Kantorovich type). Finally, some numerical experiments are presented in Section 4. These applications include quadratic (Riccati) equations and integral (Hammerstein) equations. In all these problems the proposed methods seem more efficient than second-order methods.

2. A Family of High-Order Iterative Methods

As was indicated in the introduction, we are interested in the study of the family of iterative methods as follows 𝑦_𝑛=𝑥_𝑛−𝑓𝑥_𝑛𝑓^′𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−1+𝐿_𝑓𝑥_𝑛+𝐿_𝑓𝑥_𝑛²𝑔_𝑓𝑥_𝑛𝑓𝑦_𝑛𝑓^′𝑥_𝑛. (2.1)

Note that the method (2.1) is equivalent to iterate the function 𝑀_𝑓 given by the following 𝑀_𝑓(𝑥)𝑓(𝑥)=𝑥−𝑓^′(𝑥)−1+𝐿_𝑓(𝑥)+𝐿_𝑓(𝑥)²𝑔_𝑓(𝑥)𝑓𝑥−𝑓(𝑥)/𝑓^′(𝑥)𝑓^′(𝑥), (2.2) that is, 𝑥_𝑛+1=𝑀_𝑓𝑥_𝑛. (2.3)

Particular examples of schemes included in the family with nonsmooth functions 𝑔_𝑓(𝑥) are(1)Halley𝑥_𝑛+1=𝑥_𝑛−11+(1/2)𝐿_𝑓𝑥_𝑛𝑓𝑥_𝑛𝑓^′𝑥_𝑛, (2.4)(2)Super-Halley𝑥_𝑛+1=𝑥_𝑛−𝐿1+^𝑓𝑥_𝑛21−𝐿_𝑓𝑥_𝑛𝑓𝑥_𝑛𝑓^′𝑥_𝑛, (2.5)(3)Chebyshev𝑥_𝑛+1=𝑥_𝑛−11+2𝐿_𝑓𝑥_𝑛𝑓𝑥_𝑛𝑓^′𝑥_𝑛, (2.6)(4)Chebyshev like methods. For 0≤𝛼≤2, we consider the following 𝛼-methods 𝑥_𝑛+1=𝑥_𝑛−11+2𝐿_𝑓𝑥_𝑛+𝛼𝐿_𝑓𝑥_𝑛²𝑓𝑥_𝑛𝑓^′𝑥_𝑛, (2.7)(5)Two-step 𝑦_𝑛=𝑥_𝑛−𝑓𝑥_𝑛𝑓^′𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−𝑓𝑦_𝑛𝑓^′𝑥_𝑛. (2.8)

These methods have order of convergence three that is small than the estimate 4=min{2+2,2⋅2} in Traub’s theorem (since 𝑔_𝑓(𝑥) is nonsmooth). For instance the above two-step method admits 𝑔_𝑓(𝑥)=−1/𝐿_𝑓(𝑥). Indeed, all these methods have the function 𝑓 in the denominator.

On the other hand, considering different smooth functions 𝑔_𝑓(𝑥), the following schemes are also particular examples in the family.(1) The two-step method (𝑔_𝑓(𝑥)=0) 𝑦𝑀4∶^𝑛=𝑥_𝑛−𝑓𝑥_𝑛𝑓^′𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−1+𝐿_𝑓𝑥_𝑛𝑓𝑦_𝑛𝑓^′𝑥_𝑛, (2.9)  has order four.(2) The two-step method (𝑔_𝑓(𝑥)=(1/2)((5/2)−𝐿_𝑓′(𝑥))) 𝑦𝑀5∶^𝑛=𝑥_𝑛−𝑓𝑥_𝑛𝑓′𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−1+𝐿_𝑓𝑥_𝑛+1252−𝐿_𝑓′𝑥_𝑛𝐿_𝑓𝑥_𝑛²𝑓𝑦_𝑛𝑓′𝑥_𝑛. (2.10) has order five.(3) We should start with other iterative functions Φ(𝑥) and develop a similar analysis. For instance, starting with Chebyshev’s method we can consider the method (𝑔_𝑓(𝑥)=(1/2)(3−𝐿_𝑓′(𝑥_𝑛))) 𝑦𝑀6∶^𝑛=𝑥_𝑛−11+2𝐿_𝑓𝑥_𝑛𝑓𝑥_𝑛𝑓′𝑥_𝑛,𝑥_𝑛+1=𝑦_𝑛−1+𝐿_𝑓𝑥_𝑛+123−𝐿_𝑓′𝑥_𝑛𝐿_𝑓𝑥_𝑛²𝑓𝑦_𝑛𝑓′𝑥_𝑛 (2.11) that has order six [3]. We use this scheme only in the numerical section.

3. Semilocal Convergence

Several techniques are usually considered to study the convergence of iterative methods, as we can see in the following papers [4, 17–20]. Among these, the two most common are the based on the majorant principle and on recurrence relations.

In this section, we analyze the semilocal convergence of the introduced family (1.2) under a generalization of Kantorovich conditions.

Namely, we assume that:(C1) Let 𝑥₀∈Ω such that Γ₀=𝐹^′(𝑥₀)⁻¹ exists and ‖Γ₀‖≤𝛽.(C2)‖Γ₀𝐹(𝑥₀)‖≤𝜂.(C3)‖𝐹^′′(𝑥)‖≤𝑀 for all 𝑥∈Ω.(C4)‖𝐹^′′(𝑥)−𝐹^′′(𝑦)‖≤𝐾‖𝑥−𝑦‖, 𝐾>0, 𝑥,𝑦∈Ω.

Under these hypotheses it is possible to find a cubic polynomial in an interval [𝑎,𝑏] such that 𝑝(𝑎)>0>𝑝(𝑏), 𝑝^′(𝑡)<0, 𝑝^′′(𝑡)>0 and 𝑝^′′′(𝑡)>0 in [𝑎,𝑡^∗], with 𝑡^∗ the unique simple solution of 𝑝(𝑡)=0, and verifying the following hypotheses:

For 𝑡₀∈[𝑎,𝑏] and 𝑝(𝑡₀)>0.(H1)‖Γ₀‖≤−1/𝑝^′(𝑡₀),(H2)‖Γ₀𝐹(𝑥₀)‖≤−𝑝(𝑡₀)/𝑝^′(𝑡₀),(H3)‖𝐹^′′(𝑥)‖≤𝑝^′′(𝑡) for all 𝑥∈Ω, ‖𝑥−𝑥₀‖≤𝑡−𝑡₀≤𝑡^∗−𝑡₀,(H4)‖𝐹^′′(𝑥)−𝐹^′′(𝑦)‖≤|𝑝^′′(𝑢)−𝑝^′′(𝑣)|, with ‖𝑥−𝑦‖≤|𝑢−𝑣|, 𝑥,𝑦∈Ω and 𝑢,𝑣∈[𝑎,𝑡^∗].

Some immediate properties of the polynomial may be obtained from the conditions above imposed:(1)𝑝(𝑡) is decreasing in the interval [𝑎,𝑡^∗], since 𝑝′(𝑡)<0 in that interval.(2)𝑝(𝑡)>0 in [𝑎,𝑡^∗[.(3)𝑝′(𝑡) is increasing and 𝑝(𝑡) is convex in [𝑎,𝑡^∗], since we have 𝑝′′(𝑡)>0 in [𝑎,𝑡^∗].(4)𝑝^′′(𝑡) is increasing in [𝑎,𝑡^∗], since 𝑝^′′′(𝑡)>0 in that interval.

From these properties it follows the next:(a)The Newton map associate to 𝑝(𝑡), 𝑁_𝑝(𝑡)=𝑡−(𝑝(𝑡)/𝑝′(𝑡)), is increasing in [𝑎,𝑡^∗[,  𝑁_𝑝(𝑡^∗)=𝑡^∗ and 𝑁^′_𝑝(𝑡^∗)=0.(b)The function 𝐿_𝑝(𝑡)=𝑝(𝑡)𝑝^′′(𝑡)/𝑝^′(𝑡)²>0 in [𝑎,𝑡^∗[, since 𝑝(𝑡) and 𝑝^′′(𝑡) are strictly positive in that interval. Furthermore, 𝐿_𝑝(𝑡^∗)=0, since 𝑝(𝑡^∗)=0 and 𝑝′(𝑡^∗)≠0.

In this paper, as in [21, page 43], we consider as the function 𝑝(𝑡) the following polynomial: 𝑝(𝑡)𝐾∶=6𝑡³+𝑀2𝑡²−1𝛽𝜂𝑡+𝛽, (3.1) assuming 𝜂≤4𝐾+𝑀²𝛽−𝑀𝛽𝑀²+2𝐾𝛽3𝛽𝐾𝑀+𝑀².+2𝐾𝛽 (3.2)

If this last condition holds, then the cubic polynomial 𝑝(𝑡) has two roots 𝑡^∗ and 𝑡^∗∗ (𝑡^∗≤𝑡^∗∗). We can choose 𝑎 and 𝑏 such that 0<𝑎<𝑡^∗ and 𝑏>2/(𝑀𝛽+(𝑀²𝛽²+2𝐾𝛽)).

Moreover, we need some extra conditions associated to the operator 𝐺_𝐹 and the function 𝑔_𝑝. We assume:(Hg1)‖𝐿_𝐹(𝑥)²𝐺_𝐹(𝑥)‖≤𝐿_𝑝(𝑡)²𝐺_𝑝(𝑡), for ‖𝑥−𝑥₀‖≤𝑡−𝑡₀≤𝑡^∗−𝑡₀,(Hg2)1+𝐿_𝑝(𝑡)+𝐿_𝑝(𝑡)²𝑔_𝑝(𝑡)≥0,(Hg3)𝑚^′(𝑡)>0 in [𝑎,𝑡^∗[, where 𝑚(𝑡)𝑝(𝑡)=𝑡−𝑝^′(𝑡)−1+𝐿_𝑝(𝑡)+𝐿_𝑝(𝑡)²𝑔_𝑝(𝑡)𝑝𝑡−𝑝(𝑡)/𝑝^′(𝑡)𝑝^′(𝑡). (3.3)

All the methods considered in the above section have associated functions 𝑔_𝑝 that verify the three last conditions. With the two last hypotheses on 𝑔_𝑝 and the definition of 𝑝, following [21, Corollary 2.2.2 in page 31], the next result holds.

Proposition 3.1. The sequence 𝑠_𝑛=𝑡_𝑛−𝑝𝑡_𝑛𝑝^′𝑡_𝑛,𝑡_𝑛+1=𝑠_𝑛−1+𝐿_𝑝𝑡_𝑛+𝐿_𝑝𝑡_𝑛²𝑔_𝑝𝑡_𝑛𝑝𝑠_𝑛𝑝^′𝑡_𝑛, (3.4) starting from the above 𝑡₀ converges monotonically to 𝑡^∗ the real simple solution of 𝑝(𝑡)=0 in [𝑎,𝑏].

We are now ready to prove the desired semilocal convergence.

Theorem 3.2. Let us assume 𝑥₀∈Ω and 𝑡₀∈[𝑎,𝑡^∗] verifying the hypotheses (H1)–(H4) and (Hg1)–(Hg3) with 𝜂≤4𝐾+𝑀²𝛽−𝑀𝛽𝑀²+2𝐾𝛽3𝛽𝐾𝑀+𝑀².+2𝐾𝛽 (3.5)
If 𝐵(𝑥₀,𝑡^∗−𝑡₀)⊂Ω then the sequence (1.2) is well defined and converges to 𝑥^∗ the unique solution of 𝐹(𝑥)=0 in 𝐵(𝑥₀,𝑡^∗−𝑡₀).
Moreover, 𝑥^∗−𝑥_𝑛≤𝑡^∗−𝑡_𝑛,𝑛≥0, (3.6) where 𝑠_𝑛=𝑡_𝑛−𝑝𝑡_𝑛𝑝^′𝑡_𝑛,𝑡_𝑛+1=𝑠_𝑛−1+𝐿_𝑝𝑡_𝑛+𝐿_𝑝𝑡_𝑛²𝑔_𝑝𝑡_𝑛𝑝𝑠_𝑛𝑝^′𝑡_𝑛. (3.7)

Proof. By an induction process, it is possible to verify that(i)‖𝐹^′(𝑥_𝑛)⁻¹‖≤−1/𝑝^′(𝑡_𝑛),(ii)‖𝐹(𝑥_𝑛)‖≤𝑝(𝑡_𝑛), and then,(iii)‖𝐿_𝐹(𝑥_𝑛)‖≤𝐿_𝑝(𝑡_𝑛), and(iv)‖𝑥_𝑛+1−𝑥_𝑛‖≤𝑡_𝑛+1−𝑡_𝑛.
The case 𝑛=0 follows from the initial conditions on 𝑥₀ and 𝑡₀.
We now assume that the conditions are valid for 𝑛 and we check them for 𝑛+1.(i) 𝐹^′𝑥_𝑛+1=𝐹^′𝑥_𝑛𝐼−𝐹^′𝑥_𝑛⁻¹𝐹^′𝑥_𝑛−𝐹^′𝑥_𝑛+1. (3.8) Applying Taylor’s theorem: 𝐹^′𝑥_𝑛⁻¹𝐹^′𝑥_𝑛−𝐹^′𝑥_𝑛+1≤𝐹^′𝑥_𝑛⁻¹𝐹^′′𝑥_𝑛+12𝐾𝑥_𝑛−𝑥_𝑛+1⋅𝑥_𝑛−𝑥_𝑛+11≤−𝑝^′𝑡_𝑛𝑝^′′𝑡_𝑛+12𝐾𝑡_𝑛+1−𝑡_𝑛⋅𝑡_𝑛+1−𝑡_𝑛1=−𝑝^′𝑡_𝑛𝑝^′𝑡_𝑛+1−𝑝^′𝑡_𝑛𝑝=1−^′𝑡_𝑛+1𝑝^′𝑡_𝑛<1, (3.9)because 𝑝^′(𝑡) is increasing.  By applying the general invertibility criterion, 𝐹^′(𝑥_𝑛+1) is invertible, and 𝐹^′𝑥_𝑛+1⁻¹≤𝐼−𝐹^′𝑥_𝑛⁻¹𝐹^′𝑥_𝑛−𝐹^′𝑥_𝑛+1⁻¹𝐹^′𝑥_𝑛⁻¹≤𝐹^′𝑥_𝑛⁻¹1−𝐹^′𝑥_𝑛⁻¹𝐹^′𝑥_𝑛−𝐹^′𝑥_𝑛+11≤−𝑝^′𝑡_𝑛1−1−𝑝^′𝑡_𝑛+1/𝑝^′𝑡_𝑛1=−𝑝^′𝑡_𝑛+1. (3.10)(ii) Using the following Taylor expansion 𝐹𝑦_𝑛=𝐹𝑥_𝑛+𝐹^′𝑥_𝑛𝑦_𝑛−𝑥_𝑛+1𝐹2!^′′𝑥_𝑛𝑦_𝑛−𝑥_𝑛²+^𝑦𝑛_𝑥𝑛𝐹^′′(𝑥)−𝐹^′′𝑥_𝑛𝑦_𝑛−𝑥𝑑𝑥, (3.11)and by the definition of the method 𝐹^′𝑥_𝑛𝑦_𝑛−𝑥_𝑛=−𝐹^′𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛, (3.12) we obtain that 𝐹𝑦_𝑛=1𝐹2!^′′𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹(𝑥_𝑛)²+^𝑦𝑛_𝑥𝑛𝐹^′′(𝑥)−𝐹^′′𝑥_𝑛𝑦_𝑛−𝑥𝑑𝑥, (3.13) and since 𝑝𝑠_𝑛=1𝑝2!^′′𝑡_𝑛𝑝^′𝑡_𝑛⁻¹𝑝𝑡_𝑛²+^𝑠𝑛_𝑡𝑛𝑝^′′(𝑥)−𝑝^′′𝑡_𝑛𝑠_𝑛−𝑡𝑑𝑡, (3.14) we conclude that 𝐹𝑦_𝑛≤𝑝𝑠_𝑛. (3.15) Similarly from the following expansion 𝐹𝑥_𝑛+1=𝐹𝑥_𝑛+𝐹^′𝑥_𝑛𝑥_𝑛+1−𝑥_𝑛+1𝐹2!^′′𝑥_𝑛𝑥_𝑛+1−𝑥_𝑛²+^𝑥𝑛+1_𝑥𝑛𝐹^′′(𝑥)−𝐹^′′𝑥_𝑛𝑥_𝑛+1−𝑥𝑑𝑥, (3.16) the definition of the method, the main hypotheses on 𝐺_𝐹 and the induction process, we obtain, using that ‖𝐹(𝑦_𝑛)‖≤𝑝(𝑠_𝑛) and that 𝐹^′′𝑥_𝑛𝑥_𝑛+1−𝑥_𝑛²=𝐹^′′𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛+𝐹^′′𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛𝐼+𝐿_𝐹𝑥_𝑛+𝐿_𝐹𝑥_𝑛²𝐺_𝐹𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑦_𝑛+𝐹^′′𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑥_𝑛𝐼+𝐿_𝐹𝑥_𝑛+𝐿_𝐹𝑥_𝑛²𝐺_𝐹𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑦_𝑛+𝐹^′′𝑥_𝑛𝐼+𝐿_𝐹𝑥_𝑛+𝐿_𝐹𝑥_𝑛²𝐺_𝐹𝑥_𝑛𝐹^′𝑥_𝑛⁻¹𝐹𝑦_𝑛², (3.17) the desired inequality: 𝐹𝑥_𝑛+1≤𝑝𝑡_𝑛+1. (3.18)
In this situation, the theorem holds by applying the previous estimates directly to the formulas that describe the methods, we refer [21, page 41-42] for more details.

The estimates given in the present paper are optimal in the sense that the sequence associated to 𝑝 verifies the inequalities with equalities.

4. Numerical Experiments

We consider several problems where the presented high-order methods can be considered as a good alternative to second-order methods.

4.1. Approximation of Riccati’s Equations

In this first example, we consider quadratic equations, therefore the second Fréchet derivative is constant. Particular cases of this type of equations, that appear in many applications, are Riccati’s equations [22–24]. For instance, if we consider the problem of calculating feedback controls for systems modeled by partial differential or delay differential equations, a classical controller design objective will be to find a control 𝑢(𝑡) for the state 𝑥(𝑡) such that the following objective function ^∞₀⟨𝐶𝑥(𝑡),𝐶𝑥(𝑡)⟩+𝑢^∗𝑅𝑢(𝑡)𝑑𝑡 (4.1) is minimized, where 𝑅 is a positive defined matrix and the observation 𝐶∈ℒ(𝑋,ℝ^𝑑). In practice, the control is calculated through approximation. This leads to solving an algebraic Riccati equation 𝐴^∗𝑃+𝑃𝐴−𝑃𝐵𝑅⁻¹𝐵^∗𝑃=−𝐶^∗𝐶 (4.2) for a feedback operator 𝐾=−𝑅⁻¹𝐵^∗𝑃, (4.3) see [25, 26] for more details.

In the general case, an algebraic Riccati’s equation is given by [27] 𝑅(𝑋)=𝑋𝐷𝑋−𝑋𝐴−𝐴^𝑇𝑋−𝐶=0, (4.4) where 𝐷,𝐴,𝐶∈ℝ^𝑛×𝑛 are given matrix, 𝐷 symmetric and 𝑋∈ℝ^𝑛×𝑛 is the unknown.

In this case, 𝑅^′(𝑋)𝑌=𝑋𝐷−𝐴^𝑇(),𝑅𝑌+𝑌𝐷𝑋−𝐴_′′(𝑋)𝑌𝑍=𝑌𝐷𝑍+𝑍𝐷𝑌. (4.5)

In particular, the second derivative is constant. In this case, the Kantorovich conditions for Newton’s methods have the compact form 𝑅⁻¹𝑋₀𝑅^′′𝑋₀𝑅⁻¹𝑋₀𝑅𝑋₀≤12. (4.6) Moreover, this hypothesis also gives the convergence for the high-order methods [22].

Then, using a matricial norm 𝑅^′′(𝑋)‖𝐷𝑌𝑍‖,𝑅𝑌𝑍≤2‖‖‖‖_′′(𝑋)≤2‖𝐷‖. (4.7)

Given a symmetric initial guess 𝑋₀∈ℝ^𝑛×𝑛, to obtain 𝑅^′(𝑋₀)⁻¹ we solve the equation 𝑅^′𝑋₀𝑌=𝑋₀𝐷−𝐴^𝑇𝑌+𝑌𝐷𝑋₀−𝐴=𝑍. (4.8) This equation has solution if 𝐷𝑋₀−𝐴 is stable [27], that is, all its eigenvalues have negative real part. In this following case 𝑌=𝑅^′𝑋₀⁻¹𝑍=−^∞₀exp𝐷𝑋₀−𝐴^𝑇𝑡𝑍exp𝐷𝑋₀−𝐴𝑡𝑑𝑡. (4.9)

Next, to illustrate the previous results, we consider the following algebraic Riccati equation (4.4) with matrix 𝐷=100010001=𝐶,𝐴=001010100, (4.10) and the starting point 𝑋₀=−3250100010−4. (4.11) In this case, the algebraic Riccati equation has exact solution 𝑋^∗=−√√20101−√2010−2. (4.12) Besides, from the aforesaid starting point it follows that 𝐷𝑋₀−𝐴 is a stable matrix.

Now, considering the stopping criterion ‖𝑋_𝑛−𝑋^∗‖<10⁻⁵⁰ in Table 1, we obtain the errors ‖𝑋_𝑛−𝑋^∗‖. If we now analyze the following computational order of convergence [28]: 𝑋𝜌≈ln^𝑛+1−𝑋^∗𝑋_𝑛−𝑋^∗𝑋/ln^𝑛−𝑋^∗𝑋_𝑛−1−𝑋^∗,𝑛∈ℕ, (4.13) we observe that method M6 has computationally the order of convergence at least six. See Table 2, where 𝜌_𝑁, 𝜌_CH and 𝜌_M6 denote, respectively, the computational order of convergence of the three last methods.

Table 1 

Errors for the Newton, Chebyshev and M6 methods.

Table 2 

The computational order of convergence for the Newton, Chebyshev and M6 methods.

In comparison with the classical Newton’s method, the extra computational cost per iteration of method M6, is only two new evaluations of the operator 𝐹, and two extra matrix-vector multiplications. Moreover, the same as Newton’s method only a 𝐿𝑈 decomposition is necessary. Thus, M6 is more efficient.

See [29] for more details.

4.2. Approximation of Hammerstein Equations

We consider an important special case of integral equation, the following Hammerstein equation 𝑢(𝑠)=𝜓(𝑠)+¹₀𝐻(𝑠,𝑡)𝑓(𝑡,𝑢(𝑡))𝑑𝑡. (4.14) These equations are related with boundary value problems for differential equations. For some of them, high-order methods using second derivatives are useful for their effective (discretized) solution.

The discrete version of (4.14) is 𝑥^𝑖=𝜓𝑡_𝑖+^𝑚_𝑗=0𝛾_𝑗𝐻𝑡_𝑖,𝑡_𝑗𝑓𝑡_𝑗,𝑥^𝑗,𝑖=0,1,…,𝑚, (4.15) where 0≤𝑡₀<𝑡₁<⋯<𝑡_𝑚≤1 are the grid points of some quadrature formula ∫¹₀𝑓(𝑡)𝑑𝑡≈∑^𝑚_𝑗=0𝛾_𝑗𝑓(𝑡_𝑗), and 𝑥^𝑖=𝑥(𝑡_𝑖).

The second Fréchet derivative of the associated discrete system is diagonal by blocks.

Let the following Hammerstein equation 𝑥(𝑠)1=1−4¹₀𝑠1𝑡+𝑠𝑥(𝑡)𝑑𝑡,𝑠∈[0,1]. (4.16)

The discretization of this equation verifies the Lipschitz condition of our Kantorovich theorem [4].

We consider 𝑚=20 in the quadrature trapezoidal formula and as exact solution the obtained numerically by Newton method. In Table 3, we summarize the numerical results for different methods in the family: Newton, Halley, and M4. We consider as initial guess 𝑥₀(𝑠)=1.5.

Table 3 

Errors for the Newton, Halley, and M4 methods.

Since the second derivative is diagonal by blocks, its application has a computational cost of order 𝑂(𝑚²). Thus, the computational cost in each iteration of the three schemes is, for 𝑚 sufficiently big, of the same order (𝑂(𝑚³) due to the 𝐿𝑈 decomposition). Note that we only have to do a factorization in each iteration of the three schemes. As conclusion, the scheme M4 (order four) is the most efficient for 𝑚 sufficiently big.

See [30] for other-related problems.

5. Conclusions

Summing up, in this paper we have studied a family of high-order iterative methods. Mainly, the theoretical analysis we did allows to ensure convergence conditions for all these schemes. We established priori error bounds for them and consequently their order. We have presented different applications where we may add that in these cases the analyzed high-order methods are more efficient than simpler second-order methods.

Acknowledgment

S. Amat, C. Bermúdez, and S. Busquier were supported in part by MINECO-FEDER MTM2010-17508 and 08662/PI/08. S. Plaza was supported in part by Fondecyt (Grant no. 1095025).