# A Simple Alternating Direction Method for the Conic Trust Region Subproblem

Authors: Honglan Zhu; Qin Ni
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

A simple alternating direction method is used to solve the conic trust region subproblem of unconstrained optimization. By use of the new method, the subproblem is solved by two steps in a descent direction and its orthogonal direction, the original conic trust domain subproblem into a one-dimensional subproblem and a low-dimensional quadratic model subproblem, both of which are very easy to solve. Then the global convergence of the method under some reasonable conditions is established. Numerical experiment shows that the new method seems simple and effective.


---

## Body

## 1. Introduction

In this paper, we consider the unconstrained optimization problem(1)minx∈Rnfx,where f(x) is continuously differentiable. The trust region method is a very effective method for the unconstrained optimization problem (1) (see [1–6]). Traditional trust region methods are based on a quadratic model and the corresponding quadratic program subproblem is, at the kth iteration,(2)mins∈Rnϱks=gkTs+12sTBks,(3)s.t.s≤Δk,where xk is the current iterate point, gk=∇f(xk), Bk is symmetric and an approximation to the Hessian of f(x), · refers to the Euclidean norm, and Δk is the trust region radius at the kth iteration. There are many methods that can be used to solve the subproblems (2)-(3). The simple, low cost, and effective methods are dogleg methods (see [7–11]). Now, we recall the simple dogleg algorithm for solving trust region subproblem with the quadratic model as in [12].Algorithm 1.
⁡
S tep 0. Input the data of the kth iteration i.e., gk, Bk and Δk.
S tep 1. Compute skN=-Bk-1gk. If skN≤Δk, then s∗=skN, and stop.
S tep 2. Compute skc=-(gkTgk/gkTBkgk)gk. If skc≥Δk, then s∗=-Δkgk/gk, and stop. Otherwise, go to Step 3.
S tep 3. Compute(4)d=skN-skc2,e=skN-skcTskc,f=skc2-Δk2.then s∗=skc+λ(skN-skc), where λ=(-e+e2-df)/d.In 1980, Davidon first proposed the conic model (see [13]). It is an alternative model to substitute the quadratic model. For optimization problems; if the objective function has a strong nonquadratic or its curvature changes severely, then the conic model is better than the quadratic model both the effect of data fitting and the result of numerical calculations. In addition, the conic model can supply enough freedom to make best use of both information of gradients and function values in iterate points. In view of these good properties of conic model, it has attracted wide attention of many scholars [14–28]. In [28], Ni proposed a new trust region subproblem and gave the optimality conditions for the trust region subproblems of a conic model. That is, at the kth iteration, the trial step sk is computed by solving the following conic model trust region subproblem(5)mins∈Rnϕks=gkTs1-akTs+sTBks21-akTs2,(6)s.t.s∈S,where(7)S=s∣s≤Δk,1-akTs≥ε0,horizon vector ak∈Rn, Bk is symmetric and positive semidefinite, and ε0 (0<ε0<1) is a sufficiently small positive number. We note that the conic model ϕk(s) has a denominator and the shape of the trust region S is irregular; therefore the conic trust region subproblems (5)-(7) are not easy to search for the descent point and difficult to solve. The trust region method often does not require the exact solution of trust region subproblem but only requires an approximate solution. The dogleg method for solving the trust region subproblem based on conic model is an approximate solution method; however its calculation is relatively complicated.In this paper, we continue to study the subproblems (5)-(7). In order to find an easy way to solve, inspired by the alternating direction method of multipliers (ADMM), we consider obtaining the approximate solution by the two-step search in two orthogonal directions in the trust region S. ADMM is an algorithm that solves convex optimization problems by breaking them into smaller pieces, each of which are then easier to handle. Because of its significant efficiency and easy implementation, it has recently found wide application in a number of areas (see [29–45]).In the following, we use the alternating orthogonal direction search method to find the approximate solution of the subproblems (5)-(7). The rest of this paper is organized as follows. In the next section, the motivation and description of the simple alternating direction search method are presented. In Section 3, we give the quasi-Newton method based on the conic model for solving unconstrained optimization problems and prove its global convergence properties. The numerical results are provided in Section 4.


