Coursehero >> Wisconsin >> Wisconsin >> CS 726

Course Hero has millions of student submitted documents similar to the one below including study guides, homework solutions, papers, exam answer keys and textbook solutions.

Journal IMA of Numerical Analysis (2006) 26, 604 627 doi:10.1093/imanum/drl006 Advance Access publication on March 24, 2006 The cyclic Barzilai Borwein method for unconstrained optimization YU-HONG DAI State Key Laboratory of Scienti c and Engineering Computing, Institute of Computational Mathematics and Scienti c/Engineering computing, Academy of Mathematics and Systems Science, Chinese Academy of Sciences, PO Box 2719, Beijing 100080, People s Republic of China WILLIAM W. HAGER Department of Mathematics, University of Florida, Gainesville, FL 32611, USA KLAUS SCHITTKOWSKI Department of Computer Science, University of Bayreuth, 95440 Bayreuth, Germany AND HONGCHAO ZHANG Department of Mathematics, University of Florida, Gainesville, FL 32611, USA [Received on 2 June 2005; revised on 11 January 2006] In the cyclic Barzilai Borwein (CBB) method, the same Barzilai Borwein (BB) stepsize is reused for m consecutive iterations. It is proved that CBB is locally linearly convergent at a local minimizer with positive de nite Hessian. Numerical evidence indicates that when m > n/2 3, where n is the problem dimension, CBB is locally superlinearly convergent. In the special case m = 3 and n = 2, it is proved that the convergence rate is no better than linear, in general. An implementation of the CBB method, called adaptive cyclic Barzilai Borwein (ACBB), combines a non-monotone line search and an adaptive choice for the cycle length m. In numerical experiments using the CUTEr test problem library, ACBB performs better than the existing BB gradient algorithm, while it is competitive with the well-known PRP+ conjugate gradient algorithm. Keywords: unconstrained optimization; gradient method; convex quadratic programming; non-monotone line search. 1. Introduction In this paper, we develop a cyclic version of the Barzilai Borwein (BB) gradient type method (Barzilai & Borwein, 1988) for solving an unconstrained optimization problem min f (x), Email: dyh@lsec.cc.ac.cn Email: hager@math.u .edu Email: klaus.schittkowski@uni-bayreuth.de Email: hzhang@math.u .edu x Rn , (1.1) c The author 2006. Published by Oxford University Press on behalf of the Institute of Mathematics and its Applications. All rights reserved. THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 605 where f is continuously differentiable. Gradient methods start from an initial point x0 and generate new iterates by the rule xk+1 = xk k gk , (1.2) k 0, where gk = f (xk )T is the gradient, viewed as a column vector, and k is a stepsize computed by some line search algorithm. In the steepest descent (SD) method, which can be traced back to Cauchy (1847), the exact steplength is given by k arg min f (xk gk ). R (1.3) It is well known that SD can be very slow when the Hessian of f is ill-conditioned at a local minimum (see Akaike, 1959; Forsythe, 1968). In this case, the iterates slowly approach the minimum in a zigzag fashion. On the other hand, it has been shown that if the exact SD step is reused in a cyclic fashion, then the convergence is accelerated. Given an integer m 1, which we call the cycle length, cyclic SD can be expressed as m +i SD = m +1 for i = 1, . . . , m, (1.4) SD = 0, 1, . . ., where k is the exact steplength given by (1.3). Formula (1.4) was rst proposed in Friedlander et al. (1999), while the particular choice m = 2 was also investigated in Dai (2003) and Raydan & Svaiter (2002). The analysis in Dai & Fletcher (2005a) shows that if m > n , cyclic SD is 2 likely R-superlinearly convergent. Hence, SD is accelerated when the stepsize is repeated. Let BB denote the original method of Barzilai & Borwein (1988). In this paper, we develop a cyclic Barzilai Borwein (CBB) method. The basic idea of Barzilai and Borwein is to regard the ma1 trix D( k ) = k I as an approximation of the Hessian 2 f (xk ) and impose a quasi-Newton property on D( k ): k = arg min D( )sk 1 yk 1 2 , R (1.5) where sk 1 = xk xk 1 , yk 1 = gk gk 1 and k BB k = 2. The proposed stepsize, obtained from (1.5), is T sk 1 sk 1 T sk 1 yk 1 . (1.6) Other possible choices for the stepsize k include Dai (2003), Dai & Fletcher (2006), Dai & Yang (2001, 2003), Friedlander et al. (1999), Grippo & Sciandrone (2002), Raydan & Svaiter (2002) and Sera ni et al. (2005). In this paper, we refer to (1.6) as the BB formula. The gradient method (1.2) corresponding to the BB stepsize (1.6) is called the BB method. Due to their simplicity, ef ciency and low memory requirements, BB-like methods have been used in many applications. Glunt et al. (1993) present a direct application of the BB method in chemistry. Birgin et al. (1999) use a globalized BB method to estimate the optical constants and the thickness of thin lms, while in Birgin et al. (2000) further extensions are given, leading to more ef cient projected gradient methods. Liu & Dai (2001) provide a powerful scheme for solving noisy unconstrained optimization problems by combining the BB method and a stochastic approximation method. The projected BBlike method turns out to be very useful in machine learning for training support vector machines (see Sera ni et al., 2005; Dai & Fletcher, 2006). Empirically, good performance is observed on a wide variety of classi cation problems. 606 Y.-H. DAI ET AL. The superior performance of cyclic SD, compared to the ordinary SD, as shown in Dai & Fletcher (2005a), leads us to consider the CBB method: m +i BB = m +1 for i = 1, . . . , m, (1.7) where m 1 is again the cycle length. An advantage of the CBB method is that for general non-linear functions, the stepsize is given by the simple formula (1.5) in contrast to the non-trivial optimization problem associated with the SD step (1.3). In Friedlander et al. (1999), the authors obtain global convergence of CBB when f is a strongly convex quadratic. Dai (2003) establishes the R-linear convergence of CBB for a strongly convex quadratic. In Section 2, we prove the local R-linear convergence for the CBB method at a local minimizer of a general non-linear function. In Section 3, numerical evidence for strongly convex quadratic functions indicates that the convergence is superlinear if m > n/2 3. In the special case m = 3 and n = 2, we prove that the convergence is at best linear, in general. In Section 4, we propose an adaptive method for computing an appropriate cycle length, and we obtain a globally convergent non-monotone scheme by using a modi ed version of the line search developed in Dai & Zhang (2001). This new line search, an adaptive analogue of Toint s (1997) scheme for trust region methods, accepts the original BB stepsize more often than does Raydan s (Raydan, 1997) strategy for globalizing the BB method. We refer to Raydan s globalized BB implementation as the GBB method. Numerical comparisons with the PRP+ algorithm and the SPG2 algorithm (Birgin et al., 2000) (one version of the GBB method) are given in Section 4 using the CUTEr test problem library (Bongartz et al., 1995). Throughout this paper, we use the following notations. is the Euclidean norm of a vector. The subscript k is often associated with the iteration number in an algorithm. The letters i, j, k, , m and n, either lower or upper case, designate integers. The gradient f (x) is a row vector, while g(x) = f (x)T is a column vector; here T denotes transpose. The gradient at the iterate xk is gk = g(xk ). We let 2 f (x) denote the Hessian of f at x. The ball with centre x and radius is denoted B (x). 2. Local linear convergence In this section, we prove R-linear convergence for the CBB method. In Liu & Dai (2001), it is proposed that R-linear convergence for the BB method applied to a general non-linear function could be obtained from the R-linear convergence results for a quadratic by comparing the iterates associated with a quadratic approximation to the general non-linear iterates. In our R-linear convergence result for the CBB method, we make such a comparison. The CBB iteration can be expressed as xk+1 = xk k gk , where k = siT si siT yi , i = (k), and (k) = m (k 1)/m , (2.2) (2.1) k 1. For r R, r denotes the largest integer j such that j r . We assume that f is two times Lipschitz continuously differentiable in a neighbourhood of a local minimizer x where the Hessian H = 2 f (x ) is positive de nite. The second-order Taylor approximation f to f around x is given by 1 f (x) = f (x ) + (x x )T H (x x ). 2 (2.3) THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 607 We will compare an iterate xk+ j generated by (2.1) to a CBB iterate xk, j associated with f and the starting point x k,0 = xk . More precisely, we de ne xk,0 = xk , xk, j+1 = xk, j k, j gk, j , where k, j = k si si , T j 0, (2.4) if (k + j) = (k) i = (k + j), otherwise. siT yi Here sk+ j = xk, j+1 xk, j , gk, j = H (xk, j x ) and yk+ j = gk, j+1 gk, j . We exploit the following result established in Dai (2003, Theorem 3.2): x LEMMA 2.1 Let { k, j : j 0} be the CBB iterates associated with the starting point xk,0 = xk and the in (2.3), where H is positive de nite. Given two arbitrary constants C2 > C1 > 0, there quadratic f exists a positive integer N with the following property: For any k 1 and k,0 [C1 , C2 ], xk,N x 1 xk,0 x . 2 (2.5) In our next lemma, we estimate the distance between xk, j and xk+ j . Let B (x) denote the ball with centre x and radius . Since f is two times Lipschitz continuously differentiable and 2 f (x ) is positive de nite, there exist positive constants , and 2 > 1 such that f (x) H (x x ) and 1 y T 2 f (x)y yT y 2 for all y Rn and x B (x ). (2.7) x x 2 for all x B (x ) (2.6) Note that if xi and xi+1 B (x ), then the fundamental theorem of calculus applied to yi = gi+1 gi yields 1 2 siT si siT yi 1 . 1 (2.8) Hence, when the CBB iterates lie in B (x ), condition (2.5) of Lemma 2.1 is satis ed with C1 = 1/ 2 and C2 = 1/ 1 . If we de ne g(x) = f (x)T , then the fundamental theorem of calculus can also be used to deduce that g(x) = g(x) g(x ) for all x B (x ). 2 x x (2.9) 608 Y.-H. DAI ET AL. LEMMA 2.2 Let {x j : j k} be a sequence generated by the CBB method applied to a function f with a local minimizer x , and assume that the Hessian H = 2 f (x ) is positive de nite with (2.7) satis ed. Then for any xed positive integer N , there exist positive constants and with the following property: For any x k B (x ), k [ 1 , 1 ], [0, N ] with 2 1 xk, j x we have xk+ j B (x ) for all j [0, ]. Proof. Throughout the proof, we let c denote a generic positive constant, which depends on xed constants, such as N , 1 , 2 or , but not on k, the choice of xk B (x ) or the choice of k [ 1 , 1 ]. To facilitate the proof, we also show that 2 1 g(xk+ j ) g(xk, j ) sk+ j | k+ j k, j | c xk x 2 , c xk x 1 xk,0 x 2 and for all j [0, max{0, 1}], xk x (2.10) xk+ j xk, j 2 (2.11) (2.12) (2.13) c xk x , , (2.14) for all j [0, ], where g(x) = f (x)T = H (x x ). The proof of (2.11) (2.14) is by induction on . For since xk,0 = xk . By (2.6), we have = 0, we take = . Relation (2.11) is trivial xk x 2 , g(xk ) g(xk,0 ) = g(xk ) g(xk ) which gives (2.12). Since = and xk B (x ), it follows from (2.9) that sk = k gk 2 xk x , 1 which gives (2.13). Relation (2.14) is trivial since k,0 = k . Now, proceeding by induction, suppose that there exist L [1, N ) and > 0 with the property that if (2.10) holds for any [0, L 1], then (2.11) (2.14) are satis ed for all j [0, ]. We wish to show that for a smaller choice of > 0, we can replace L by L + 1. Hence, we suppose that (2.10) holds for all j [0, L]. Since (2.10) holds for all j [0, L 1], it follows from the induction hypothesis and (2.13) that xk+L+1 x xk x + L sk+i i=0 c xk x . (2.15) Consequently, by choosing smaller, if necessary, we have xk+L+1 B (x ) when xk B (x ). By the triangle inequality, xk+L+1 xk,L+1 = xk+L k+L g(xk+L ) [xk,L k,L g(xk,L )] xk+L xk,L + | k,L | g(xk+L ) g(xk,L ) + | k+L k,L | g(xk+L ) . (2.16) THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 609 We now analyse each of the terms in (2.16). By the induction hypothesis, the bound (2.11) with j = L holds, which gives xk+L xk,L c xk x 2 . (2.17) By the de nition of , either k,L = k [ 1 , 1 ] or 2 1 k,L = In this latter case, 1 2 siT si siT H si = siT si siT yi 1 . 1 siT si siT yi , i = (k + L). Hence, in either case, k,L [ 1 , 1 ]. It follows from (2.12) with j = L that 2 1 | k,L | g(xk+L ) g(xk,L ) 1 g(xk+L ) g(xk,L ) 1 c xk x 2 . (2.18) Also, by (2.14) with j = L and (2.9), we have | k+L k,L | g(xk+L ) Utilizing (2.15) (with L replaced by L 1) gives | k+L k,L | g(xk+L ) c xk x 2 . (2.19) c xk x xk+L x . We combine (2.16) (2.19) to obtain (2.11) for j = L + 1. Note that in establishing (2.11), we exploited (2.12) (2.14). Consequently, to complete the induction step, each of these estimates should be proved for j = L + 1. Focusing on (2.12) for j = L + 1, we have g(xk+L+1 ) g(xk,L+1 ) g(xk+L+1 ) g(xk+L+1 ) + g(xk+L+1 ) g(xk,L+1 ) = g(x k+L+1 ) g(xk+L+1 ) + H (xk+L+1 xk,L+1 ) g(xk+L+1 ) H (xk+L+1 x ) + 2 xk+L+1 xk,L+1 g(xk+L+1 ) H (xk+L+1 x ) + c xk x 2 , since H 2 by (2.7). The last inequality is due to (2.11) for j = L + 1, which was just established. Since we chose small enough that xk+L+1 B (x ) (see (2.15)), (2.6) implies that g(xk+L+1 ) H (xk+L+1 x ) Hence, g(xk+L+1 ) g(xk,L+1 ) xk+L+1 x 2 c xk x 2 . c xk x 2 , which establishes (2.12) for j = L + 1. 610 Y.-H. DAI ET AL. Observe that k+L+1 equals either k [ 1 , 1 ] or (siT si )/(siT yi ), where k + L i= 2 1 (k + L + 1) > k. In this latter case, since xk+ j B (x ) for 0 j L + 1, it follows from (2.8) that 1 . k+L+1 1 Combining this with (2.9), (2.15) and the bound (2.13) for j sk+L+1 = k+L+1 g(xk+L+1 ) L, we obtain c xk x . 2 xk+L+1 x 1 Hence, (2.13) is established for j = L + 1. Finally, we focus on (2.14) for j = L + 1. If (k + L + 1) = (k), then k,L+1 = k+L+1 = k , so we are done. Otherwise, (k + L + 1) > (k), and there exists an index i (0, L] such that k+L+1 = By (2.11) and the fact that i T sk+i sk+i T sk+i yk+i and k,L+1 = sk+i sk+i T sk+i yk+i T . L, we have sk+i sk+i c xk x 2 . Combining this with (2.13), and choosing smaller, if necessary, gives T T |sk+i sk+i sk+i sk+i | = |2sk+i (sk+i sk+i ) sk+i sk+i T 2 | c xk x 3 . (2.20) Since k,i [ 1 , 1 ], we have 2 1 sk+i = k,i gk,i Furthermore, by (2.10), it follows that sk+i 1 1 xk,0 x = xk x . 2 2 2 2 (2.21) 1 H (xk,i x ) 2 1 xk,i x . 2 Hence, combining (2.20) and (2.21) gives 1 T sk+i sk+i sk+i sk+i T = T |sk+i sk+i sk+i sk+i | T sk+i sk+i T c xk x . (2.22) Now let us consider the denominators of k+i and k,i . Observe that T T T sk+i yk+i sk+i yk+i = sk+i (yk+i yk+i ) + (sk+i sk+i )T yk+i T = sk+i (yk+i yk+i ) + (sk+i sk+i )T H sk+i . (2.23) By (2.11) and (2.13), we have |(sk+i sk+i )T H sk+i | = |(sk+i sk+i )T H sk+i (sk+i sk+i )T H (sk+i sk+i )| c xk x 3 (2.24) THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 611 for suf ciently small. Also, by (2.12) and (2.13), we have T |sk+i (yk+i yk+i )| sk+i ( gk+i+1 gk,i+1 + gk+i gk,i ) c xk x 3 . (2.25) Combining (2.23) (2.25) yields T |sk+i yk+i sk+i yk+i | T c xk x 3 . (2.26) Since xk+i and xk+i+1 B (x ), it follows from (2.7) that T T sk+i yk+i = sk+i (gk+i+1 gk+i ) 1 sk+i 2 = 1 | k+i |2 gk+i 2 . (2.27) By (2.8) and (2.7), we have | k+i |2 gk+i Finally, (2.10) gives xk+i x Combining (2.27) (2.29) yields T sk+i yk+i 2 2 1 gk+i 2 2 2 = 1 g(xk+i ) g(x ) 2 2 2 2 1 2 2 xk+i x 2 . (2.28) 1 xk x 2 . 4 (2.29) 3 1 4 2 2 xk x 2 . (2.30) Combining (2.26) and (2.30) gives 1 Observe that | k+L+1 k,L+1 | = T sk+i sk+i T sk+i yk+i sk+i yk+i T T sk+i yk+i = T |sk+i yk+i sk+i yk+i | T T sk+i yk+i c xk x . (2.31) sk+i sk+i T sk+i yk+i T = k,L+1 1 sk+i yk+i T T sk+i yk+i T sk+i sk+i sk+i yk+i T T sk+i yk+i sk+i sk+i T 1 1 1 = where a =1 T sk+i sk+i sk+i sk+i T 1 |a(1 b) + b| 1 1 (|a| + |b| + |ab|), 1 (2.32) T sk+i sk+i sk+i sk+i T and b = 1 sk+i yk+i T T sk+i yk+i . 612 Together, (2.22), (2.31) and (2.32) yield Y.-H. DAI ET AL. | k+L+1 k,L+1 | c xk x for suf ciently small. This completes the proof of (2.11) (2.14). THEOREM 2.3 Let x be a local minimizer of f , and assume that the Hessian 2 f (x ) is positive de nite. Then there exist positive constants and and a positive constant c < 1 with the property that for all starting points x0 , x1 B (x ), x0 = x1 , the CBB iterates generated by (2.1) and (2.2) satisfy xk x ck x 1 x . Proof. Let N > 0 be the integer given in Lemma 2.1, corresponding to C1 = 1 and C2 = 1 , 1 2 and let 1 and 1 denote the constants and given in Lemma 2.2 Let 2 denote the constant c in 1 , (2.13). In other words, these constant 1 , 1 and 2 have the property that whenever xk x k [ 1 , 1 ] and 2 1 xk, j x we have x k+ j xk, j sk+ j 1 x k x 2 , 2 x k x , (2.33) (2.34) (2.35) 1 xk,0 x 2 for 0 j 1 < N, xk+ j B (x ), for all j [0, ]. Moreover, by the triangle inequality and (2.34), it follows that x k+ j x for all j [0, ]. We de ne = min{ 1 , , (4 1 ) 1 }. (N 2 + 1) xk x = 3 xk x , 3 = (N 2 + 1), (2.36) (2.37) For any x0 and x1 B (x ), we de ne a sequence 1 = k1 < k2 < in the following way: Starting with the index k1 = 1, let j1 > 0 be the smallest integer with the property that 1 1 xk1 ,0 x = x1 x . 2 2 Since x0 and x1 B (x ) B (x ), it follows from (2.8) that xk1 , j1 x 1,0 = 1 = By Lemma 2.1, j1 T s0 s0 T s0 y0 [ 1 , 1 ]. 2 1 N . De ne k2 = k1 + j1 > k1 . By (2.33) and (2.37), we have xk2 x = xk1 + j1 x 1 x k 1 x 2 xk1 + j1 xk1 , j1 + xk1 , j1 x + Since x1 x 1 xk1 ,0 x 2 1 3 2 = 1 x k 1 x + x k1 x x k1 x . 2 4 , it follows that xk2 B (x ). By (2.35), x j B (x ) for 1 j (2.38) k1 . THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 613 Now, proceed by induction. Assume that ki has been determined with xki B (x ) and x j B (x ) for 1 j ki . Let ji > 0 be the smallest integer with the property that xki , ji x 1 1 xki ,0 x = x ki x . 2 2 3 x ki x . 4 N (i 1) + 1. Hence, i x k1 x k/N . Set ki+1 = ki + ji > ki . Exactly as in (2.38), we have xki+1 x Again, xki+1 B (x ) and x j B (x ) for j [1, ki+1 ]. For any k [ki , ki+1 ), we have k ki + N 1 N i, since ki Also, (2.36) gives xk x 3 x k i x 3 3 4 (k/N ) 1 3 3 4 i 1 x1 x = ck x 1 x , where = This completes the proof. 3. The CBB method for convex quadratic programming In this section, we give numerical evidence which indicates that when m is suf ciently large, the CBB method is superlinearly convergent for a quadratic function f (x) = 1T x Ax bT x, 2 (3.1) 4 3 3 and c = 3 4 1/N < 1. where A Rn n is symmetric and positive de nite and b Rn . Since CBB is invariant under an orthogonal transformation and since gradient components corresponding to identical eigenvalues can be combined (see, e.g. Dai & Fletcher, 2005b), we assume without loss of generality that A is diagonal: A = diag( 1 , 2 , . . . , n ) with 0 < 1 < 2 < < n . (3.2) In the following subsections, we give an overview of the experimental convergence results; we then show in the special case m = 2 and n = 3 that the convergence rate is no better than linear, in general. Finally, we show that the convergence rate for CBB is strictly faster than that of SD. We obtain some further insights by applying our techniques to cyclic SD. 3.1 Asymptotic behaviour and cycle number gk+1 = (I k A)gk . In the quadratic case, it follows from (1.2) and (3.1) that (3.3) 614 Y.-H. DAI ET AL. TABLE 1 Transition to superlinear convergence n Superlinear m Linear m (i) 2 1 3 3 2 4 2 1 5 4 3 6 4 3 8 5 4 10 6 5 12 7 6 14 8 7 If gk denotes the ith component of the gradient gk , then by (3.3) and (3.2), we have gk+1 = (1 k i )gk , (i) (i) (i) i = 1, 2, . . . , n. (i) (3.4) We assume that gk = 0 for all suf ciently large k. If gk = 0, then by (3.4), component i remains zero during all subsequent iterations; hence, it can be discarded. In the BB method, starting values are needed for x0 and x1 in order to compute 1 . In our study of CBB, we treat 1 as a free parameter. In our numerical experiments, 1 is the exact stepsize (1.3). For different choices of the diagonal matrix (3.2) and the starting point x1 , we have evaluated the convergence rate of CBB. By the analysis given in Friedlander et al. (1999) for positive de nite quadratics or by the result given in Theorem 2.3 for general non-linear functions, the convergence rate of the iterates is at least linear. On the other hand, for m suf ciently large, we observe experimentally that the convergence rate is superlinear. For xed dimension n, the value of m where the convergence rate makes a transition between linear and non-linear is shown in Table 1. More precisely, for each value of n, the convergence rate is superlinear when m is greater than or equal to the integer given in the second row of Table 1. The convergence is linear when m is less than or equal to the integer given in the third row of Table 1. The limiting integers appearing in Table 1 are computed in the following way: For each dimension, we randomly generate 30 problems, with eigenvalues uniformly distributed on (0, n], and 50 starting points a total of 1500 problems. For each test problem, we perform 1000n CBB iterations, and we plot log(log( gk )) versus the iteration number. We t the data with a least squares line, and we compute the correlation coef cient to determine how well the linear regression model ts the data. If the correlation coef cient is 1 (or 1), then the linear t is perfect, while a correlation coef cient of 0 means that the data are uncorrelated. A linear t in a plot of log(log( gk )) versus the iteration number indicates superlinear convergence. For m large enough, the correlation coef cients are between 1.0 and 0.98, indicating superlinear convergence. As we decrease m, the correlation coef cient abruptly jumps to the order of 0.8. The integers shown in Table 1 re ect the values of m where the correlation coef cient jumps. Based on Table 1, the convergence rate is conjectured to be superlinear for m > n/2 3. For n < 6, the relationship between m and n at the transition between linear and superlinear convergence is more complicated, as seen in Table 1. Graphs illustrating the convergence appear in Fig. 1. The horizontal axis in these gures is the iteration number, while the vertical axis gives log(log( gk )). Here represents the sup norm. In this case, straight lines correspond to superlinear convergence the slope of the line re ects the convergence order. In Fig. 1, the bottom two graphs correspond to superlinear convergence, while the top two graphs correspond to linear convergence for these top two examples, a plot of log( gk ) versus the iteration number is linear. 3.2 Analysis for the case m = 2 and n = 3 The theoretical veri cation of the experimental results given in Table 1 is not easy. We have the following partial result in connection with the column m = 2. THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 615 FIG. 1. Graphs of log(log( gk )) versus k; (a) 3 n 6 and m = 3; (b) 6 n 9 and m = 4. THEOREM 3.1 For n = 3, there exists a choice for the diagonal matrix (3.2) and a starting guess x1 with the property that k+8 = k for each k, and the convergence rate of CBB with m = 2 is at most linear. Proof. To begin, we treat the initial stepsize 1 as a variable. For each k, we de ne the vector u k by uk = (i) (i) (gk )2 , gk 2 i = 1, . . . , n. (3.5) The above de nition is important and is used for some other gradient methods (see Forsythe, 1968; Dai & Yang, 2001). For the case m = 2, we can obtain by (2.1), (2.2) and (3.4) and the de nition of u k that u 2k+1 = for all k (i) (1 2k 1 i )4 u 2k 1 n =1 (1 2k 1 (i) )4 u 2k 1 () (3.6) 1 and i = 1, . . . , n. In the same fashion, we have 2k+1 = n 2 (i) i=1 (1 2k 1 i ) u 2k 1 . n 2 (i) i=1 i (1 2k 1 i ) u 2k 1 (3.7) We want to force our examples to satisfy u9 = u1 and 9 = 1 . (3.8) For k 1, a subsequent iteration of the method is uniquely determined by u 2k 1 and 2k 1 . It follows from (3.8) that u 8k+1 = u 1 and 8k+1 = 1 for all k 1, and hence a cycle occurs. For any i and j, let bi j be de ned by bi j = 1 2i 1 j . (3.9) Henceforth, we focus on the case n = 3 speci ed in the statement of the Theorem 3.1. To satisfy relation (3.8), we impose the following condition on the stepsizes { 1 , 3 , 5 , 7 }: 4 bi j = , i=1 j = 1, 2, 3, (3.10) 616 Y.-H. DAI ET AL. where > 0 is a positive number. By (3.6) and (3.10), we know that the rst equation of (3.8) is satis ed. On the other hand, (3.6), (3.7), 9 = 1 and the de nition of (3.9) imply the following system of linear equations for u 1 : 2 2 2 b11 b21 b12 b22 b13 b23 u (1) 1 42 42 42 b12 b22 b32 b13 b23 b33 (2) b11 b21 b31 u (3.11) T u1 = 4 4 2 1 = 0. b b b b41 b4 b4 b2 b42 b4 b4 b2 b43 (3) 11 21 31 12 22 32 13 23 33 u1 5 b4 b4 b2 5 b4 b4 b2 5 b4 b4 b2 b11 21 31 41 b12 22 32 42 b13 23 33 43 2 1 The above system has three variables and four equations. Multiplying the jth column by b1 j b2 j b4 j for j = 1, 2, 3 and using condition (3.10), it follows that the rank of the coef cient matrix T is the same as the rank of the 4 3 matrix B with entries bi j . By the de nition of bi j , the rank of T is at most 2; hence, the linear system (3.11) has a non-zero solution u 1 . To complete the construction, u 1 should satisfy the constraints u 1 > 0, and (1) (2) (i) i = 1, 2, 3, (3.12) u 1 + u 1 + u 1 = 1. The above conditions are ful lled if we look for a solution { 1 , 3 , 5 , 7 } of (3.10) such that 1 1 1 , 3 ( 1 , 2 ) 1 1 and 5 , 7 ( 2 , 3 ). (3) (3.13) (3.14) In this case, we may choose 2 1 u 1 = t b11 b21 b12 b13 b43 b42 2 1 , b12 b22 b13 b11 b41 b43 2 1 , b13 b23 b11 b12 b42 b41 T , (3.15) where t > 0 is a scaling factor such that (3.13) holds. Therefore, if we choose { 1 , 3 , 5 , 7 } satisfying (3.10) and (3.14) and furthermore u 1 from (3.15), relation (3.8) holds. Hence, we have that u 8+i = u i and 8+i = i for all i 1. Now we discuss a possible choice of 0 > in (3.10). Speci cally, we are interested in the maximal value of such that (3.10) and (3.14) hold. By continuity assumption, we know that suitable solutions exist for any (0, ). This leads to the maximization problem 4 max : i=1 1 1 1 1 bi j = ( j = 1, 2, 3); 1 , 3 ( 1 , 2 ), 5 , 7 ( 2 , 3 ) . (3.16) To solve (3.16), we consider the Lagrangian function 3 4 L( , 1 , 3 , 5 , 7 , 1 , 2 , 3 ) = + j=1 j i=1 (1 2i 1 j ) , (3.17) THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 617 where { j } are the multipliers corresponding to equality constraints. Since at a KKT point of (3.16) the partial derivatives of L are zero, we require { i } to satisfy relation (3.10), 1 + 2 + 3 = 1 and 3 4 j j j=1 =1 =i (1 2 1 j ) = 0 (i = 1, 2, 3, 4). (3.18) Dividing each relation in (3.18) by and using (3.10), we obtain the following linear equations for = ( 1 , 2 , 3 )T : H = 0, where H R4 3 with h i j = j bi 1 . j (3.19) To guarantee that system (3.19) has a non-zero solution , the rank of the coef cient matrix H must be at most 2. Let H3,3 denote the submatrix formed by the rst three rows of H . By direct calculation, we obtain det(H3,3 ) = 1 2 3 ( 1 2 )( 2 3 )( 3 1 )( 1 3 )( 3 5 )( 5 1 ) . i, j {1,2,3} bi j (3.20) Thus, det(H3,3 ) = 0 and inequality constraints (3.14) lead to 1 = 3 . Similarly, we can get 5 = 7 . From (3.10), we know that (3.16) achieves its maximum = at 1 = 3 = ( ) 1 , 5 = 7 = ( + ) 1 , ( 1 2 )2 ( 2 3 )2 ( 1 2 + 2 3 + 3 1 2 )2 2 (3.21) (3.22) where = 1 + 3 , = 1 + 2 , 1 = ( 1 )2 and 2 = ( 2 )2 . From the continuity argument, we 2 2 know that there exist cyclic examples of the CBB method with m = 2 for any (0, ). For example, we may consider the following symmetric subfamily of examples with (0, 1 ]: 2 1 , 5 = [ 1 + (1 ) 2 ] 1 , 3 , 7 = [ (1 ) 1 + 2 ] 1 . (3.23) It is easy to check that the above { i } satis es (3.10) and (3.14). When moves from 0 to 1 , we can see 2 that the value moves from 0 to . Now we present some numerical examples. Suppose that 1 = 1, 2 = 5 and 3 = 8. Because of 9 (3.22), we choose 1 = 3 = 1 and 5 = 7 = 1 from where the maximizer = 49 is found. From 2 7 972 28 1 (3.11), we get u 1 = ( 1001 , 1001 , 1001 )T . By the de nition of u 1 , the previous discussions and choosing T with any t > 0 and = 1 , the CBB method with m = 2 produces g1 = t ( 18 3, 2 7, 1) 1 2 cycling of the sequence given by {u i } and { i }. By assuming that the Hessian matrix is A = diag(1, 5, 8), we also compute the sequences {u 2k 1 } and { 2k 1 } generated by (3.6) and (3.7), respectively. Initial values for u 1 and 1 are obtained by an SD step at u 0 , i.e. 1 = 0 = uTu0 0 u T Au 0 0 , u1 = (i) (1 0 i )2 (u 0 )2 (1 0 )2 (u 0 )2 () (i) (i = 1, 2, 3). 618 Y.-H. DAI ET AL. For different u 0 , we see that different cycles are obtained, which are numerically stable. In Table 2, the 1 1 index k can be different for each vector u 0 so that k+1 , k+3 ( 1 , 2 ). 3.3 Comparison with SD The analysis in Section 3.2 shows that CBB with m = 2 is at best linearly convergent. By (3.4) and (3.10), we obtain gk+8 2 = gk 2 for all k 1, (3.24) where is the parameter in (3.10). The above relation implies that the convergence rate of the method only depends on the value . Furthermore, Table 2 tells us that this value of is related to the starting point. It may be very small or relatively large. The maximal possible value of is the in (3.21). In the 3D case, we get gk+1 2 3 1 gk 3 + 1 2 (3.25) for the SD method (see Akaike, 1959). It is not dif cult to show that < 3 1 3 + 1 4 . (3.26) Thus, we see that CBB with m = 2 is faster than the SD method if n = 3. This result could be extended to the arbitrary dimensions since we observe that CBB with m = 2 generates similar cycles for higherdimensional quadratics. The examples provided in Section 3.2 for CBB with m = 2 are helpful in understanding and analysing the behaviour of other non-monotone gradient methods. For example, we can also use the same technique to construct cyclic examples for the alternate step (AS) gradient method, at least theoretically. The AS method corresponds to the cyclic SD method (1.4) with m = 2. In fact, if we de ne u k as in (3.5), we obtain for all k 1 2k 1 = u 2k 1 u 2k 1 () () , u 2k+1 = (i) (1 2k 1 i )4 u 2k 1 (1 2k 1 )4 u 2k 1 () (i) (3.27) TABLE 2 Different choices of u 0 generate different cycles uT 0 (1, 2, 3) (1, 3, 2) (2, 1, 3) (2, 3, 1) (3, 1, 2) (3, 2, 1) 1 k+1 1 k+3 1 k+5 1 k+7 4.2186 10 6 1.2890 10 1 1.5024 10 2 1.3706 10 1 1.6018 10 3 9.4127 10 4 4.9103 3.2088 1.1099 1.5797 4.9846 1.0015 1.0000 1.3409 1.2764 2.0807 1.0026 4.9912 8.0000 6.9100 5.0197 5.7248 7.9086 7.8776 5.0008 7.2058 7.9938 7.7683 7.7458 7.8866 THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 619 for i = 1, . . . , n. For any n with u 2n+1 = u 1 and 2n+1 = 1 , we require the stepsizes { 2k 1 : k = 1, . . . , n 1} to satisfy n 1 bi j = , i=1 j = 1, . . . , n, (3.28) where bi j is given by (3.9). At the same time, we obtain the following linear equations for u 1 : i 1 T u 1 = 0, where T R(n 1) n with Ti j = bi j =1 b 4j . (3.29) The above system (3.29) has n variables, but n 1 equations. If there is a positive solution u 1 , then we () may scale the vector and obtain another positive solution u 1 = cu 1 with u 1 = 1, which completes the construction of a cyclic example. Here we present a 5D example. We rst x 1 = 1, 3 = 0.1, 5 = 0.2 and 7 = 0.0625, and then choose = (0.73477, 1.3452, 4.2721, 10.554, 16.154) which are ve roots of the equation 4 k=1 (1 2k 1 w) = 0.2. Therefore, we get the matrix 9.5537 15.154 461.26 32451 . 0.08696 16870 0.04056 1800.5 0.26523 0.34515 3.2721 0.00458 0.01228 65.659 T = 0.00311 0.00582 1.7964 0.00184 0.00208 0.00406 The system T u 1 = 0 has the positive solution u 1 = (5.6163 105 , 3.3397 105 , 7.3848 103 , 9.9533 102 , 1.0)T which leads to u 1 = (6.2128 10 1 , 3.6945 10 1 , 8.1693 10 3 , 1.1011 10 3 , 1.1062 10 6 )T . Therefore, if we choose the above initial vector u 1 , we get u 10k+1 = u 1 and 10k+1 = 1 for all k 1, and hence, the AS method falls into a cycle. Unlike CBB with m = 2, we have not found any cyclic example for the AS method which are numerically stable. 4. An adaptive cyclic Barzilai Borwein method In this section, we examine the convergence speed of CBB for different values of m [1, 7], using quadratic programming problems of the form: f (x) = 1T x Ax, 2 A = diag( 1 , . . . , n ). (4.1) We will see that the choice for m has a signi cant impact on performance. This leads us to propose an adaptive choice for m. The BB algorithm with this adaptive choice for m and a non-monotone line search is called adaptive cyclic Barzilai Borwein (ACBB). Numerical comparisons with SPG2 and with conjugate gradient codes using the CUTEr test problem library are given later in Section 4. 620 4.1 A numerical investigation of CBB Y.-H. DAI ET AL. We consider the test problem (4.1) with four different condition numbers C for the diagonal matrix, C = 102 , C = 103 , C = 104 and C = 105 , and with three different dimensions n = 102 , n = 103 and n = 104 . We let 1 = 1, n = C, the condition number. The other diagonal elements i , 2 i n 1, (i) are randomly generated on the interval (1, n ). The starting points x1 , i = 1, . . . , n, are randomly generated on the interval [ 5, 5]. The stopping condition is gk 2 10 8 . For each case, 10 runs are made and the average number of iterations required by each algorithm is listed in Table 3 (under the columns labelled BB and CBB). The upper bound for the number of iterations is 9999. If this upper bound is exceeded, then the corresponding entry in Table 3 is F. In Table 3, we see that m = 2 gives the worst numerical results in Section 3, we saw that as m increases, convergence became superlinear. For each case, a suitably chosen m drastically improves the ef ciency of the BB method. For example, in case of n = 102 and cond = 105 , CBB with m = 7 only requires one- fth of the iterations of the BB method. The optimal choice of m varies from one test case to another. If the problem condition is relatively small (cond = 102 , 103 ), a smaller value of m (3 or 4) is preferred. If the problem condition is relatively large (cond = 104 , 105 ), a larger value of m is more ef cient. This observation is the motivation for introducing an adaptive choice for m in the CBB method. Our adaptive idea arises from the following considerations. If a stepsize is used in nitely often in the gradient method, namely, k , then under the assumption that the function Hessian A has no T T multiple eigenvalues, the gradient gk must approximate an eigenvector of A, and gk Agk /gk gk tends to the corresponding eigenvalue of A (see Dai, 2003). Thus, it is reasonable to assume that repeated use of a BB stepsize leads to good approximations of eigenvectors of A. First, we de ne k = T gk Agk . gk Agk (4.2) TABLE 3 Comparing CBB(m) method with an ACBB method CBB n 102 Cond 102 103 104 105 102 103 104 105 102 103 104 105 BB 147 505 1509 5412 147 505 1609 5699 156 539 1634 6362 m=2 m=3 m=4 m=5 m=6 m=7 219 156 145 150 160 166 2715 468 364 376 395 412 F 1425 814 852 776 628 F 5415 3074 1670 1672 1157 274 160 158 162 166 181 1756 548 504 493 550 540 F 1862 1533 1377 1578 1447 F 6760 4755 3506 3516 2957 227 162 166 167 170 187 3200 515 551 539 536 573 F 1823 1701 1782 1747 1893 F 6779 5194 4965 4349 4736 Adaptive M =5 136 367 878 2607 150 481 1470 4412 156 497 1587 4687 M = 10 134 349 771 1915 145 460 1378 3187 156 505 1517 4743 103 104 THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 621 If gk is exactly an eigenvector of A, we know that k = 1. If k 1, then gk can be regarded as an SD BB approximation of an eigenvector of A and k k . In this case, it is worthwhile to calculate a new BB so that the method accepts a step close to the SD step. Therefore, we test the condition BB stepsize k k , (4.3) where (0, 1) is constant. If the above condition holds, we calculate a new BB stepsize. We also introduce a parameter M, and if the number of cycles m > M, we calculate a new BB stepsize. Numerical results for this ACBB with = 0.95 are listed under the column Adaptive of Table 3, where two values of M = 5, 10 are tested. From Table 3, we see that the adaptive strategy makes sense. The performance with M = 5 or M = 10 is often better than that of the BB method. This motivates the use of a similar strategy for designing an ef cient gradient algorithms for unconstrained optimization. 4.2 Non-monotone line search and cycle number As mentioned in Section 1, the choice of the stepsize k is very important for the performance of a gradient method. For the BB method, function values do not decrease monotonically. Hence, when implementing BB or CBB, it is important to use a non-monotone line search. T Assuming that dk is a descent direction at the kth iteration (gk dk < 0), a common termination condition for the steplength algorithm is f (xk + k dk ) T fr + k gk dk , (4.4) where fr is the so-called reference function value and (0, 1) a constant. If fr = f (xk ), then the line search is monotone since f (xk+1 ) < f (xk ). The non-monotone line search proposed in Grippo et al. (1986) chooses fr to be the maximum function value for the M most recent iterates. That is, at the kth iteration, we have fr = f max = 0 i min{k,M 1} max f (xk i ). (4.5) This non-monotone line search is used by Raydan (1997) to obtain GBB. Dai & Schittkowski (2005) extended the same idea to a sequential quadratic programming method for general constrained nonlinear optimization. An even more adaptive way of choosing fr is proposed by Toint (1997) for trust region algorithms and then extended by Dai & Zhang (2001). Compared with (4.5), the new adaptive way of choosing fr allows big jumps in function values, and is therefore very suitable for the BB algorithm (see Dai & Fletcher, 2005b, 2006; Dai & Zhang, 2001). The numerical results which we report in this section are based on the non-monotone line search algorithm given in Dai & Zhang (2001). The line search in our paper differs from the line search in Dai & Zhang (2001) in the initialization of the stepsize. Here, the starting guess for the stepsize coincides with the prior BB step until the cycle length has been reached; at which point, we recompute the step using the BB formula. In each subsequent subiteration, after computing a new BB step, we replace (4.4) with f (xk + k dk ) min{ f max , fr } + k gk dk , T where fr is the reference value given in Dai & Zhang (2001) and k is the initial trial stepsize (the previous BB step). It is proved in Dai & Zhang (2001, Theorem 3.2) that the criteria given for choosing 622 Y.-H. DAI ET AL. the non-monotone stepsize ensures convergence in the sense that lim inf gk = 0. k We now explain how we decided to terminate the current cycle, and recompute the stepsize using the BB formula. Note that the re-initialization of the stepsize has no effect on convergence, it only affects the initial stepsize used in the line search. Loosely, we would like to compute a new BB step in any of the following cases: R1. The number of times m the current BB stepsize has been reused is suf ciently large: m where M is a constant. R2. The following non-quadratic analogue of (4.3) is satis ed: T sk yk sk 2 yk M, 2 , (4.6) where < 1 is near 1. We feel that Condition (4.6) should only be used in a neighbourhood of a local minimizer, where f is approximately quadratic. Hence, we only use Condition (4.6) when the stepsize is suf ciently small: sk 2 < min c1 f k+1 ,1 , gk+1 (4.7) where c1 is a constant. R3. The current step sk is suf ciently large: sk 2 max c2 f k+1 ,1 , gk+1 (4.8) where c2 is a constant. R4. In the previous iteration, the BB step was truncated in the line search. That is, the BB step had to be modi ed by the non-monotone line search routine to ensure convergence. Nominally, we recompute the BB stepsize in any of the cases R1 R4. One case where we prefer to retain the current stepsize is the case where the iterates lie in a region where f is not strongly convex. T Note that if sk yk < 0, then there exists a point between xk and xk+1 where the Hessian of f has negative eigenvalues. In detail, our rules for terminating the current cycle and re-initializing the BB stepsize are the following conditions. 4.2.1 Cycle termination/stepsize initialization. T T1. If any of the conditions R1 through R4 are satis ed and sk yk > 0, then the current cycle is terminated and the initial stepsize for the next cycle is given by T sk sk T sk yk k+1 = max min , min where min < max are xed constants. , max , THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 623 T2. If the length m of the current cycle satis es m 1.5M, then the current cycle is terminated and the initial stepsize for the next cycle is given by k+1 = max{1/ gk+1 , k }. T Condition T2 is a safeguard for the situation where sk yk < 0 in a series of iterations. 4.3 Numerical results In this subsection, we compare the performance of our ACBB stepsize algorithm, denoted ACBB, with the SPG2 algorithm of Birgin et al. (2000, 2001), with the PRP+ conjugate gradient code developed by Gilbert & Nocedal (1992) and with the CG DESCENT code of Hager & Zhang (2005b, to appear). The SPG2 algorithm is an extension of Raydan s (1997) GBB algorithm which was downloaded from the TANGO web page maintained by Ernesto Birgin. In our tests, we set the bounds in SPG2 to in nity. The PRP+ code is available at http://www.ece.northwestern.edu/ nocedal/software.html. The CG DESCENT code is found at http://www.math.u .edu/ hager/papers/CG. The line search in the PRP+ code is a modi cation of subroutine CSRCH of Mor & Thuente (1994), which employs e various polynomial interpolation schemes and safeguards in satisfying the strong Wolfe conditions. CG DESCENT employs an approximate Wolfe line search. All codes are written in Fortran and compiled with f77 under the default compiler settings on a Sun workstation. The parameters used by CG DESCENT are the default parameter values given in Hager & Zhang (2006) for version 1.1 of the code. For SPG2, we use parameter values recommended on the TANGO web page. In particular, the steplength was restricted to the interval [10 30 , 1030 ], while the memory in the non-monotone line search was 10. The parameters of the ACBB algorithm are min = 10 30 , max = 1030 , c1 = c2 = 0.1 and M = 4. For the initial iteration, the starting stepsize for the line search was 1 = 1/ g1 . The parameter values for the non-monotone line search routine from Dai & Zhang (2001) were = 10 4 , 1 = 0.1, 2 = 0.9, = 0.975, L = 3, M = 8 and P = 40. Our numerical experiments are based on the entire set of 160 unconstrained optimization problem available from CUTEr in the fall of 2004. As explained in Hager & Zhang (2006), we deleted problems that were small or problems where different solvers converged to different local minimizers. After the deletion process, we were left with 111 test problems with dimension ranging from 50 to 104 . Nominally, our stopping criterion was the following: f (xk ) max{10 6 , 10 12 f (x0 ) }. (4.9) In a few cases, this criterion was too lenient. For example, with the test problem PENALTY1, the computed cost still differs from the optimal cost by a factor of 105 when Criterion (4.9) is satis ed. As a result, different solvers obtain completely different values for the cost, and the test problem would be discarded. By changing the convergence criterion to f (xk ) 10 6 , the computed costs all agreed to six digits. The problems for which the convergence criterion was strengthened were DQRTIC, PENALTY1, POWER, QUARTC and VARDIM. The CPU time in seconds and the number of iterations, function evaluations and gradient evaluations for each of the methods are posted at the following web site: http://www.math.u .edu/ hager/papers/CG. Here we analyse the performance data using the pro les of Dolan & Mor (2002). That is, we plot the e fraction p of problems for which any given method is within a factor of the best time. In a plot of 624 Y.-H. DAI ET AL. FIG. 2. Performance based on CPU time. TABLE 4 Number of times each method was fastest (time metric, stopping criterion (4.9)) Method CG DESCENT ACBB PRP+ SPG2 Fastest 70 36 9 9 performance pro les, the top curve is the method that solved the most problems in a time that was within a factor of the best time. The percentage of the test problems for which a method is the fastest is given on the left axis of the plot. The right side of the plot gives the percentage of the test problems that were successfully solved by each of the methods. In essence, the right side is a measure of an algorithm s robustness. In Fig. 2, we use CPU time to compare the performance of the four codes ACBB, SPG2, PRP+ and CG DESCENT. Note that the horizontal axis in Fig. 2 is scaled proportional to log2 ( ). The best performance, relative to the CPU time metric, was obtained by CG DESCENT, the top curve in Fig. 2, followed by ACBB. The horizontal axis in the gure stops at = 16 since the plots are essentially at for larger values of . For this collection of methods, the number of times any method achieved the best time is shown in Table 4. The column total in Table 4 exceeds 111 due to ties for some test problems. The results of Fig. 2 indicate that ACBB is much more ef cient than SPG2, while it performed better than PRP+, but not as well as CG DESCENT. From the experience in Raydan (1997), the GBB algorithm, with a traditional non-monotone line search (Grippo et al., 1986), may be affected signi cantly by nearly singular Hessians at the solution. We observe that nearly singular Hessians do not affect ACBB signi cantly. In fact, Table 3 also indicates that ACBB becomes more ef cient as the problem becomes more singular. Furthermore, since ACBB does not need to calculate the BB stepsize at every iteration, CPU time is saved, which can be signi cant when the problem dimension is large. For this test set, we found that the average cycle length for ACBB was 2.59. In other words, the BB step is re-evaluated after two or three iterations, on average. This memory length is smaller than the memory length that works THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 625 TABLE 5 CPU times for selected problems Problem FLETCHER FLETCHER BDQRTIC VARDIM VARDIM Dimension 5000 1000 1000 10000 5000 ACBB 9.14 1.32 0.37 0.05 0.02 CG DESCENT 989.55 27.27 3.40 2.13 0.92 well for quadratic function. When the iterates are far from a local minimizer of a general non-linear function, the iterates may not behave like the iterates of a quadratic. In this case, better numerical results are obtained when the BB stepsize is updated more frequently. Even though ACBB did not perform as well as CG DESCENT for the complete set of test problems, there were some cases where it performed exceptionally well (see Table 5). One important advantage of the ACBB scheme over conjugate gradient routines such as PRP+ or CG DESCENT is that in many cases, the stepsize for ACBB is either the previous stepsize or the BB stepsize (1.5). In contrast, with conjugate gradient routines, each iteration requires a line search. Due to the simplicity of the ACBB stepsize, it can be more ef cient when the iterates are in a regime where the function is irregular and the asymptotic convergence properties of the conjugate gradient method are not in effect. One such application is bound-constrained optimization problems as components of x reach the bounds, these components are often held xed, and the associated partial derivative change discontinuously. In Hager & Zhang (2005a) ACBB is combined with CG DESCENT to obtain a very ef cient active set algorithm for box-constrained optimization problems. 5. Conclusion and discussion In this paper, we analyse the CBB method. For general non-linear functions, we prove linear convergence. For convex quadratic functions, our numerical results indicate that when m > n/2 3, CBB is likely to be R-superlinear. For the special case n = 3 and m = 2, the convergence rate, in general, is no better than linear. By utilizing non-monotone line search techniques, we develop an ACBB stepsize algorithm for general non-linear unconstrained optimization problems. The test results in Fig. 2 indicate that ACBB is signi cantly faster than SPG2. Since the mathematical foundations of ACBB and the conjugate gradient algorithms are completely different, the performance seems to depend on the problem. Roughly speaking, if the objective function is close to quadratic, the conjugate gradient routines seem to be more ef cient; if the objective function is highly non-linear, then ACBB is comparable to or even better than conjugate gradient algorithms. Acknowledgements Constructive and detailed comments by the referees are gratefully acknowledged and appreciated. Y-HD was supported by the Alexander von Humboldt Foundation under grant CHN/1112740 STP and Chinese National Science Foundation grants 10171104 and 40233029. WWH and HZ were supported by US National Science Foundation grant no. 0203270. REFERENCES AKAIKE, H. (1959) On a successive transformation of probability distribution and its application to the analysis of the optimum gradient method. Ann. Inst. Stat. Math. Tokyo, 11, 1 17. 626 Y.-H. DAI ET AL. BARZILAI, J. & BORWEIN, J. M. (1988) Two point step size gradient methods. IMA J. Numer. Anal., 8, 141 148. BIRGIN, E. G., CHAMBOULEYRON, I. & MART NEZ, J. M. (1999) Estimation of the optical constants and the I thickness of thin lms using unconstrained optimization. J. Comput. Phys., 151, 862 880. BIRGIN, E. G., MART NEZ, J. M. & RAYDAN, M. (2000) Nonmonotone spectral projected gradient methods for I convex sets. SIAM J. Optim., 10, 1196 1211. BIRGIN, E. G., MART NEZ, J. M. & RAYDAN, M. (2001) Algorithm 813: SPG software for convex-constrained I optimization. ACM Trans. Math. Softw., 27, 340 349. BONGARTZ, I., CONN, A. R., GOULD, N. I. M. & TOINT, P. L. (1995) CUTE: constrained and unconstrained testing environments. ACM Trans. Math. Softw., 21, 123 160. CAUCHY, A. (1847) M thode g n rale pour la r solution des syst` ms d equations simultan es. Comp. Rend. Sci. e ee e e e Paris, 25, 46 89. DAI, Y. H. (2003) Alternate stepsize gradient method. Optimization, 52, 395 415. DAI, Y. H. & FLETCHER, R. (2005a) On the asymptotic behaviour of some new gradient methods. Math. Prog., 103, 541 559. DAI, Y. H. & FLETCHER, R. (2005b) Projected Barzilai-Borwein methods for large-scale box-constrained quadratic programming. Numer. Math., 100, 21 47. DAI, Y. H. & FLETCHER, R. (2006) New algorithms for singly linearly constrained quadratic programs subject to lower and upper bounds. Math. Prog., 106, 403 421. DAI, Y. H. & SCHITTKOWSKI, K. (2005) A sequential quadratic programming algorithm with non-monotone line search. Technical Report. Department of Mathematics, University of Bayreuth (submitted). DAI, Y. H. & YANG, X. Q. (2006) A new gradient method with an optimal stepsize property. Comput. Optim. Appl., 33, 73 88. DAI, Y. H. & YUAN, Y. (2003) Alternate minimization gradient method. IMA J. Numer. Anal., 23, 377 393. DAI, Y. H. & ZHANG, H. (2001) An adaptive two-point stepsize gradient algorithm. Numer. Algorithms, 27, 377 385. DOLAN, E. D. & MOR E , J. J. (2002) Benchmarking optimization software with performance pro les. Math. Program., 91, 201 213. FORSYTHE, G. E. (1968) On the asymptotic directions of the s-dimensional optimum gradient method. Numer. Math., 11, 57 76. FRIEDLANDER, A., MART NEZ, J. M., MOLINA, B. & RAYDAN, M. (1999) Gradient method with retards and I generalizations. SIAM J. Numer. Anal., 36, 275 289. GILBERT, J. C. & NOCEDAL, J. (1992) Global convergence properties of conjugate gradient methods for optimization. SIAM J. Optim., 2, 21 42. GLUNT, W., HAYDEN, T. L. & RAYDAN, M. (1993) Molecular conformations from distance matrices. J. Comput. Chem., 14, 114 120. GRIPPO, L., LAMPARIELLO, F. & LUCIDI, S. (1986) A nonmonotone line search technique for Newton s method. SIAM J. Numer. Anal., 23, 707 716. GRIPPO, L. & SCIANDRONE, M. (2002) Nonmonotone globalization techniques for the Barzilai-Borwein gradient method. Comput. Optim. Appl., 23, 143 169. HAGER, W. W. & ZHANG, H. (2005a) A new active set algorithm for box constrained optimization. SIAM J. Optim. (to appear). HAGER, W. W. & ZHANG, H. (2005b) A new conjugate gradient method with guaranteed descent and an ef cient line search. SIAM J. Optim., 16, 170 192. HAGER, W. W. & ZHANG, H. (2006) Algorithm 851: CG DESCENT, a conjugate gradient method with guaranteed descent. ACM Trans. Math. Softw., 32. LIU, W. B. & DAI, Y. H. (2001) Minimization algorithms based on supervisor and searcher cooperation. J. Optim. Theory Appl., 111, 359 379. THE CBB METHOD FOR UNCONSTRAINED OPTIMIZATION 627 MOR E , J. J. & THUENTE, D. J. (1994) Line search algorithms with guaranteed suf cient decrease. ACM Trans. Math. Softw., 20, 286 307. RAYDAN, M. (1997) The Barzilai and Borwein gradient method for the large scale unconstrained minimization problem. SIAM J. Optim., 7, 26 33. RAYDAN, M. & SVAITER, B. F. (2002) Relaxed steepest descent and Cauchy-Barzilai-Borwein method. Comput. Optim. Appl., 21, 155 167. SERAFINI, T. ZANGHIRATI, G. & ZANNI, L. (2005) Gradient projection methods for quadratic programs and applications in training support vector machines. Optim. Methods Softw., 20, 353 378. TOINT, P. L. (1997) A non-monotone trust region algorithm for nonlinear optimization subject to convex constraints. Math. Prog., 77, 69 94.

Find millions of documents here - Study Guides, Homework Solutions, Papers, Exam Answer Keys and more. Course Hero has millions of course related materials that will enable you to learn better, faster and get an A in all your courses.
Below is a small sample set of documents: