A Modified Hestenes-Stiefel-Type Derivative-Free Method for Large-Scale Nonlinear Monotone Equations

Abstract

The goal of this paper is to extend the modified Hestenes-Stiefel method to solve large-scale nonlinear monotone equations. The method is presented by combining the hyperplane projection method (Solodov, M.V.; Svaiter, B.F. A globally convergent inexact Newton method for systems of monotone equations, in: M. Fukushima, L. Qi (Eds.) Reformulation: Nonsmooth, Piecewise Smooth, Semismooth and Smoothing Methods, Kluwer Academic Publishers. 1998, 355–369) and the modified Hestenes-Stiefel method in Dai and Wen (Dai, Z.; Wen, F. Global convergence of a modified Hestenes-Stiefel nonlinear conjugate gradient method with Armijo line search. Numer Algor. 2012, 59, 79–93). In addition, we propose a new line search for the derivative-free method. Global convergence of the proposed method is established if the system of nonlinear equations are Lipschitz continuous and monotone. Preliminary numerical results are given to test the effectiveness of the proposed method.

1. Introduction

In this paper, we consider the problem of finding numerical solutions of the following large-scale nonlinear equations

$$F\left(x\right)=0,$$

(1)

where the function $F:R^n⟶R^n$ is monotone and continuous. If $F\left(x\right)$ is monotone, it implies that the following inequality holds

$$\langle F\left(x\right)-F\left(y\right),x-y\rangle \ge 0,\forall x,y\in R^n.$$

(2)

Nonlinear monotone equations can be applied in different fields, for example, they are used as subproblems in the generalized proximal algorithms with Bregman distances [1]. Some monotone variational inequality problems can be converted into nonlinear monotone equations [2]. Monotone systems of equations can also be applied in $L_1$-norm regularization sparse optimization problems (see [3,4]) and discrete mathematics such as graph theory (see [5,6]).

Being aware of the important applications of nonlinear monotone equations, in recent years, many scholars have paid attention to propose efficient algorithms for solving problem (1). These algorithms are mainly divided into the following categories.

Each of the Newton-type method, Levenberg-Marquardt method, and quasi-Newton method enjoy fast local convergence property, and are attractive (see [7,8,9,10,11,12,13]). However, for large-scale problems, a drawback of these methods is that at each iteration these algorithms require computing a large-scale linear system of equations by using approximate systems of equations or a Jacobian matrix. The large demand of storage for matrix results in improper handling of large-scale nonlinear monotone systems.

In recent years, gradient-type algorithms have attracted the attention of many scholars. The main reasons are low storage requirements, easy implementation, and global convergence under mild conditions. For example, the spectral gradient method [14] only needs to use gradient information which is easy but effective for an optimization problem. From a different perspective, the spectral gradient method [14] is extended to solve nonlinear monotone equations by combining it with the projection method (see [15,16]).

In addition, for large-scale unconstrained optimization problems, the conjugate gradient method (CG) is another easy but effective method, due to the two attractive features: one is the low memory requirement, the other is strong global convergence properties. In recent years, the conjugate gradient method has achieved rich results (see [17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36]) from the perspective of the sufficient descent property, qausi-Newton direction and conjugacy condition. Inspired by the extension of spectral gradient method to nonlinear monotone equations, CG methods have been applied in solving the nonlinear monotone equations problem. For example, PRP-type method ([37,38,39,40,41,42,43]); Perry conjugate gradient method ([44]); Liu-Storey type method [45], among others.

In this paper, we will focus on extending the Hestenes-Stiefel (HS) CG method to solve large-scale nonlinear monotone equations. To the best of our knowledge, the HS CG method [46] is generally considered as the most efficient CG method for computing performance. However, the HS CG method does not enjoy sufficient descent property. Based on the modified secant equation [10], Dai and Wen [47] propose a modified HS conjugate gradient method that can generate sufficient descent directions (i.e., $c>0$ exists such that $g_k^T{d}_k<-c{\parallel g_k\parallel }^2$). Global convergence results under the Armijo line search are obtained in Dai and Wen [47]. Hence, we aim to present a derivative-free method to solve the nonlinear monotone Equations (1). The proposed method can be seem as a further study of the modified HS CG method in Dai and Wen [47] for unconstrained optimization problems.

Our paper makes two contributions to large-scale nonlinear monotone equations. Firstly, a new line search is proposed for the derivative-free method. A significant advantage of this line search is that it is easier to obtain the search stepsize. Secondly, we propose a derivative-free method for solving large-scale nonlinear monotone equations which combines the modified Hestenes-Stiefel method in Dai and Wen [47] for unconstrained optimization problems and the hyperplane projection method [13]. A good property of the proposed method is that it is suitable to solve large-scale nonlinear monotone equations due to its lower storage requirement.

The rest of the article is organized as follows. In Section 2, we give the algorithm and prove the sufficient descent property. The global convergence is proved in Section 3. We report the numerical results In Section 4. The last Section gives the conclusion.

2. Algorithm and the Sufficient Descent Property

In this section, we will present the derivative-free method for solving problem (1), that is a combination of the modified Hestenes-Stiefel method [47] and the hyperplane projection method [13]. Different from the traditional conjugate gradient method, the iteration sequence $\left\{x_{k+1}\right\}$ is obtained in two steps at each iteration.

