# A Cauchy Point Direction Trust Region Algorithm for Nonlinear Equations

Authors: Zhou Sheng; Dan Luo
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

In this paper, a Cauchy point direction trust region algorithm is presented to solve nonlinear equations. The search direction is an optimal convex combination of the trust region direction and the Cauchy point direction with the sufficiently descent property and the automatic trust region property. The global convergence of the proposed algorithm is proven under some conditions. The preliminary numerical results demonstrate that the proposed algorithm is promising and has better convergence behaviors than the other two existing algorithms for solving large-scale nonlinear equations.


---

## Body

## 1. Introduction

The system of nonlinear equations is one of the most important mathematical models, with extensive applications in chemical equilibrium problems [1], problem of hydraulic circuits of power plants [2], L-shaped beam structures [3], and image restoration problems [4]. Some ℓ1-norm regularized optimization problems in compressive sensing [5] are also obtained by reformulating the systems of nonlinear equations. In this paper, consider the following nonlinear system of equation:(1)Fx=0,where F:ℜn⟶ℜn is continuously differentiable. Let(2)fx=12Fx2,where f:ℜn⟶ℜ. It is clear that the nonlinear system of equations (1) is equivalent to the unconstrained optimization problem as follows:(3)minx∈ℜnfx.Many literature studies have studied a variety of numerical algorithms to solve the nonlinear system of equations, such as the conjugate gradient algorithm [6–9], Levenberg–Marquardt algorithm [10–13], and trust region algorithm [14–21]. The trust region method plays an important role in the area of nonlinear optimization. It was proposed by Powell in [22] to solve the nonlinear optimization problem by using an iterative structure. Recently, some trust region-based methods have been developed to solve nonlinear equations and have shown very promising results. For instance, Fan and Pan [16] discussed a regularized trust region subproblem with an adaptive trust region radius, which possesses quadratic convergence under the local error bound condition. Yuan et al. [17] proposed a trust region-based algorithm, in which the Jacobian is updated by the BFGS formula, for solving nonlinear equations. The algorithm showed global convergence at a quadratic convergence rate. This algorithm does not compute the Jacobian matrix in each interrogation, which significantly reduces the computation burden. Esmaeili and Kimiaei [19] studied an adaptive trust region algorithm with the moderate radius size by employing the nonmonotone technique, and the global and quadratic local convergence was given. However, the maximal dimension of the test problems was only 300. Qi et al. [23] proposed an asymptotic Newton search direction method which is an adequate combination of the projected gradient direction and the projected trust region direction for solving box-constrained nonsmooth equations and proved the global and quadratic local convergence. Motivated by the search direction by Qi et al., Yuan et al. [24] used the idea of the search direction to solve nonsmooth unconstrained optimization problems and also showed global convergence.Inspired by the fact that the search direction by Qi et al. [23] speeds up the local convergence rate compared to the gradient direction, we propose an efficient Cauchy point direction trust region algorithm to solve (1) in this paper. The search direction is an optimal convex combination of the trust region direction and the Cauchy point direction with the automatic trust region property and the sufficiently descent property, which can improve the numerical performance. We also show the global convergence of the proposed algorithm. Numerical results indicate that this algorithm has better behavior convergence compared with that of the classical trust region algorithm and the adaptive regularized trust region algorithm. However, the local convergence is not given.The remainder of this paper is organized as follows. In Section2, we introduce our motivation and a new search direction. In Section 3, we give a Cauchy point direction trust region algorithm for solving problem (1). The global convergence of the proposed algorithm is proven in Section 4. Preliminary numerical results are reported in Section 5. We conclude our paper in Section 6.Throughout this paper, unless otherwise specified,⋅ denotes the Euclidean norm of vectors or matrices.


## 2. Motivation and Search Direction

In this section, motivated by the definition of the Cauchy point which is parallel to the steepest descent direction, we present a search direction with the Cauchy point. LetΔk be the trust region radius and Δk>0.


### 2.1. Cauchy Point Direction

Noting thatJkTJk is a semidefinite matrix, we define the Cauchy point at the current point xk:(4)dkCΔk≔−ckΔkJkTFkJkTFk,where(5)ck=1,ifFkTJkJkT2Fk=0,min1,JkTFk3ΔkFkTJkJkT2Fk,ifFkTJkJkT2Fk>0,Fk=Fxk, and Jk=∇Fxk is the Jacobian of Fx or its approximation.The following lemma gives a nice property of the Cauchy point directiondkC.Proposition 1 (automatic trust region property).
Suppose thatdkC is defined by (4), then dkC≤Δk.Proof.
It is easy to get the result of this proposition from0<ck≤1. This completes the proof.
From the definition ofdkC and the above proposition, we can determine the descent property and the automatic trust region property of the Cauchy point direction dkC.


### 2.2. Trust Region Direction

Let the trial stepdkTRΔk be a solution of the following trust region subproblem:(6)mind∈ℜnmkd≔12Fk+Jkd2,s.t.d≤Δk,where Δk is the same as the above definition.From the famous result of Powell [22], we give the following lemma and omit its proof.Lemma 1.
LetdkTRΔk be a solution of the trust region subproblem (6), then(7)fk−mkdkTRΔk=12Fk2−12Fk+JkdkTRΔk2≥12JkTFkminΔk,JkTFkJkTJk.As mentioned by Qi et al. [23] and Yuan et al. [24], when xk is far from the optimization solution, the trust region direction dkTRΔk may not possess the descent property. However, it follows from the definition of the Cauchy point direction that dkCΔk always possesses the descent property. On the contrary, if the search direction is defined by using a combination of dkC and dkTR, it will be a descent direction. Hence, we give a new search direction by using a convex combination form in the next section.


### 2.3. Search Direction

Letλk be a solution of the one-dimensional minimization problem as follows:(8)minλ∈0,1hkλ≔12Fk+JkλdkCΔk+1−λdkTRΔk2,and the search direction is defined as(9)dkΔk=λkdkCΔk+1−λkdkTRΔk.For the above one-dimensional minimization problem, we can use the gold section method, Fibonacci method, parabolic method [25], etc. The idea of the search direction follows from Qi et al. [23]. However, we use the Cauchy point direction which possesses the automatic trust region property to get optimal λk by using the gold section method for problem (6) in this paper.


## 2.1. Cauchy Point Direction

Noting thatJkTJk is a semidefinite matrix, we define the Cauchy point at the current point xk:(4)dkCΔk≔−ckΔkJkTFkJkTFk,where(5)ck=1,ifFkTJkJkT2Fk=0,min1,JkTFk3ΔkFkTJkJkT2Fk,ifFkTJkJkT2Fk>0,Fk=Fxk, and Jk=∇Fxk is the Jacobian of Fx or its approximation.The following lemma gives a nice property of the Cauchy point directiondkC.Proposition 1 (automatic trust region property).
Suppose thatdkC is defined by (4), then dkC≤Δk.Proof.
It is easy to get the result of this proposition from0<ck≤1. This completes the proof.
From the definition ofdkC and the above proposition, we can determine the descent property and the automatic trust region property of the Cauchy point direction dkC.


## 2.2. Trust Region Direction

Let the trial stepdkTRΔk be a solution of the following trust region subproblem:(6)mind∈ℜnmkd≔12Fk+Jkd2,s.t.d≤Δk,where Δk is the same as the above definition.From the famous result of Powell [22], we give the following lemma and omit its proof.Lemma 1.
LetdkTRΔk be a solution of the trust region subproblem (6), then(7)fk−mkdkTRΔk=12Fk2−12Fk+JkdkTRΔk2≥12JkTFkminΔk,JkTFkJkTJk.As mentioned by Qi et al. [23] and Yuan et al. [24], when xk is far from the optimization solution, the trust region direction dkTRΔk may not possess the descent property. However, it follows from the definition of the Cauchy point direction that dkCΔk always possesses the descent property. On the contrary, if the search direction is defined by using a combination of dkC and dkTR, it will be a descent direction. Hence, we give a new search direction by using a convex combination form in the next section.


## 2.3. Search Direction

Letλk be a solution of the one-dimensional minimization problem as follows:(8)minλ∈0,1hkλ≔12Fk+JkλdkCΔk+1−λdkTRΔk2,and the search direction is defined as(9)dkΔk=λkdkCΔk+1−λkdkTRΔk.For the above one-dimensional minimization problem, we can use the gold section method, Fibonacci method, parabolic method [25], etc. The idea of the search direction follows from Qi et al. [23]. However, we use the Cauchy point direction which possesses the automatic trust region property to get optimal λk by using the gold section method for problem (6) in this paper.


