A Trust-Region Method for Nonsmooth Nonconvex Optimization

Different types of nonsmooth trust-region methods have already been proposed and analyzed for the general optimization problem (1.2) throughout the last two decades. Several of these nonsmooth trust-region methods utilize abstract model functions on a theoretical level which means that the model function is typically not specified. In [26], a nonsmooth trust-region method is proposed for (1.2) under the assumption that $\psi$ is regular. A nonsmooth trust-region algorithm for general problems is investigated in [64]. In this work, an abstract first-order model is considered that is not necessarily based on subgradient information or directional derivatives. Extending the results in [64], Grohs and Hosseini propose a Riemannian trust-region algorithm, see [34]. Here, the objective function is defined on a complete Riemannian manifold. All mentioned methods derive global convergence under an assumption similar to the concept of a “strict model” stated in [59]. Using this concept, a nonsmooth bundle trust-region algorithm with global convergence is constructed in [4].

In [20], a hybrid approach is presented using simpler and more tractable quadratic model functions. The method switches to a complicated second model if the quadratic model is not accurate enough and if it is strictly necessary. In [3], a quadratic model function is analyzed where the first-order term is derived from a suitable approximation of the steepest descent direction and the second-order term is updated utilizing a BFGS scheme. The authors apply an algorithmic approach proposed in [49] to compute the approximation of the $\epsilon$-subdifferential and steepest descent direction. Another class of methods employs smoothing techniques. In [32], the authors first present a smooth trust-region method without using derivatives, and then, in the nonsmooth case, use this methodology after smoothing the objective function. Furthermore, trust-region algorithms for nonsmooth problems can be developed based on smooth merit functions for the problem. In [66], a nonsmooth convex optimization is investigated and the Moreau envelope is considered as a smooth merit function. A smooth trust-region method is performed on the smooth merit function, where the second-order term of the model function is again updated by the BFGS formula.

Bundle methods are an important and related class of methods for nonsmooth problems [44, 52, 51, 39, 45, 42, 40]. The ubiquitous cutting-plane model in bundle methods is polyhedral, i.e., the supremum of a finite affine family. This model builds approximations of convex functions based on the subgradient inequality. In [67], an efficient bundle technique for convex optimization has been proposed; in [24], a convex bundle method is derived to deal with additional noise, i.e., the case when the objective function and the subgradient can not be evaluated exactly. Different modifications of the bundle ideas for nonconvex problems have been established in [59, 67]. In [35] and [60], the authors consider bundle methods for nonsmooth nonconvex optimization when the function values and the subgradients of $\psi$ can only be evaluated inexactly.

Local convergence properties and rates for nonsmooth problems are typically studied utilizing additional and more subtle structures. In this regard, some fundamental and helpful concepts are the idea of an “active manifold” and the family of “partly smooth” functions introduced by [46]. In particular, the problem (1.1) has been investigated when the nonsmooth term $\varphi$ is partly smooth relative to a smooth active manifold. The so-called finite activity identification is established for forward-backward splitting methods by [48] and, more recently, for SAGA/Prox-SVRG by [63]. After the identification, those algorithms enter a locally linear convergence regime. In [36], the authors use partial smoothness and prox-regularity to identify the active constraints after finitely many iterations, which is an extension of other works on finite constraint identification, see [17, 15, 74]. After identifying the active manifold, the nonsmooth problem may become a smooth optimization on a Riemannian manifold. Some algorithms and analysis for Riemannian trust-region methods were studied in the literature; see e.g., [2, 1, 38, 8, 6]. There are numerous applications of Riemannian trust-region methods, such as eigenproblems [6, 7], low-rank matrix completion [13], and tensor problems [14].

We note that for composite programs, nonsmooth trust-region methods have also been studied in the literature, such as [28, 29, 79, 76, 33] for $\min\ g(x)+h(f(x))$ where $f$ and $g$ are smooth and $h$ is convex, [78, 9] for $\min\ h(f(x))$ where $f$ is smooth and $h$ is convex, [25] for $\min\ h(f(x))$ where $f$ is locally Lipschitz and $h$ is smooth and convex and [16] for $\min\ g(x)+h(f(x))$ where $f$ is smooth and $g$ and $h$ are convex. For the problem (1.1), there are also other efficient methods, such as gradient-type methods [30, 57], semismooth Newton methods, [55, 47, 75], proximal Newton methods [61, 43, 72], or forward-backward envelope-based (quasi-)Newton methods [62, 69].