In the first step, the algorithm produces an iterative sequence $\{z_k=x_k+{\alpha }_k{d}_k\}$, where $d_k$ is the search direction, and ${\alpha }_k>0$ is the steplength obtained by a suitable line search. For most iterative algorithms of optimization problems, the line search plays an important role in convergence analysis and numerical calculation. Zhang and Zhou [15] obtained the steplength ${\alpha }_k>0$ by the following Armijio-type line search: calculating the search steplength ${\alpha }_k=max\{\beta {\rho }^i:i=0,1,\dots ,\}$ such that

$$-F{(x_k+{\alpha }_k{d}_k)}^T{d}_k\ge \sigma {\alpha }_k{\parallel d_k\parallel }^2,$$

(3)

where $\beta$ is some initial attempt for ${\alpha }_k$, $\beta >0$ and $\rho \in (0,1)$.

In addition, Li and Li [39] introduced an alternative line search, that is, computing the search steplength ${\alpha }_k=max\{\beta {\rho }^i:i=0,1,\dots ,\}$ such that

$$-F{(x_k+{\alpha }_k{d}_k)}^T{d}_k\ge \sigma {\alpha }_k\parallel F(x_k+{\alpha }_k{d}_k)\parallel ·\parallel d_k{\parallel }^2,$$

(4)

From the above introduction, we can see that the steplength ${\alpha }_k$ is obtained by calculating ${\alpha }_k=max\{\beta {\rho }^i:i=0,1,\dots ,\}$ such that (3) or (4) is satisfied. If the point $x_k$ is far from the solution, the obtained steplength ${\alpha }_k$ may be very small. Taking this into account, we present the following line search rule where the steplength ${\alpha }_k$ is obtained by computing ${\alpha }_k=max\{\beta {\rho }^i:i=0,1,\dots ,\}$ such that

$$-F{(x_k+{\alpha }_k{d}_k)}^T{d}_k\ge \sigma {\alpha }_k\min \{\parallel d_k{\parallel }^2,\parallel F(x_k+{\alpha }_k{d}_k)\parallel \parallel d_k{\parallel }^2,-F{\left(x_k\right)}^T{d}_k\}.$$

(5)

In the second step, $\left\{x_{k+1}\right\}$ can be determined by using $x_k,z_k,F\left(z_k\right)$ via the hyperplane projection method [13]. Now, let’s introduce how to generate $\left\{x_{k+1}\right\}$ via the hyperplane projection method [13]. Along the search direction $d_k>0$, we can generate a point $z_k=x_k+{\alpha }_k{d}_k$ by a suitable line search such that

$$F{\left(z_k\right)}^T(x_k-z_k)>0.$$

(6)

On the other hand, the monotonicity of F implies that for any solution $x^{\ast }(F\left(x^{\ast }\right)=0)$, the following inequality holds

$$F{\left(z_k\right)}^T(x^{\ast }-z_k)=-{(F\left(z_k\right)-F\left(x^{\ast }\right))}^T(z_k-x^{\ast })\le 0.$$

(7)

From (6) and (7), we can see that $\{F{\left(z_k\right)}^T(x_k-z_k)>0\}$ holds for any $x_k$ and $\{F{\left(z_k\right)}^T(x^{\ast }-z_k)\le 0\}$ holds for the solution $x^{\ast }$. Therefore, from (6) and (7), there is a hyperplane

$${\mathit{H}}_k=\{x\in R^n|F{\left(z_k\right)}^T(x-z_k)=0\},$$

(8)

which can strictly separate the current point $x_k$ from the $x^{\ast }$ (zero point) of equation in (1).

Following Solodov and Svaiter [13] and Zhang and Zhou [15], taking the projection of $x_k$ onto the hyperplane (8) as the next iterate $x_{k+1}$ is a reasonable choice. In detail, the next iterate point $x_{k+1}$ can be computed by

$$x_{k+1}=x_k-\frac{F{\left(z_k\right)}^T(x_k-z_k)}{\parallel F\left(z_k\right){\parallel }^2}F\left(z_k\right).$$

(9)

In what following, we pay our attention to the search direction which plays a crucial role in an iterative algorithm. Our main starting point is to extend the search direction of Dai and Wen [47] to nonlinear monotone equations problem (1). Similar to Dai and Wen [47] for unconstrained optimization, we give the search direction as

$$d_k=\left\{\begin{array}{ccc}-F_0,\hfill & \mathrm{if}\hfill & k=0,\hfill \\ -F_k+{\beta }_k^{\mathrm{NHZ}}{d}_{k-1},\hfill & \mathrm{if}\hfill & k\ge 1,\hfill \end{array}\right.$$

(10)

where

$${\beta }_k^{\mathrm{NHZ}}=\frac{F_k^T{y}_{k-1}}{d_{k-1}^T{w}_{k-1}}-\mu \frac{\parallel y_{k-1}{\parallel }^2}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}{F}_k^T{d}_{k-1},\mu >\frac{1}{4},$$

(11)

$$w_{k-1}=y_{k-1}+\gamma {\overline{s}}_{k-1},\gamma >0,y_{k-1}=F\left(x_k\right)-F\left(x_{k-1}\right),{\overline{s}}_{k-1}=z_{k-1}-x_{k-1}={\alpha }_{k-1}{d}_{k-1}.$$

(12)

For simplicity, we refer to (10) and (11) as NHZ method hereafter.