## 2. Range ofε0 and the Approximate Solution of the Subproblem

In this section, we will modify the range ofε0 and give the motivation and description of the algorithm. We note that the conic model ϕk(s) has one more parameter ak than ϱk(s). Therefore, ϕk(s) can make full use of the existing function information to satisfy more interpolation conditions by using of the function values and the gradient values. All of these will improve the effectiveness of the algorithm. In general, ak chooses a descent direction, such as gk-1, gk, or sk-1 (see [13–17]). For convenience, we omit the index k of ak, gk, and Bk in this section. Therefore, in this paper we assume that(8)aTg<0and B is positive (abbreviated as B≻0).Let(9)Bag=aTBa-aTaaTg,(10)ε¯0=min⁡1,aTBaBag.From (8), we have a≠0 and(11)Bag>0.Although in principle we are seeking the optimal solution of the subproblems (5)-(7), it is enough to find an approximate solution sk in the feasible region and guarantee a sufficient reduction in the model and the global convergence. Therefore, in order to simplify this algorithm we choose ε0 in (7) such that it satisfies(12)0<ε0<ε¯0.In the following, we consider the alternating direction search method to solve the subproblems (5)-(7) by making full use of the parameters a. The new method is divided into two steps. First, we search in the direction of a and then search in the direction y which is perpendicular to a.Let(13)s=τa+y,where τ≥0,a≠0,y∈Rn, and aTy=0. Then, the solving process of subproblems (5)-(7) is divided into the following two stages.In the first stage, sety=0 and then s=τa. Substituting it into (5)-(7), we have(14)minτ∈Rρτ=τaTg1-τaTa+τ2aTBa21-τaTa2,(15)s.t.τ∈Ω,where(16)Ω=τ∣0≤τ≤τΔ∩τ∣τ≤τdorτ≥τu,(17)τΔ=Δa,τd=1-ε0a2,τu=1+ε0a2.By the direct computation, we have that the derivative ofρ(τ) is(18)ρ′τ=Bagτ+aTg-a6τ-τm3,where(19)τm=1a2.From (17), we know that 0<τd<τm<τu and then τm∉Ω. Since Bag>0, then ρ(τ) has only one stationary point(20)τcp=-aTgBag.By simple calculation, the following lemmas can be easily obtained.Lemma 2.
SupposeaTg<0,B≻0 and 0<ε0<ε¯0. Then τcp<τd.Proof.
Combining with the definition ofε¯0 as in (10), we easily know that(21)0<ε0<aTBaBag.From (17) and (20), we have(22)τcp-τd=ε0-aTBa/Baga2<0.Lemma 3.
Under the same conditions as Lemma2, then ρ(τ) is monotonically increasing in the trust region (τcp,τm); ρ(τ) is monotonically decreasing for 0<τ<τcp and τ>τm.Proof.
From (19) and (20), we know that(23)τm-τcp=aTBaBaga2>0.Then combining with (18) we can obtain that the lemma obviously holds.Letτ∗ be the solution of the subproblems (14)-(15). Then, we can get the following theorems by analysis.Theorem 4.
Under the same conditions as Lemma2, then the solution of the subproblems (14)-(15) is(24)τ∗=min⁡τcp,τΔ,ifτΔ≤τd,τcp,ifτΔ>τd.Proof.
(1) IfτΔ≤τd then from (16) we know that Ω=[0,τΔ]. From Lemmas 2 and 3, we can obtain that if τcp≤τΔ then τ∗=τcp and if τΔ<τcp then τ∗=τΔ. Therefore, τ∗=min⁡τcp,τΔ.
(2) Ifτd<τΔ<τu, then Ω=[0,τd]. From Lemmas 2 and 3, we can similarly get τ∗=τcp.
(3) IfτΔ≥τu, then Ω=[0,τd]∪[τu,τΔ]. From Lemmas 2 and 3, we know that (25)τ∗=arg⁡min⁡ρτcp,ρτΔ.Because of τΔ≥τu, then from (17) we have 1-Δa≤-ε0<0. Since aTg<0, then(26)ρτΔ=ΔaTga1-Δa+Δ2aTBa2a21-Δa2>0.Then τ∗=τcp. The theorem is proved.It is worth noting that ifτ∗=τΔ then from (17) we have τ∗a=Δ. Therefore, for this case we set s∗=τ∗a and exit the calculation of subproblem. Otherwise, we know that τ∗a is inside the trust region S. Then, we should carry out the calculation of the second stage below.We sets=τ∗a+y and substitute it into ϕk(s). Then the subproblems (5)-(7) become(27)miny∈Rnψy=gTτ∗a+y1-τ∗aTa+τ∗a+yTBτ∗a+y21-τ∗aTa2,(28)s.t.y≤Δ~,aTy=0,where(29)Δ~=Δ2-τ∗2a2.In order to remove the equality constraint in (28), we use the null space technique. That is, for a≠0 then there exist n-1 mutually orthogonal unit vectors q,q,…,qn-1 orthogonal to the parameter vector a. Set Q=[q,q,…,qn-1] and y=Qu, where u∈Rn-1. Then (27)-(28) can be simplified as the following subproblem:(30)minu∈Rn-1ψ~u=g~Tu+12uTB~u,(31)s.t.u≤Δ~,where(32)g~=QTg1-τ∗aTa+τ∗QTBa1-τ∗aTa2,B~=QTBQ1-τ∗aTa2Set gk=g~, Bk=B~, and Δk=Δ~. By Algorithm 1, we can obtain the solution u∗ of the subproblems (30)-(31). Then y∗=Qu∗ and s∗=τ∗a+y∗. Thus, the subproblems (5)-(7) are solved approximately.Now we could give the alternating direction search method for solving the conic trust region subproblems (5)-(7) as follows.Algorithm 5.
Givena,g,B, and Δ,
S tep 1. If aTg≥0, then set a=0. Then solve the subproblems (5)-(7) by Algorithm 1 to get sk, and stop.
S tep 2. Compute Bag and ε-0 by (9) and (10). let ε0=0.9ε-0.
S tep 3. Compute τΔ, τd and τcp by (17) and (20).
S tep 4. Compute τ∗ by (24).
S tep 5. If τ∗=τΔ, then sk=τΔa, and stop; otherwise, compute Q, Δ~, g~, and B~ by (29) and (32).
S tep 6. Set gk=g~, Bk=B~ and Δk=Δ~. Then solve the subproblems (30)-(31) by Algorithm 1 to get u∗.
S tep 7. Set y∗=Qu∗ and sk=τ∗a+y∗, and stop.In order to discuss the lower bound of predicted reduction in each iteration, we define the following predicted reduction:(33)preds=ϕ0-ϕs,pred1τ=ρ0-ρτ,(34)pred2y=ψ0-ψy,pred3u=ψ~0-ψ~uNow we should prove the following theorem to guarantee the global convergence of the algorithm proposed in the next section.Theorem 6.
Under the same conditions as Lemma2, if sk is obtained from the above Algorithm 5, then there exists a positive constant c3 such that(35)predsk≥12c3gmin⁡Δ,gB.Proof.
(1) Ifsk is obtained by Step 1 of Algorithm 5, then from Nocedal and Wright [46] we have(36)predsk≥12c1gmin⁡Δ,gB,where c1∈(0,1].
(2) Ifsk is obtained by Step 5 of Algorithm 5, then sk=τΔa and τ∗=τΔ, where τ∗ as defined in (24). By computation, we have(37)predsk=pred1τΔ=-ρτΔ=-ΔΔBag+2aaTg-Δa2aTg2a21-Δa2.From (24), we know that if τ∗=τΔ then τΔ≤τcp and τΔ≤τd. And then from (17) and (20), we can obtain(38)ΔBag+aaTg≤0and(39)1-Δa≥ε0>0.Combining with (37)-(39), we know that(40)predsk≥-ΔaTg2a1-Δa≥ϵΔg2,where(41)ϵ=-aTgag.For aTg<0, then 0<ϵ≤1 holds obviously.
(3) Ifsk is obtained by Step 7 of Algorithm 5, then sk=τ∗a+Qu∗, where τ∗≠τΔ. From (24), we know that τ∗=τcp. Combining with (33) and (34), we have (42)predsk=pred1τcp+pred3u∗.Because u∗ is obtained by Algorithm 1, then from [46] we have(43)pred3u∗≥12c2g~min⁡Δ~,g~B~where c2∈(0,1] and Δ~, g~, and B~ are defined by (29) and (32). Thus,(44)predsk≥pred1τcp=-τcpaTg1-τcpa2-τcp2aTBa21-τcpa22=aTg22aTBa≥ϵ2g22B,where the second equality is from (20) and the last inequality is from (41).
Therefore, the theorem follows from (36), (40), and (44) with(45)c3=min⁡c1,ϵ2.