## 3. New Algorithm

In this section, letdkΔk be the optimal solution of (3). Then, the actual reduction is defined as(10)Aredk≔fk−fxk+dkΔk,and the predicted reduction is(11)Predk≔mk0−mkdkΔk=12Fk2−12Fk+JkdkΔk2.Thus,rk is defined by the ratio between Aredk and Predk:(12)rk≔AredkPredk.We now list the detailed steps of the Cauchy point direction trust region algorithm in Algorithm1.Algorithm 1.
Cauchy point direction trust region algorithm for solving (3).Step 1. Given initial point x0∈ℜn and parameters 0≤μ1<μ2<1, 0<η1<1<η2, 0<δ<1, 0≤ϵ≪1, and Δ0>0. Let k≔0.Step 2. If the termination condition JkTFk≤ϵ is satisfied at the iteration point xk, then stop. Otherwise, go to Step 3.Step 3. Solve the trust region subproblem (6) to obtain dkTRΔk and compute ck and dkCΔk from (4) and (5), respectively.Step 4. To solve the one-dimensional minimization problem (8) giving λk by using Algorithm 2, calculate the search direction from (6).Step 5. Compute Aredk,Predk,andrk, if(13)fk−12Fk+Jkdk2≥−δFkJkTdkCΔk,rk≥μ1,are satisfied; letxk+1≔xk+dkΔk, and update Δk by(14)Δk+1≔η1Δk,ifrk<μ1,Δk,ifμ1≤rk<μ2,η2Δk,ifrk≥μ2.Step 6. If rk<μ1, return to Step 3, and let xk+1≔xk. Otherwise, return to Step 2. Let k≔k+1.
In Step 4 of Algorithm1, we use the gold section algorithm as given in Algorithm 2.Algorithm 2.
Gold section algorithm [25] for solving problem (8).Step 1. For the given initial parameters a0=0,b0=1,ϵ′>0, compute(15)p0≔a0+0.382b0−a0,q0≔a0+0.618b0−a0,andhkp0,hkq0. Let k≔0.Step 2. If hkpk≤hkqk, go to Step 3. Otherwise, go to Step 4.Step 3. If the termination condition qk−ak≤ϵ′ is satisfied, stop. Otherwise, let(16)ak+1≔ak,bk+1≔qk,hkqk+1≔hkpk,qk+1≔pk,pk+1≔ak+1+0.382bk+1−ak+1,and computehkpk+1 and go to Step 5.Step 4. If the termination condition bk−pk≤ϵ′ is satisfied, stop. Otherwise, let(17)ak+1≔pk,bk+1≔bk,hkpk+1≔hkqk,pk+1≔qk,qk+1≔ak+1+0.618bk+1−ak+1,and computehkqk+1 and go to Step 5.Step 5. Let k≔k+1; return to Step 2.


## 4. Global Convergence

