Hybrid Modification of Accelerated Double Direction Method

Abstract

We present a hybridization of the accelerated gradient method with two vector directions. This hybridization is based on the usage of a chosen three-term hybrid model. Derived hybrid accelerated double direction model keeps preferable properties of both included methods. Convergence analysis demonstrates at least linear convergence of the proposed iterative scheme on the set of uniformly convex and strictly convex quadratic functions. The results of numerical experiments confirm better performance profile in favor of derived hybrid accelerated double direction model when compared to its forerunners.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

The main goal herein is to derive an efficient optimization method for minimization of an objective function $f:{\mathbb{R}}^n\to \mathbb{R}$ . Therewith, we assume the function $f$ is uniformly convex and twice continuously differentiable. Furthermore, for the gradient and the Hessian of the function $f$ at the $k$ -th iterative point we use the next notation: $$\begin{array}{c}\text{(1)}& g_k\left(x\right)=▽f\left(x_k\right),G_k\left(x\right)={▽}^{\mathrm{2}}f\left(x_k\right).\end{array}$$ The general form of the iterations for finding the extreme values of the objective function $f$ is given by the next expression: $$\begin{array}{}\text{(2)}& x_{k+\mathrm{1}}=x_k+t_k{d}_k,\end{array}$$ where $x_k$ is the current, $x_{k+\mathrm{1}}$ is the next iterative point, $t_k$ is iterative step length value, and $d_k$ is an iterative vector direction which leads us to the solution of the problem. Certainly, $t_k$ and $d_k$ are the most important issues of an optimization model (2) and they generate the efficiency of a relevant method. For that reason, the way of defining these two crucial elements is of great importance for each minimization scheme.

In one of the first algorithms for solving unconstrained optimization problems, denoted as the steepest descent gradient method which is exposed by Cauchy, the iteration is defined as $$\begin{array}{}\text{(3)}& x_{k+\mathrm{1}}=x_k-t_k{g}_k,\end{array}$$ and here, the descent direction is simply presented as the negative gradient vector, while iterative step size value is calculated by the exact line search formula: $$\begin{array}{}\text{(4)}& t_k=\arg \underset{t>\mathrm{0}}{\min }f\left(x_k+t{d}_k\right).\end{array}$$

Furthermore, in the general Newton method $$\begin{array}{}\text{(5)}& x_{k+\mathrm{1}}=x_k-G_k^{-\mathrm{1}}{g}_k{t}_k,\end{array}$$ the vector direction is calculated as the product of the inverse Hessian and the gradient of the objective function. Defining the vector direction in this way guarantees fast convergence properties, but still, practical computing of the function Hessian and its inverse can be difficult. And so, many modified Newton, Newton-conjugate, and quasi Newton schemes were developed in which the calculation of the Hessian and its inverse is, somehow, avoided.

In the quasi-Newton methods the Hessian of the goal function or its inverse is approximated by the adequately defined matrix. Using this type of methods we generally reduce the time of computations since we avoid the complicated calculations in deriving the Hessian of the objective function. Nevertheless, the methods of quasi Newton type preserve good properties of the Newton method. For these reasons, in this paper, we propose the method of quasi Newton type where the value of the iterative step size parameter $t_k$ is obtained by the inexact Backtracking line search procedure.

In the second section we give an overview of some accelerated gradient methods and hybrid iterations. We elaborate the deriving of the hybrid accelerated double direction method and restate the algorithm in the third section of this paper. In the fourth section we give a convergence analysis regarding the proposed iteration. Numerical experiments and comparison are presented in the last section of this paper.

2. Preliminaries: Accelerated Gradient Methods and Hybrid Iterations

The authors in [1] rightfully detected a class of accelerated gradient descent methods, defined by the general iterative scheme $$\begin{array}{}\text{(6)}& x_{k+\mathrm{1}}=x_k-{\gamma }_k^{-\mathrm{1}}{t}_k{g}_k.\end{array}$$ In the previous expression, ${\gamma }_k$ presents an iterative acceleration parameter which improves performance of the relevant method. A common way to determine this parameter is through the features of the second-order Taylor’s series taken on appropriate scheme (6). Acceleration parameters that were computed in such way are applied in the methods described in [1–5]. According to the iteration form (6), we can conclude that the accelerated gradient methods are of the quasi-Newton type in which the approximation of the Hessian, i.e., its inverse, is obtained by the scalar matrix ${\gamma }_{k}I$ , where $I$ is appropriate identity matrix and ${\gamma }_k=\gamma (x_k,x_{k-\mathrm{1}})$ is the matching acceleration parameter. Here are several expressions for defining the acceleration parameters of some accelerated gradient schemes: $$\begin{array}{}\text{(7)}& {\gamma }_{k+\mathrm{1}}^{SM}=\mathrm{2}{\gamma }_k\frac{{\gamma }_k\left[f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)\right]+t_k{‖g_k‖}^{\mathrm{2}}}{t_k^{\mathrm{2}}{‖g_k‖}^{\mathrm{2}}}\text{.}\end{array}$$ (SM method [1]) $$\begin{array}{}\text{(8)}& {\gamma }_{k+\mathrm{1}}^{ADD}=\mathrm{2}\frac{f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)-{\alpha }_k{g}_k^{\mathrm{T}}\left({\alpha }_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}{{\left({\alpha }_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}\left({\alpha }_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)},\end{array}$$ (ADD method [4]) $$\begin{array}{}\text{(9)}& {\gamma }_{k+\mathrm{1}}^{ADSS}=\mathrm{2}\frac{f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)+\left({\alpha }_k{{\gamma }_k}^{-\mathrm{1}}+{\beta }_k\right){‖g_k‖}^{\mathrm{2}}}{{\left({\alpha }_k{{\gamma }_k}^{-\mathrm{1}}+{\beta }_k\right)}^{\mathrm{2}}{‖g_k‖}^{\mathrm{2}}},\end{array}$$ (ADSS method [2]) $$\begin{array}{}\text{(10)}& {\gamma }_{k+\mathrm{1}}^{TADSS}=\mathrm{2}\frac{f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)+{\psi }_k{‖g_k‖}^{\mathrm{2}}}{{\psi }_k^{\mathrm{2}}{‖g_k‖}^{\mathrm{2}}},\end{array}$$ (TADSS method [5]) $$\begin{array}{}\text{(11)}& {\gamma }_{k+\mathrm{1}}^{HSM}=\mathrm{2}{\gamma }_k\frac{{\gamma }_k\left[f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)\right]+\left({\alpha }_k+\mathrm{1}\right)t_k{‖g_k‖}^{\mathrm{2}}}{{\left({\alpha }_k+\mathrm{1}\right)}^{\mathrm{2}}{t}_k^{\mathrm{2}}{‖g_k‖}^{\mathrm{2}}},\end{array}$$ (HSM method [3]).