More developments about trust-region methods may be found in review papers such as [80].

In this work, we propose and investigate a trust-region method for nonsmooth composite programs. The approach utilizes quadratic model functions to approximate the underlying nonsmooth objective function. This methodology leads to classical and tractable trust-region subproblems that can be solved efficiently by standard optimization methods if the second-order information is symmetric. We also discuss an efficient subproblem solver in the case that the second-order information does not stem from a symmetric matrix. The linear part of our proposed quadratic model can be based on the steepest descent direction or other directions such as proximal gradient-type descent directions. Our algorithm contains the following steps: computation of a model function and (approximate) solution of the associated trust-region subproblem; acceptance test of the calculated step; determination of a suitable stepsize by a cheap method followed by some stepsize safeguards and a second acceptance test (if the first test is not successful); update of the trust-region radius and a modification step via a novel truncation mechanism.

In order to control the approximation error between the quadratic model function and the nonsmooth objective function, we define a stepsize safeguard strategy that tries to avoid points along a specific direction at which the directional derivative is not continuous. Specifically, this strategy tries to guarantee that the objective function is directionally differentiable along a specific direction. Since a direct implementation of such a strategy can yield arbitrarily small stepsizes, we consider functions which can be truncated and propose an additional truncation step that allows to enlarge the stepsize. This modification is an essential and new part in our global convergence theory. We verify that the family of functions that can be truncated is rich and contains many important examples, such as the $\ell_{1}$-norm, $\ell_{\infty}$-norm or group sparse-type penalty terms. Moreover, we provide a detailed global convergence analysis of the proposed trust-region framework. In particular, we show that every accumulation point of a sequence generated by our algorithm is a stationary point. Global convergence of nonsmooth trust-region methods typically requires a certain uniform accuracy assumption on the model which coincides with the concept of the already mentioned strict model proposed by [59]. Our assumptions are similar to these standard requirements and can be verified for a large family containing polyhedral problems and group lasso. Furthermore, we also show how a strict model – aside from utilizing the original objective function – can be constructed.

We analyze the local properties of the nonsmooth trust-region method for (1.1) when $\varphi$ is a partly smooth function. In particular, it is possible to establish quadratic convergence of our approach in this case. We assume that the underlying manifold is an affine subspace and that a strict complementary condition holds. After the finite activity identification, we transfer our problem to a smooth problem in the affine subspace by proving that an appropriate choice of the first-order and the second-order model coincides with the Riemannian framework. Results from Riemannian trust-region theory can then be applied to derive local quadratic convergence. Additionally, if the nonsmooth term is polyhedral, it can be shown that the Riemannian Hessian can be computed without knowing the underlying manifold.

The rest of this paper is organized as follows. In Section 2, descent directions and several properties of $\psi$ and $\varphi$ used in our algorithm are discussed. In Section 3, we present the nonsmooth trust-region framework. In Section 4, the global convergence of our method is established. In Section 5, we show fast local convergence by studying the nonsmooth composite program for partly smooth $\varphi$. Some preliminary numerical experiments are presented in Section 6.

In this work, the expression $\|\cdot\|$ denotes the $\ell_{2}$-norm and the Frobenius norm for vectors and matrices, respectively. For $x\in\mathbb{R}^{n}$ and $r>0$, $B_{r}(x):=\{y\in\mathbb{R}^{n}:\|y-x\|<r\}$ denotes the open ball with radius $r$ around $x$. Let $\Lambda$ be a given symmetric positive definite matrix. The proximal operator is defined via

$$\mathrm{prox}_{\varphi}^{\Lambda}(z)=\mathop{\operatorname*{argmin}}_{y\in\mathbb{R}^{n}}\varphi(y)+\frac{1}{2}\|y-z\|_{\Lambda}^{2},$$

where $\|x\|_{\Lambda}^{2}:=x^{T}\Lambda x$. We also slightly abuse the notation and write $\mathrm{prox}_{\varphi}^{\lambda}=\mathrm{prox}_{\varphi}^{\Lambda}$ for $\Lambda=\lambda I$.

The directional derivative of a function $h:\mathbb{R}^{n}\rightarrow\mathbb{R}$ at $x$ along $d$ is denoted by

$$h^{\prime}(x;d):=\lim_{t\rightarrow 0^{+}}\frac{h(x+td)-h(x)}{t},$$