Further, in this paper, the function F is assumed to satisfy the following assumptions witch are often utilized in convergence analysis for nonlinear monotone equations (see, [37,38,39,40,41,42,43,44,45,48]).

Assumption A1.

($A_1$)

The F is a monotone function:

$${(F\left(x\right)-F\left(y\right))}^T(x-y)\rangle \ge 0,\forall x,y\in R^n.$$

(13)

($A_2$)

The F is Lipschitz continuous function, namely, there exists a $L>0$ such that

$$\parallel F\left(x\right)-F\left(y\right)\parallel \le L\parallel x-y\parallel ,\forall x,y\in R^n.$$

(14)

In what following, we will describe the proposed Algorithm 1.

Algorithm 1: NHZ derivative-free method.

Step 0: Given $x_0\in R^n$ as an initial point, and the constants $\epsilon >0,\beta >0,\sigma >0,\rho \in (0,1)$. Then set k := 0.  Step 1: Calculate $F\left(x_k\right)$. If $\parallel F\left(x_k\right)\parallel \le \epsilon$, stop the algorithm. Otherwise, go to step 2.  Step 2: Determine the search direction $d_k$ by (10), (11) and (12).   Step 3: Calculate the search steplength ${\alpha }_k$ by (5). Let $z_k=x_k+{\alpha }_k{d}_k$.  Step 4: Calculate $F\left(z_k\right)$. If $\parallel F\left(z_k\right)\parallel \le \epsilon$, stop the algorithm. Otherwise, calculate $x_{k+1}$ by using the projection (9). Set k := k + 1 and go to Step 1.

In what follow, we will show that the proposed NHZ derivative-free method enjoys the sufficient descent property which plays an important role in proof of convergence. From now on, we use $F_k$ to denote $F\left(x_k\right)$.

Theorem 1.

The search direction $d_k$ generated by (10), (11) and (12) is sufficient descent direction. That is, if $d_k^T{w}_k\ne 0$, then we have

$$F_k^T{d}_k\le -(1-\frac{1}{4\mu }){\parallel F_k\parallel }^2,\mu >\frac{1}{4}.$$

(15)

Proof. 

When $k=0$, we have

$$F_0^T{d}_0=-\parallel F_0{\parallel }^2\le -(1-\frac{1}{4\mu }){\parallel F_0\parallel }^2.$$

It is obvious that (15) is satisfied for $k=0$. □

Now we will show that the sufficient descent condition (15) holds for $k\ge 1$. We can obtain from (10) and (11) that

$$\begin{array}{ccc}\hfill F_k^T{d}_k& =& -\parallel F_k{\parallel }^2+{\beta }_k^{}{F}_k^T{d}_{k-1}\hfill \\ & =& -\parallel F_k{\parallel }^2+\left\{\frac{F_k^T{y}_{k-1}}{d_{k-1}^T{w}_{k-1}}-\mu \frac{\parallel y_{k-1}{\parallel }^2}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}{F}_k^T{d}_{k-1}\right\}{F}_k^T{d}_{k-1}\hfill \\ & =& \frac{F_k^T{y}_{k-1}\left(d_{k-1}^T{w}_{k-1}\right)\left(F_k^T{d}_{k-1}\right)-\parallel F_k{\parallel }^2{\left(d_{k-1}^T{w}_{k-1}\right)}^2-\mu {\parallel y_{k-1}\parallel }^2{\left(F_k^T{d}_{k-1}\right)}^2}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}.\hfill \end{array}$$

Define

$$u_k=\frac{1}{\sqrt{2\mu }}\left(d_{k-1}^T{w}_{k-1}\right)F_k,v_k=\sqrt{2\mu }\left(F_k^T{d}_{k-1}\right)y_{k-1}.$$

(16)

By using the Equation (16) and the inequality $u_k^T{v}_k\le 1/2(\parallel u_k{\parallel }^2+\parallel v_k{\parallel }^2)$, we have

$$\begin{array}{ccc}\hfill F_k^T{d}_k& =& \frac{u_k^T{v}_k-1/2(\parallel u_k{\parallel }^2+\parallel v_k{\parallel }^2)}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}-(1-\frac{1}{4\mu })\frac{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}{\parallel F_k\parallel }^2\hfill \\ & \le & -(1-\frac{1}{4\mu }){\parallel F_k\parallel }^2.\hfill \end{array}$$

Thus (15) holds for $k\ge 1$.

3. Global Convergence Analysis

Now, we will investigate the global convergence of Algorithm 1. Firstly, we give the following Lemma which shows the line search strategy (5) is well-defined if the search directions $\left\{d_k\right\}$ satisfy the sufficient descent property.

Lemma 1.

If the iterative sequences $\left\{x_k\right\}$ and $\left\{z_k\right\}$ are generated by the Algorithm 1. Then, there always exists a steplength ${\alpha }_k$ satisfying the line search (5).

Proof. 

Assume that, for any nonnegative integer i for $\beta {\rho }^i$, the line search strategy (5) does not hold in the $k_0$-th iterate, then we have

$$-F{(x_{k_0}+\beta {\rho }^i{d}_{k_0})}^T{d}_{k_0}<\sigma \beta {\rho }^i\min \{\parallel d_{k_0}{\parallel }^2,\parallel F(x_{k_0}+\beta {\rho }^i{d}_{k_0})\parallel \parallel d_{k_0}{\parallel }^2,-F{\left(x_{k_0}\right)}^T{d}_{k_0}\}.$$