An interesting concept of merging iterations through the hybrid expression was suggested in some research articles (see [6–8]). Some of representations are given by the next set of iterations: $$\begin{array}{c}\text{(12)}& u_{\mathrm{1}}=u\in \mathbb{C},u_{k+\mathrm{1}}=T{u}_k,\hspace{1em}k\in \mathbb{N},\\ \text{(13)}& v_{\mathrm{1}}=v\in \mathbb{C},v_{k+\mathrm{1}}=\left(\mathrm{1}-{\alpha }_k\right)v_k+{\alpha }_{k}T{v}_k,\hspace{1em}k\in \mathbb{N},\\ \text{(14)}& z_{\mathrm{1}}=z\in \mathbb{C},z_{k+\mathrm{1}}=\left(\mathrm{1}-{\alpha }_k\right)z_k+{\alpha }_{k}T{y}_k,\\ y_k=\left(\mathrm{1}-{\beta }_k\right)z_k+{\beta }_{k}T{z}_k,\hspace{1em}k\in \mathbb{N},\end{array}$$ where $T:\mathbb{C}\to \mathbb{C}$ is a mapping defined on nonempty convex subset $\mathbb{C}$ of a normed space $\mathbb{E}$ , $v_k,z_k$ , and $y_k$ present the sequences defined by proposed iterations, and $\{{\alpha }_k\},\{{\beta }_k\}\in (\mathrm{0,1})$ .

In [9] it was proved that the hybrid method $$\begin{array}{c}\text{(15)}& x_{\mathrm{1}}=x\in \mathbb{R},x_{k+\mathrm{1}}=T{y}_k,\\ y_k=\left(\mathrm{1}-{\alpha }_k\right)x_k+{\alpha }_{k}T{x}_k,\hspace{1em}k\in \mathbb{N},\end{array}$$ proposed by Picard, Mann, and Ishikawa, upgrades the hybrid models mentioned above. The authors of [3] used the advantages of the hybrid model (15) and derived a hybrid version of the accelerated gradient SM method from [1], termed the HSM method and defined by $$\begin{array}{}\text{(16)}& x_{k+\mathrm{1}}=x_k-\left({\alpha }_k+\mathrm{1}\right)t_k{\gamma }_k^{-\mathrm{1}}{g}_k.\end{array}$$ Numerical tests from [3] confirmed that the hybrid model (16) upgrades its forerunner SM iterative rule.

3. HADD Algorithm

We are motivated by the confirmed advantages which were approved in [3] when the scheme (15) was applied on the SM method. As a result, the hybrid SM model (called HSM) was defined and tested in [3]. Herein, we apply the same hybridization strategy to the accelerated double direction method (ADD method, shortly), introduced in [4]. Derived scheme will be based on the hybrid scheme (15) and with that it keeps the accelerated features of the ADD iterations.

In order to complete the presentation, we start from the ADD iteration: $$\begin{array}{}\text{(17)}& x_{k+\mathrm{1}}=x_k-t_k{{\gamma }_k}^{-\mathrm{1}}{g}_k+t_k^{\mathrm{2}}{d}_k,\end{array}$$ where $t_k$ is appropriately defined step size, the first direction vector is given by $-{{\gamma }_k}^{-\mathrm{1}}{g}_k=-{\left({\gamma }_k^{ADD}\right)}^{-\mathrm{1}}{g}_k$ , and the second one, $d_k$ , is determined based on the next procedure Second direction. That procedure was introduced in [4], which was derived as a practical appearance of the more general procedure considered in [10]. The procedure Second direction is restated in Algorithm 1.

Algorithm 1: Procedure Second direction (calculation of the second direction vector $d_k$ ).

Require: Objective function $f(x)$ , gradient $g_k$ and stepsize $t_k$ .

1: Compute