if it exists. In the composite case $\psi=f+\varphi$ with smooth $f$ and real-valued convex $\varphi$, the directional derivative $\psi^{\prime}(x;d)$ is well-defined for all $x,d\in\mathbb{R}^{n}$ and the subdifferential of $\psi$ is defined via

$$\partial\psi(x)=\{\nabla f(x)\}+\partial\varphi(x)=\{a\in\mathbb{R}^{n}\mid\langle a,d\rangle\leq\psi^{\prime}(x;d),\ \forall\ d\in\mathbb{R}^{n}\},$$

where $\partial\varphi(x)$ is the usual subdifferential of a convex function. The steepest descent direction of $\psi$ is defined as

$$d_{s}(x):=\begin{cases}\mathop{\operatorname*{argmin}}_{\|d\|\leq 1}\psi^{\prime}(x;d)&\text{if }0\notin\partial\psi(x),\\ 0&\text{if }0\in\partial\psi(x).\end{cases}$$

In this paper, we will repeatedly work with the following normalization condition for a direction $d(x)$:

$$\|d(x)\|=\begin{cases}1&\text{if }0\notin\partial\psi(x),\\ 0&\text{if }0\in\partial\psi(x).\end{cases}$$

(2.1)

We can see that $d_{s}(x)$ satisfies the property (2.1). We say that a point $x^{*}$ is a stationary point of problem (1.2) if $\psi^{\prime}(x^{*};d_{s}(x^{*}))=0$, i.e., if and only if, $0\in\partial\psi(x^{*})$.

If the objective function $\psi$ is smooth, then the first-order information is carried in the gradient $\nabla\psi$. We notice that the gradient is directly connected to the steepest descent direction $d_{s}(x)=-\nabla\psi(x)/\|\nabla\psi(x)\|$ and thus, we can use the following representation $-\nabla\psi(x)=\psi^{\prime}(x;d_{s}(x))d_{s}(x)$. Motivated by this observation, a natural choice of the first-order information and extension in the nonsmooth case is $g(x)=\psi^{\prime}(x;d_{s}(x))d_{s}(x)$. Since in some cases it might be hard or expensive to directly compute $g(x)=\psi^{\prime}(x;d_{s}(x))d_{s}(x)$, we can also utilize a general descent direction $d(x)$ satisfying (2.1) instead. In our model function, we will work with directions of the form $g(x)=u(x)d(x)$ where $d(x)$ is a descent direction, $\psi^{\prime}(x;d(x))<0$, satisfying (2.1) and $u(x)$ is an upper bound of $\psi^{\prime}(x,d(x))$ with

$$u(x)\begin{cases}\in[\psi^{\prime}(x,d(x)),0)&\text{if }0\notin\partial\psi(x),\\ =0&\text{if }0\in\partial\psi(x).\end{cases}$$

(2.2)

This implies $g(x)=0$ if and only if $0\in\partial\psi(x)$, i.e., if $x$ is a stationary point of (1.2). The direction $g(x)$ plays a similar role as the gradient in the smooth case. We would call $g(x)$ as pseudo-gradient. Our aim in the rest of this subsection is to propose several strategies in the settings of composite programs (1.1) for computing and choosing the functions $u(x)$ and $d(x)$.

If we utilize a smooth model $m_{k}$, once we have selected the descent direction $d$, (directionally) noncontinuous points of $\psi^{\prime}(\cdot,d)$ will contribute to the inaccuracy of $m_{k}$. Hence, we should keep the stepsize relatively small to avoid those points. For any $x,d\in\mathbb{R}^{n}$ with $\|d\|=1$, we set the stepsize safeguard $\Gamma(x,d)$ to guarantee that $t\mapsto\psi^{\prime}(x+td;d)$ is continuous on $(0,\Gamma(x,d))$, which is equivalent to saying that $t\mapsto\psi\left(x+td\right)$ is continuously differentiable on $t\in(0,\Gamma(x,d))$.

We prefer to choose the largest possible value of $\Gamma(x,d)$:

$$\Gamma(x,d)=\Gamma_{\max}(x,d):=\sup\left\{T>0:\tilde{\psi}^{\prime}_{x,d}(t):=\psi^{\prime}(x+td;d)\in C^{1}(0,T)\right\},$$

(2.11)

since it can intuitively lead to faster convergence. We will see that this choice works well for polyhedral problems, where $\varphi$ is the supremum of several affine functions, such as in $\ell_{1}$- and $\ell_{\infty}$-optimization.