(17)

Letting $i\mapsto \infty$, we have from the continuity of F and $\rho \in (0,1)$ that

$$-F{(x_{k_0}+\beta {\rho }^i{d}_{k_0})}^T{d}_{k_0}<0,$$

(18)

which contradicts (15). The proof is completed. □

The next Lemma indicates that the line search strategy (5) provides a lower bound for steplength ${\alpha }_k$.

Lemma 2.

If the iterative sequences $\left\{x_k\right\}$ and $\left\{z_k\right\}$ are generated by the Algorithm 1. Then, we can obtain that

$${\alpha }_k\ge \mathit{min}\left\{\beta ,\frac{\delta \rho }{(L+\sigma )\parallel }\frac{\parallel F_k{\parallel }^2}{\parallel d_k{\parallel }^2}\right\},$$

(19)

where $\delta =1-\frac{1}{4\mu }$.

Proof. 

Tf ${\alpha }_k=\beta$, it is obviously that (19) holds. Suppose that ${\alpha }_k\ne \beta$, then we can obtain that ${\alpha }_k^{{}^{\prime }}={\rho }^{-1}{\alpha }_k$ does not satisfy the line search process (5). That is,

$$-F{\left(z_k^{{}^{\prime }}\right)}^T{d}_k<\sigma {\alpha }_k^{{}^{\prime }}\min \{\parallel d_k{\parallel }^2,\parallel F\left(z_k^{{}^{\prime }}\right)\parallel \parallel d_k{\parallel }^2,-F{\left(x_k^{{}^{\prime }}\right)}^T{d}_k\}\le \sigma {\alpha }_k^{{}^{\prime }}{\parallel d_k\parallel }^2,$$

(20)

where $z_k^{{}^{\prime }}=x_k+{\alpha }_k^{{}^{\prime }}{d}_k.$ □

From the sufficient descent condition (15), we have that

$$(1-\frac{1}{4\mu })\parallel F_k{\parallel }^2\doteq \delta {\parallel F_k\parallel }^2\le -F_k^T{d}_k.$$

(21)

From the Lipschitz continuity of F (14), (20) and (21), we can obtain that

$$\begin{array}{ccc}\hfill -F_k^T{d}_k& =& {(F\left(z_k^{{}^{\prime }}\right)-F\left(x_k\right))}^T{d}_k-F{\left(z_k^{{}^{\prime }}\right)}^T{d}_k\hfill \\ & \le & \parallel F\left(z_k^{{}^{\prime }}\right)-F\left(x_k\right)\parallel \parallel d_k\parallel +\sigma {\alpha }_k^{{}^{\prime }}{\parallel d_k\parallel }^2\hfill \\ & \le & L\parallel z_k^{{}^{\prime }}-x_k\parallel \parallel d_k\parallel +\sigma {\alpha }_k^{{}^{\prime }}{\parallel d_k\parallel }^2\hfill \\ & =& L{\alpha }_k^{{}^{\prime }}\parallel d_k{\parallel }^2+\sigma {\alpha }_k^{{}^{\prime }}{\parallel d_k\parallel }^2\hfill \\ & =& {\rho }^{-1}{\alpha }_k(L+\sigma ){\parallel d_k\parallel }^2.\hfill \end{array}$$

Therefore, the above inequalities and (21) imply

$${\alpha }_k\ge \frac{\delta \rho }{(L+\sigma )\parallel }\frac{\parallel F_k{\parallel }^2}{\parallel d_k{\parallel }^2}.$$

This shows that Lemma about the search steplength ${\alpha }_k$ holds.

The next Lemma is proved by Solodov and Svaiter (see Lemma 2.1 in [13]), which can also hold for Algorithm 1. Now, we give this lemma without proof, because its proof is similar in Solodov & Svaiter [13].

Lemma 3.

Assume the function F is monotone and the Lipschitz continuous condition (14) holds. If the iterative sequences $\left\{x_k\right\}$ is generated by the Algorithm 1, then for any $x^{\ast }$, such that $F\left(x^{\ast }\right)=0$, we can obtain that

$$\parallel x_{k+1}-x^{\ast }{\parallel }^2\le \parallel x_k-x^{\ast }{\parallel }^2-{\parallel x_{k+1}-x_k\parallel }^2.$$

In particular, the iterative sequence $\left\{x_k\right\}$ is bounded and

$$\sum _{k=0}^{\infty }{\parallel x_{k+1}-x_k\parallel }^2<\infty .$$

(22)

Remark 1.

The above Lemma 2 confirms that the sequence $\{\parallel x_k-x^{\ast }\parallel \}$ decreases with k. In addition, (22) implies that

$$\underset{k\to \infty }{lim}{\parallel x_{k+1}-x_k\parallel }^{}=0.$$

(23)

Theorem 2.

If the iterative sequences $\left\{x_k\right\}$ is generated by the Algorithm 1, then, we can have that

$$\underset{k\to \infty }{lim}{\alpha }_k\parallel d_k\parallel =0.$$

(24)

Proof. 

We can obtain from (5) and (9) that, for any k,

$$\parallel x_{k+1}-x_k\parallel =\frac{|F{\left(z_k\right)}^T(x_k-z_k)|}{\parallel F\left(z_k\right)\parallel }=\frac{-{\alpha }_{k}F{\left(z_k\right)}^T{d}_k}{\parallel F\left(z_k\right)\parallel }\ge \sigma {\alpha }_k^2{\parallel d_k\parallel }^2.$$