## 3. The Algorithm and Its Convergence

In this section, we propose a quasi-Newton method with a conic model for unconstrained minimization and prove its convergence under some reasonable conditions. In order to solve the problem (1), we approximate f(x) with a conic model of the form(46)mks=fk+gkTs1-akTs+12sTBks1-akTs2,where fk=f(xk), gk=∇f(xk), Bk∈Rn×n, and ak∈Rn are parameter vectors.Now we give the simple alternating direction trust region algorithm based on conic model (46).Algorithm 7.
⁡
S tep 0. Choose parameters ϵ,ε,∈(0,1), 0<η1<η2<1, 0<δ1<1<δ2, and Δ->0; give a starting point x0∈Rn, B0∈Rn×n, a0∈Rn, and an initial trust region radius Δ0∈(0,Δ-]; set k=0.
S tep 1. Compute fk and gk. If gk<ε, and then stop with xk as the approximate optimal solution; otherwise go to Step 2.
S tep 2. Set a=ak, g=gk, B=Bk, and Δ=Δk. Then solve the subproblem (5)-(7) by Algorithm 5 to get one of the approximate solution sk.
S tep 3. Compute the ratio of predicted reduction and the actual reduction (47)rk=aredskpredsk,where(48)aredsk=fxk-fxk+sk,(49)predsk=-gkTsk1-akTsk-12skTBksk1-akTsk2.
S tep 4. If rk≤η1, then set Δk=δ1Δk and go to Step 2. If rk>η1, then set xk+1=xk+sk and choose the new trust region bound satisfying(50)Δk=min⁡δ2Δk,Δ-,ifrk≥η2,sk=Δk,Δk,otherwise.
S tep 5. Generate ak+1 and Bk+1; set k=k+1, and go to Step 1.The choice of parameterak+1 in the cone model method is crucial. In general, ak+1 and Bk+1 are chosen to satisfy certain interpolation conditions, which means that the conic model function interpolates both the function values and the gradient values of the objective function at xk and xk+1. The choice of the parameters ak+1 and Bk+1 can refer to [13–17] and [47–49], respectively. In this paper, we are not prepared to study the specific iterative formulas of ak+1 and Bk+1 in depth and directly adopt the choice of ak+1 in [16] and the choice of Bk+1 in [49].From (49) and Theorem 6, we have(51)predsk≥12c3gkmin⁡Δk,gkBk,where c3 as defined by (45). The following theorem guarantees that the Algorithm 7 is globally convergent.Theorem 8.
Under the same conditions as Lemma2, suppose that the level set(52)Lx0=x∣fx≤fx0and the sequence ak, gk, and Bk is all uniformly bounded, Bk is symmetric and positive definite, and f is twice continuously differentiable in L(x0). Then for any ε>0, Algorithm 7 terminates in finite number of iterations, that is,(53)limk→∞⁡gk=0.Proof.
We give the proof by contradiction. Suppose that there isε1>0 such that(54)gk≥ε1,∀k.From the hypothesis, we have(55)ak≤a¯,gk≤g¯,0<Bk≤B¯.Combining with (51)-(55), we have(56)predsk≥12c3min⁡Δk,ε1B¯≥12ζΔk,where the first inequality follows from(57)minp,q≥pqp+q,∀p,q>0,and the second inequality is from Δk≤Δ¯ and(58)ζ=c3ε1ε1+B-Δ-.From Step 4 of Algorithm 7 and (56), we obtain that for all k(59)fk-fk+1≥η1predsk≥12η1ζΔk.Since f(x) is bounded from below and fk+1<fk, then we have(60)∞>∑kfk-fk+1≥∑k12η1ζΔk,which implies that(61)∑k=1∞Δk<∞And then(62)limk→∞⁡Δk=0,limk→∞⁡sk=0.
On the other hand, from (55) and (62) we can get(63)11-akTsk=1+akTsk+osk,(64)skTBksk21-akTsk2=12skTBksk+osk2.Then from (55) to (63), we have (65)fk-fxk+sk-predsk=fk-fxk+sk+1+akTskgkTsk+12skTBksk+osk2=-12skT∇2fxk+ϑksksk+akTskgkTsk+12skTBksk+osk2≤12M1+2a-g-+B¯+O1sk2≤12Q+O1Δk2,where ϑk∈(0,1) and Q=M1+2a-g-+B¯. Combining with (56) and (65), we can get that(66)rk-1=fk-fxk+skpredsk-1≤Q+O1ζΔk.From (62) we have rk→1. Hence, there is a sufficiently large positive number K such that ∀k≥K and(67)rk>η1holds. From Step 4 of Algorithm 7, it follows that (68)Δk+1≥Δk,∀k≥K,which is a contradiction to (62). The theorem is proved.