However, in some other cases, we may need to set $\Gamma(x,d)$ more carefully. For example, for the group lasso problem $\min_{X\in\mathbb{R}^{n_{1}\times n_{2}}}f(X)+\varphi(X)$, where $f$ is smooth and $\varphi$ is given by $\varphi(X)=\sum_{i=1}^{n_{2}}\left\|X_{i}\right\|$, an appropriate choice for $X=(X_{1},X_{2},\cdots,X_{n_{2}})\in\mathbb{R}^{n_{1}\times n_{2}}$ and $D=(D_{1},D_{2},\cdots,D_{n_{2}})\in\mathbb{R}^{n_{1}\times n_{2}}$ with $\|D\|=1$ is:

$$\Gamma(X,D)=\min\left\{\Gamma_{\max}(X,D),\min_{X_{i}\neq 0}\frac{\|X_{i}\|^{1+\sigma}}{1-\theta_{i}^{2}},\min_{X_{i}\neq 0}\frac{\|X_{i}\|}{\max\{-2\theta_{i},0\}}\right\},\quad\sigma>0.$$

(2.12)

Here, $\theta_{i}$ is given by $\theta_{i}:=\langle X_{i},D_{i}\rangle/(\|X_{i}\|\cdot\|D_{i}\|)$ and we use $c/0:=+\infty$ if $c>0$. The term $\Gamma_{\max}(X,D)$ is defined as in (2.11). This $\Gamma(X,D)$ is specifically designed to overcome some technical difficulties; see, e.g., Lemma 4.4.

Next, let us define the function $\Gamma:\mathbb{R}^{n}\rightarrow\mathbb{R}^{+}$ via

$$\Gamma(x):=\inf_{d\in\mathbb{R}^{n},\ \|d\|=1}\Gamma_{\max}(x,d).$$

(2.13)

The scalar $\Gamma(x)$ is important in our convergence analysis as it provides a lower bound for the stepsize safeguard. For the composite program (1.1) with $\varphi(x)=\|x\|_{1}$, $\Gamma$ can be simply calculated as follows $\Gamma(x)=\min\{|x_{i}|:x_{i}\neq 0\}$ where $\min\emptyset:=+\infty$. Further examples for $\Gamma$ can be found in Appendix A.2.

Since in our algorithmic design we utilize simple, linear-quadratic models to approximate the nonsmooth function $\psi$, we need to introduce stepsize safeguards that allow to intrinsically control the accuracy of the model. However, if the “safeguard” $\Gamma(x)$ is very small, the resulting step might be close to the old iterate and the algorithm can start to stagnate. In order to prevent such an undesirable behavior, we discuss an additional modification step that allows to increase $\Gamma(x)$.

Specifically, given a point $x$, first we want to find a point $x^{\prime}$ near $x$ such that $\Gamma(x^{\prime})$ is relatively large. Let us consider the simplest case where $\psi=f+\varphi$ and $\varphi(x)=\|x\|_{1}$. If $x$ has a nonzero component with small absolute values, then $\Gamma(x)$ is also small. So we can replace those components with 0 and get a new point $x^{\prime}$ satisfying $\Gamma(x^{\prime})>\Gamma(x)$. Since only some components with small absolute values are truncated to 0, the point $x^{\prime}$ is close to $x$. In more general cases, we define a class of functions that allow similar operations:

In Definition 2.5, $\Gamma(T(x,a))\geq a$ means that we can make the value of $\Gamma(\cdot)$ larger by performing truncation and $\|T(x,a)-x\|\leq\kappa a$ implies that the change caused by $T(\cdot,a)$ can be controlled. Example 2.6 shows that $\varphi(x)=\|x\|_{1}$ can be truncated and we present more examples ($\ell_{\infty}$-optimization and group lasso) in the appendix.

Let us mention that, for a smooth regularizer $\varphi$, all properties discussed above are satisfied, since the stepsize safeguard can be chosen as $+\infty$ and no truncation is needed.

Recall that in the classical trust-region method for a smooth optimization problem, $\min_{x\in\mathbb{R}^{n}}\psi(x)$, the model function is $m_{k}(s)=\psi(x^{k})+\langle\nabla\psi(x^{k}),s\rangle+\frac{1}{2}\langle s,B^{k}s\rangle$. As mentioned at the beginning of Section 2.2, a natural extension is