$$d_k(t)=\left\{\begin{array}{cc}{d}_k^{\ast }\hfill & \mathrm{i}\mathrm{f}\hspace{0.17em}\hspace{0.17em}k\le m-\mathrm{1}\hfill \\ {\displaystyle \sum }_{i=\mathrm{2}}^m‍t^{i-\mathrm{1}}{d}_{k-i+\mathrm{1}}^{\ast }\hfill & \mathrm{i}\mathrm{f}\hspace{0.17em}\hspace{0.17em}k\ge m,\hfill \end{array}\right.$$

where $t=t_k$ is the step size, $d_k^{\ast }$ is the solution of the problem ${\mathrm{m}\mathrm{i}\mathrm{n}}_{x\in R}{\Phi }_k(d)$ , and

$${\Phi }_k(d)=\nabla f{\left(x_k\right)}^{\mathrm{T}}d+\frac{\mathrm{1}}{\mathrm{2}}{\gamma }_{k+\mathrm{1}}I=g_k^{\mathrm{T}}d+\frac{\mathrm{1}}{\mathrm{2}}{\gamma }_{k+\mathrm{1}}I.$$

2: Return $d_k=d_k^{\ast }$ .

Remark 1.

For further investigation within this topic, the second direction $d_k$ in the ADD iteration can be defined differently. For example, in [11] the authors proposed directional k-step Newton methods for solving a single nonlinear equation in n-variables. Accordingly, they established the semi-local convergence analysis for these models, based on two different approaches. The first one is based on recurrent relations, while the other, more preferable, is established using recurrent functions. Using one (or both) approaches from [11] in determining the second direction $d_k$ in the ADD method as well as in its hybrid version can be an interesting topic in further research.

Applying the hybrid scheme (15) on the iterative rule (17), we get the hybrid iterative scheme $$\begin{array}{c}\text{(18)}& x_{\mathrm{1}}=x\in \mathbb{R},x_{k+\mathrm{1}}=T{y}_k=y_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k+t_k^{\mathrm{2}}{d}_k,\\ y_k=\left(\mathrm{1}-{\alpha }_k\right)x_k+{\alpha }_{k}T{x}_k=\left(\mathrm{1}-{\alpha }_k\right)x_k+{\alpha }_k\left(x_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k+t_k^{\mathrm{2}}{d}_k\right)\\ =x_k-{\alpha }_k{t}_k{\gamma }_k^{-\mathrm{1}}{g}_k+{\alpha }_k{t}_k^{\mathrm{2}}{d}_k,\hspace{1em}k\in \mathbb{N}.\end{array}$$ After replacing the third expression from the set of equations (18) into the second one, the next iterative rule follows: $$\begin{array}{}\text{(19)}& x_{k+\mathrm{1}}=x_k-\left({\alpha }_k+\mathrm{1}\right)t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\left({\alpha }_k+\mathrm{1}\right)t_k^{\mathrm{2}}{d}_k.\end{array}$$

To simplify further calculation, we will use a constant value for the parameter ${\alpha }_k\in (\mathrm{0,1})$ in (19), just like the authors did in [1, 9]. So, in (19), instead of ${\alpha }_k+\mathrm{1}\in (\mathrm{1,2})$ we simply take $\alpha \in (\mathrm{1,2})$ . Now we can restate a hybrid ADD method, or the HADD iterative scheme, as follows: $$\begin{array}{}\text{(20)}& x_{k+\mathrm{1}}=x_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\alpha t_k^{\mathrm{2}}{d}_k,\hspace{1em}\alpha \in \left(\mathrm{1,2}\right).\end{array}$$

Yet, we need to determine the iterative value of the accelerated parameter ${\gamma }_k={\gamma }_k^{HADD}.$ As we mentioned previously, this parameter can be appropriately defined using Taylor’s expansion of the proposed iteration (20) in two successive iterative points: $$\begin{array}{}\text{(21)}& f\left(x_{k+\mathrm{1}}\right)\approx f\left(x_k\right)+g_k^{\mathrm{T}}\alpha \left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}{\nabla }^{\mathrm{2}}f\left(\xi \right)\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right).\end{array}$$ The parameter $\xi$ in the previous expansion fulfills the condition $$\begin{array}{c}\text{(22)}& \xi \in \left[x_k,x_{k+\mathrm{1}}\right],\xi =x_k+\kappa \left(x_{k+\mathrm{1}}-x_k\right)=x_k+\kappa \left(\alpha t_k^{\mathrm{2}}{d}_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right),\hspace{1em}\mathrm{0}\le \kappa \le \mathrm{1}.\end{array}$$ In the next relation we substitute the value ${▽}^{\mathrm{2}}f(\xi )$ from (21) by the scalar matrix ${\gamma }_{k+\mathrm{1}}I$ , which leads to $$\begin{array}{}\text{(23)}& f\left(x_{k+\mathrm{1}}\right)=f\left(x_k\right)+g_k^{\mathrm{T}}\alpha \left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{\gamma }_{k+\mathrm{1}}{\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right).\end{array}$$ From (23), it is possible to derive the approximation factor ${\gamma }_{k+\mathrm{1}}={\gamma }_{k+\mathrm{1}}^{HADD}$ of the HADD scheme: $$\begin{array}{}\text{(24)}& {\gamma }_{k+\mathrm{1}}^{HADD}={\gamma }_{k+\mathrm{1}}=\mathrm{2}\frac{f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)-\alpha g_k^{\mathrm{T}}\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)}{{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}.\end{array}$$

With the aim of preserving the Second-Order Necessary Condition and Second-Order Sufficient Condition, we assume positivity of the acceleration parameter: ${\gamma }_{k+\mathrm{1}}>\mathrm{0}.$ In practical computation, it is possible that (24) generates negative value for ${\gamma }_{k+\mathrm{1}}.$ We resolve this situation by taking ${\gamma }_{k+\mathrm{1}}=\mathrm{1}$ in such cases. As a consequence, then the first vector direction becomes the negative gradient vector $-g_k.$ In this special case, the next iterative point of the iteration (20) becomes $$\begin{array}{}\text{(25)}& x_{k+\mathrm{2}}=x_{k+\mathrm{1}}-\alpha t_{k+\mathrm{1}}{g}_{k+\mathrm{1}}+\alpha t_{k+\mathrm{1}}^{\mathrm{2}}{d}_{k+\mathrm{1}}.\end{array}$$

In order to present the main HADD algorithm, we need two additional auxiliary procedures. The first one is previously displayed Algorithm 1, by which we calculate the second vector direction, $d_k$ . The second procedure is the Backtracking line search algorithm for calculating the iterative step size value.

Algorithm 3 describes the main algorithm, termed the HADD algorithm.

Algorithm 2: The Backtracking line search procedure.

Require: Objective function $f(x)$ , the direction $d_k$ of the search at the point $x_k$ , and numbers $\mathrm{0}<\sigma <\mathrm{0.5}$ and $\beta \in (\mathrm{0,1})$ .

1: $t=\mathrm{1}$ .

2: While $f(x_k+t{d}_k)>f(x_k)+\sigma t{g}_k^T{d}_k$ , take $t≔t\beta$ .

3: Return $t_k=t$ .

Algorithm 3: HADD method.

Require: $\mathrm{0}<\rho <\mathrm{1}$ , $\mathrm{0}<\tau <\mathrm{1}$ , $x_{\mathrm{0}}\in {\mathbb{R}}^n$ and $\alpha \in (\mathrm{1,2})$ .

1: Set $k=\mathrm{0}.$ For given $x_{\mathrm{0}}$ , take ${\gamma }_{\mathrm{0}}=\mathrm{1}$ and calculate $f(x_{\mathrm{0}})$ and $g(x_{\mathrm{0}})$ .

2: If $‖g_k‖<ϵ$ , then return $x_k$ , $f(x_k)$ and stop the algorithm, else continue by the next step.

3: Compute the iterative step length, $t_k$ , by Algorithm 2.

4: Compute second vector direction, $d_k$ , by the rule Second direction described in Algorithm 1.

5: Compute $x_{k+\mathrm{1}}$ using (20), then calculate $f(x_{k+\mathrm{1}})$ and $g_{k+\mathrm{1}}.$

6: Calculate the approximation parameter ${\gamma }_{k+\mathrm{1}}$ using (24).

7: If ${\gamma }_{k+\mathrm{1}}<\mathrm{0}$ take ${\gamma }_{k+\mathrm{1}}=\mathrm{1}$ .

8: Set $k≔k+\mathrm{1}$ and go to Step 2.

9: Return $x_{k+\mathrm{1}}$ and $f(x_{k+\mathrm{1}}).$

4. Convergence of the HADD Method

The convergence properties of the established HADD iterative method are considered on the set of uniformly convex and strictly convex quadratic functions. In the case of uniformly convex functions the statements are the same as exposed in [1, 4]. For that reason, we just restate the following lemma, in which decreasing of the objective function in two successive points is estimated with respect to the HADD scheme. Thereupon, the upcoming theorem confirms linear convergence of our hybrid accelerated model.

Lemma 2.

Suppose the function $f$ is twice continuously differentiable and uniformly convex on ${\mathbb{R}}^n$ . With that, let the sequence $\{x_k\}$ be generated by Algorithm 3. Then the next estimation is true $$\begin{array}{}\text{(26)}& f\left(x_k\right)-f\left(x_{k+\mathrm{1}}\right)\ge \mu {‖g_k‖}^{\mathrm{2}},\end{array}$$ where $$\begin{array}{}\text{(27)}& \mu =\left\{\frac{\sigma }{M},\frac{\sigma \left(\mathrm{1}-\sigma \right)}{L}\beta \right\}.\end{array}$$

Theorem 3.

For the twice continuously differentiable and uniformly convex function $f$ on ${\mathbb{R}}^n$ and the sequence $\{x_k\}$ generated by Algorithm 3, the following holds: $$\begin{array}{}\text{(28)}& \underset{k\to \infty }{\mathrm{l}\mathrm{i}\mathrm{m}}‖g_k‖=\mathrm{0}.\end{array}$$ Therewith, the sequence $\{x_k\}$ converges to the optimal solution at least linearly.

We show now that the iteration (20) is convergent regarding the set of strictly convex quadratic functions $$\begin{array}{}\text{(29)}& f\left(x\right)=\frac{\mathrm{1}}{\mathrm{2}}{x}^{\mathrm{T}}Ax-b^{\mathrm{T}}x.\end{array}$$ In (29), it is assumed that $A$ is a real $n\times n$ symmetric positive definite matrix and that vector $b\in {\mathbb{R}}^n$ is given. The smallest and the largest eigenvalues of the matrix $A$ , respectively, are denoted by ${\lambda }_{\mathrm{1}}$ and ${\lambda }_n$ .

Lemma 4.

Let $f$ be the strictly convex quadratic function defined by (29), where $A\in {\mathbb{R}}^{n\times n}$ is a symmetric positive definite matrix. Let ${\lambda }_{\mathrm{1}}$ and ${\lambda }_n$ be the smallest and the largest eigenvalues of $A$ . Then, the following inequalities are valid for the hybrid accelerated gradient model (20): $$\begin{array}{}\text{(30)}& {\lambda }_{\mathrm{1}}\le \frac{{\gamma }_{k+\mathrm{1}}}{t_{k+\mathrm{1}}}\le \frac{\mathrm{4}{\lambda }_n}{\sigma },\hspace{1em}k\in \mathbb{N},\end{array}$$

Proof.

Let us calculate the difference in two successive iterative points of the goal function (29): $$\begin{array}{}\text{(31)}& f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)=\frac{\mathrm{1}}{\mathrm{2}}{x}_{k+\mathrm{1}}^{\mathrm{T}}A{x}_{k+\mathrm{1}}-b^{\mathrm{T}}{x}_{k+\mathrm{1}}-\frac{\mathrm{1}}{\mathrm{2}}{x}_k^{\mathrm{T}}A{x}_k+b^{\mathrm{T}}{x}_k.\end{array}$$

