# An Enhanced Matrix-Free Secant Method via Predictor-Corrector Modified Line Search Strategies for Solving Systems of Nonlinear Equations

Authors: M. Y. Waziri; Z. A. Majid
Journal: International Journal of Mathematics and Mathematical Sciences
 (2013)

---

## Abstract

Diagonal updating scheme is among the cheapest Newton-like methods for solving system of nonlinear equations. Nevertheless, the method has some shortcomings. In this paper, we proposed an improved matrix-free secant updating scheme via line search strategies, by using the steps of backtracking in the Armijo-type line search as a step length predictor and Wolfe-Like condition as corrector. Our approach aims at improving the overall performance of diagonal secant updating scheme. Under mild assumptions, the global convergence results have been presented. Numerical experiments verify that the proposed approach is very promising.


---

## Body

## 1. Introduction

Consider the problem(1)F(x)=0,
with F:Rn→Rn.The mappingF is assumed to satisfy the following standard assumptions:(1)
Fis continuously differentiable in an open convex set Φ;(2)
there exists a solutionx* of (1) in Φ such that F(x*)=0 and F′(x*)≠0;(3)
the JacobianF′(x) is local Lipschitz continuous at x*.The famous scheme for finding the solution to (1) is the Newton method. The method is simple to implement, and it produces an iterative sequence {xk} from any given initial guess x0 in the neighborhood of x* via
(2)xk+1=xk-(F′(xk))-1F(xk),
where k=0,1,2…n.The attractive features of this method are that it is easy to implement and converges rapidly [1]. However, the major difficulty of Newton-type method is matrix storage requirements especially when handling large systems of nonlinear equations [1]. To overcome such difficulty, the simple modification on the Newton method is the fixed Newton method. Fixed Newton method for the determination of solution x* is given by
(3)xk+1=xk-(F′(x0))-1F(xk),k=0,1,2,….
This method avoids computation and storing the Jacobian in each iterations (except at k=0). However, it still requires to solve systems of n linear equations and consumes more CPU time as the system’s dimension increases [2, 3].Quasi-Newton method is another variant of Newton-type methods. It replaces the Jacobian or its inverse with an approximation which can be updated at each iteration and is given as(4)xk+1=xk-Bk-1F(xk),
where the matrix Bk is the approximation of the Jacobian at xk. The main idea behind quasi-Newton’s method is to eliminate the evaluation cost of the Jacobian matrix [4, 5]. However, the most critical part of such quasi-Newton method is storing the full matrix of the approximate Jacobian, which can be a very expensive task as the dimension of systems increases [6]. In contrast, this paper presents an improved variant of quasi-Newton update via the steps of backtracking in the Armijo-type line search and Wolfe–like condition and incorporating restarting strategy whenever the updating matrix is singular or nearly singular. The anticipation has been to improve the overall performance of the diagonal quasi-Newton method.We organized the rest of this paper as follows. In the next section, we present the proposed method. Section3 presents convergence results. Numerical experiments are reported in Section 4, and finally, conclusion is given in Section 5.


## 2. The Improved Secant Diagonal Updating