(25)

In particular, it follows from (23) and (25) that

$$\underset{k\to \infty }{lim}{\alpha }_k\parallel d_k\parallel =0.$$

(26)

 □

Lemma 4.

If the iterative sequences $\left\{x_k\right\}$ is generated by the Algorithm 1, and $x^{\ast }$ satisfies $F\left(x^{\ast }\right)=0,z_k^{{}^{\prime }}=x_k+{\alpha }_k^{{}^{\prime }}{d}_k,{\alpha }_k^{{}^{\prime }}={\rho }^{-1}{\alpha }_k$. Then, $\{\parallel F\left(z_k^{{}^{\prime }}\right)\parallel \}$ and $\{\parallel F_k\parallel \}$ are bounded, i.e, there is a constant $M\ge 0$, such that

$$\parallel F\left(z_k^{{}^{\prime }}\right)\parallel \le M,\parallel F_k\parallel \le M.$$

(27)

Proof. 

By the Lemma 2, we have

$$\parallel x_k-x^{\ast }\parallel \le \parallel x_0-x^{\ast }\parallel .$$

From (26), we have that there is a constant $M_1>0$ such that ${\alpha }_k\parallel d_k\parallel \le M_1$. Hence

$$\parallel z_k^{{}^{\prime }}-x^{\ast }\parallel \le \parallel x_k-x^{\ast }\parallel +{\alpha }_k^{{}^{\prime }}\parallel d_k\parallel \le \parallel x_0-x^{\ast }\parallel +{\rho }^{-1}{\alpha }_k\parallel d_k\parallel \le \parallel x_0-x^{\ast }\parallel +M_1.$$

(28)

Since the function $F\left(x\right)$ is Lipschitz continuous, we can easily obtain the following two inequalities