Including the iteration (20) we continue computations: $$\begin{array}{}\text{(32)}& f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)=\frac{\mathrm{1}}{\mathrm{2}}{\left(x_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\alpha t_k^{\mathrm{2}}{d}_k\right)}^{\mathrm{T}}A\left(x_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\alpha t_k^{\mathrm{2}}{d}_k\right)-b^{\mathrm{T}}\left(x_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\alpha t_k^{\mathrm{2}}{d}_k\right)-\frac{\mathrm{1}}{\mathrm{2}}{x}_k^{\mathrm{T}}A{x}_k+b^{\mathrm{T}}{x}_k=\frac{\mathrm{1}}{\mathrm{2}}{x}_k^{\mathrm{T}}A{x}_k-\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A{x}_k+\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k^{\mathrm{2}}{d}_k^{\mathrm{T}}A{x}_k-\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A{x}_k+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{2}}{g}_k^{\mathrm{T}}A{g}_k-\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k^{\mathrm{3}}{\gamma }_k^{-\mathrm{1}}{d}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k^{\mathrm{2}}{d}_k^{\mathrm{T}}A{x}_k-\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{\gamma }_k^{-\mathrm{1}}{d}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{4}}{d}_k^{\mathrm{T}}A{d}_k-b^{\mathrm{T}}{x}_k+\alpha t_k{\gamma }_k^{-\mathrm{1}}{b}^{\mathrm{T}}{g}_k-\alpha t_k^{\mathrm{2}}{b}^{\mathrm{T}}{d}_k-\frac{\mathrm{1}}{\mathrm{2}}{x}_k^{\mathrm{T}}A{x}_k+b^{\mathrm{T}}{x}_k.\end{array}$$ Applying the equality $g_k=A{x}_k-b$ and the symmetry property of $A$ , we get $$\begin{array}{}\text{(33)}& f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)=-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A{x}_k+\alpha t_k^{\mathrm{2}}{d}_k^{\mathrm{T}}A{x}_k-{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{\gamma }_k^{-\mathrm{1}}{d}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{2}}{g}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{4}}{d}_k^{\mathrm{T}}A{d}_k+\alpha t_k{\gamma }_k^{-\mathrm{1}}{b}^{\mathrm{T}}{g}_k-\alpha t_k^{\mathrm{2}}{b}^{\mathrm{T}}{d}_k.\end{array}$$