In this section, the global convergence of Algorithm1 is proven. Firstly, it is easy to have the following lemma from Lemma 1.Lemma 2.
LetdkΔk be defined by (9), where λk is a solution of problem (8). Then,(18)Predk≥12JkTFkminΔk,JkTFkJkTJk.Proof.
From the definition ofdkΔk, we have mkdkΔk≤mkdkTRΔk. It follows from the definition of Predk and Lemma 1 that Predk≥fk−mkdkTRΔk holds. Thus, (16) holds, and therefore, the proof is completed.
The next lemma indicates an important property of the Cauchy point directiondkC, which is the sufficiently descent property.Lemma 3.
JkTFkTdkCΔk=−ckΔkJkTFk holds, for any Δk>0.Proof.
For anyΔk, from (9), we have(19)JkTFkTdkCΔk=FkTJkdkCΔk=−ckΔkFkTJkJkTFkJkTFk=−ckΔkJkTFk,and note that 0<ck≤1; thus,(20)JkTFkTdkCΔk=−ckΔkJkTFk<0.
This completes the proof.
To derive the global convergence of Algorithm1, the following assumption is considered.Assumption 1.
The JacobianJx of Fx is bounded; i.e., there exists a positive constant M>0 such that Jx≤M, for any x∈ℜn.
In the sequel, we show that “Steps 3–6” of Algorithm1 do not cycle in the inner loop infinitely.Lemma 4.
Suppose that Assumption1 holds, let xk be generated by Algorithm 1. If xk is not a stationary point of problem (3), then Algorithm 1 is well defined; i.e., the inner cycle of Algorithm 1 can be left in a finite number of internal iterates.Proof.
It follows fromxk which is not a stationary point that there exists a positive constant ϵ0>0 such that(21)JkTFk≥ϵ0,∀k.
Since0<ck≤1 and Assumption 1 holds, there exists a positive constant α>0 such that(22)JkdkC≤ckΔkJk≤αΔk.
Letting0<Δk≤2/α21−δϵ0ck, then from (4) and (22), we obtain(23)fk−12Fk+JkdkΔk2≥fk−12Fk+JkdkCΔk2=−FkTJkdkCΔk−12JkdkCΔk2=−δFkTJkdkCΔk−1−δFkTJkdkCΔk−12JkdkCΔk2≥−δFkTJkdkCΔk+1−δckΔkJkTFk−12α2Δk2,and noting (21), we have 1−δckΔkJkTFk−1/2α2Δk2≥1−δckΔk−1/2α2Δk2≥0, which implies(24)fk−12Fk+JkdkΔk2≥−δFkTJkdkCΔk,and therefore, (13) holds. Next, we prove that (14) is also satisfied. By contradiction, we assume that “Step 3–Step 6” of Algorithm 1 do cycle in the inner loop infinitely; then, Δk⟶0.
From the definition ofdkCΔk, we get(25)12Fk+JkdkCΔk2=12Fk2+FkTJkdkCΔk+12JkdkCΔk2≤12Fk2−ckΔkJkTFk+12α2Δk2.
It follows from (13) and Lemma 3 that there exists a positive constant β>0 such that(26)fk−12Fk+JkdkΔk2≥−δFkTJkdkCΔk=δckΔkJkTFk≥βΔk.
By the definitiondkTRΔk, we have dkTR≤Δk, and noting Proposition 1 and 0<ck≤1, further, we get(27)dkΔk≤Δk.
Thus, from (26) and (27), we have(28)rk−1=Aredk−PredkPredk=1/2Fk+JkdkΔk2−1/2Fk+JkdkΔk+OdkΔk22Predk=Fk+JkdkΔkOdkΔk2+OdkΔk4fk−1/2Fk+Jkdk2≤OdkΔk2βΔk≤OΔk2βΔk⟶0,and therefore rk⟶1, which means that rk≥μ2 when Δk is sufficiently small. However, the definition of Δk in Step 5 of Algorithm 1 contradicts the fact that Δk⟶0. Hence, the cycling between Step 3 and Step 6 terminates finitely. This completes the proof.
We can get the global convergence of Algorithm1 based on the above lemmas.Theorem 1.
Suppose that Assumption1 holds and xk be generated by Algorithm 1, then(29)liminfk⟶∞JkTFk=0.Proof.
Suppose that the conclusion is not true, then there exists a positive constantϵ1>0 such that for any k,(30)JkTFk≥ϵ1.
From Lemma2, Assumption 1, and the above inequation, we have(31)Predk≥12JkTFkminΔk,JkTFkJkTJk≥12ϵ1minΔk,ϵ1M2.
Noting the updating rule ofΔk, rk≥μ1, (31), and Lemma 4, we deduce that(32)+∞≥∑k=1∞12Fk2−12Fk+12,≥∑k=1∞μ1Predk,≥∑k=1∞12μ1ϵ1minΔk,ϵ1M2,≥∑k=1∞12μ1ϵ1minΔ,ϵ1M2,=+∞,which yields a contradiction that this means (30) is not true; i.e., the result (29) of the theorem holds, and therefore, the proof is completed.


## 5. Numerical Results

In this section, the numerical results of the proposed algorithm for solving problem (1) are reported. We denote the proposed Algorithm 1 by CTR and compare it with some existing algorithms. The adaptive trust region algorithm by Fan and Pan [16] and the classical trust region algorithm [26] are denoted by ATR and TTR, respectively. In numerical experiments, all parameters of Algorithm 1 are chosen as follows: μ1=0.1, μ2=0.9, η1=0.25, η2=3, ϵ=10−5, and c0=Δ0=1δ=0.9. We choose ϵ′=10−6 in Algorithm 2. Steihaug method [27] is employed to obtain the approximative trial step dk for solving the trust region subproblem (6). We will terminate algorithms if the number of iterations is larger than 2000.All codes were implemented in MATLAB R2010b. The numerical experiments were performed on a PC with an Intel Pentium (R) Dual-Core CPU at 3.20 GHz and 2.00 GB of RAM and using the Windows 7 operating system.We list the test problems in Table1 that can also be found in [28]. In Tables 1–3, “No.” denotes the number of the test problem, “Problem name” denotes the name of the test problem, “x0” denotes the initial point, “n” denotes the dimension of the problem, “NI” is the total number of iterations, “NF” is the number of the function evaluations, “Fk” refers to the optimum function values, and “λmean” denotes the mean value of λ. It is worth mentioning that “∗” represents the stopping of the algorithms in situations that the number of iterations exceeds the maximal iterations, but the terminal conditions cannot be satisfied yet.Table 1
Test problems.
No.Problem namex01Exponential function 1n/n−1,…,n/n−12Exponential function 21/n2,…,1/n23Extended Rosenbrock function5,1,…,5,14Trigonometric function101/100n,…,101/100n5Singular function1,…,16Logarithmic function1,…,17Broyden tridiagonal function−1,…,−18Trigexp function0,…,09Function 15−1,…,−110Strictly convex function 11/n,…,1/n11Strictly convex function 21,…,112Zero Jacobian function100n−100/n,n−1000n−500/360n2,…,n−1000n−500/360n213Linear function-full rank100,…,10014Penalty function1/3,…,1/315Brown almost-linear function1/2,…,1/216Extended Powell singular function1.5×10−4,…,1.5×10−417Tridiagonal system12,…,1218Extended Freudenstein and Roth function6,3,…,6,319Discrete boundary value problem−n/n+12,…,i−n−1/n+12,…,−1/n+1220Troesch problem1/2,…,1/2Table 2
Numerical results.
No.nCTRATRTTRNI/NFFkλmeanNI/NFFkNI/NFFk15006/73.18E − 063.32E − 077/86.73E − 066/73.18E − 0610005/68.95E − 063.32E − 077/84.74E − 065/68.95E − 0620005/66.32E − 063.32E − 077/83.34E − 065/66.32E − 0625002/33.50E − 061.57E − 062/33.38E − 062/33.38E − 0610002/36.67E − 063.32E − 072/36.69E − 062/36.69E − 0620003/49.13E − 063.32E − 073/88.39E − 063/58.44E − 06350012/162.16E − 073.31E − 0112/203.89E − 0616/219.44E − 10100011/155.92E − 072.96E − 0112/205.51E − 0616/211.39E − 09200011/151.78E − 073.50E − 0112/207.82E − 0616/211.84E − 0941009/117.33E − 065.28E − 012001/22933.03E−03∗2001/22703.48E−03∗550022/254.39E − 062.96E − 0118/195.84E − 0619/208.32E − 06100025/276.37E − 064.35E − 0123/585.68E − 0623/244.33E − 06200024/285.24E − 064.50E − 0119/205.53E − 0621/235.23E − 0665005/62.12E − 069.81E − 045/64.91E − 087/82.74E − 0810005/61.44E − 065.05E − 015/61.91E − 077/81.58E − 0720006/77.69E − 076.08E − 015/67.56E − 087/85.08E − 0875006/72.57E − 091.59E − 016/71.51E − 096/74.28E − 1110006/72.63E − 061.63E − 016/71.03E − 086/73.83E − 0620006/75.50E − 062.26E − 016/72.60E − 086/75.43E − 06810014/164.31E − 082.92E − 012001/22813.76E+00∗2001/22683.76E+00∗9100014/185.42E − 061.46E − 0214/202.76E − 0716/221.36E − 09105006/75.54E − 084.14E − 016/71.33E − 076/71.09E − 0710006/76.04E − 083.71E − 016/72.49E − 076/79.91E − 0720006/71.16E − 074.47E − 016/74.33E − 077/82.32E − 07111008/93.22E − 072.21E − 018/93.38E − 088/91.15E − 06125005/63.62E − 128.36E − 0116/175.81E − 0618/195.96E − 0610005/64.15E − 118.82E − 0116/177.38E − 0618/199.31E − 0620006/83.42E − 078.29E − 0116/178.23E − 0619/202.82E − 06135009/102.41E − 077.19E − 014/57.67E − 069/101.04E − 0810009/103.93E − 077.62E − 014/58.32E − 069/103.37E − 08200010/117.01E − 067.30E − 015/61.00E − 1010/111.35E − 07145004/56.96E − 061.09E − 027/105.22E − 065/65.34E − 06100012/131.79E − 071.77E − 018/102.39E − 066/82.65E − 06200010/115.39E − 061.15E − 0111/161.70E − 0610/115.32E − 07155003/42.85E − 093.61E − 013/42.48E − 064/52.49E − 0610004/51.32E − 075.56E − 013/46.34E − 076/78.80E − 0820004/51.40E − 085.04E − 014/53.74E − 076/73.63E − 07165001/21.10E − 063.32E − 071/21.10E − 061/21.10E − 0610001/21.55E − 063.32E − 071/21.55E − 061/21.55E − 0620001/22.20E − 063.32E − 071/22.20E − 061/22.20E − 0617100034/449.22E − 093.58E − 0190/2196.85E − 0963/754.89E − 10Table 3
Numerical results (continued).
No.nCTRATRTTRNI/NFFkλmeanNI/NFFkNI/NFFk185006/83.57E − 074.13E − 017/85.84E − 087/86.12E − 0910008/101.14E − 063.17E − 017/87.87E − 088/97.10E − 0820008/103.08E − 063.17E − 017/81.07E − 078/92.29E − 08195002/39.77E − 093.32E − 073/46.77E − 092/39.77E − 0910002/39.10E − 103.32E − 073/46.17E − 102/39.09E − 1020002/34.47E − 103.32E − 073/43.02E − 102/34.47E − 102010004/51.90E − 073.32E − 078/98.94E − 084/51.87E − 07It is shown from Tables2 and 3 that the proposed algorithm is efficient to solve these test problems and is competitive with the other two algorithms in NI and NF for most problems. For problems 1, 2, 16, 19, and 20, λmean is close to zero, which implies that the trust region direction was used many times. Moreover, the number of iteration is the same as that of algorithm TTR, which means that our algorithm speeds up the local convergence of the iterates. It can be seen that λmean is not too large for other test problems. However, for problems 4 and 8, algorithms ATR and TTR are not efficient to get the optimal solution when NI>2000. We used the performance profiles by Dolan and Moré [29] to compare these algorithms with respect to NI and NF again, as shown in Figures 1 and 2. Therefore, it is noted from Tables 2 and 3 and Figures 1 and 2 that the proposed algorithm is well promising and has better performance in test problems in comparison with two existing algorithms. In order to investigate the behavior of the convergence rate of three algorithms, we plot the curve 0−k−Lgxk−x∗ on problem 10 with n=1000 and n=2000 in Figure 3, where x∗ is the optimal point. Figure 3 indicates that Algorithm 1 converges faster than others.Figure 1
Performance profiles of the algorithms for the test problems (NI).Figure 2
Performance profiles of the algorithms for the test problems (NF).Figure 3
The behavior of the convergence: (a)n = 1000 and (b) n = 2000.
(a)(b)


## 6. Conclusions

In this paper, we presented a Cauchy point direction trust region algorithm for solving the nonlinear system of equations. This algorithm combines the trust region direction and the Cauchy point direction with the descent property and the automatic trust region. The optimal convex combination parameters were determined by using the gold section algorithm which is an exact one-dimensional linear search algorithm. The new algorithm was proven to have the descent property and the trust region property. We also showed the global convergence of the proposed algorithm under suitable condition. More importantly, the numerical results demonstrated the faster convergence of the proposed algorithm over two existing methods, which is very promising to solve large-scale nonlinear equations. However, the convergence rate of the proposed algorithm remains unclear. Therefore, the demonstration of the theoretical convergence rate of the proposed algorithm, application in image restoration problems, and compressive sensing will be the topics among others in our future study.



---

*Source: 4072651-2020-08-24.xml*