$$m_{k}(s)=\psi(x^{k})+\langle\psi^{\prime}(x^{k};d_{s}(x^{k}))d_{s}(x^{k}),s\rangle+\frac{1}{2}\langle s,B^{k}s\rangle.$$

(3.1)

This model function is still quadratic and fits the objective function well along the steepest descent direction, which means that they have the same directional derivative in this direction. Though approximating the nonsmooth function $\psi$ with a quadratic function might not lead to good trust-region models in general, we can design specific quadratic models that fit $\psi$ well along certain directions.

One may wish to use different descent directions. Or it might be expensive to compute the steepest descent direction $d_{s}(x)$ and its Clarke’s generalized directional derivative $\psi^{\prime}(x;d_{s}(x))$. Therefore we use a descent direction $d(x)$ satisfying (2.1) and $u(x)$ satisfying (2.2) instead. We can now define our model function

$$m_{k}(s)=\psi_{k}+\langle g^{k},s\rangle+\frac{1}{2}\langle s,B^{k}s\rangle,$$

(3.2)

where $\psi_{k}=\psi(x^{k})$ and $g^{k}=g(x^{k})=u(x^{k})d(x^{k})$, and the associated trust-region subproblem is given by

$$\min_{s}\ m_{k}(s)=\psi_{k}+\langle g^{k},s\rangle+\frac{1}{2}\langle s,B^{k}s\rangle\quad\text{s.t.}\quad\|s\|\leq\Delta_{k}.$$

(3.3)

This subproblem is quadratic and coincides with the classical approaches if $B^{k}$ is symmetric. The matrices $B^{k}$ in the model (3.2) are typically chosen to capture or approximate the second-order information of the objective function $\psi$. However, such a careful choice is only required when certain local convergence properties should be guaranteed. Similar to classical trust-region methods, global convergence of the method will generally not be affected by $\{B^{k}\}$ and flexible choices of $B^{k}$ are possible under some mild boundedness conditions that will be specified later.

We remark that the major difference between our algorithm and other nonsmooth trust-region methods in the literature is the first-order information in the model function. The methods in [64, 20] employ first-order terms that tend to be more complicated than a simple linear function, in order to approximate the nonsmooth objective function well and to satisfy certain accuracy assumptions. Although the model function in [3] is quadratic, their first-order term needs to be built using a steepest descent direction based on the $\epsilon$-subdifferential of $\psi$.

An important concept for solving (3.3) is the so-called Cauchy point, which is defined via

$$s^{k}_{C}:=-\alpha_{k}^{C}g^{k}\ \text{and}\ \alpha_{k}^{C}:=\mathop{\operatorname*{argmin}}_{0\leq t\leq\Delta_{k}/\|g^{k}\|}m_{k}(-tg^{k}).$$

The Cauchy point is computational inexpensive [58, Algorithm 4.2] and it leads to sufficient reduction of the model function (Cauchy decrease condition):

$$m_{k}(0)-m_{k}(s^{k}_{C})\geq\frac{1}{2}\|g^{k}\|\min\left\{\Delta_{k},\frac{\|g^{k}\|}{\|B^{k}\|}\right\},$$

(3.4)

see, e.g., in [58, Lemma 4.3]. Furthermore, it can be shown that

$$m_{k}(0)-m_{k}(\alpha_{k}\bar{s}^{k}_{C})\geq\frac{\alpha_{k}}{\|s^{k}_{C}\|}\cdot\frac{1}{2}\|g^{k}\|\min\left\{\Delta_{k},\frac{\|g^{k}\|}{\|B^{k}\|}\right\},$$

(3.5)

where $\bar{s}^{k}_{C}=s^{k}_{C}/\|s^{k}_{C}\|$ and $0<\alpha_{k}\leq\|s^{k}_{C}\|$.

In our algorithm, we need to generate an approximate solution of (3.3) that achieves a similar model descent compared to the Cauchy descent condition (3.4) in some sense. More precisely, we need to recover a solution $s^{k}$ satisfying

$$m_{k}(0)-m_{k}(s^{k})\geq\frac{\gamma_{1}}{2}\|g^{k}\|\min\left\{\Delta_{k},\gamma_{2}\|g^{k}\|\right\},$$

(3.6)

where $\gamma_{1},\gamma_{2}>0$ are constants which do not depend on $k$, and

$$m_{k}(0)-m_{k}(s^{k})\geq(1-\ell(\|s^{k}\|))[m_{k}(0)-m_{k}(s^{k}_{C})],$$