The right hand side of the previous expression can be further transformed as follows: $$\begin{array}{}\text{(34)}& f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)=\alpha t_k{\gamma }_k^{-\mathrm{1}}\left(b^{\mathrm{T}}{g}_k-x_k^{\mathrm{T}}A{g}_k\right)+\alpha t_k^{\mathrm{2}}\left(x_k^{\mathrm{T}}A{d}_k-b^{\mathrm{T}}{d}_k\right)-{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{\gamma }_k^{-\mathrm{1}}{d}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}\alpha t_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{2}}{g}_k^{\mathrm{T}}A{g}_k+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{4}}{d}_k^{\mathrm{T}}A{d}_k=\alpha t_k{g}_k^{\mathrm{T}}\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)-\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{d}_k^{\mathrm{T}}A\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)=\left(\alpha t_k{g}_k^{\mathrm{T}}-\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A+\frac{\mathrm{1}}{\mathrm{2}}{\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{d}_k^{T}A\right)\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right).\end{array}$$ The replacement of $f(x_{k+\mathrm{1}})-f(x_k)$ by the right hand side of (34) into (24) leads us to $$\begin{array}{}\text{(35)}& {\gamma }_{k+\mathrm{1}}=\mathrm{2}\frac{\left(\alpha t_k{g}_k^{\mathrm{T}}-\left(\mathrm{1}/\mathrm{2}\right){\alpha }^{\mathrm{2}}{t}_k^{\mathrm{2}}{\gamma }_k^{-\mathrm{1}}{g}_k^{\mathrm{T}}A+\left(\mathrm{1}/\mathrm{2}\right){\alpha }^{\mathrm{2}}{t}_k^{\mathrm{3}}{d}_k^{T}A\right)\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)-\alpha t_k{g}_{\mathrm{k}}^T\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}{{\alpha }^{\mathrm{2}}{\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}\left(t_k^{\mathrm{2}}{d}_k-t_k{\gamma }_k^{-\mathrm{1}}{g}_k\right)}.\end{array}$$