## 4. Numerical Tests

In this section, Algorithm7 is tested with some standard test problems from [16, 50]. The names of the 16 test problems are listed in Table 1. All the computations are carried out in Matlab R2016b on a microcomputer in double precision arithmetic. These tests use the same stopping criterion gk≤10-5. The columns in the Tables have the following meanings: No. denotes the numbers of the test problems; n is the dimension of the test problems; Iter is the number of iterations; nf is the number of function evaluations performed; ng is the number of gradient evaluations; fk is the final objective function value; g is the Euclidean norm of the final gradient; CPU(s) denotes the total iteration time of the algorithm in seconds.Table 1
Test functions.
No. Problem No. Problem 1 Cube 2 Penalty-I 3 Beale 4 Conic 5 Extended powell 6 Variably Dimensioned 7 Rosenbrock 8 Extended Trigonometric 9 Tridiagonal Exponential 10 Brent 11 Troesch 12 Cragg and Levy 13 Broyden Tridiagonal 14 Brown 15 Discrete Boundary Value 16 Extended TrigonometricThe parameters in these algorithms are(69)a0=0,B0=I,ε0=ϵ=10-5,Δ0=1,Δ-=10,η1=0.01,η2=0.75,δ1=0.5,δ2=2.In order to analyze the effectiveness of our new algorithm, we compare Algorithm7 with the alternating direction trust region method based on conic model (abbreviated as ADCTR) in [12]. The numerical results of ADCTR and Algorithm 7 are listed in Table 2. We note that the optimal value of these test problems is f∗=0. From Table 2, we can see that the performance of Algorithm 7 is feasible and effective. For the above 16 problems, Algorithm 7 is better than the ADCTR for 13 tests and is somewhat bad for 4 tests, and the two algorithms are same in efficiency for the other 1 tests. Therefore, it seems that Algorithm 7 is better than algorithm ADCTR in [12].Table 2
Numerical results.
No. n Starting point Algorithm Iter nf/ng f k g CPU (s) 1 2 (2,1.5) ADCTR 52 53/43 1.0377e-15 1.9477e-06 0.064681 Algorithm7 56 57/37 7.7243e-13 5.0126e-06 0.049091 2 2 (1,2) ADCTR 10 11/11 9.0831e-06 8.9419e-06 0.048493 Algorithm7 8 9/9 9.0831e-06 8.9461e-06 0.033654 3 2 (2,-2) ADCTR 18 19/18 9.0379e-15 9.3925e-07 0.053525 Algorithm7 16 17/17 7.3180e-16 4.3875e-08 0.039889 4 2 (-1,1) ADCTR 16 17/13 1.1407e-12 2.1360e-06 0.050445 Algorithm7 9 10/10 9.8953e-13 1.9895e-06 0.038502 5 4 (3,-1,0,1) ADCTR 41 42/34 4.8648e-09 4.5887e-06 0.062011 Algorithm7 107 108/83 9.8794e-12 3.9079e-06 0.054692 6 4 0.25∗(3,2,1,0) ADCTR 32 33/29 2.3856e-14 3.0965e-07 0.066287 Algorithm7 27 28/27 3.6994e-13 5.7059e-06 0.059686 7 2 (-1.2,1) ADCTR 50 51/49 1.5486e-14 5.4101e-06 0.064227 Algorithm7 47 48/40 2.4987e-19 2.1924e-08 0.048847 8 4 0.25∗(1,1,1,1) ADCTR 47 48/34 7.9158e-15 4.1153e-07 0.076865 Algorithm7 28 29/22 3.0282e-04 8.2932e-06 0.085045 9 4 1.5∗(1,1,1,1) ADCTR 7 8/8 8.1577e-12 4.5905e-06 0.058505 Algorithm7 6 7/7 5.2383e-13 2.5151e-06 0.054092 10 4 10∗(1,1,1,1) ADCTR 81 82/58 5.8024e-18 4.6604e-07 0.089702 Algorithm7 56 57/40 7.1113e-18 2.9344e-07 0.077486 11 4 (-1,-1,-1,-1) ADCTR 64 65/52 9.9974e-14 4.4986e-06 0.082793 Algorithm7 51 52/41 2.0655e-13 5.6223e-06 0.060827 12 4 (-1,-1,-1,-1) ADCTR 48 49/43 1.1247e-08 5.2578e-06 0.068215 Algorithm7 128 129/125 4.3892e-08 6.6956e-06 0.054621 13 4 (-1,-1,-1,-1) ADCTR 35 36/19 1.4498e-11 5.0442e-06 0.063276 Algorithm7 38 39/19 1.4483e-11 6.4239e-06 0.047787 14 4 (0,1,0,1) ADCTR 91 92/52 0.1998e-06 2.5916e-07 0.089294 Algorithm7 63 64/41 1.3996e-05 5.7261e-07 0.052852 15 4 -0.08∗(2, 3,3,2) ADCTR 23 24/15 2.0042e-12 8.2898e-06 0.061544 Algorithm7 19 20/12 2.3574e-12 6.3180e-06 0.043328 16 4 0.25∗(1,1,1,1) ADCTR 14 15/15 3.0282e-04 4.9068e-06 0.048488 Algorithm7 14 15/11 3.6700e-04 8.1833e-06 0.049998


## 5. Conclusions

The algorithm ADCTR and Algorithm7 are similar; that is, the idea of alternating direction method is used to solve the conic trust region subproblem. However, Algorithm 7 in this paper takes into account the special property that the parameter vector a is generally taken as the descending direction. Thus, under the assumption of aTg<0, the calculation of Algorithm 7 is simpler to calculate and has shorter CPU time, better calculation effect, and also global convergence.However, there are still many aspects worthy of further study, for example, weakening the positive definite condition ofBk, using algorithms to solve large-scale problems, calculation of convergence rate, and so on.



---

*Source: 5358191-2018-12-18.xml*