(3.7)

where $\ell:\mathbb{R}^{+}\rightarrow[0,1]$ is chosen as a monotonically decreasing function with $\lim_{\Delta\rightarrow 0^{+}}\ell(\Delta)=0$. The classical choice $\ell(\Delta)\equiv 0$ is also allowed.

In the trust-region framework, we will work with the parameters $0<\eta\leq\eta_{1}<\eta_{2}<1$, $0<r_{1}<1<r_{2}$, and $\Delta_{\max}>0$. Let $s^{k}$ denote the generated solution of (3.3). Similar to the classical trust-region method, we define the ratio between actual reduction and predicted reduction as

$$\rho^{1}_{k}=\frac{\psi(x^{k})-\psi(x^{k}+s^{k})}{m_{k}(0)-m_{k}(s^{k})}.$$

(3.8)

If the proposed step $x^{k}+s^{k}$ is “successful”, i.e., $\rho^{1}_{k}\geq\eta_{1}$, we accept the step, i.e., $\tilde{x}^{k}=x^{k}+s^{k}$, and update the trust-region radius $\Delta_{k}$ as

$$\Delta_{k+1}=\begin{cases}\min\{{\Delta}_{\max},r_{2}\Delta_{k}\}&\text{if }\rho^{1}_{k}>\eta_{2},\\ \Delta_{k}&\text{otherwise}.\end{cases}$$

(3.9)

If $\rho^{1}_{k}<\eta_{1}$, we can introduce an additional stepsize strategy to refine the step. Specifically, we now consider the normalized descent direction $\bar{s}^{k}:={s^{k}}/{\|s^{k}\|}$. In the following, we will always use the notation $\bar{x}:={x}/{\|x\|}$. Instead of setting the stepsize as $\|s^{k}\|$ and working with $s^{k}$ directly, we calculate $\alpha_{k}$ via

$$\alpha_{k}=\min\left\{\Gamma(x^{k},\bar{s}^{k}),\|s^{k}\|\right\},$$

(3.10)

where $\Gamma(x^{k},\bar{s}^{k})$ is the stepsize safeguard. If

$$m_{k}(0)-m_{k}(\alpha_{k}\bar{s}^{k})\geq\frac{\alpha_{k}}{2\|s^{k}\|}(m_{k}(0)-m_{k}(s^{k})),$$

(3.11)

which means that the modified step yields sufficient descent, we use the direction $\bar{s}^{k}$ and the stepsize $\alpha_{k}$; otherwise we set

$$s^{k}=s^{k}_{C},\ \bar{s}^{k}=\bar{s}^{k}_{C},\text{ and }\alpha_{k}=\min\left\{\Gamma(x^{k},\bar{s}^{k}_{C}),\|s^{k}_{C}\|\right\}.$$

(3.12)

If $m_{k}$ is convex (which, e.g., can be ensured when the matrix $B^{k}$ is chosen to be positive semidefinite) then (3.11) holds automatically and the latter case will not occur. In (3.12) we utilize the Cauchy point and the corresponding stepsize as a simpler gradient-based step. As we have seen such a step can always guarantee certain descent properties. Next, we perform a second ratio test

$$\rho^{2}_{k}=\frac{\psi(x^{k})-\psi(x^{k}+\alpha_{k}\bar{s}^{k})}{m_{k}(0)-m_{k}(\alpha_{k}\bar{s}^{k})}.$$

(3.13)

According to this ratio, we update the trust-region radius as

$$\Delta_{k+1}=\begin{cases}r_{1}\Delta_{k}&\text{if }\rho^{2}_{k}<\eta_{1},\\ \min\{{\Delta}_{\max},r_{2}\Delta_{k}\}&\text{if }\rho^{2}_{k}>\eta_{2},\\ \Delta_{k}&\text{otherwise},\end{cases}$$

(3.14)

and decide whether to accept the proposed step

$$\tilde{x}^{k}=\begin{cases}x^{k}+\alpha_{k}\bar{s}^{k}&\text{if }\rho_{k}^{2}\geq\eta,\\ x^{k}&\text{if }\rho_{k}^{2}<\eta.\end{cases}$$

(3.15)

We declare the step as “subsuccessful” if $\rho^{1}_{k}<\eta_{1}$ while $\rho^{2}_{k}\geq\eta$, i.e., even if the original step is unsuccessful, the refined version can still provide some descent which is essential to guarantee convergence.