After some calculations, we obtain $$\begin{array}{}\text{(36)}& {\gamma }_{k+\mathrm{1}}=\frac{{\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}A\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}{{\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}^{\mathrm{T}}\left(t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k\right)}.\end{array}$$ Previous expression confirms that ${\gamma }_{k+\mathrm{1}}$ is the Rayleigh quotient of the real symmetric matrix $A$ at the vector $t_k{d}_k-{\gamma }_k^{-\mathrm{1}}{g}_k$ , which leads us to the conclusion $$\begin{array}{}\text{(37)}& {\lambda }_{\mathrm{1}}\le {\gamma }_{k+\mathrm{1}}\le {\lambda }_n,\hspace{1em}k\in \mathbb{N}.\end{array}$$ The left hand side in inequalities (30) arrives from the fact $\mathrm{0}\le t_{k+\mathrm{1}}\le \mathrm{1}$ . To prove the right hand side of (30), we use the estimation [[3], eq. (3.8)]: $$\begin{array}{}\text{(38)}& t_k>\frac{\beta \left(\mathrm{1}-\sigma \right){\gamma }_k}{\alpha L}.\end{array}$$ Previous inequality implies $$\begin{array}{}\text{(39)}& \frac{{\gamma }_{k+\mathrm{1}}}{t_{k+\mathrm{1}}}<\frac{L\alpha }{\beta \left(\mathrm{1}-\sigma \right)}.\end{array}$$ We can approximate the Lipschitz constant $L$ by the largest eigenvalue ${\lambda }_n$ and use the fact that $\alpha \le \mathrm{2},\hspace{0.17em}\hspace{0.17em}\mathrm{0}<\sigma <\mathrm{0.5}$ and $\beta \in (\sigma ,\mathrm{1})$ . Then (39) is restated to $$\begin{array}{}\text{(40)}& \frac{{\gamma }_{k+\mathrm{1}}}{t_{k+\mathrm{1}}}<\frac{L\alpha }{\beta \left(\mathrm{1}-\sigma \right)}<\frac{\mathrm{2}{\lambda }_n}{\sigma ·\mathrm{0.5}}=\frac{\mathrm{4}{\lambda }_n}{\sigma }.\end{array}$$ Estimation of the Lipschitz constant $L$ by the largest eigenvalue ${\lambda }_n$ is certainly valid since the matrix $A$ is symmetric and $g(x)=Ax-b.$ From these two facts we conclude that $$\begin{array}{}\text{(41)}& ‖g\left(x\right)-g\left(y\right)‖=‖Ax-Ay‖=‖A\left(x-y\right)‖\le ‖A‖‖x-y‖={\lambda }_n‖x-y‖,\end{array}$$ which completes the proof.

Theorem 5.

Let the iterations (19) be applied on strictly convex quadratic function $f$ given by the expression (29). Suppose that the condition $$\begin{array}{}\text{(42)}& {\lambda }_n<\frac{\mathrm{2}{\lambda }_{\mathrm{1}}}{\alpha }\end{array}$$ holds for the largest and the smallest eigenvalues of symmetric positive definite matrix $A$ . Then, the following estimations are true: $$\begin{array}{c}\text{(43)}& {\left(p_i^{k+\mathrm{1}}\right)}^{\mathrm{2}}\le {\delta }^{\mathrm{2}}{\left(p_i^k\right)}^{\mathrm{2}},{\left(q_i^{k+\mathrm{1}}\right)}^{\mathrm{2}}\le {\lambda }_n^{\mathrm{2}}{\left(q_i^k\right)}^{\mathrm{2}},\end{array}$$ where $$\begin{array}{}\text{(44)}& \delta =\max \left\{\mathrm{1}-\frac{\sigma {\lambda }_{\mathrm{1}}}{\mathrm{2}{\lambda }_n},\frac{{\lambda }_n}{{\lambda }_{\mathrm{1}}}-\mathrm{1}\right\}\end{array}$$ and $p_i^k,q_i^k\in \mathbb{R},k,i,n\in \mathbb{N}.$ Therewith $$\begin{array}{}\text{(45)}& \underset{k\to \infty }{\mathrm{l}\mathrm{i}\mathrm{m}}‖g_k‖=\mathrm{0}.\end{array}$$

Proof.

Let us consider the orthonormal system of eigenvectors $\{v_{\mathrm{1}},v_{\mathrm{2}},\dots ,v_n\}$ of the matrix $A$ . Thereon, we construct the sequence of values $\{x_k\}$ by applying Algorithm 3 on strictly convex quadratic function $f$ defined by (29). Then, for some $k\in \mathbb{N}$ and for some constants $p_{\mathrm{1}}^k,p_{\mathrm{2}}^k,\dots ,p_n^k,q_{\mathrm{1}}^k,q_{\mathrm{2}}^k,\dots ,q_n^k\in \mathbb{R}$ it follows that $$\begin{array}{c}\text{(46)}& g_k=A{x}_k-b={\displaystyle \sum }_{i=\mathrm{1}}^n{p}_i^k{v}_i,d_k={\displaystyle \sum }_{i=\mathrm{1}}^n{q}_i^k{v}_i.\end{array}$$

Applying (20), further we conclude that $$\begin{array}{}\text{(47)}& g_{k+\mathrm{1}}=A\left(x_k-\alpha t_k{\gamma }_k^{-\mathrm{1}}{g}_k+\alpha t_k^{\mathrm{2}}{d}_k\right)-b=\left(I-\alpha t_k{\gamma }_k^{-\mathrm{1}}A\right)g_k+\alpha t_k^{\mathrm{2}}A{d}_k.\end{array}$$ Having in mind the representation (46), one can verify $$\begin{array}{}\text{(48)}& g_{k+\mathrm{1}}={\displaystyle \sum }_{i=\mathrm{1}}^n{p}_i^{k+\mathrm{1}}{v}_i={\displaystyle \sum }_{i=\mathrm{1}}^n\left(\left(\mathrm{1}-\alpha t_k{\gamma }_k^{-\mathrm{1}}{\lambda }_i\right)p_i^k+\alpha t_k^{\mathrm{2}}{\lambda }_i{q}_i^k\right)v_i={\displaystyle \sum }_{i=\mathrm{1}}^n\left(\mathrm{1}-\alpha t_k{\gamma }_k^{-\mathrm{1}}{\lambda }_i\right)p_i^k{v}_i+\alpha t_k^{\mathrm{2}}{\displaystyle \sum }_{i=\mathrm{1}}^n{\lambda }_i{q}_i^k{v}_i.\end{array}$$ To prove $$\begin{array}{}\text{(49)}& {\left(p_i^{k+\mathrm{1}}\right)}^{\mathrm{2}}\le {\delta }^{\mathrm{2}}{\left(p_i^k\right)}^{\mathrm{2}},{\left(q_i^{k+\mathrm{1}}\right)}^{\mathrm{2}}\le {\lambda }_n^{\mathrm{2}}{\left(q_i^k\right)}^{\mathrm{2}},\end{array}$$ i.e., the inequalities (43), we only need to verify that $$\begin{array}{}\text{(50)}& \left|\mathrm{1}-\alpha \frac{{\lambda }_i}{{\gamma }_k{t}_k^{-\mathrm{1}}}\right|\le \delta \end{array}$$ since $\left|{\lambda }_i\right|\le {\lambda }_n$ for all $i\in \{\mathrm{1,2},\dots ,n\}$ . There are two possibilities:

(1) $\alpha {\lambda }_i\le t_k^{-\mathrm{1}}{\gamma }_k.$

This case implies the next set of inequalities: $$\begin{array}{}\text{(51)}& \mathrm{1}>\alpha \frac{{\lambda }_i}{{\gamma }_k{t}_k^{-\mathrm{1}}}=\frac{\alpha {\lambda }_i}{{\gamma }_k}·t_k\ge \frac{\alpha {\lambda }_{\mathrm{1}}}{{\gamma }_k}\frac{\beta \left(\mathrm{1}-\sigma \right){\gamma }_k}{L\alpha }=\frac{{\lambda }_{\mathrm{1}}\beta \left(\mathrm{1}-\sigma \right)}{L}>\frac{{\lambda }_{\mathrm{1}}\sigma ·\mathrm{0.5}}{{\lambda }_n}=\frac{{\lambda }_{\mathrm{1}}}{\mathrm{2}{\lambda }_n}\sigma .\end{array}$$

As a consequence, we can conclude that $$\begin{array}{}\text{(52)}& \mathrm{1}-\alpha \frac{{\lambda }_i}{{\gamma }_k{t}_k^{-\mathrm{1}}}\le \mathrm{1}-\frac{{\lambda }_{\mathrm{1}}}{\mathrm{2}{\lambda }_n}\sigma \le \delta .\end{array}$$

(2) $\alpha {\lambda }_i\ge t_k^{-\mathrm{1}}{\gamma }_k.$

In this case, one can verify the following estimations: $$\begin{array}{}\text{(53)}& \mathrm{1}<\frac{\alpha {\lambda }_i}{{\gamma }_k{t}_k^{-\mathrm{1}}}\le \frac{\alpha {\lambda }_n}{{\lambda }_{\mathrm{1}}}\Rightarrow \left|\mathrm{1}-\alpha \frac{{\lambda }_i}{{\gamma }_k{t}_k^{-\mathrm{1}}}\right|\le \frac{\alpha {\lambda }_n}{{\lambda }_{\mathrm{1}}}\le \delta .\end{array}$$ The representation (46) and the fact that $\{v_{\mathrm{1}},v_{\mathrm{2}},\dots ,v_n\}$ is an orthonormal system of eigenvectors lead to the next conclusion $$\begin{array}{}\text{(54)}& {‖g_k‖}^{\mathrm{2}}={\displaystyle \sum }_{i=\mathrm{1}}^n{\left(d_i^k\right)}^{\mathrm{2}}.\end{array}$$ Now, knowing that the parameter $\delta$ under condition ${\lambda }_n<\mathrm{2}{\lambda }_{\mathrm{1}}/\alpha$ satisfies $\mathrm{0}<\delta <\mathrm{1}$ , we confirm that the final statement is true. $$\begin{array}{}\text{(55)}& \underset{k\to \infty }{\mathrm{l}\mathrm{i}\mathrm{m}}‖g_k‖=\mathrm{0}\end{array}$$

Remark 6.

The assumption (42) used in the previous theorem is required in order to prove that the HADD process is convergent for the strictly convex quadratics. Therewith, knowing that the hybrid parameter $\alpha \in (\mathrm{1,2})$ implies ${\lambda }_n/{\lambda }_{\mathrm{1}}$ $\in (\mathrm{1,2})$ points to the conclusion that Theorem 5 is applicable to very few cases. However, this is not entirely so since we choose only one particular value $\alpha \in (\mathrm{1,2})$ for the practical computations. Regarding this matter, the authors in [3] numerically confirmed that the optimal value of the hybrid parameter $\alpha$ is the one close to the left limit of the interval $(\mathrm{1,2})$ , i.e., the value which is very close to 1. Therefore, we choose $\alpha =\mathrm{1.1}$ for numerical tests displayed in the next section. Choosing the similar values for hybrid parameter $\alpha$ , the condition (42) becomes very close to the condition ${\lambda }_n<\mathrm{2}{\lambda }_{\mathrm{1}}$ , used in [12], under which Q-linear convergence rate of the preconditioned BB method was established.

5. Computational Tests and Comparisons

The performance of the C $++$ implementation of derived HADD model is investigated on a set of 630 test unconstrained optimization problems picked from [13]. We conduct the testings on a Workstation Intel Celeron $\mathrm{1.6}$ GHz. The following stopping criteria are used: $$\begin{array}{c}\text{(56)}& ‖g_k‖\le \mathrm{1}{\mathrm{0}}^{-\mathrm{6}},\frac{\left|f\left(x_{k+\mathrm{1}}\right)-f\left(x_k\right)\right|}{\mathrm{1}+\left|f\left(x_k\right)\right|}\le \mathrm{1}{\mathrm{0}}^{-\mathrm{16}}.\end{array}$$ The values of the Backtracking parameters are set up as follows: $\sigma =\mathrm{0.000}\mathrm{1}$ and $\beta =\mathrm{0.8}$ .