$$\parallel F\left(z_k^{{}^{\prime }}\right)\parallel \le \parallel F\left(z_k^{{}^{\prime }}\right)-F\left(x^{\ast }\right)\parallel \le L\parallel z_k^{{}^{\prime }}-x^{\ast }\parallel \le L(\parallel x_0-x^{\ast }\parallel +{\rho }^{-1}{M}_1,$$

and

$$\parallel F_k\parallel \le \parallel F\left(x_k\right)-F\left(x^{\ast }\right)\parallel \le L\parallel x_k-x^{\ast }\parallel \le L\parallel x_0-x^{\ast }\parallel .$$

Let $M=max\{L\parallel x_0-x^{\ast }\parallel ,L(\parallel x_0-x^{\ast }\parallel +{\rho }^{-1}{M}_1\}$. We can obtain (27). □

In what following, we will give the global convergence theorem for our proposed method.

Theorem 3.

If the iterative sequences $\left\{x_k\right\}$ is generated by the Algorithm 1. We can obtain that

$$\underset{k\to \infty }{liminf}\parallel F_k\parallel =0.$$

(29)

Proof. 

We will prove this Theorem by contradiction. Assume that (29) is not true. Then it implies there is a constant $\epsilon >0$, s.t. $\parallel F_k\parallel >\epsilon$.

Since $F_k\ne 0$, we have from (15) that $d_k\ne 0$. Hence, the monotonicity of F together with (10) implies

$$\begin{array}{ccc}\hfill {\overline{s}}_{k-1}^T{w}_{k-1}& =& \langle F\left(z_{k-1}\right)-F\left(x_{k-1}\right),z_{k-1}-x_{k-1}\rangle +\gamma {\overline{s}}_{k-1}^T{\overline{s}}_{k-1}\hfill \\ & \ge & \gamma {\overline{s}}_{k-1}^T{\overline{s}}_{k-1}.\hfill \end{array}$$

This together with the definition of ${\overline{s}}_{k-1}$ implies

$$d_{k-1}^T{w}_{k-1}\ge \gamma {\alpha }_{k-1}{\parallel d_{k-1}\parallel }^2.$$

(30)

We have from (11), (12) and (30) that

$$\begin{array}{ccc}\hfill |{\beta }_k^{\mathrm{NHZ}}|& =& |\frac{F_k^T{y}_{k-1}}{d_{k-1}^T{w}_{k-1}}-\mu \frac{\parallel y_{k-1}{\parallel }^2}{{\left(d_{k-1}^T{w}_{k-1}\right)}^2}{F}_k^T{d}_{k-1}|\hfill \\ & \le & \frac{L{\alpha }_{k-1}\parallel d_{k-1}\parallel \parallel F_k\parallel }{\gamma {\alpha }_{k-1}{\parallel d_{k-1}\parallel }^2}+\mu \frac{L^2{\alpha }_{k-1}^2{\parallel d_{k-1}\parallel }^2}{{\gamma }^2{\alpha }_{k-1}^2{\parallel d_{k-1}\parallel }^4}\parallel F_k\parallel \parallel d_{k-1}\parallel \hfill \\ & \le & (\frac{L}{\gamma }+\mu \frac{L^2}{{\gamma }^2})\frac{\parallel F_k\parallel }{\parallel d_{k-1}\parallel }.\hfill \end{array}$$

Therefore, from (10) and (30), we can obtain

$$\begin{array}{ccc}\hfill \parallel d_k\parallel & \le & \parallel F_k\parallel +|{\beta }_k^{\mathrm{NHZ}}|\parallel d_{k-1}\parallel \hfill \\ & \le & \parallel F_k\parallel +(\frac{L}{\gamma }+\mu \frac{L^2}{\gamma })\parallel F_k\parallel \hfill \\ & \le & (1+\frac{L}{\gamma }+\mu \frac{L^2}{{\gamma }^2})M.\hfill \end{array}$$

Define $C=(1+\frac{L}{\gamma }+\mu \frac{L^2}{{\gamma }^2})M$. Then, we can have $\parallel d_k\parallel \le C$. □

It follows from Lemmas 2, Lemmas 3, $\parallel F_k\parallel \ge \epsilon$ and $\parallel d_k\parallel \ge \epsilon$ that for all k sufficiently large,

$$\begin{array}{ccc}\hfill {\alpha }_k\parallel d_k\parallel & \ge & \min \left\{\beta ,\frac{\delta \rho }{(L+\sigma )\parallel }\frac{\parallel F_k{\parallel }^2}{\parallel d_k{\parallel }^2}\right\}\parallel d_k\parallel \hfill \\ & \ge & \min \left\{\beta \epsilon ,\frac{\delta \rho {\epsilon }^2}{(L+\sigma )\parallel C}\right\}>0.\hfill \end{array}$$

It is obvious that the above inequality contradicts with (24). That is, (29) holds. And we complete this proof.

4. Numerical Experiments

Now, we will give some numerical experiments to test the numerical performance of our proposed method. We try to test the NHZ Algorithm 1 and compare it’s performance with the spectral gradient (SG) method [15] and the MPRP method in [39]. In the testing experiments, all codes were written in Matlab R2018a, and run on a Lenovo PC with 4 GB RAM memory.

To obtain better numerical performance, we select the following initial steplength as in [39] and [43]

$$\beta =|\frac{F{\left(x_k\right)}^T{d}_k}{d_k^T(F(x_k+ϵd_k)-F\left(x_k\right))/ϵ}|\approx |\frac{F{\left(x_k\right)}^T{d}_k}{d_k^T\nabla F\left(x_k\right)d_k}|,ϵ={10}^{-8}.$$

(31)

We set $\rho =0.5,\sigma =2$. Further, let $\beta =1$ if $\beta <{10}^{-4}$.

Following the MPRP method in [39], we terminate the iterative process if the following condition

$$\min \{\parallel F\left(x_k\right)\parallel ,\parallel F\left(z_k\right)\parallel \}\le atol+rtol\parallel F\left(x_0\right)\parallel$$

is satisfied, and $rtol=atol={10}^{-4}$.

The numerical performance of SG, MPRP, and NHZ methods are tested by using the following five nonlinear monotone equations problem with different various sizes and initial points.

Problem 1

([49]). The specific expression of the function $F\left(x\right)$ is defined as

$$F_i\left(x\right)=2x_i-sin\left(\right|x_i\left|\right),i=1,\dots ,n.$$

Problem 2

([49]). The specific expression of the function $F\left(x\right)$ is defined as

$$\begin{array}{ccc}\hfill F_1\left(x\right)& =& 2x_1+sin\left(x_1\right)-1,\hfill \\ \hfill F_i\left(x\right)& =& -2x_{i-1}+2x_i+sin\left(x_i\right)-1,i=2,\dots ,n-1,\hfill \\ \hfill F_n\left(x\right)& =& 2x_n+sin\left(x_n\right)-1.\hfill \end{array}$$

Problem 3

([50]). The specific expression of the function $F\left(x\right)$ is defined as

$$F_i\left(x\right)=x_i-sin\left(x_i\right),i=1,\dots ,n.$$

Problem 4

([50]). The specific expression of the function $F\left(x\right)$ is defined as

$$F\left(x\right)=Ax+g\left(x\right)$$

where $g\left(x\right)={(e^{x_1}-1,e^{x_2}-1,\dots ,e^{x_n}-1)}^T$

$$A=\left(\begin{array}{ccccc}2& -1& & & \\ -1& 2& -1& & \\ & \ddots & \ddots & \ddots & \\ & & \ddots & \ddots & -1\\ & & & -1& 2\end{array}\right)$$

Problem 5.

The specific expression of the function $F\left(x\right)$ is defined as

$$F\left(x\right)=Ax+\left|X\right|-B$$

where $\left|X\right|=\left(\right|x_1|,|x_2|,\dots ,|x_n{\left|\right)}^T,B={(1,1,\dots ,1)}^T$, and

$$A=\left(\begin{array}{ccccc}2& -1& & & \\ -1& 2& -1& & \\ & \ddots & \ddots & \ddots & \\ & & \ddots & \ddots & -1\\ & & & -1& 2\end{array}\right).$$

The numerical results for 5 tested problems are reported in

Table 1. Numerical results for the tested Problem 1 with various sizes and given initial points.

Initial	Dim.		SG			MPRP			NHZ
		Time	Iter	Feval	Time	Iter	Feval	Time	Iter	Feval
$x_1$	1000	1.16	16	37	0.81	12	39	0.65	10	29
$x_2$	1000	1.16	16	37	0.81	12	39	0.65	10	29
$x_3$	1000	0.77	13	26	0.55	8	27	0.48	7	22
$x_4$	1000	1.24	15	43	0.85	9	42	0.78	7	35
$x_5$	1000	1.15	14	44	0.53	6	28	0.48	5	24
$x_6$	1000	1.15	14	44	0.53	6	28	0.48	5	24
$x_1$	5000	6.96	18	37	5.32	13	42	4.55	11	32
$x_2$	5000	6.96	17	37	5.01	12	39	3.85	11	29
$x_3$	5000	4.16	11	22	3.01	7	24	2.43	6	20
$x_4$	5000	10.45	18	62	7.38	12	60	5.99	11	48
$x_5$	5000	6.88	14	44	3.40	6	28	2.84	5	24
$x_6$	5000	6.96	15	45	3.40	6	28	2.84	5	24
$x_1$	10,000	27.77	18	38	21.15	13	42	15.72	11	32
$x_2$	10,000	27.62	17	36	19.65	12	39	14.36	11	27
$x_3$	10,000	14.90	10	20	11.92	7	24	8.96	6	20
$x_4$	10,000	44.65	20	69	35.15	14	72	30.98	12	68
$x_5$	10,000	29.55	15	45	13.86	6	28	11.22	5	24
$x_6$	10,000	29.55	15	45	13.86	6	28	11.22	5	24

,

Table 2. Numerical results for the tested Problem 2 with various sizes and given initial points.

Initial	Dim.		SG			MPRP			NHZ
		Time	Iter	Feval	Time	Iter	Feval	Time	Iter	Feval
$x_1$	1000	0.67	326	646	0.17	196	591	0.12	155	568
$x_2$	1000	0.64	364	714	0.11	203	612	0.10	180	584
$x_3$	1000	0.64	343	691	0.19	195	588	0.17	174	555
$x_4$	1000	1.18	468	941	0.15	243	741	0.14	220	709
$x_5$	1000	0.67	479	988	0.17	194	585	0.15	162	549
$x_6$	1000	0.65	431	857	0.16	189	571	0.14	165	529
$x_1$	5000	10.38	723	2047	8.19	487	1465	8.10	468	1448
$x_2$	5000	16.87	753	2622	12.85	769	2310	12.63	744	2205
$x_3$	5000	9.01	824	2849	7.98	469	1411	7.65	442	1324
$x_4$	5000	18.74	1023	3261	15.96	929	2840	15.51	901	2781
$x_5$	5000	10.05	1226	3053	7.59	453	1363	7.50	442	1320
$x_6$	5000	11.32	836	2476	6.32	371	1122	6.20	338	1103
$x_1$	10,000	47.90	960	2002	34.59	539	1622	10.01	460	1508
$x_2$	10,000	82.62	1334	4066	65.06	1023	3072	60.79	1001	3003
$x_3$	10,000	50.99	833	2469	34.30	516	1553	32.32	501	1502
$x_4$	10,000	56.82	2042	6668	39.65	1668	5075	36.31	1602	5003
$x_5$	10,000	49.63	832	2268	31.57	497	1497	30.25	436	1405
$x_6$	10,000	45.70	850	1706	25.36	396	1202	22.24	375	1106

,

Table 3. Numerical results for the tested Problem 3 with various sizes and given initial points.

Initial	Dim.		SG			MPRP			NHZ
		Time	Iter	Feval	Time	Iter	Feval	Time	Iter	Feval
$x_1$	1000	1.16	16	37	0.81	12	39	0.62	10	28
$x_2$	1000	1.17	17	36	0.83	12	39	0.71	11	28
$x_3$	1000	0.77	11	24	0.57	8	27	0.49	7	28
$x_4$	1000	1.25	14	44	0.88	9	42	0.75	7	32
$x_5$	1000	1.16	13	42	0.56	6	28	0.48	5	22
$x_6$	1000	1.16	13	42	0.57	6	28	0.48	5	22
$x_1$	5000	6.98	17	36	5.42	13	42	4.63	11	32
$x_2$	5000	6.98	17	36	5.11	12	39	3.95	11	30
$x_3$	5000	4.29	10	22	3.12	7	24	2.34	6	20
$x_4$	5000	10.57	19	64	7.46	12	60	6.25	11	52
$x_5$	5000	6.99	13	42	3.52	6	28	3.92	5	24
$x_6$	5000	6.99	13	42	3.52	6	28	3.92	5	24
$x_1$	10,000	27.78	17	36	21.35	13	42	15.97	11	32
$x_2$	10,000	27.79	17	36	19.75	12	39	15.86	11	30
$x_3$	10,000	15.65	9	26	11.99	7	24	9.98	6	19
$x_4$	10,000	44.85	20	69	35.36	14	72	29.98	12	60
$x_5$	10,000	29.89	14	45	13.98	6	28	12.56	5	24
$x_6$	10,000	29.89	14	45	13.98	6	28	13.59	6	24

,

Table 4. Numerical results for the tested Problem 4 with various sizes and given initial points.

Initial	Dim.		SG			MPRP			NHZ
		Time	Iter	Feval	Time	Iter	Feval	Time	Iter	Feval
$x_1$	1000	0.20	219	431	0.06	50	216	0.05	38	168
$x_2$	1000	0.28	261	463	0.06	56	252	0.05	46	185
$x_3$	1000	0.28	224	329	0.05	34	152	0.03	32	137
$x_4$	1000	0.22	263	529	0.07	100	421	0.06	96	399
$x_5$	1000	0.28	183	403	0.06	42	187	0.05	40	177
$x_6$	1000	0.28	212	424	0.06	60	261	0.05	48	218
$x_1$	5000	2.15	263	456	1.11	48	209	1.05	47	183
$x_2$	5000	2.45	225	378	1.19	46	224	0.92	38	169
$x_3$	5000	1.65	122	267	0.62	27	117	0.65	29	128
$x_4$	5000	3.41	265	558	2.59	109	483	2.49	104	455
$x_5$	5000	2.86	290	467	1.21	53	231	1.07	44	189
$x_6$	5000	2.97	231	477	1.20	54	234	1.16	48	213
$x_1$	10,000	5.30	278	502	3.96	45	195	3.83	42	185
$x_2$	10,000	6.26	237	574	4.27	41	210	3.84	38	158
$x_3$	10,000	5.62	275	585	1.93	42	96	2.62	35	142
$x_4$	10,000	18.15	333	596	10.86	117	533	9.45	109	498
$x_5$	10,000	13.52	341	595	4.34	49	212	3.85	44	186
$x_6$	10,000	13.55	336	553	4.89	56	246	3.78	48	195

and

Table 5. Numerical results for the tested Problem 5 with various sizes and given initial points.

Initial	Dim.		SG			MPRP			NHZ
		Time	Iter	Feval	Time	Iter	Feval	Time	Iter	Feval
$x_1$	1000	0.89	119	289	0.66	47	199	0.44	38	168
$x_2$	1000	0.78	122	263	0.45	22	105	0.44	24	98
$x_3$	1000	0.69	130	235	0.35	48	209	0.28	38	120
$x_4$	1000	0.85	190	249	0.47	37	165	0.34	35	98
$x_5$	1000	0.75	194	248	0.55	94	237	0.45	66	192
$x_6$	1000	1.22	225	462	0.79	174	396	0.75	142	372
$x_1$	5000	2.32	113	260	1.22	51	221	0.98	42	172
$x_2$	5000	2.92	128	270	0.56	22	105	0.58	28	96
$x_3$	5000	3.80	228	412	1.11	47	200	0.79	44	144
$x_4$	5000	3.50	216	424	1.20	48	206	0.79	44	142
$x_5$	5000	3.00	226	443	1.17	47	206	0.81	44	122
$x_6$	5000	6.57	461	881	5.06	308	707	4.25	262	628
$x_1$	10,000	5.92	66	209	3.90	44	191	3.42	38	184
$x_2$	10,000	6.86	68	218	51.1	22	105	49.1	21	98
$x_3$	10,000	5.76	60	181	4.24	47	204	3.23	38	132
$x_4$	10,000	11.84	69	227	10.5	48	209	8.52	44	148
$x_5$	10,000	10.55	68	221	4.02	45	196	3.82	42	168
$x_6$	10,000	12.46	89	326	10.1	74	262	7.83	68	232

respectively, where the given initial points are $x_1={(0.1,\dots ,0.1)}^T,x_2={(1,\dots ,1)}^T,x_3=(1,\frac{1}{2},$$\dots ,\frac{1}{n}{)}^T,x_4={(-10,\dots ,-10)}^T,x_5={(-0.1,\dots ,-0.1)}^T,x_6={(-1,\dots ,-1)}^T.$ Meanwhile, the numerical results are listed in the Table 1, Table 2, Table 3, Table 4 and Table 5 where (Time) is the CPU time (in seconds), (Iter) is the number of iterations and (Feval) is the number of function evaluations.

Table 1, Table 2, Table 3, Table 4 and Table 5 report the numerical results of the proposed algorithm, the spectral gradient (SG) method [15] and the MPRP method in [39] with five tested problems where the test indicators of Time, Iter and Feval are used to evaluate numerical performance.

Comparing the CPU time within the tested algorithms, we note that our proposed algorithm needs lower time consuming than both the spectral gradient (SG) method [15] and the MPRP method in [39] for all tested problems, and the difference is substantial and significant especially for large-scale problems. In addition, the MPRP method in [39] needs lower CPU time than the spectral gradient (SG) method [15]. Assessing the number of iterations within the tested algorithms, we find that the NHZ method requires fewer iterations than the spectral gradient (SG) method [15] and the MPRP method in [39] for all tested problems. We also note that our proposed algorithm requires fewer number of function evaluations for all tested problems, and the difference is substantial and significant.

In sum, from the numerical results in Table 1, Table 2, Table 3, Table 4 and Table 5, it isn’t difficult to see that the proposed algorithm performs better than the spectral gradient (SG) method [15] and the MPRP method in [39] for three test indicators which implies that the modified Hestenes-Stiefel-based derivative-free method is computationally efficient for nonlinear monotone equations.

5. Conclusions

This paper aims to present a modified Hestenes-Stiefel method to solve the nonlinear monotone equations which combines the hyperplane projection method [13] and the modified Hestenes-Stiefel method in Dai and Wen [47]. In the proposed method, the search direction satisfies sufficient descent conditions. A new line search is proposed for the derivative-free method. Under appropriate conditions, the proposed method converges globally. The given numerical results show the presented method is more efficient compared to the methods proposed by the spectral gradient method Zhang and Zhou [15] and the MPRP method in Li and Li [39].

In addition, we also expect that our proposed method and its further modifications could produce new applications for problems in relevant areas of symmetric equations [51], image processing [52], and finance [53,54,55].