An exact algorithm for graph partitioning

AN EXA CT ALGORITHM F OR GRAPH P AR TITIONING ∗ WILLIAM W. HAG ER † , DZUNG T. PHAN ‡ , AND HONGCHAO ZHANG § Abstract. An exact algorithm is presen ted f or solving edge we ighted graph partitioning prob- lems. The algorithm is based on a branch and b ound method applied to a cont inuous quadratic programming formulation of the problem. Low er b ounds are obtained b y decomposing the ob jective function into con v ex and conc av e parts an d replac ing the c onca ve part by an affine underestimate. It is sho wn that the b est affine underestimate can b e expressed in terms of the cen ter and the radius of the small est sphere con taining the feasible set. The concav e term is obtained either by a constan t diagonal shift asso ciated with the smallest eigenv alue of the ob jectiv e function Hessian, or by a di- agonal shift obtained b y solving a s emidefinite programming problem. Num er ical results sho w that the proposed algorithm is competitive with state-of-the-art graph partitioning codes. AMS sub ject cla ssifications. 90C35, 90C20, 90C27, 90C46 Key w ords. graph partitioning, min-cut, quadratic programming, branch and b ound, affine underestimate 1. In tro duction. Given a graph with edge weights, the g raph partitioning prob- lem is to partition the vertices in to t wo sets satisfying sp ecified size constraints, while minimizing the s um of the w eights of the edges that connect the vertices in the tw o sets. Graph pa rtitioning problems a rise in many are a s including VLSI design, data mining, parallel computing, and spa rse matrix factor izations [14, 22, 27, 36]. The graph partitioning problem is NP-hard [11]. There are tw o general classes of metho ds for the graph partitioning pro blem, exact metho ds whic h compute the optimal partition, and heuristic metho ds which try to quickly compute an a pproximate so lution. Heuristic metho ds include sp ectral metho ds [18], g e o metric methods [12], multilev el schemes [19], optimiza tio n-based metho ds [8], and methods that employ randomizatio n techniques such as genetic algorithms [34]. Softw a re which implements heuristic methods includes Metis ([24, 25, 26]), Chaco [17], Part y [32], PaT oH [4], SCOTCH [31], J ostle [3 7], Zoltan [6 ], and HUND [13]. This pap er develops an exact algo rithm for the graph partitioning pro blem. In earlier work, Brunetta, Conforti, a nd Rinaldi [3] propo s e a branch-and-cut scheme based on a linear prog ramming relax a tion and subseque nt cuts based on separatio n techn iques. A column g e ne r ation approa ch is developed by Johnson, Mehrotra, and Nemhauser [22], while Mitchell [28] develops a po lyhedral appro ach. Karisch, Rendl, and Clausen [23] de velop a bra nch-and-bound metho d utilizing a semidefinite pro- gramming relaxation to obtain a low er b ound. Sensen [33] develops a br anch-and- bo und metho d based on a lower b ound obtained by solving a mult icommo dity flow problem. ∗ Nov ember 09, 2009. Thi s material is based upon wo rk supported b y the National Science F oundation under Gran t 0620286. † hager@ma th.ufl.edu , http ://www.math.ufl.edu/ ∼ hager , PO Bo x 11810 5, Departmen t of Mathe- matics, Univ ersity of Fl orida, Gaine sville, FL 32611-8105. Phone (352) 392-0281. F ax (352) 392-8357. ‡ dphan@ma th.ufl.edu , h ttp://www.math.ufl.edu/ ∼ dphan , PO Box 118105, Department of Math- ematics, Universit y of Flor i da, Gainesville, FL 32611-8105. Phone (352) 392-0281. F ax (352) 392- 8357. § hozhang@ math.lsu.ed u , ht tp://www.math.lsu.edu/ ∼ hozhang, Departmen t of Mathemat ics, 140 Lock ett Hall , Center f or Computation and T ec hnology , Louisiana State Universit y , Baton Rouge, LA 70803-4918. Phone (225) 578-1982. F ax (225) 578-4276. 1 2 W. W. HAGER, D. T. PHAN, H. ZHANG In this paper , we develop a branch-and-bound algorithm based o n a quadratic progra mming (QP) for m ulation of the graph pa rtitioning problem. The ob jective function o f the QP is express ed as the sum of a conv ex and a co ncav e function. W e cons ide r t wo different tec hniques for ma king this decomposition, one based on eigenv a lues and the other based o n se midefinite progr amming. In e a ch case, we give an a ffine under estimate for the conca ve function, which leads to a tra ctable lo wer bo und in the bra nch and b ound algo rithm. The pap er is or ganized as follows. In Section 2 we r eview the c ontin uo us q uadratic progra mming formulation of the graph partitioning problem dev elop ed in [14] and w e explain ho w to asso ciate a solution of the contin uous problem with the solution to the discrete pro blem. In Section 3 we discuss approaches for decompo sing the o b jectiv e function for the QP int o the sum of co nv e x and a concave functions , a nd in ea ch c a se, we show how to gener a te an affine low er bo und for the concav e pa rt. Section 4 gives the bra nch-and-bound alg o rithm, while Section 5 provides necessary and sufficien t conditions for a lo cal minimizer. Section 6 compa res the p erfor mance of the new branch-and-bo und algo rithm to earlier results given in [23] and [33]. Notation. Throughout the pap er , k · k de no tes the Euclidian norm. 1 is the vector whose entries are a ll 1. The dimension will be clea r from context. If A ∈ R n × n , A  0 means that A is p o sitive semidefinite. W e let e i denote the i -th co lumn of the ident ity matrix ; a gain, the dimension will be clear from context. If S is a set, then |S | is the num b er of elements in S . The gradient ∇ f ( x ) is a row vector. 2. Con tin uous quadratic programming form ulation. L e t G b e a graph with n vertices V = { 1 , 2 , · · · , n } , and let a ij be a w eight asso cia ted with the e dg e ( i , j ). When there is no edg e b etw een i a nd j , we set a ij = 0. F or each i and j , we a ssume that a ii = 0 a nd a ij = a j i ; in other w ords, w e co nsider an undirected gra ph without self lo ops (a simple, undirec ted graph). The sign of the w eights is not re stricted, and in fact, a ij could b e negative, as it would be in the max-cut problem. Giv en integers l and u such that 0 ≤ l ≤ u ≤ n , we wish to partition the vertices into tw o disjoint sets, with b etw een l a nd u vertices in one set, while minimizing the sum of the weigh ts as so ciated with e dg es co nnecting vertices in different s ets. The edges connecting the tw o sets in the pa rtition are r eferred to a s the cut edges, and the optimal partition minimizes the sum of the weigh ts of the cut edg es. Hence, the gr aph par titioning pr oblem is also called the min- cut problem. In [1 4] we show that for a suitable choice of the diago na l matrix D , the graph partitioning problem is eq uiv alent to the following contin uous q uadratic prog ramming problem: minimize f ( x ) := ( 1 − x ) T ( A + D ) x sub ject to 0 ≤ x ≤ 1 , l ≤ 1 T x ≤ u, (2.1) where A is the matrix with elements a ij . Supp ose x is binary a nd let us define the sets V 0 = { i : x i = 0 } and V 1 = { i : x i = 1 } . (2.2) It can b e c heck ed that f ( x ) is the sum of the weights of the cut edg es asso cia ted with the partition (2.2). Hence, if we add the restrictio n tha t x is binary , then (2 .1) GRAPH P AR TITIONING 3 is exactly equiv alent to finding the par tition which minimizes the weigh t of the cut edges. Note, thoug h, that there are no binary constraints in (2.1). T he equiv a lence betw een (2.1) a nd the graph partitioning pro blem is a s follows (see [14, Thm. 2.1]): Theorem 2.1. If t he diagonal matrix D is chosen so that d ii + d j j ≥ 2 a ij and d ii ≥ 0 (2.3) for e ach i and j , then ( 2.1 ) has a binary solut ion x and the p artition given by ( 2.2 ) is a min-cut. The g eneralizatio n of this r esult to m ultiset pa r titioning is given in [15]. The condition (2.3) is sa tis fied, for ex ample, by the choice d j j = ma x { 0 , a 1 j , a 2 j , . . . , a nj } for each j . The pro of of Theorem 2 .1 was based on showing that any solution to (2.1) could be tr a nsformed to a binary solutio n witho ut c hanging the ob jective function v alue. With a mo dification of this idea, an y fea sible p oint can b e tr ansformed to a binary feas ible p oint without increa sing the o b jective function v alue. W e now g ive a constructive pro of of this re sult, w hich is used when we solve (2.1). Corollar y 2.2. If x is fe asible in ( 2.1 ) and the dia gonal matrix D satisfies ( 2.3 ) , then ther e exists a binary y with f ( y ) ≤ f ( x ) and y i = x i whenever x i is binary. Pr o of . W e first show how to find z w ith the pr op erty that z is feasible in (2.1 ), f ( z ) ≤ f ( x ), 1 T z is int eger, and the only comp onents of z and x which differ ar e the fractional compo nents o f x . If 1 T x = u o r 1 T x = l , then w e are do ne since l and u ar e int egers ; hence, we assume that l < 1 T x < u . If all comp onents of x are binary , then we are done, s o supp o se that there ex ists a no nbinary compo nent x i . Since a ii = 0, a T aylor expa nsion o f f gives f ( x + α e i ) = f ( x ) + α ∇ f ( x ) i − α 2 d ii , where e i is the i -th column of the ident ity matrix. The q ua dratic term in the expansio n is no npo sitive since d ii ≥ 0. If the first der iv ative term is negative, then increa se α ab ov e 0 until either x i + α b ecomes 1 or 1 T x + α is an integer. Since the first deriv ative term is neg ative and α > 0 , f ( x + α e i ) < f ( x ). If 1 T x + α b eco mes an integer, then we a r e done. If x i + α b eco mes 1, then w e reach a po int x 1 with o ne more binary comp onent and with an ob jective function v alue no larger than f ( x ). If the first deriv ative term is nonnegative, then decrea se α b elow 0 un til either x i + α b eco mes 0 or 1 T x + α is a n integer. Since the first der iv ative term is no nnegative and α < 0, f ( x + α e i ) ≤ f ( x ). If 1 T x + α b eco mes an integer, then we a re done. If x i + α b ecomes 0, then we r each a point x 1 with one more binar y comp onent and w ith a smaller v alue for the cost function. In this latter ca se, we cho o se ano ther nonbinary co mpo nent of x 1 and rep ea t the pr o cess. Hence, there is no loss of generality in assuming that 1 T x is an integer. Suppo se that x is not binar y . Since 1 T x is an integer, x m ust hav e a t least tw o nonbinary comp onents, say x i and x j . Aga in, expanding f is a T aylor series gives f ( x + α ( e i − e j )) = f ( x ) + α ( ∇ f ( x ) i − ∇ f ( x ) j ) + α 2 (2 a ij − d ii − d j j ) . By (2.3), the quadra tic term is no npo sitive for any choice of α . If the first deriv a tive term is negative, then we incr ease α ab ov e 0 until either x i + α reaches 1 or x j − α reach 4 W. W. HAGER, D. T. PHAN, H. ZHANG 0. Since the first deriv a tive term is negative and α > 0, we hav e f ( x + α ( e i − e j )) < f ( x ). If the first deriv ative ter m is no nneg ative, then we decr e a se α b elow 0 until either x i + α r eaches 0 or x j − α r each 1. Since the first deriv a tive term is nonnegative and α < 0, it follows that f ( x + α ( e i − e j )) ≤ f ( x ). In either ca se, the v alue of the cost function do es not increa s e, and we reach a feasible p oint x 1 with 1 T x 1 int eger and with a t least o ne more binar y comp onent. If x 1 is not binar y , then x 1 m ust have at lea st tw o nonbinary comp onents; hence, the adjustment pro cess can b e contin ued un til all the comp o nents of x are bina ry . These adjustments to x do no t increa s e the v alue o f the cost function and we only alter the fractio nal comp onents o f x . This completes the pro o f. 3. Con v ex low er b ounds for the ob jectiv e function. W e compute an exa ct solution to the contin uous formu lation (2.1) of gr aph partitioning problem using a branch and b ound algorithm. The b o unding pr o cess r equires a low er b ound for the ob jective function when restricted to the intersection of a b ox and t wo ha lf spaces. This lower b o und is obtained by writing the ob jective function as the sum of a conv ex and a conc ave function and by replacing the co ncav e part by the b est affine underes- timate. Two different strateg ies a re g iven for decomp osing the ob jective function. 3.1. Lo w er bound based on minim um eige nv alue. Let us decomp ose the ob jective function f ( x ) = ( 1 − x ) T ( A + D ) x in the fo llowing wa y: f ( x ) = ( f ( x ) + σ k x k 2 ) − σ k x k 2 , where σ is the maximum o f 0 and the larges t eigenv alue o f A + D . This represents a DC (difference conv ex) decomp osition (see [2 0]) since f ( x ) + σ k x k 2 and σ k x k 2 are bo th conv ex. The concav e ter m −k x k 2 is under estimated by an affine function ℓ to obtain a conv ex underestimate f L of f g iven by f L ( x ) =  f ( x ) + σ k x k 2  + σ ℓ ( x ) . (3.1) W e now conside r the problem of finding the b est affine underestimate ℓ for the con- cav e function −k x k 2 ov er a given co mpa ct, conv ex set denoted C . The set of affine underestimators for −k x k 2 is given by S 1 = { ℓ : R n → R such that ℓ is a ffine and − k x k 2 ≥ ℓ ( x ) for a ll x ∈ C } . The b est affine underes timate is a solution of the problem min ℓ ∈S 1 max x ∈C −  k x k 2 + ℓ ( x )  . (3.2) The follo wing result generalizes Theor em 3 .1 in [16] where w e determine the best affine underestimate for − k x k 2 ov er an ellipso id. Theorem 3. 1. L et C ⊂ R n b e a c omp act, c onvex set and let c b e the c enter and r b e the ra dius of the smal lest spher e c ontaining C . Th is smal lest spher e is un ique and a solution of ( 3.2 ) is ℓ ∗ ( x ) = − 2 c T x + k c k 2 − r 2 . F u rthermor e, min ℓ ∈S 1 max x ∈C −  k x k 2 + ℓ ∗ ( x )  = r 2 . GRAPH P AR TITIONING 5 Pr o of . T o b egin, we will show tha t the minimization in (3.2) can be restricted to a compact set. Clea rly , when carr y ing out the minimization in (3.2), we should restr ic t our attention to thos e ℓ which touch the function h ( x ) := −k x k 2 at some p oint in C . Let y ∈ C denote the po int o f contact. Since h ( x ) ≥ ℓ ( x ) a nd h ( y ) = ℓ ( y ), a low er bo und for the err or h ( x ) − ℓ ( x ) over x ∈ C is h ( x ) − ℓ ( x ) ≥ | ℓ ( x ) − ℓ ( y ) | − | h ( x ) − h ( y ) | . If M is the difference betw een the max imum and minim um v a lue of h ov er C , then we hav e h ( x ) − ℓ ( x ) ≥ | ℓ ( x ) − ℓ ( y ) | − M . (3.3) An upp er b ound for the minimum in (3.2) is obtained by the linear function ℓ 0 which is co nstant on C , with v alue equal to the minimum o f h ( x ) over x ∈ C . If w is a po int where h attains its minimum over C , then we hav e max x ∈C h ( x ) − ℓ 0 ( x ) = max x ∈C h ( x ) − h ( w ) = M . Let us restrict o ur attention to the linear functions ℓ which achiev e an ob jective function v alue in (3.2) which is at least as small as that of ℓ 0 . F or these ℓ a nd for x ∈ C , we hav e h ( x ) − ℓ ( x ) ≤ max x ∈C h ( x ) − ℓ ( x ) ≤ max x ∈C h ( x ) − ℓ 0 ( x ) = M . (3.4) Combining (3.3) and (3.4) g ives | ℓ ( x ) − ℓ ( y ) | ≤ 2 M . (3.5) Thu s, when we c a rry out the minimization in (3.2), we should r estrict our a ttent ion to linear functions which touch h a t some p oint y ∈ C and with the c hange in ℓ acr oss C satisfying the b o und (3.5) for all x ∈ C . This tells us that the minimization in (3.2) can be restric ted to a compact set, and that a minimizer must exist. Suppo se that ℓ attains the minimum in (3 .2). Let z b e a p oint in C where h ( x ) − ℓ ( x ) achieves its maximum. A T aylor e xpansion around x = z gives h ( x ) − ℓ ( x ) = h ( z ) − ℓ ( z ) + ( ∇ h ( z ) − ∇ ℓ )( x − z ) − k x − z k 2 . Since ℓ ∈ S 1 , h ( x ) − ℓ ( x ) ≥ 0 for all x ∈ C . It follows that h ( z ) − ℓ ( z ) ≥ − ( ∇ h ( z ) − ∇ ℓ )( x − z ) + k x − z k 2 . (3.6) Since C is conv ex, the first-o rder optimality conditions for z g ive ( ∇ h ( z ) − ∇ ℓ )( x − z ) ≤ 0 for all x ∈ C . It follows from (3.6) that h ( z ) − ℓ ( z ) ≥ k x − z k 2 (3.7) for all x ∈ C . There exists x ∈ C suc h that k x − z k ≥ r o r else z would b e the center of a smaller sphere containing C . Hence, (3.7 ) implies that h ( z ) − ℓ ( z ) ≥ r 2 . 6 W. W. HAGER, D. T. PHAN, H. ZHANG ||x−c’|| < ||x − c|| C H x c c’ Fig. 3.1 . Supp ose c 6∈ C It follows that max x ∈C h ( x ) − ℓ ( x ) ≥ h ( z ) − ℓ ( z ) ≥ r 2 . (3.8) W e now observe that for the sp ecific linear function ℓ ∗ given in the sta temen t of the theorem, (3.8 ) beco mes an equality , which implies the optimality of ℓ ∗ in (3.2). Expand h in a T aylor series ar ound x = c to obtain h ( x ) = −k c k 2 − 2 c T ( x − c ) − k x − c k 2 = − 2 c T x + k c k 2 − k x − c k 2 . Subtract ℓ ∗ ( x ) = − 2 c T x + k c k 2 − r 2 from bo th sides to o btain h ( x ) − ℓ ∗ ( x ) = r 2 − k x − c k 2 . (3.9) If c ∈ C , then the maximum in (3.9) ov er x ∈ C is a ttained by x = c for which h ( c ) − ℓ ∗ ( c ) = r 2 . Consequently , (3.8) b eco mes a n eq ua lity fo r ℓ = ℓ ∗ , which implies the optimality of ℓ ∗ in (3.2). W e can show that c ∈ C as follo ws: Supp ose c 6∈ C . Since C is compa ct and conv ex, there exists a hyperplane H strictly separating c and C – see Figur e 3.1 If c ′ is the pro jection of c o nto H , then k x − c ′ k < k x − c k for all x ∈ C . (3.10) Let x ′ ∈ C b e the p oint which is far thest from c ′ and let x ∈ C b e the p oint far thest from c . Hence, k x − c k = r . By (3.10), we have k x ′ − c ′ k < k x − c k = r ; it fo llows that the spher e with center c ′ and ra dius k x ′ − c ′ k co ntains C and ha s ra dius sma lle r than r . This contradicts the assumption that r w as the sphere of smallest radius containing C . GRAPH P AR TITIONING 7 The uniqueness of the smallest sphere containing C is as follo ws: Supp o se that there exist t wo different smallest s pheres S 1 and S 2 containing C . Let S 3 be the smallest spher e co nt aining S 1 ∩ S 2 . Since the diameter o f the intersection is strictly less than the diameter of S 1 or S 2 , we co ntradict the ass umption that S 1 and S 2 were spheres of smallest ra dius co ntaining C . Remark 1. A lthough the smal lest spher e c ontaining C in The or em 3.1 is unique, the b est line ar u n der estimator of h ( x ) = −k x k 2 is not unique. F or example, supp ose a and b ∈ R n and C is the line se gment C = { x ∈ R n : x = α a + (1 − α ) b , α ∈ [0 , 1] } . Alo ng this line se gment, h is a c onc ave quadr atic in one variable. The b est affine under est imate along the line se gment c orr esp onds to the line c onne cting the ends of the quadr atic r est r icte d to the line se gment. Henc e, in R n +1 , any hyp erplane which c ontains t he p oints ( h ( a ) , a ) and ( h ( b ) , b ) le ads to a b est affine under est imate. Remark 2. L et C b e the b ox B = { x ∈ R n : p ≤ x ≤ q } . The diameter of B , the distanc e b etwe en the p oints in B with gr e atest sep ar ation, is k p − q k . Henc e, the s m al lest spher e c ontaining B ha s r adius at le ast k p − q k / 2 . If x ∈ B , then | x i − ( p i + q i ) / 2 | ≤ ( q i − p i ) / 2 for every i . Conse qu ently, k x − ( p + q ) / 2 k ≤ k p − q k / 2 and the spher e with c enter c = ( p + q ) / 2 and ra dius r = k p − q k / 2 c ontains B . It fol lows that t his is the smal lest spher e c ontaining B sinc e any other spher e must have r adius at le ast k p − q k / 2 . Remark 3. Finding t he smal lest spher e c ontaining C may not b e e asy. However, the c enter and r adius of any spher e c ontaining C yields an affine under est imate for k x k 2 over C . That is, if S is a spher e with C ⊂ S , then t he b est affine u nder estimate for −k x k 2 over S is also an affine u nder estimate for −k x k 2 over C . 3.2. Lo w er b ound based on s emidefini te programming. A differen t DC decomp osition of f ( x ) = ( 1 − x ) T ( A + D ) x is the following: f ( x ) = ( f ( x ) + x T Λx ) − x T Λx , where Λ is a diago nal matrix with i -th diag onal element λ i ≥ 0. W e would like to ma ke the s e cond term x T Λx a s small as p oss ible while keeping the first ter m f ( x ) + x T Λx conv ex. This sugges ts the following semidefinite prog ramming pro blem minimize P n i =1 λ i sub ject to Λ − ( A + D )  0 , Λ  0 , (3.11) where λ is the diagonal of Λ . If the diagonal of A + D is no nnegative, then the inequality Λ  0 can b e dro pped since it is implied b y the ine q uality Λ − ( A + D )  0 . As b efo r e, we se e k the b est linear underestimate of the concave function − x T Λx ov er a compact, conv ex se t C . If any o f the λ i v anish, then reorder the comp onents of x so that x = ( y , z ) where z co rresp onds to the co mpo nents o f λ i that v anish. Let Λ + be the principal s ubmatrix o f Λ corr esp onding to the p o sitive diagonal elements, and define the set C + = { y : ( y , z ) ∈ C for some z } . 8 W. W. HAGER, D. T. PHAN, H. ZHANG The proble m o f finding the b est linear underestima te fo r − x T Λx ov er C is essentially equiv alent to finding the b est linear underestimate for − y T Λ + y over the C + . Hence, there is no loss of g enerality in assuming that the diago nal of Λ is strictly p ositive. As a consequence of Theo rem 3.1, we hav e Corollar y 3.2. Supp ose the diago nal of Λ is st rictly p ositive and let c b e the c enter and r the ra dius of the u nique smal lest spher e c ontaining the set Λ 1 / 2 C := { Λ 1 / 2 x : x ∈ C } . The b est line ar un der estimate of − x T Λx over the c omp act, c onvex set C is ℓ ∗ ( x ) = − 2 c T Λ 1 / 2 x + k c k 2 − r 2 . F u rthermor e, min ℓ ∈S 2 max x ∈C −  x T Λx + ℓ ∗ ( x )  = r 2 , wher e S 2 = { ℓ : R n → R such that ℓ is affine and − x T Λx ≥ ℓ ( x ) for al l x ∈ C } . Pr o of . With the change of v a riables y = Λ 1 / 2 x , an affine function in x is trans- formed to a n affine function in y and co nv e r sely , an affine function in y is transformed to a n a ffine function in x . Hence, the problem of finding the b est affine underestimate for − x T Λx ov er C is equiv alent to the pro blem o f finding the b est a ffine under estimate for −k y k 2 ov er Λ 1 / 2 C . Apply The o rem 3.1 to the transformed pro blem in y , and then transform back to x . Remark 4. If C is the b ox { x ∈ R n : 0 ≤ x ≤ 1 } , then Λ 1 / 2 C is also a b ox to which we c an apply t he observation in R emark 2. In p articular, we have c = 1 2 Λ 1 / 2 1 = 1 2 λ 1 / 2 and r = k Λ 1 / 2 1 k / 2 = k λ 1 / 2 k / 2 . (3.12) Henc e, k c k 2 − r 2 = 0 and we have ℓ ∗ ( x ) = − λ T x . Remark 5. L et us c onsider the set C = { x ∈ R n : 0 ≤ x ≤ 1 , 1 T x = b } , wher e 0 < b < n . Determining the smal lest spher e c ontaining Λ 1 / 2 C may not b e e asy. However, as indic ate d in Rema rk 3, any spher e c ontaining Λ 1 / 2 C yield s an under est imate for x T Λx . Observe t hat Λ 1 / 2 C = { y ∈ R n : 0 ≤ y ≤ λ 1 / 2 , y T λ − 1 / 2 = b } . As observe d in R emark 4 , the c enter c and r adius r of the smal lest spher e S c ontaining the set { y ∈ R n : 0 ≤ y ≤ λ 1 / 2 } ar e given in (3.12). The interse ction of this spher e with t he hyp erplane y T λ − 1 / 2 = b is a lower dimensional spher e S ′ whose c enter c ′ is the pr oje ction of c onto the GRAPH P AR TITIONING 9 hyp erplane. S ′ c ontains C sinc e C is c ontaine d in b oth the original spher e S and the hyp erplane. With a little algebr a, we obtain c ′ = 1 2 λ 1 / 2 +  b − . 5 n P n i =1 λ − 1 i  λ − 1 / 2 . By t he Pythagor e an The or em, the r adius r ′ of t he lower dimensional spher e S ′ is r ′ = v u u t . 25 n X i =1 λ i ! − ( b − . 5 n ) 2 P n i =1 λ − 1 i . Henc e, by Cor ol lary 3.2, an un der estimate of − x T Λx is given by ℓ ( x ) = − λ T x +  n − 2 b P n i =1 λ − 1 i  1 T x + k c ′ k 2 − ( r ′ ) 2 . Sinc e 1 T x = b when x ∈ C , it c an b e shown, after some algebr a, t hat ℓ ( x ) = − λ T x (al l the c onstants in the affine fun ction c anc el). Henc e, the affine under estimate ℓ ∗ c ompute d in R emark 4 for the unit b ox and the affine under est imate ℓ c ompute d in this r emark for the unit b ox interse ct the hyp erplane 1 T x = b ar e t he same. 4. Branc h and b ound algorithm. Since the contin uo us quadratic pro gram (2.1) has a binary solution, the branching pro cess in the branch and b ound algo rithm is based on setting v aria ble s to 0 or 1 and reducing the problem dimensio n (we do not employ bisections of the feasible region as in [1 6]). W e b egin by constr ucting a linear o rdering of the vertices of the graph according to an estimate for the difficulty in placing the vertex in the partition. F o r the numerical exp eriments, the order was based on the total weigh t o f the edges connecting a vertex to the adjacent vertices. If tw o vertices v 1 and v 2 hav e weights w 1 and w 2 resp ectively , then v 1 precedes v 2 if w 1 > w 2 . Let v 1 , v 2 , . . . , v n denote the or dered vertices. Level i in the br anch and b ound tree cor resp onds to s e tting the v i -th comp onent o f x to the v alues 0 or 1 . Each leaf at level i repr esents a sp ecific selection of 0 and 1 v alues for the v 1 through v i -th comp onents of x . Hence, a leaf at level i has a lab el of the for m τ = ( b 1 , b 2 , . . . , b i ) , b j = 0 or 1 for 1 ≤ j ≤ i. (4.1) Corresp o nding to this lea f, the v a lue of the v j -th comp onent of x is b j for 1 ≤ j ≤ i . Let T k denote the br anch and b ound tree at iteration k a nd le t E ( T k ) deno te the leav es in the tree. Supp ose τ ∈ E ( T k ) lies at level i in T k as in (4 .1 ). Let x τ denote the vector gotten by r emoving comp onents v j , 1 ≤ j ≤ i , from x . The v j -th comp onent of x has the pre-assigned bina r y v alue b j for 1 ≤ j ≤ i . After taking into account these assigned binary v alues, the quadratic pr oblem reduces to a low er dimensional problem in the v ariable x τ of the form minimize f τ ( x τ ) sub ject to 0 ≤ x τ ≤ 1 , l τ ≤ 1 T x τ ≤ u τ , where u τ = u − i X j =1 b j and l τ = l − i X j =1 b j . 10 W. W. HAGER, D. T. PHAN, H. ZHANG ( v ) 2 ( v ) 3 ( v ) 1 Level 1 Level 2 Level 3 (0) (1) (0,1) (0,0) (1,1) (1,0) (1,0,1) (1,0,0) Fig. 4.1 . Br anch and b ound tr e e Using the techniques develop ed in Section 3 , we replace f τ by a conv ex low er bo und denoted f L τ and we consider the co nv ex problem minimize f L τ ( x τ ) sub ject to 0 ≤ x τ ≤ 1 , l τ ≤ 1 T x τ ≤ u τ . (4.2) Let M ( τ ) denote the optimal ob jective function v alue for (4.2). A t iteration k , the leaf τ ∈ E ( T k ) for which M ( τ ) is sma llest is used to branch to the next level. If τ has the form (4.1), then the branching pro ces ses g enerates the tw o new leav es ( b 1 , b 2 , . . . , b i , 0) and ( b 1 , b 2 , . . . , b i , 1) . (4.3) An illustration inv olving a 3-level bra nch and bo und tree app ear s in Fig ure 4.1. During the br anch and bound pro cess, w e m ust also co mpute an upp er b ound for the minimal ob jective function v a lue in (2.1). This upper b o und is obtained using a heuristic techn ique based on the gradient pro jection a lg orithm and sphere approximations to the feasible se t. These heur istics for generating an upp er b o und will be describ ed in a s eparate pap er. As p ointed o ut earlie r , many heuristic techniques are av aila ble (for example, Metis ([24, 25, 26]), Chaco [17], a nd Party [32]). An a dv antage of o ur quadr atic pro g ramming based heuris tic is that we s tart at the so lutio n to the low e r bo unding problem, a solution which typically has fractional entries and which is a feas ible star ting p oint for (2.1). Conse q uent ly , the upp er b ound is no larg e r than the ob jective function v alue assoc iated with the optimal point in the lo wer-bound problem. Con v ex quadratic branc h and b ound (CQB) 1. Initialize T 0 = ∅ a nd k = 0. E v aluate b oth a low er bo und for the solution to (2.1) and an upp er deno ted U 0 . 2. Cho ose τ k ∈ E ( T k ) such that M ( τ k ) = min { M ( τ ) : τ ∈ E ( T k ) } . If M ( τ k ) = U k , then stop, an optimal solution of (2.1) has bee n found. 3. Assuming that τ k has the for m (4.1), let T k +1 be the tree obtained by bra nch- ing at τ k and adding tw o new leav es as in (4 .3); also see Figure 4 .1. Ev a luate low e r b ounds for the quadr a tic progra mming pro blems (4.2) ass o ciated with the tw o new leaves, and ev alua te a n impr ov ed upp er b ound, denoted U k +1 , by using solutio ns to the low er b ound problems as starting g uesses in a descent metho d applied to (2.1 ). 4. Replace k by k + 1 and re turn to step 2. Conv ergence is assured since there a re a finite num b er of binary v alues for the comp onents of x . In the w orst case, the branch a nd b ound a lg orithm will build all 2 n +1 − 1 no de s of the tree. GRAPH P AR TITIONING 11 5. Necessary and s ufficient optimali t y conditi ons. W e use the gradient pro jection algorithm to obtain an upper b ound for a solution to (2.1). Since the gradient pro jection algor ithm can terminate a t a stationary point, we need to b e able to dis ting uish b etw e en a stationar y p o int a nd a lo cal minimizer, and at a stationa r y po int which is no t a lo ca l minimizer, we need a fast way to compute a descent direction. W e begin b y s ta ting the first-order optimality conditions. Given a sca lar λ , define the vector µ ( x , λ ) = ( A + D ) 1 − 2( A + D ) x + λ 1 , and the set-v a lued maps N : R → 2 R and M : R → 2 R N ( ν ) =    R if ν = 0 { 1 } if ν < 0 { 0 } if ν > 0 , M ( ν ) =    R if ν = 0 { u } if ν > 0 { l } if ν < 0 . F or any v ector µ , N ( µ ) is a vector of sets whos e i -co mp o nent is the set N ( µ i ). The first-order optimality (Ka rush-Kuhn-T uck e r) conditions asso c iated with a lo c al minimizer x of (2.1) ca n b e wr itten in the following way: F or so me scalar λ , we ha ve 0 ≤ x ≤ 1 , x ∈ N ( µ ( x , λ )) , l ≤ 1 T x ≤ u , and 1 T x ∈ M ( λ ) . (5.1) The first and third conditions in (5.1) are the cons traints in (2.1 ), while the re main- ing t wo conditions corres p o nd to complementary slackness and statio narity of the Lagra ng ian. In [14] w e give a necessar y and sufficient optimality conditions for (2.1), which we now review. Giv en any x that is fea s ible in (2.1), let us define the sets U ( x ) = { i : x i = 1 } , L ( x ) = { i : x i = 0 } , and F ( x ) = { i : 0 < x i < 1 } . W e also introduce subsets U 0 and L 0 defined by U 0 ( x , λ ) = { i ∈ U ( x ) : µ i ( x , λ ) = 0 } and L 0 ( x , λ ) = { i ∈ L ( x ) : µ i ( x , λ ) = 0 } . Theorem 5.1. Supp ose that l = u and D is chosen so that d ii + d j j ≥ 2 a ij . (5.2) for al l i and j . A n e c essary and sufficient c ondition for x to b e a lo c al m inimizer in ( 2.1 ) is that the fol lowing al l hold: (P1) F or some λ , the first-or der c onditions ( 5.1 ) ar e satisfie d at x . (P2) F or e ach i and j ∈ F ( x ) , we have d ii + d j j = 2 a ij . (P3) Consider t he thr e e sets U 0 ( x , λ ) , L 0 ( x , λ ) , and F ( x ) . F or e ach i and j in two differ ent sets , we have d ii + d j j = 2 a ij . In trea ting the situation l < u , a n additional condition concer ning the dual multi- pliers λ and µ in the first-o rder optimality c onditions (5 .1 ) en ters int o the statement of the res ult: (P4) If λ = µ i ( x , 0) = 0 fo r some i , then d ii = 0 in any of the fol lowing thr e e c ases: (a) l < 1 T x < u . (b) x i > 0 and 1 T x = u . 12 W. W. HAGER, D. T. PHAN, H. ZHANG (c) x i < 1 and 1 T x = l . Corollar y 5.2. Supp ose that l < u and D is chosen so that d ii + d j j ≥ 2 a ij and d ii ≥ 0 (5.3) for al l i and j . A n e c essary and sufficient c ondition for x to b e a lo c al m inimizer in ( 2.1 ) is that (P1) – (P4) al l hold. Based on The o rem 5.1 and Cor ollary 5 .2, we can easily chec k whether a given stationary p oint is a local minimizer. This is in contrast to the g eneral quadr a tic progra mming problem for which deciding whether a given p oint is a lo cal minimizer is NP-har d (see [29, 30]). W e now obser ve that when x is a statio nary p oint and when any of the conditions (P2)–(P4) are violated, then a des cent dir ection is readily av aila ble. Proposition 5 .3. Supp ose t hat x is a stationary p oint for ( 2.1 ) and ( 5.3 ) holds. If either (P 2) or (P3) is violate d, then d = e i − e j , with an appr opriate choic e of sign, is a desc ent dir e ction. If l < u , λ = 0 = µ i ( x , 0) , and d ii > 0 , then d = e i , with an appr opriate choic e of sign, is a desc ent dir e ction in any of t he c ases (a)–(c) of (P4) . Pr o of . The Lagra ngian L ass o ciated with (2.1) has the fo rm L ( x ) = f ( x ) + λ ( 1 T x − b ) − X i ∈L µ i x i − X i ∈U µ i ( x i − 1 ) , (5.4) where b = u if λ > 0, b = l if λ < 0, and µ stands fo r µ ( x , λ ). The sets L and U denote L ( x ) and U ( x ) resp ectively . By the first-order optimality conditions (5.1), we hav e L ( x ) = f ( x ) and ∇ L ( x ) = 0 . E xpanding the Lag rangia n around x gives L ( x + y ) = L ( x ) + ∇ L ( x ) y + 1 2 y T ∇ 2 L ( x ) y = f ( x ) − y T ( A + D ) y . W e substitute for L using (5.4) to obtain f ( x + y ) = L ( x + y ) − λ ( 1 T ( x + y ) − b ) + X i ∈L µ i ( x i + y i ) + X i ∈U µ i ( x i + y i − 1 ) = f ( x ) − λ 1 T y − y T ( A + D ) y + X i ∈L µ i y i + X i ∈U µ i y i . (5.5) If (P2) is v iolated, then there a re indices i a nd j ∈ F ( x ) such that d ii + d j j > 2 a ij . W e inser t y = α ( e i − e j ) in (5.5) to o btain f ( x + α ( e i − e j )) = f ( x ) + α 2 (2 a ij − d ii − d j j ) . (5.6) Since the co efficient o f α 2 is negative, d = e i − e j is a descent direc tio n for the ob jective function. Since 0 < x i < 1 and 0 < x j < 1, feas ibilit y is preserved for α sufficiently small. In a similar manner , if (P3) is vio la ted b y indices i and j , then (5.6) again holds and d = ± ( e i − e j ) is aga in a des cent direction where the sign is chosen appropria tely to prese r ve feasibility . F or example, if i ∈ L 0 ( x ) and j ∈ U 0 ( x ), then x i = 0 and x j = 1. Consequen tly , x + α ( e i − e j ) is feasible if α > 0 is sufficient ly small. GRAPH P AR TITIONING 13 Finally , supp ose that l < u , λ = 0 = µ i ( x , 0), and d ii > 0. Substituting y = α e i in (5.5) yields f ( x + α e i ) = f ( x ) − α 2 d ii . Since the co efficient d ii of α 2 is p ositive, d = ± e i is a descent direction. Moreover, in any of the ca ses (a)–(c) o f (P4), x + α d is feasible for s ome α > 0 with an appropria te choice of the sign of d . W e now give a necessa ry a nd sufficient conditio n for a lo c a l minimizer to b e strict. When a lo ca l minimizer is not s trict, it may b e p o ssible to move to a neig hboring po int which has the sa me ob jective function v alue but which is not a lo ca l minimizer . Corollar y 5.4. If x is a lo c al minimizer for ( 2.1 ) and ( 5.3 ) holds, then x is a strict lo c al minimizer if and only if the fol lowing c onditions hold: (C1) F ( x ) is empty. (C2) ∇ f ( x ) i > ∇ f ( x ) j for every i ∈ L ( x ) and j ∈ U ( x ) . (C3) If l < u , the first-or der optimality c onditions ( 5 .1 ) h old fo r λ = 0 , and Z := { i : ∇ f ( x ) i = 0 } 6 = ∅ , t hen either (a) 1 T x = u and x i = 0 for al l i ∈ Z or (b) 1 T x = l and x i = 1 for al l i ∈ Z . Pr o of . Throug hout the pro o f, we let µ , F , L and U denote µ ( x , λ ), F ( x ), L ( x ), and U ( x ) r esp ectively , wher e x is a lo cal minimizer for (2.1) and the pa ir ( x , λ ) satisfies the first-or der optimality conditions (5.1). T o b eg in, supp os e that x is a strict lo ca l minimizer o f (2.1). T ha t is, f ( y ) > f ( x ) w hen y is a feasible p oint near x . If F has at lea st tw o elements, then by (P2) o f Theo rem 5 .1, d ii + d j j = 2 a ij for each i and j ∈ F . Since the fir st-order o ptimality conditions (5.1) ho ld at x , it follows from (5.6) that f ( x + α ( e i − e j )) = f ( x ) (5.7) for all α . Since this vio lates the assumption that x is a strict lo cal minimizer, w e conclude that |F | ≤ 1. If 1 T x = u or 1 T x = l , then since u a nd l a r e int egers , it is not p os sible for x to hav e just o ne fr a ctional co mpo nent. Consequently , F is empty . If l < 1 T x < u , then by complementary slackness, λ = 0. Supp os e that |F | = 1 and i ∈ F . B y (P4) of Coro llary 5 .2, d ii = 0 . Ag ain, by (5.5) it follows that f ( x + α e i ) = f ( x ) for all α . This violates the ass umption that x is a strict lo cal minimizer of (2.1). Hence, F is empty . By the first-or der co nditions (5.1), there exis ts λ such that µ i ( x , λ ) ≥ 0 ≥ µ j ( x , λ ) (5.8) for all i ∈ L and j ∈ U . If this inequality b ecomes an equa lity for so me i ∈ L and j ∈ U , then µ i = 0 = µ j , and b y (P3) o f Co rollar y 5.2, we hav e d ii + d j j = 2 a ij . Again, (5.7) violates the assumption that x is a strict lo cal minimizer. Hence, one of the inequalities in (5.8) is strict. The λ o n each side of (5.8 ) is cancelled to obtain (C2). Suppo se that l < u , λ = 0, and Z := { i : ∇ f ( x ) i = 0 } 6 = ∅ . When λ = 0, we hav e µ ( x , 0) = ∇ f ( x ). Hence, Z = { i : µ i ( x , 0) = 0 } 6 = ∅ . It follows from (P4) that in any of the cases (a)–(c), we have d ii = 0. In particular , if l < 1 T x < u , then by 14 W. W. HAGER, D. T. PHAN, H. ZHANG (5.5), we have f ( x + α e i ) = f ( x ) for all α . Again, this violates the assumption that x is a s trict lo cal minimum. Similarly , if for some i ∈ Z , either x i > 0 and 1 T x = u or x i < 1 and 1 T x = l , the ide ntit y f ( x + α e i ) = f ( x ) implies that we violate the strict lo cal optimality of x . This establishes (C3). Conv ersely , suppo se that x is a lo cal minim izer and (C1)–(C3) hold. W e will show that ∇ f ( x ) y > 0 whenever y 6 = 0 a nd x + y feasible in (2.1) . (5.9) As a result, by the mean v alue theo r em, f ( x + y ) > f ( x ) when y is sufficiently small. Hence, x is a s trict lo cal minimizer. When x + y is feas ible in (2.1), we have y i ≥ 0 for all i ∈ L and y i ≤ 0 for all i ∈ U . (5.10) By the first-or der optimality co ndition (5.1 ), µ i ≥ 0 for a ll i ∈ L and µ i ≤ 0 for a ll i ∈ U . Hence, we hav e ( ∇ f ( x ) + λ 1 T ) y = µ T y = X i ∈L µ i y i + X i ∈U µ i y i ≥ 0 . (5.11) W e now co nsider three cases . First, s uppo se that 1 T y = 0 a nd y 6 = 0 . By (C1) F is empty and hence, by (5.1 0), y i > 0 for some i ∈ L and y j < 0 for some j ∈ U . After adding λ to ea ch side in the inequality in (C2), it follows that either min i ∈L µ i ≥ 0 > max j ∈U µ j (5.12) or min i ∈L µ i > 0 ≥ max j ∈U µ j . (5.13) Combining (5.1 1), (5.12), and (5.13) gives ∇ f ( x ) y ≥ µ i y i − µ j y j > 0 since either µ i > 0 or µ j < 0 , and y i > 0 > y j . Second, supp os e that 1 T y 6 = 0 and λ 6 = 0. T o be s pe cific, suppo se that λ > 0. By complementary slackness, 1 T x = u . Since x + y is feasible in (2.1) and 1 T y 6 = 0, we m ust have 1 T y < 0 . Hence, by (5.1 1), ∇ f ( x ) y > 0. The case λ < 0 is similar. Finally , cons ider the case 1 T y 6 = 0 and λ = 0. In this case, we must hav e l < u . If the set Z in (C3) is empty , then ∇ f ( x ) i = µ i 6 = 0 for a ll i , and by (5.11), ∇ f ( x ) y > 0. If Z 6 = ∅ , then by (C3), either 1 T x = u and x i = 0 for all i ∈ Z o r 1 T x = l and x i = 1 for all i ∈ Z . T o b e spec ific, supp ose that 1 T x = u a nd x i = 0 for all i ∈ Z . Aga in, since x + y is feasible in (2.1) and 1 T y 6 = 0, we hav e 1 T y < 0. If U = ∅ , then x = 0 since F = ∅ . Since 1 T y < 0, w e contradict the feasibilit y of x + y . Hence, U 6 = ∅ . Since 1 T y < 0 , there exists j ∈ U such that y j < 0. Since Z ⊂ L , it follows from (5.12) tha t µ j < 0. B y (5.11) ∇ f ( x ) y ≥ µ j y j > 0. The ca se 1 T x = l and x i = 1 for all i ∈ Z is similar . This completes the pro of of (5.9 ), and the c o rollar y has b een established. 6. Numerical results. W e inv estigate the per formance o f the br a nch and bound algorithm based o n the lower b ounds in Section 3 using a series of test problems. The co des w ere written in C a nd the exp e riments were conducted on an In tel Xeon Quad- Core X5355 2.66 GHz computer using the Lin ux o pe r ating system. Only o ne of the 4 GRAPH P AR TITIONING 15 pro cesso r s was used in the ex pe riments. T o ev a luate the low er b ound, we solve (4.2) by the gradient pro jection method with an ex act linesea rch and Barzilai-Bor wein steplength [1]. The stopping cr iterion in our exp eriments was k P ( x k − g k ) − x k k ≤ 10 − 4 , where P deno tes the pro jection on to the feasible set and g k is the gradient of the ob jective function at x k . The s olution of the semidefinite progr amming problem (3.11) was obtained using V ersion 6.0 .1 of the CSDP co de [2] av ailable at ht tps://pr o jects.coin-or.o rg/Cs dp/ W e compare the p erformance of our a lgorithm with results r ep orted by K a risch, Rendl, and Cla usen in [23] and b y Sensen in [33]. Since these earlier results were obtained on different computer s, we obtained estimates for the cor resp onding running time o n our computer using the LINP ACK benchmarks [7]. Since our co mputer is r o ughly 30 times fas ter than the HP 9000 /735 used in [23] and it is r oughly 7 times fas ter than the Sun UltrSP AR C-I I 4 00Mhz machine used in [3 3], the e a rlier CPU times were divided by 3 0 and 7 resp ectively to obtain the estimated running time o n our co mputer. Note that the same int erior -p oint a lgorithm that we use, which is the ma in ro utine in the CSDP co de, was use d to solve the semidefinite r e laxation in [23]. The test pr o blems w ere ba sed on the gr aph bisection problem where l = u = n/ 2. Two different data sets w ere us ed for the A matrice s in the numerical exp eriments. Most of the test pr oblems came fro m the librar y of Brunetta, Conforti, and Rinaldi [3] which is av ailable at ftp://ftp.math.unipd.it/pub/Mis c/equicut . Some of the test matrices w ere from the UF Sparse Matrix Libra ry ma intained by Timothy Davis: ht tp://www.cise .ufl.edu/ resear ch/sparse/matrices/ Since this se c ond s et of ma tr ices is not directly connected with gr aph partitioning, we create an A for gra ph partitioning as follows: If the matrix S from the library was symmetric, then A was the adjacency matrix defined as fo llows: the diagonal o f A is zero, a ij = 1 if s ij 6 = 0, and a ij = 0 otherwise. If S was no t symmetr ic, then A was the adjacency matrix o f S T S . 6.1. Lo w er b ound comparison. Our num erica l study b egins with a compar i- son of the lo wer b ound o f Section 3.1 based on the minim um eigenv alue of A + D and the b est affine under e s timate, and the lower bo und of Section 3.2 based o n se midef- inite prog r amming. W e la be l these tw o low er bounds LB 1 and LB 2 resp ectively . In T a ble 6.1, the first 5 g raphs cor resp ond to matrices from the UF Sparse Matrix Library , while the next 5 gr aphs were from the test set of B runetta, Conforti, and Rinaldi. The column lab eled “ O pt” is the minimum cut and while n is the pr o b- lem dimension. The numerical results indica te that the low er bound LB 2 based on semidefinite programming is generally b etter (larg er) than LB 1 . In T able 6.1 the best low e r b o und is highlighted in bold. B ased on these res ults, we use the semidefinite progra mming-based lower bo und in the numerical exp er iments which follow. 6.2. Algorithm performance. Unless stated otherwis e , the remaining test problems ca me fro m the library o f Br unetta, Conforti, and Rinaldi [3]. T able 6.2 gives r esults fo r matrices asso ciated with the finite element metho d [35]. The thr ee 16 W. W. HAGER, D. T. PHAN, H. ZHANG T ab le 6.1 Comp arison of two lower bo unds Graph n LB 1 LB 2 Opt Tina Disc a l 11 0.31 0.86 12 jg1009 9 1.55 1.72 1 6 jg1011 11 1.48 0.94 24 Stranke94 10 1.76 1.77 24 Hamrle1 32 -1.93 1.12 17 4x5t 20 -21.7 1 5.43 28 8x5t 40 -16.1 6 2.91 33 t050 30 0.90 18.54 397 2x17m 34 1.33 1.27 316 s090 6 0 -9.84 13. 10 23 8 metho ds are lab eled CQ B (our conv ex quadratic branch and b ound algor ithm), KR C (algorithm of Ka risch, Rendl, and Cla us en [23]), and SEN (algorithm of Sensen [33]). “ n ” is the problem dimens ion, “%” is the p ercent of nonzeros in the matrix, and “# no des” is the num b er of no des in the br a nch and b ound tree. The CPU time is given in s econds. The best time is highlighted in bo ld. As can b e se en in T able 6.2, CQB was fastest in 6 out of the 10 pr oblems even though the num b e r of no des in the branch and bound tre e was muc h larger . Thus b oth KRC and SEN provided muc h tighter relaxa tio ns, how ever, the time to solv e their relaxed pr oblems w as muc h lar ger than the time to optimize our conv ex q uadratics. T able 6.3 gives results for compiler design problems [1 0, 21]. F or this test set, KRC was fastest in 3 out of 5 test problems. Note though that the times for CQ B were comp etitive with KRC. T able 6 .4 gives results for binary de B r uijn gr aphs which arise in a pplications related to par allel computer a rchitecture [5, 9]. These gra phs ar e constructed by the following pro cedure. W e fir st build a dir ected graph using the Mathematica command: A = TableF orm[To AdjacencyMatrix[DeBruijnGraph[2, n]]] T o obtain the graph par titio ning test pr oblem, we a dd the Mathematica generated matrix to its transp ose and set the diagonal to 0 . F or this test set, SEN had b y far the b est p erfor mance. T able 6.5 gives results for to roidal grid graphs. These graphs are connected with an h × k grid, the num b er of vertices in the gr aph is n = hk and there are 2 hk edges whose weights a re chosen from a uniform distribution on the in terv al [1 , 10]. Since Sensen did not s olve e ither this test set, or the remaining test sets, we now compare betw een CQB and KRC. W e se e in T able 6.5 that CQB was faster than KRC in 9 of the 10 toro idal gr id cases. T able 6 .6 gives res ults for mixe d grid graphs. These are complete gr aphs as s o ci- ated with an planar h × k planar g rid; the edg es in the planar gr id received integer weigh ts uniformly dr awn from [1,10 0], while all the o ther edges needed to complete the graph received integer weigh ts uniformly drawn fr om [1,10 ]. F or these graphs, K RC was muc h faster than CQB. Notice that the graphs in this test set are completely dense. One tre nd that is seen in these n umerical exper iments is that as the graph density increases , the p erformance of CQB r elative to the other metho ds deg rades. Results for pla nar gr id graph ar e given in T able 6.7. These graphs are a s so ciated GRAPH P AR TITIONING 17 with an h × k grid. There ar e hk vertices and 2 hk − h − k edges whose weigh ts ar e int egers uniformly drawn fr om [1,10]. F or this relatively sparse test set, CQB was faster in 7 out o f 1 0 problems. T able 6.8 gives r esults for randomly generated gra phs. F or these graphs, the density is first fixed and then the edge s are ass igned in teger weight s uniformly drawn from [1 ,10]. F or this test set, CQB is fastest in 11 of 2 0 cas es. Again, observe that the relative p er formance of CQB degrades as the densit y increa ses, mainly due to the large num b er of no des in the bra nch and bo und tree. 7. Conclusions . An exact algorithm is presented for solving the gr a ph parti- tioning problem with upper and low er b ounds o n the size of each set in the partition. The a lgorithm is based on a contin uo us quadratic programming for m ulation of the discrete partitioning problem. W e show how to transform a feas ible x for the gra ph partitioning QP (2.1) to a binary feasible p oint y with an ob jective function v a lue which sa tisfies f ( y ) ≤ f ( x ). The binar y feasible po int corres p o nds to a partition of the gr aph vertices and f ( y ) is the weigh t of the cut edg es. At any stationar y p o int of (2.1 ) which is not a lo ca l minimizer, P rop osition 5.3 provides a descent direction that can b e used to strictly improv e the ob jective function v a lue. In the bra nch and bound algorithm, the ob jective function is decompos e d into the sum o f a conv e x and a concav e par t. A low er bo und for the ob jectiv e function is a chiev e d by repla cing the concav e pa rt by a n a ffine underestimate. Two different decomp ositions were cons idered, one based on the minimum eigenv alue of the matr ix in the o b jectiv e function, and the other based o n the solution to a semidefinite pro- gramming problem. The semidefinite prog ramming a pproach generally led to mu ch tighter low er b o unds. In a ser ies of numerical exp er iment s, the new algo rithm CQ B (conv ex qua dratic br a nch and b o und) was comp etitive with s tate-of-the-ar t pa rtition- ing metho ds; the rela tive p e rformance of CQB was better for spars e gra phs than for dense graphs. REFERENCES [1] J. Barzilai and J. M. Bor wein , Two p oint step size g r adient metho ds , IM A J. Numer. A nal., 8 (19 88), pp. 141–1 48. [2] B. Borchers , CSDP, A C lib r ary for semidefinite pr o gr amming , Optimization Metho ds and Soft ware, 11(1) (1999), pp. 613–623. [3] L. Brunett a, M. Confor ti, and G. Rinaldi , A br anch-and-cut algorithm for the equicut pr oblem , Math. Program., 78 (1997), pp. 243–263. [4] U. V. Ca t al y ¨ urek and C. A ykana t , Hy p er gr aph-p artitioning b ase d de c omp osition for p ar- al lel sp arse-matrix ve ct or multiplic ation , IEEE T r ansaction on Parallel and Distributed Systems, 10 (1999), pp. 673–693. [5] O. Collins, S. Dolinar, R. McEliece, and F. Pollara , A VLSI de c omp osition of the de Bruijn gr aph , J. ACM, 39 (1992), pp. 931–948. [6] K. Dev ine, E. Boman, R. Heaphy , R. Bisseling, and U. Ca t al yurek , Par al lel hyp e r gr aph p artitioning for scientific c omputing , i n Pro c. of 20th In ternational Parallel and Di stributed Pro cessing Symposium (IPDPS’06), IEEE, 2006. [7] J. J. Dongarra , Performanc e of various co mputers using standar d line ar e quations softwar e , T ech. Report cs-89-85, University of T ennessee, Kno xville, TN, USA, 2008. [8] J. F alkner, F. Rendl, and H. Wol kowic z , A c omputational study of gr aph p artitioning , Math. Programming, 66 (1994), pp. 211–240. [9] R. Feldmann, B. Monien , P. M ysliwietz, and S. T sch ¨ oke , A b etter upp e r b ound on the bise c tion width of de Bruijn networks , i n ST ACS 97, v ol. 1200 of Lecture Notes in Computer Science, Springer, Berli n, 1997, pp. 511–522. [10] C. E. Ferreira, A. Mar tin, C. C. de S ouza, R. Weisman tel, and L. A. W olsey , The no de c ap acitate d gr aph p artitioning pr oblem: A c omputational study , Math. Programming, 81 (1998), pp. 229–256. 18 W. W. HAGER, D. T. PHAN, H. ZHANG [11] M. R. Garey , D. S. J ohnson, a nd L. S tockmeyer , Some simplifie d N P -co mplete gr aph pr oblems , Theoretical Computer Scien ce, 1 (1976), pp. 237–26 7. [12] J. R. Gilber t, G. L. M iller, and S. H. Teng , Ge ometric mesh p artitioning: Implementation and exp eriments , SIAM J. Sci. Comput., 19 (1998 ), pp. 2091–211 0. [13] L. G rigori, E. Boman, S. Donf ack, and T. A. Da vis , Hyp er gra ph-b ase d unsymmetric neste d disse cti on or dering for sp arse LU factorization , SIAM J. Sci. Comput., (2008 ). under submission. [14] W. W. Hager a nd Y. Kr y l yuk , Gr aph p artitioning and c ontinuous quadr atic pr o gr amming , SIAM J. Disc. Math., 12 (1999), pp. 500–523. [15] , M ultiset gr aph p artit i oning , Math. Meth. Oper. Res., 55 (2002), pp. 1–1 0. [16] W. W. Hager and D. T. Phan , An el lipsoidal br anch and b ound algorithm for glob al opti- mization , SIAM J. Optim., 20 (2009), pp. 740–758. [17] B. Hendrickson and R. Leland , The Chac o user’s g uide - V ersion 2.0, Sandia National L ab or atories, T e chnic al R ep ort SA N D94-2692 , 1994. [18] , An impr ove d sp e c tr al g r aph p artitioning algorithm for mapping p ar al lel c omputations , SIAM J. Sci. Comput., 16 (1995), pp. 452–469. [19] , A multilevel algorithm f or p artiti oning gra phs , in Pro c. Sup ercomputing ’95, A CM, Nov . 1995. [20] R. Horst, P. M. P ardalos, an d N. V. Thoai , Intr o duction to Glob al Optimization , Kluw er Academic Publishers, Dordrech t, Holland, 1995. [21] E. L. John son, A. Mehrotra, and G. L. Nemh auser , Min-cut clustering , M ath. Program., 62 (1993) , pp. 133– 151. [22] T. Johnson , A c oncurr ent dynamic task gr aph , in Pr oc. 1993 In tl. Conf. on Parallel Pro cessing, 1993. (TR-93-011, CISE Dept., Univ. of Florida). [23] S. E. Karisch, F. Rendl, a nd J. Clausen , Solving gr aph b ise c t ion pr oblems with semidefinite pr o gr amming , INF ORMS Journal on Computing, 12 (2000), pp. 177–1 91. [24] G. Kar ypis and V. Kumar , A fast and high quality multilevel sc heme for p artitioning irr e gular gr aphs , SIAM J. Sci. Comput., 20 (1998), pp. 359 –392. [25] , Par al lel multilev el k- way p artitioning scheme for irre g ular gra phs , SIAM Review, 41 (1999), pp. 278–300. [26] , M ultilevel k-way hyp e r gr aph p artitioning , VLSI Design, 11 (2000), pp. 285–300. [27] T. Lengauer , Combinatorial Algorithms for Inte gr ated Cir cuit La yout , John Wil ey , Chichest er, 1990. [28] J. Mitchell , Br anch-and-cut for the k -way e quip artit i on pr oblem , tec h. r eport, Department of Mathematica l Sciences, Rensselaer Polytec hnic Institut e, 2001. [29] K. G. Mur ty and S . N. Kabadi , Some NP-c omplete pr oblems in q uadr atic and line ar pr o- gr amming , Math. Program., 39 (1987), pp. 117–129. [30] P. M. P ardalos a nd S . A. V a v asis , Q uadr atic pr o gr amming with one negative eigenvalue is NP-har d , J. Global Optim. , 1 (1991 ), pp. 15– 22. [31] F. Pellegrini, J. Roman, an d P. R. Amestoy , Hy bridizing neste d disse ction and halo ap- pr oximate minimum de gr e e for efficient sp arse matrix or dering , Concurrency: Pract. Exp., 12 (2000) , pp. 68–8 4. [32] R. Preis and R. Diekmann , The P AR TY p arti t ioning libr ary user guide - version 1.1 , 1996. [33] N. Sensen , L ower b ounds and exact algorithms for the gr aph p artitioning pr oblem using multi- c ommo dity fl ows , in Lecture Notes in Computer Science, vo l. 2161, Spri nger- V erlag, 2001, pp. 391–4 03. [34] A. J. S oper, C. W alsha w, and M. Cross , A c ombine d multilevel se ar ch and multilevel opti- mization appr o ach to gr aph-p artition , J. Global O ptim., 29 (2004), pp. 225–241. [35] C. Souza, R. Keunings, L. Wolsey, and O. Zone , A new appr o ach to minimizing the fr on- twidth in finit e element ca lculations , Computer Methods in Applied Mechan ics and Engi- neering, 111 (1994), pp. 323–334. [36] S.-H. Teng , Pr ovably go o d p artitioning and lo ad b alancing algorithms for p ar al lel adaptive N-b o dy simulation , SIAM J. Sci. Comput., 19 (1998), pp. 635–6 56. [37] C. W alsha w, M. Cross, and M. Everett , Par al lel dynamic gr aph p artitioning for adaptive unstructur e d meshs , J. Parallel Distri b. Comput., 47 (1997), pp. 102–10 8. GRAPH P AR TITIONING 19 T ab le 6.2 Mesh Instanc es CQB KRC SEN graph n % # no des time # no des time # no des time m4 32 1 0 22 0.0 5 1 0.03 1 0 .14 ma 54 5 8 0.16 1 0.10 1 0.28 me 60 5 13 0.20 1 0.13 1 0 .28 m6 70 5 205 0.4 7 1 1.23 1 1.43 m b 74 4 95 0.43 1 0.98 1 1.14 mc 74 5 412 0.5 2 1 1.53 1 1.43 md 80 4 101 0 .55 1 0.96 1 1.28 mf 90 4 99 0.79 1 0 .8 0 1 1.85 m1 100 3 200 1.04 15 36.50 1 3.00 m8 148 2 35 16 6.62 1 10.70 1 4.14 T ab le 6.3 Compiler Design CQB KRC SEN graph n % # no des time # no des time # no des time cd30 3 0 13 11 0.05 1 0.0 3 1 0.00 cd45 4 5 10 35 0.27 1 0.23 1 0 .57 cd47a 47 9 4 5 0.34 1 0.33 7 1.00 cd47b 4 7 9 67 0.29 35 3.73 3 1.4 3 cd61 6 1 10 95 0.86 1 0.67 6 6 .00 T ab le 6.4 de Bruijn Networks CQB KRC SEN graph n % # no des time # no des time # no des time debr5 32 1 2 57 0.11 3 0.2 0 1 0.00 debr6 64 6 7327 2.25 55 1 5.63 1 1. 00 debr7 128 3 161 4094 5 1:22:4 5 711 46 :36 1 10. 28 T ab le 6.5 T or oidal Grid: a weig hte d h × k grid with hk vertic e s and 2 hk e dges that re c eiv e d inte ger weig hts uniformly dr awn fr om [1,10] CQB KRC graph n % # no des time # no des time 4x5t 20 21 13 0.01 1 0.03 6x5t 30 14 46 0.05 1 0.10 8x5t 40 10 141 0.1 6 1 0.20 21x2t 4 2 10 18 0.0 2 1 0 .17 23x2t 4 6 9 78 0.15 33 4.16 4x12t 4 8 9 69 0.17 3 0.56 5x10t 5 0 8 129 0.24 1 0.2 0 6x10t 6 0 7 992 0.54 43 11.66 7x10t 7 0 6 844 0.68 47 19.06 10x8t 8 0 5 420 0.91 45 31.46 20 W. W. HAGER, D. T. PHAN, H. ZHANG T ab le 6.6 Mixe d Grid Gr aphs CQB KR C graph n % # no des time # no des time 2x10m 20 100 150 0.03 1 0.03 6x5m 30 1 00 24 7 6 0.20 1 0.03 2x17m 34 100 42 410 2.12 21 0.96 10x4m 40 100 51 713 3.74 2 0.06 5x10m 50 100 3 5887 97 296.19 1 0. 06 T ab le 6.7 Planar Grid CQB KRC graph n % # no des time # no des time 10x2g 20 15 10 0.01 1 0.03 5x6g 3 0 11 44 0.0 5 1 0 .10 2x16g 32 9 23 0.0 6 1 0 .13 18x2g 36 8 19 0.08 1 0. 06 2x19g 38 8 53 0.2 9 49 1.8 3 5x8g 4 0 9 24 0.08 1 0.06 3x14g 42 8 31 0.1 4 5 0 .60 5x10g 50 7 178 0.34 1 0.30 6x10g 60 6 224 0. 35 57 1 0.63 7x10g 70 5 271 0. 63 61 1 8.56 T ab le 6.8 R andomly Gener ate d Gr aphs CQB KRC graph n % # no des time # no des time v090 20 1 0 12 0.01 1 0.03 v000 20 1 00 952 0.02 1 0.03 t090 30 10 10 0.05 1 0.03 t050 30 50 5081 0.32 17 0.73 t000 30 1 00 1226 70 3.79 3 0.20 q090 40 1 0 89 0 .14 1 0. 13 q080 40 2 0 914 0.2 4 3 1 2 .30 q030 40 7 0 5 54652 32.2 3 23 2 . 06 q020 40 8 0 1364 517 72.58 7 0.8 3 q010 40 9 0 4344 123 217 .16 1 3 1.36 q000 40 1 00 81 8698 4 3 80.72 1 0. 13 c090 50 10 397 0. 2 9 1 0.33 c080 50 20 14290 2.20 45 6.13 c070 50 30 1 3629 0 15.70 49 8 .06 c030 50 70 228587 29 2756.26 51 5.46 c290 52 10 340 0. 3 4 1 0.40 c490 54 10 14 4 3 0.54 1 5 3.3 0 c690 56 10 34 0 5 0.82 3 1 .00 c890 58 10 13385 2.66 71 17.53 s090 60 1 0 828 3 2.01 37 9.9 0