It is well known that it is not always feasible to compute the full elements of the Jacobian matrix of the given nonlinear function or it may be very expensive. We often have to approximate the Jacobian matrix by some other approach, and the famous method of doing so is quasi-Newton's method [4]. The basic idea underlining this approach has been to reduce the evaluation cost of the Jacobian matrix. This new approach generates a sequence of points {xk} via
(5)xk+1=xk-αkBkF(xk),
where αk is a step length, and Bk is a diagonal approximation of the inverse Jacobian matrix which can be updated in each iterations. To achieve this, we incorporate some new line search strategies, via the steps of backtracking in the Armijo-type line search as predictor and then improved via Wolfe-like condition as a corrector. The Armijo rule is among the inexact line search methods which guarantees a sufficient degree of accuracy to ensure the algorithm convergence. Nevertheless, the scheme needs high floating points operations and function call. In this paper, we present a simple line search strategy, which is less computational cost, floating points operations and CPU time consumptions respectively compared to classical Armijor line search. The approach is as follows.Givenϵ∈(0,1) and σ<1, the proposed approach finds the appropriate α, such that
(6)∥F(xk+αkdk)∥≤σ∥F(xk)∥.
In addition, the new strategy is implemented in an iterative way using a fixed initial value of α as follows.Algorithm 1 (Armijo-Like).Step  0. Set K=0, α0>0, and σ<1.Step 1.
∥
F
(
x
k
+
1
)
∥
≤
σ
∥
F
(
x
k
)
∥. Choose αk as the step size; stop. Otherwise αk+1=αk/2, k=k+1.Waziri et al. [7] have set the step length (αk=1, for all  k); this approach is mostly used in many Newton-like methods. Here, we continue in the spirit of diagonal updating, using a new line search strategy to obtain a good step length (αk) in every iteration, anticipating to produce a more accurate approximation of the Jacobian inverse matrix and then employing restating strategy whenever the updating matrix is undefined. To this end, Bk+1 would be obtained almost similar to the diagonal updating scheme presented in [7] in which sk=-αkBkF(xk) instead of sk=-BkF(xk). Now, let the deviation between Bk+1 and Bk denoted as Φk=Bk+1-Bk be minimized under some norms; the optimal solution is given as
(7)Φk=(ykTsk-ykTBkyk)tr(Ψk2)Ψk,
where Ψk=diag((yk(1))2,(yk(2))2,…,(yk(n))2), ∑i=1n‍(yk(i))4=tr(Ψk2), yk=F(xk+1)-F(xk), sk=xk+1-xk, and tr is the trace operation, respectively. The updated formula for the proposed matrix Bk is given as [7]
(8)Bk+1=Bk+(ykTsk-ykTBkyk)tr(Ψk2)Ψk.
To safeguard on the possibilities of generating undefined Bk+1, we let  Bk+1=B0, whenever this situation happens
(9)Bk+1={Bk+(ykTsk-ykTBkyk)tr(Ψk2)Ψk,tr(Ψk2)≠0,B0,otherwise.Now, we can describe the algorithm for our proposed approach as follows.Algorithm 2 (EMFM).Step 1.
Choose an initial guessx0, σ∈(0,1), B0=In, α0>0 and let k∶=0.Step 2.
ComputeF(xk), and if ∥F(xk)∥≤10-4, stop.Step 3.
Computed=-F(xk)Bk.Step 4.
If∥F(xk+αkdk)∥≤σ∥F(xk)∥, retain αk and go to 5. Otherwise, set αk+1=αk/2 and repeat 4.Step 5.
If∥F(xk+αkdk)-F(xk)∥≥∥F(xk+αkdk)∥-∥F(xk)∥, retain αk and go to 6. Otherwise, set αk+1=αk×1.1 and repeat 5.Step 6.
Letxk+1=xk+αkdk.Step 7.
If∥xk+1-xk∥2+∥F(xk)∥2≤10-4, stop. Otherwise go to Step 8.Step 8.
If∥ΔFk∥2≥ϵ1, where ϵ1=10-4, compute Bk+1, and if not, Bk+1=B0.Step 9.
Setk∶=k+1 and go to 2.


## 3. Convergence Analysis

We present the convergence result of EMFM method by proving the existence of the step lengthαk>0. We will make the following assumptions on nonlinear system F.Assumption 3.
(i)F is differentiable in an open convex set E in ℜn.(ii)
There existsx*∈E such that (x*)=0, and F′(x) is continuous for all x.(iii)
F(x) satisfies Lipschitz condition of order one; that is, there exists a positive constant μ such that
(10)∥F(x)-F(y)∥≤μ∥x-y∥,
for all x,y∈ℜn.(iv)
There exist constantsc1≤c2 such that c1∥ω∥2≤ωTF′(x)ω≤c2∥ω∥2 for all x∈E and ω∈ℜn.To this end, we proceed by given the following result for the step length generated by the proposed strategies.Theorem 4.
Assume thatF is a strictly convex function. Suppose that the new strategies are employed with dk≠0 and positive αk exist, for all sufficiently large k. Then, the iterates {xk} generated by the line search algorithm have the property that
(11)lim∥F(xk+αkdk)∥=0,ask⟶∞.Proof.
From condition (iii) of Assumption3 and the fact that sk=xk+1-xk, we have
(12)∥F(xk)-F(xk+αkdk)∥F≤μ∥sk∥F=μ∥αdk∥F≤μ|α|∥dk∥F.
Recall that
(13)∥F(x)-F(xk+αkdk)∥F≤∥F(x)∥F-∥F(xk+αkdk)∥F.
Equation (13) gives
(14)∥F(xk)∥F-∥F(xk+αkdk)∥F≤αμ∥dk∥F=μ∥αdk∥F≤μ|α|∥dk∥F.
Since ∥dk∥≠0 and ∥F(xk+αkdk)∥≤σ∥F(xk)∥, hence, it follows that
(15)∥F(xk)∥F-σ∥F(xk)∥F≤αμ∥dk∥F,(1-σ)∥F(xk)∥F≤αμ∥dk∥F.
After little simplifications, we obtain
(16)α≥(1-σ)∥F(xk)∥Fμ∥dk∥F>0.
Hence, αk exist and are positive.
We continue to show the convergence of the iterates{xk} by recalling that
(17)∥F(xk+αkdk)∥≤σ∥F(xk)∥.
Then, we have
(18)∥F(xk+αkdk)∥≤σ2∥F(xk)∥≤σ3∥F(xk)∥
for a finite k, and (18) yields
(19)∥F(xk+αkdk)∥≤σk∥F(xk)∥.
Since the algorithm terminates at ∥F(xk)∥=0, it follows from (19) that
(20)∥F(xk+αkdk)∥∥F(xk)∥≤σk.
We have
(21)lim∥F(xk+αkdk)∥∥F(xk)∥≤limσk,ask⟶∞.
Therefore, due to σ<1, it implies that
(22)limσk⟶0,ask⟶∞.
hence
(23)lim∥F(xk+αkdk)∥=0,ask⟶∞.To show the convergence results of EMFM method, we require to show that the updating matrixBk is bounded above and below by some positive constants. Hence, we can state the following result on the boundedness of {∥Φk∥F} by assuming that, without loss of generality, the updating matrix (8) is always used; then, we have the following.Theorem 5.
LetF satisfy Assumption 3, and let ∥yk∥≠0 for all finite k. Let {Bk} be the sequence generated by (8). If the given nonsingular B0 satisfies
(24)θ≤B0i≤ω,i=1,2,...,n,
for some constants θ and ω, then the sequence {∥Bk∥F} is bounded for all finite k.Proof.
SinceBk+1=Bk+Φk, it follows that
(25)∥Bk+1∥F≤∥Bk∥F+∥Φk∥F.
For k=0 and assuming B0=I, we have
(26)|Φ0(i)|=|y0Ts0-y0TB0y0tr(Ψ02)(y0(i))2|≤|y0Ts0-y0TB0y0|tr(Ψ02)(y0(max))2,
where (y0(max))2 is the largest element among (y0(i))2,i=1,2,…,n.
After multiplying (26) by (y0(max))2/(y0(max))2 and substituting tr(Ψ02)=∑i=1n‍(y0(i))4, we have
(27)|Φ0(i)|≤|y0Ts0-y0TB0y0|(y0(max))2∑i=1n(y0(i))4(y0(max))4.
Since (y0(max))4/∑i=1n‍(y0(i))4≤1, then (27) turns into
(28)|Φ0(i)|≤|y0TF′(x)y0-y0TB0y0|(y0(max))2.
From Assumption 3 and B0=I, (28) becomes
(29)|Φ0(i)|≤|c-1|(y0Ty0)(y0(max))2,
where c=max{|c1|,|c2|}.
Since(y0(i))2≤(y0(max))2 for i=1,…,n, it follows that
(30)|Φ0(i)|≤n|c-1|(y0(max))2(y0(max))2.
Hence, we obtain
(31)∥Φ0∥F≤n3/2|c-1|.
Suppose that η=n3/2|c-1|; then,
(32)∥Φ0∥F≤η.
From the fact that ∥B0∥F=n, it follows that
(33)∥B1∥F≤β,
where β=n+η>0.


## 4. Numerical Results

In this section, we consider some benchmark problems to illustrate the performance of the method proposed in this paper for solving large-scale systems of nonlinear equations when compared to some Newton-like methods. The computations are performed in MATLAB 7.0 using double precision computer, and the stopping rule used is(34)∥sk∥+∥F(xk)∥≤10-4.
The identity matrix has been chosen as an initial approximate Jacobian inverse. We further design the codes to terminate whenever one of the following happens:(i)
the number of iteration is at least 250, but no point ofxk that satisfies (34) is obtained;(ii)
CPU time in seconds reaches 250;(iii)
There is insufficient memory to initiate the run.The performances of these methods are be compared in terms of number of iterations and CPU time in seconds. In the following, some details on the benchmarks test problems are presented.Problem 1.
System ofn nonlinear equations is as follows:
(35)fi(x)=xi-3xi(sinxi3-0.66)+2,i=1,2,…,n,x0=(3,3,…,3).Problem 2.
Extended Trigonometric function of Spedicator [8] is as follows:
(36)fi(x)=n-∑j=1n‍cosxj+i(1-cosxi)-sinxi+exp(1-cosxi)-∑j=1n‍(xj+1),i=1,…,n,x0=(1n,1n,…,1n).Problem 3.
System ofn nonlinear equations is as follows:
(37)fi(x)=xi-0.1xi+12,fn(x)=xn-0.1x12,i=1,…,n-1,x0=(7,7,…,7).
(A1) Nonadiabatic Stirred Tank Reactors. The first application is the model of [9]. The model deals with two continuous nonadiabatic stirred tank reactors. The reactors are in a series, at steady state having a recycle component, and with an exothermic first-order irreversible reaction. By eliminating certain variables, the model results into two nonlinear systems of equations as follows:
(38)f1=(1-λ)[D10(1+β1)-x1]exp(10x11+(10x1/γ))-x1,f2=x1-(1+β2)x2+(1-λ)×[D10-β1x1-(1+β2x2)]exp(10x21+(10x2/γ)).
The dimensionless temperatures of the two reactors are represented by x1 and x2. The parameters λ, γ, β1, β2, and D are given as 1, 1000, 2, 2, and 22, respectively.
(A2) Navigation by Range Measurements. Consider two beacons determining position by measuring the distances r1 and r2, where (p1,p2) is the position of beacon 1, (q1,q2) is the position of beacon 2, and (u,v) is an unknown point (see Figure 1). This can be modeled into a two nonlinear systems of equations as follows:
(39)f1=(p1-u)2+(q1-v)2-r1=0,f1=(p2-u)2+(q2-v)2-r2=0.
The parameters (p1,q1), (p2,q2), r1, and r2 are given as (10,10), (10,-10), 14, and 16, respectively.Figure 1Table1 shows the number of iterations and CPU time for these five methods, respectively. In Table 1, the value “Dim” denotes the dimension of the systems, “NI” denotes number of iterations, while “CPU” is the CPU time in seconds respectively. We analyze the performance of each method via execution time, floating points operations, and storage locations. One can observe that EMFM has the smallest number of iterations compared to the classical diagonal updating (I-VDN) proposed by Waziri et al. [7]. This shows that the line search strategies presented in this paper have increased the convergence speed of the classical diagonal updating method.Table 1
Numerical comparison of NM, FN, BM, I-VDN, and EMFM methods.
Prob
Dim
NM
FN
BM
I-VDN
EMFM
NI
CPU
NI
CPU
NI
CPU
NI
CPU
NI
CPU
1
25
4
0.031
11
0.016
8
0.015
8
0.006
6
0.001
2
25
4
0.062
11
0.031
8
0.016
8
0.002
7
0.001
3
25
6
0.031
—
—
10
0.031
10
0.001
7
0.001
1
50
4
0.046
11
0.031
8
0.046
8
0.011
6
0.008
2
50
4
0.109
11
0.062
8
0.034
8
0.010
7
0.005
3
50
6
0.156
—
—
10
0.124
10
0.008
7
0.004
1
100
4
0.064
11
0.032
8
0.048
8
0.014
6
0.010
2
100
4
0.125
11
0.046
8
0.032
8
0.015
7
0.011
3
100
6
0.508
—
—
10
0.187
10
0.031
7
0.015
1
1000
4
1.919
12
0.7332
9
2.122
9
0.033
6
0.018
2
1000
4
2.277
12
0.998
8
1.778
8
0.049
5
0.031
3
1000
6
88.031
—
—
10
13.057
10
0.041
7
0.024If we compare the performance of all methods, in terms of CPU time, it is clear that EMFM method consumes less CPU time than the others and still keeping memory requirement and CPU time in seconds to onlyO(n). All five methods are able to obtain the solution (x*=-1.1243,1.5001) of A2, but EMFM method consumes less CPU time in second (0.001) compared to the other 4 methods. Moreover, for A1, still proposed method has shown a promising performance with less storage locations (2 locations) whereas NM, FN, and BM, respectively, required 4 locations for each.


## 5. Conclusion

From the fact that there is a rapid development in research on solving nonlinear systems, nevertheless, the dimension of the nonlinear system is most of the times so large that it requires a very costly arithmetic operations when using some other Newton-like methods; so we usually employ cheap iterative approach, and the good candidate is classical diagonal updating. Notwithstanding, the classical updating scheme has some lapses. One is that it usually needs high number of iterations. The other is that the method generally has slow convergence behavior. This paper proposes an enhanced diagonal secant updating scheme based on the steps of backtracking in the Armijo-type line search and then improved via Wolfe-like condition. Our approach aims at improving the overall performance of diagonal secant updating scheme. The algorithm only requires to store a row vector while ignoring all the off and low diagonal elements and therefore largely reduces memory locations. In addition, as it uses two line search strategies (predictor and corrector) to obtain a new iterates point, the spectral properties of the diagonal updating scheme is improved, and rapid convergence property is gained.Computational experiment suggests that it is very vital for diagonal updating scheme to use line search strategy. EMFM method has very good solving speed and the best performance among the Newton-like methods. Finally, it can be concluded that this approach would certainly be quite useful for solving large-scale systems of nonlinear equations.



---

*Source: 814587-2013-02-19.xml*