We compare the hybrid accelerated HADD method with its forerunner ADD scheme, as well as with the hybrid accelerated HSM method. The number of function evaluations is the performance profile measured in all tests. The dominance of the ADD method regarding the number of iterations among the other comparative models was confirmed in [4]. However, from that research we do not have any information about the behavior of the ADD method when the number of function evaluations is involved. With respect to this parameter, the HSM scheme upgrades the accelerated SM method as well as Nesterov’s line search algorithm; see [3]. For these reasons, our experimental goal is to numerically prove better performance feature of the HADD method, considering the number of function evaluations, when compared with the ADD and the HSM method.

In Table 1, we display the number of problems, out of 630, for which an algorithm achieved the minimum number of function evaluations. In the same table, we also display the total number of problems for which all three algorithms achieved an equal number of function evaluations. Based on the results displayed in this table, it is obvious that the HADD scheme convincingly outperforms the other two comparative models.

Table 1

Comparison between the HADD, ADD, and HSM methods regarding the minimal number of function evaluations.

Comparative methods	HADD	ADD	HSM	=
Number of function evaluations	142	0	18	50

For more clear visualization of the performance of the HADD algorithm versus the ADD and the HSM algorithms, we display in Figure 1 the Dolan-Moré’s performance profile subject to the number of function evaluations metric. As we can see, the HADD scheme is more robust and therewith more efficient than the other two methods.

[figure omitted; refer to PDF]

Obtained numerical results confirm that applied hybridization process is a good way to improve some important characteristics of chosen accelerated methods. Preferable outcomes of the HADD scheme, regarding analyzed characteristic, come from the properly chosen hybrid value $\alpha$ , together with derived accelerated parameter ${\gamma }^{HADD}$ . Good convergent properties of defined HADD process can be a reason for applying proposed hybridization on some other gradient and accelerated gradient models.

6. Conclusion

We present a hybrid accelerated double direction gradient method for solving unconstrained optimization problems. The HADD method is derived by applying good properties of the hybrid representation introduced in [3] in conjunction with the form of double direction optimization model with accelerated parameter presented in [4]. The convergence of defined optimization model is provided on the set of uniformly convex and strictly convex quadratic functions.

The HADD scheme reserves preferable features of both forerunner methods. Therewith, according to conducted numerical experiments, it outperforms the ADD and HSM methods regarding the requested number of function evaluations. We evaluated the Dolan-Moré performance profiles of comparative methods and showed that the HADD iteration is the most efficient compared to the other two algorithms.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

The first and the third authors acknowledge support from the internal research project IS01-17 supported by the Faculty of Sciences and Mathematics, University of Priština, Serbia. The second author gratefully acknowledges support from the Project supported by Ministry of Education and Science of Republic of Serbia, Grant No. 174013.

References

[1] P. S. Stanimirović, M. B. Miladinović, "Accelerated gradient descent methods with line search," Numerical Algorithms, vol. 54 no. 4, pp. 503-520, DOI: 10.1007/s11075-009-9350-8, 2010.

[2] M. J. Petrović, "An accelerated double step size model in unconstrained optimization," Applied Mathematics and Computation, vol. 250, pp. 309-319, DOI: 10.1016/j.amc.2014.10.104, 2015.

[3] M. J. Petrovic', V. Rakocevic', N. Kontrec, S. Panic', D. Ilic', "Hybridization of accelerated gradient descent method," Numer. Algor, 2018.

[4] M. J. Petrović, P. S. Stanimirović, "Accelerated Double Direction Method for Solving Unconstrained Optimization Problems," Mathematical Problems in Engineering, vol. 2014,DOI: 10.1155/2014/965104, 2014.

[5] P. S. Stanimirović, G. V. Milovanović, M. J. Petrović, "A Transformation of Accelerated Double Step Size Method for Unconstrained Optimization," Mathematical Problems in Engineering, vol. 2015,DOI: 10.1155/2015/283679, 2015.

[6] S. Ishikawa, "Fixed points by a new iteration method," Proceedings of the American Mathematical Society, vol. 44, pp. 147-150, DOI: 10.1090/S0002-9939-1974-0336469-5, 1974.

[7] W. R. Mann, "Mean value methods in iteration," Proceedings of the American Mathematical Society, vol. 4, pp. 506-510, DOI: 10.1090/S0002-9939-1953-0054846-3, 1953.

[8] E. Picard, "Memoire sur la theorie des equations aux derivees partielles et la methode des approximations successives," Journal de Mathématiques Pures et Appliquées, vol. 6, pp. 145-210, 1890.

[9] S. H. Khan, "A Picard-Mann hybrid iterative process," Fixed Point Theory and Applications, vol. 2013, article 69,DOI: 10.1186/1687-1812-2013-69, 2013.

[10] N. I. Djuranović-Miličić, M. Gardaševic-Filipović, "A multi-step curve search algorithm in nonlinear optimization: nondifferentiable convex case," Facta Universitatis, Series: Mathematics and Informatics, vol. 25, pp. 11-24, 2010.

[11] A. Kumar, D. K. Gupta, E. Martinez, S. Singh, "Directional k -step Newton methods in n variables and its semilocal convergence analysis," Mediterranean Journal of Mathematics, vol. 15 no. 2,DOI: 10.1007/s00009-018-1077-0, 2018.

[12] B. Molina, M. Raydan, "Preconditioned Barzilai-BORwein method for the numerical solution of partial differential equations," Numerical Algorithms, vol. 13 no. 1-2, pp. 45-60, DOI: 10.1007/bf02143126, 1996.

[13] N. Andrei, "An unconstrained optimization test functions collection," Advanced Modeling and Optimization, vol. 10 no. 1, pp. 147-161, 2008.