An Introduction to Semidefinite Programming: Theory, Algorithms, and Applications

Introduction to Semidefinite Programming (SDP) Robert M. Freund March, 2004 2004 Massachusetts Institute of Technology. 1 1 Introduction Semidefinite programming (SDP ) is the most exciting development in mathematical programming in the 1990’s. SDP has applications in such diverse fields as traditional convex constrained optimization, control theory, and combinatorial optimization. Because SDP is solvable via interior point methods, most of these applications can usually be solved very efficiently in practice as well as in theory. 2 Review of Linear Programming Consider the linear programming problem in standard form: LP : minimize c · x s.t. ai · x = bi , i = 1, . . . , m x ∈ n+ . Here x is a vector of n variables, and we write “c · x” for the inner-product
“ nj=1 cj xj ”, etc. n the nonnegative orthant. Also, n+ := {x ∈ n | x ≥ 0}, and we call + In fact, n+ is a closed convex cone, where K is called a closed a convex cone if K satisfies the following two conditions: • If x, w ∈ K, then αx + βw ∈ K for all nonnegative scalars α and β. • K is a closed set. 2 In words, LP is the following problem: “Minimize the linear function c · x, subject to the condition that x must solve m given equations ai · x = bi , i = 1, . . . , m, and that x must lie in the closed convex cone K = n+ .” We will write the standard linear programming dual problem as: LD : maximize m
i=1 m
s.t. i=1 yi bi yi ai + s = c s ∈ n+ . Given a feasible solution x of LP and a feasible solution (y, s) of LD, the

m duality gap is simply c·x− m i=1 yi bi = (c− i=1 yi ai ) ·x = s·x ≥ 0, because x ≥ 0 and s ≥ 0. We know from LP duality theory that so long as the primal problem LP is feasible and has bounded optimal objective value, then the primal and the dual both attain their optima with no duality gap. That is, there exists x∗ and (y ∗ , s∗ ) feasible for the primal and dual, respectively,
∗ ∗ ∗ such that c · x∗ − m i=1 yi bi = s · x = 0. 3 Facts about Matrices and the Semidefinite Cone If X is an n × n matrix, then X is a positive semidefinite (psd) matrix if v T Xv ≥ 0 for any v ∈ n . 3 If X is an n × n matrix, then X is a positive definite (pd) matrix if v T Xv > 0 for any v ∈ n , v = 0. n denote Let S n denote the set of symmetric n × n matrices, and let S+ the set of positive semidefinite (psd) n × n symmetric matrices. Similarly n denote the set of positive definite (pd) n × n symmetric matrices. let S++ Let X and Y be any symmetric matrices. We write “X  0” to denote that X is symmetric and positive semidefinite, and we write “X  Y ” to denote that X − Y  0. We write “X  0” to denote that X is symmetric and positive definite, etc. 2 n = {X ∈ S n | X  0} is a closed convex cone in n of Remark 1 S+ dimension n × (n + 1)/2. n . Pick any scalars To see why this remark is true, suppose that X, W ∈ S+ n α, β ≥ 0. For any v ∈  , we have: v T (αX + βW )v = αv T Xv + βv T W v ≥ 0, n . This shows that S n is a cone. It is also straightwhereby αX + βW ∈ S+ + n forward to show that S+ is a closed set. Recall the following properties of symmetric matrices: • If X ∈ S n , then X = QDQT for some orthonormal matrix Q and some diagonal matrix D. (Recall that Q is orthonormal means that 4 Q−1 = QT , and that D is diagonal means that the off-diagonal entries of D are all zeros.) • If X = QDQT as above, then the columns of Q form a set of n orthogonal eigenvectors of X, whose eigenvalues are the corresponding diagonal entries of D. • X  0 if and only if X = QDQT where the eigenvalues (i.e., the diagonal entries of D) are all nonnegative. • X  0 if and only if X = QDQT where the eigenvalues (i.e., the diagonal entries of D) are all positive. • If X  0 and if Xii = 0, then Xij = Xji = 0 for all j = 1, . . . , n. • Consider the matrix M defined as follows:  M= P vT v d  , where P  0, v is a vector, and d is a scalar. Then M  0 if and only if d − v T P −1 v > 0. 4 Semidefinite Programming Let X ∈ S n . We can think of X as a matrix, or equivalently, as an array of n2 components of the form (x11 , . . . , xnn ). We can also just think of X as an object (a vector) in the space S n . All three different equivalent ways of looking at X will be useful. What will a linear function of X look like? If C(X) is a linear function of X, then C(X) can be written as C • X, where 5 C • X := n  n  Cij Xij . i=1 j=1 If X is a symmetric matrix, there is no loss of generality in assuming that the matrix C is also symmetric. With this notation, we are now ready to define a semidefinite program. A semidefinite program (SDP ) is an optimization problem of the form: SDP : minimize C • X s.t. Ai • X = bi , i = 1, . . . , m, X  0, Notice that in an SDP that the variable is the matrix X, but it might be helpful to think of X as an array of n2 numbers or simply as a vector in S n . The objective function is the linear function C • X and there are m linear equations that X must satisfy, namely Ai • X = bi , i = 1, . . . , m. The variable X also must lie in the (closed convex) cone of positive semidefn . Note that the data for SDP consists of the inite symmetric matrices S+ symmetric matrix C (which is the data for the objective function) and the m symmetric matrices A1 , . . . , Am , and the m−vector b, which form the m linear equations. Let us see an example of an SDP for n = 3 and m = 2. Define the following matrices: ⎛ ⎞ 1 0 1 A1 = ⎝ 0 3 7 ⎠ , 1 7 5 ⎛ 0 A2 = ⎝ 2 8 2 6 0 6 ⎞ ⎛ ⎞ 8 1 2 3 0 ⎠ , and C = ⎝ 2 9 0 ⎠ , 4 3 0 7 and b1 = 11 and b2 = 19. Then the variable X will be the 3 × 3 symmetric matrix: ⎛ x11 X = ⎝ x21 x31 x12 x22 x32 ⎞ x13 x23 ⎠ , x33 and so, for example, C • X = x11 + 2x12 + 3x13 + 2x21 + 9x22 + 0x23 + 3x31 + 0x32 + 7x33 = x11 + 4x12 + 6x13 + 9x22 + 0x23 + 7x33 . since, in particular, X is symmetric. Therefore the SDP can be written as: SDP : minimize x11 + 4x12 + 6x13 + 9x22 + 0x23 + 7x33 s.t. x11 + 0x12 + 2x13 + 3x22 + 14x23 + 5x33 = 11 0x11 + 4x12 + 16x13 + 6x22 + 0x23 + 4x33 = 19 ⎛ x11 X = ⎝ x21 x31 x12 x22 x32 ⎞ x13 x23 ⎠  0. x33 Notice that SDP looks remarkably similar to a linear program. However, the standard LP constraint that x must lie in the nonnegative orthant is replaced by the constraint that the variable X must lie in the cone of positive semidefinite matrices. Just as “x ≥ 0” states that each of the n components 7 of x must be nonnegative, it may be helpful to think of “X  0” as stating that each of the n eigenvalues of X must be nonnegative. It is easy to see that a linear program LP is a special instance of an SDP . To see one way of doing this, suppose that (c, a1 , . . . , am , b1 , . . . , bm ) comprise the data for LP . Then define: ⎛ ai1 ⎜ 0 ⎜ Ai = ⎜ .. ⎝ . 0 0 ai2 .. . 0 ... ... .. . 0 0 .. . ⎞ ⎛ c1 ⎟ ⎜0 ⎟ ⎜ ⎟ , i = 1, . . . , m, and C = ⎜ .. ⎠ ⎝ . 0 . . . ain 0 c2 . .. ... ... .. . 0 0 .. . 0 . . . cn ⎞ ⎟ ⎟ ⎟. ⎠ Then LP can be written as: SDP : minimize C • X Ai • X = bi , i = 1, . . . , m, Xij = 0, i = 1, . . . , n, j = i + 1, . . . , n, X  0, s.t. with the association that ⎛ x1 ⎜ 0 ⎜ X = ⎜ .. ⎝ . 0 x2 .. . ... ... .. . 0 0 .. . 0 0 ... xn ⎞ ⎟ ⎟ ⎟. ⎠ Of course, in practice one would never want to convert an instance of LP into an instance of SDP . The above construction merely shows that SDP 8 includes linear programming as a special case. 5 Semidefinite Programming Duality The dual problem of SDP is defined (or derived from first principles) to be: m
SDD : maximize i=1 m
s.t. i=1 yi bi yi Ai + S = C S  0. One convenient way of thinking about this problem is as follows. Given mul
tipliers y1 , . . . , ym , the objective is to maximize the linear function m i=1 yi bi . The constraints of SDD state that the matrix S defined as S=C− m  yi Ai i=1 must be positive semidefinite. That is, C− m  yi Ai  0. i=1 We illustrate this construction with the example presented earlier. The dual problem is: 9 SDD : maximize 11y1 + 19y2 ⎛ s.t. ⎞ ⎛ ⎞ ⎛ 1 0 1 0 2 8 1 y1 ⎝ 0 3 7 ⎠ + y 2 ⎝ 2 6 0 ⎠ + S = ⎝ 2 1 7 5 8 0 4 3 2 9 0 ⎞ 3 0⎠ 7 S  0, which we can rewrite in the following form: SDD : maximize s.t. 11y1 + 19y2 ⎛ ⎞ 1 − 1y1 − 0y2 2 − 0y1 − 2y2 3 − 1y1 − 8y2 ⎜ ⎟ ⎝ 2 − 0y1 − 2y2 9 − 3y1 − 6y2 0 − 7y1 − 0y2 ⎠  0. 3 − 1y1 − 8y2 0 − 7y1 − 0y2 7 − 5y1 − 4y2 It is often easier to “see” and to work with a semidefinite program when it is presented in the format of the dual SDD, since the variables are the m multipliers y1 , . . . , ym . As in linear programming, we can switch from one format of SDP (primal or dual) to any other format with great ease, and there is no loss of generality in assuming a particular specific format for the primal or the dual. The following proposition states that weak duality must hold for the primal and dual of SDP : 10 Proposition 5.1 Given a feasible solution X of SDP and a feasible solu
tion (y, S) of SDD, the duality gap is C • X − m i=1 yi bi = S • X ≥ 0. If
m C • X − i=1 yi bi = 0, then X and (y, S) are each optimal solutions to SDP and SDD, respectively, and furthermore, SX = 0. In order to prove Proposition 5.1, it will be convenient to work with the trace of a matrix, defined below: trace(M ) = n  Mjj . j=1 Simple arithmetic can be used to establish the following two elementary identities: • trace(M N ) = trace(N M ) • A • B = trace(AT B) Proof of Proposition 5.1. For the first part of the proposition, we must show that if S  0 and X  0, then S • X ≥ 0. Let S = P DP T and X = QEQT where P, Q are orthonormal matrices and D, E are nonnegative diagonal matrices. We have: S • X = trace(S T X) = trace(SX) = trace(P DP T QEQT ) = trace(DP T QEQT P ) = n  Djj (P T QEQT P )jj ≥ 0, j=1 where the last inequality follows from the fact that all Djj ≥ 0 and the fact 11 that the diagonal of the symmetric positive semidefinite matrix P T QEQT P must be nonnegative. To prove the second part of the proposition, suppose that trace(SX) = 0. Then from the above equalities, we have n  Djj (P T QEQT P )jj = 0. j=1 However, this implies that for each j = 1, . . . , n, either Djj = 0 or the (P T QEQT P )jj = 0. Furthermore, the latter case implies that the j th row of P T QEQT P is all zeros. Therefore DP T QEQT P = 0, and so SX = P DP T QEQT = 0. q.e.d. Unlike the case of linear programming, we cannot assert that either SDP or SDD will attain their respective optima, and/or that there will be no duality gap, unless certain regularity conditions hold. One such regularity condition which ensures that strong duality will prevail is a version of the “Slater condition”, summarized in the following theorem which we will not prove: ∗ denote the optimal objective function values Theorem 5.1 Let zP∗ and zD of SDP and SDD, respectively. Suppose that there exists a feasible solution ˆ of SDP such that X ˆ  0, and that there exists a feasible solution (ˆ X y, Ŝ) of SDD such that Ŝ  0. Then both SDP and SDD attain their optimal ∗. values, and zP∗ = zD 12 6 Key Properties of Linear Programming that do not extend to SDP The following summarizes some of the more important properties of linear programming that do not extend to SDP : • There may be a finite or infinite duality gap. The primal and/or dual may or may not attain their optima. As noted above in Theorem 5.1, both programs will attain their common optimal value if both programs have feasible solutions in the interior of the semidefinite cone. • There is no finite algorithm for solving SDP . There is a simplex algorithm, but it is not a finite algorithm. There is no direct analog of a “basic feasible solution” for SDP . • Given rational data, the feasible region may have no rational solutions. The optimal solution may not have rational components or rational eigenvalues. • Given rational data whose binary encoding is size L, the norms of any L feasible and/or optimal solutions may exceed 22 (or worse). • Given rational data whose binary encoding is size L, the norms of any L feasible and/or optimal solutions may be less than 2−2 (or worse). 7 SDP in Combinatorial Optimization SDP has wide applicability in combinatorial optimization. A number of N P −hard combinatorial optimization problems have convex relaxations that are semidefinite programs. In many instances, the SDP relaxation is very tight in practice, and in certain instances in particular, the optimal solution to the SDP relaxation can be converted to a feasible solution for the original problem with provably good objective value. An example of the use of SDP in combinatorial optimization is given below. 13 7.1 An SDP Relaxation of the MAX CUT Problem Let G be an undirected graph with nodes N = {1, . . . , n}, and edge set E. Let wij = wji be the weight on edge (i, j), for (i, j) ∈ E. We assume that wij ≥ 0 for all (i, j) ∈ E. The MAX CUT problem is to determine a subset S of the nodes N for which the sum of the weights of the edges that cross from S to its complement S¯ is maximized (where (S¯ := N \ S) . We can formulate MAX CUT as an integer program as follows. Let xj = 1 for j ∈ S and xj = −1 for j ∈ S̄. Then our formulation is: M AXCU T : maximizex s.t. 1 4 n
n
i=1 j=1 wij (1 − xi xj ) xj ∈ {−1, 1}, j = 1, . . . , n. Now let Y = xxT , whereby Yij = xi xj i = 1, . . . , n, j = 1, . . . , n. Also let W be the matrix whose (i, j)th element is wij for i = 1, . . . , n and j = 1, . . . , n. Then MAX CUT can be equivalently formulated as: 14 M AXCU T : maximizeY,x 1 4 n
n
i=1 j=1 wij − W • Y xj ∈ {−1, 1}, j = 1, . . . , n s.t. Y = xxT . Notice in this problem that the first set of constraints are equivalent to Yjj = 1, j = 1, . . . , n. We therefore obtain: 1 4 M AXCU T : maximizeY,x n
n
i=1 j=1 wij − W • Y Yjj = 1, j = 1, . . . , n s.t. Y = xxT . Last of all, notice that the matrix Y = xxT is a symmetric rank-1 positive semidefinite matrix. If we relax this condition by removing the rank1 restriction, we obtain the following relaxtion of MAX CUT, which is a semidefinite program: RELAX : maximizeY 1 4 n
n
i=1 j=1 wij − W • Y Yjj = 1, j = 1, . . . , n s.t. Y  0. It is therefore easy to see that RELAX provides an upper bound on MAXCUT, i.e., 15 M AXCU T ≤ RELAX. As it turns out, one can also prove without too much effort that: 0.87856 RELAX ≤ M AXCU T ≤ RELAX. This is an impressive result, in that it states that the value of the semidefinite relaxation is guaranteed to be no more than 12% higher than the value of N P -hard problem MAX CUT. 8 SDP in Convex Optimization As stated above, SDP has very wide applications in convex optimization. The types of constraints that can be modeled in the SDP framework include: linear inequalities, convex quadratic inequalities, lower bounds on matrix norms, lower bounds on determinants of symmetric positive semidefinite matrices, lower bounds on the geometric mean of a nonnegative vector, plus many others. Using these and other constructions, the following problems (among many others) can be cast in the form of a semidefinite program: linear programming, optimizing a convex quadratic form subject to convex quadratic inequality constraints, minimizing the volume of an ellipsoid that covers a given set of points and ellipsoids, maximizing the volume of an ellipsoid that is contained in a given polytope, plus a variety of maximum eigenvalue and minimum eigenvalue problems. In the subsections below we demonstrate how some important problems in convex optimization can be re-formulated as instances of SDP . 16 8.1 SDP for Convex Quadratically Constrained Quadratic Programming, Part I A convex quadratically constrained quadratic program is a problem of the form: QCQP : minimize xT Q0 x + q0T x + c0 x s.t. xT Qi x + qiT x + ci ≤ 0 , i = 1, . . . , m, where the Q0  0 and Qi  0, i = 1, . . . , m. This problem is the same as: QCQP : minimize θ x, θ s.t. xT Q0 x + q0T x + c0 − θ ≤ 0 xT Qi x + qiT x + ci ≤ 0 , i = 1, . . . , m. We can factor each Qi into Qi = MiT Mi for some matrix Mi . Then note the equivalence: I xT MiT Mi x −ci − qiT x 0 In this way we can write QCQP as: 17 ⇔ xT Qi x + qiT x + ci ≤ 0. QCQP : minimize θ x, θ s.t. I xT M0T M0 x −c0 − q0T x + θ I xT MiT Mi x −ci − qiT x 0  0 , i = 1, . . . , m. Notice in the above formulation that the variables are θ and x and that all matrix coefficients are linear functions of θ and x. 8.2 SDP for Convex Quadratically Constrained Quadratic Programming, Part II As it turns out, there is an alternative way to formulate QCQP as a semidefinite program. We begin with the following elementary proposition. Proposition 8.1 Given a vector x ∈ k and a matrix W ∈ k×k , then 1 xT x W 0 if and only if W  xxT . Using this proposition, it is straightforward to show that QCQP is equivalent to the following semi-definite program: 18 QCQP : minimize θ x, W, θ s.t. c0 − θ 1 2 q0 ci 1 2 qi 1 T 2 q0 Q0 1 T 2 qi Qi 1 xT x W • • 1 xT x W 1 xT x W ≤0 ≤ 0 , i = 1, . . . , m 0. variables, but that Notice in this formulation that there are now (n+1)(n+2) 2 the constraints are all linear inequalities as opposed to semi-definite inequalities. 8.3 SDP for the Smallest Circumscrib ed Ellipsoid Problem A given matrix R  0 and a given point z can be used to define an ellipsoid in n : ER,z := {y | (y − z)T R(y − z) ≤ 1}. One can prove that the volume of ER,z is proportional to det(R−1 ). Suppose we are given a convex set X ∈ n described as the convex hull 19 of k points c1 , . . . , ck . We would like to find an ellipsoid circumscribing these k points that has minimum volume. Our problem can be written in the following form: M CP : minimize vol (ER,z ) R, z s.t. ci ∈ ER,z , i = 1, . . . , k, which is equivalent to: M CP : minimize − ln(det(R)) R, z s.t. (ci − z)T R(ci − z) ≤ 1, i = 1, . . . , k R  0, Now factor R = M 2 where M  0 (that is, M is a square root of R), and now M CP becomes: M CP : minimize − ln(det(M 2 )) M, z s.t. (ci − z)T M T M (ci − z) ≤ 1, i = 1, . . . , k, M  0. Next notice the equivalence: I (M ci − M z)T M ci − M z 1 0 20 ⇔ (ci − z)T M T M (ci − z) ≤ 1 In this way we can write M CP as: M CP : minimize −2 ln(det(M )) M, z I s.t. (M ci − M z)T M ci − M z 1  0, i = 1, . . . , k, M  0. Last of all, make the substitution y = M z to obtain: M CP : minimize −2 ln(det(M )) M, y I M ci − y s.t. T (M ci − y) 1  0, i = 1, . . . , k, M  0. Notice that this last program involves semidefinite constraints where all of the matrix coefficients are linear functions of the variables M and y. However, the objective function is not a linear function. It is possible to convert this problem further into a genuine instance of SDP , because there is a way to convert constraints of the form − ln(det(X)) ≤ θ to a semidefinite system. Nevertheless, this is not necessary, either from a theoretical or a practical viewpoint, because it turns out that the function f (X) = − ln(det(X)) is extremely well-behaved and is very easy to optimize (both in theory and in practice). 21 Finally, note that after solving the formulation of M CP above, we can recover the matrix R and the center z of the optimal ellipsoid by computing R = M 2 and z = M −1 y. 8.4 SDP for the Largest Inscribed Ellipsoid Problem Recall that a given matrix R  0 and a given point z can be used to define an ellipsoid in n : ER,z := {x | (x − z)T R(x − z) ≤ 1}, and that the volume of ER,z is proportional to det(R−1 ). Suppose we are given a convex set X ∈ n described as the intersection of k halfspaces {x | (ai )T x ≤ bi }, i = 1, . . . , k, that is, X = {x | Ax ≤ b} where the ith row of the matrix A consists of the entries of the vector ai , i = 1, . . . , k. We would like to find an ellipsoid inscribed in X of maximum volume. Our problem can be written in the following form: M IP : maximize vol (ER,z ) R, z s.t. ER,z ⊂ X, 22 which is equivalent to: M IP : maximize det(R−1 ) R, z s.t. ER,z ⊂ {x | (ai )T x ≤ bi }, i = 1, . . . , k R  0, which is equivalent to: M IP : maximize ln(det(R−1 )) R, z s.t. maxx {aTi x | (x − z)T R(x − z) ≤ 1} ≤ bi , i = 1, . . . , k R  0. For a given i = 1, . . . , k, the solution to the optimization problem in the ith constraint is R−1 ai x∗ = z + aTi R−1 ai with optimal objective function value aTi z + aTi R−1 ai , and so M IP can be rewritten as: 23 M IP : maximize ln(det(R−1 )) R, z aTi z + s.t. aTi R−1 ai ≤ bi , i = 1, . . . , k R  0. Now factor R−1 = M 2 where M  0 (that is, M is a square root of R−1 ), and now M IP becomes: M IP : maximize ln(det(M 2 )) M, z s.t. aTi z + aTi M T M ai ≤ bi , i = 1, . . . , k M  0, which we can re-write as: M IP : maximize 2 ln(det(M )) M, z s.t. aTi M T M ai ≤ (bi − aTi z)2 , i = 1, . . . , k bi − aTi z ≥ 0, i = 1, . . . , k M  0. Next notice the equivalence: (bi − aTi z)I (M ai )T M ai (bi − aTi z) 0 ⇔ 24 ⎧ ⎫ T T T 2 ⎪ ⎨ ai M M ai ≤ (bi − ai z) ⎪ ⎬ ⎪ ⎩ and bi − aTi z ≥ 0. ⎪ ⎭ In this way we can write M IP as: M IP : minimize −2 ln(det(M )) M, z (bi − aTi z)I M ai s.t. T (bi − aTi z) (M ai )  0, i = 1, . . . , k, M  0. Notice that this last program involves semidefinite constraints where all of the matrix coefficients are linear functions of the variables M and z. However, the objective function is not a linear function. It is possible, but not necessary in practice, to convert this problem further into a genuine instance of SDP , because there is a way to convert constraints of the form − ln(det(X)) ≤ θ to a semidefinite system. Such a conversion is not necessary, either from a theoretical or a practical viewpoint, because it turns out that the function f (X) = − ln(det(X)) is extremely well-behaved and is very easy to optimize (both in theory and in practice). Finally, note that after solving the formulation of M IP above, we can recover the matrix R of the optimal ellipsoid by computing R = M −2 . 25 8.5 SDP for Eigenvalue Optimization There are many types of eigenvalue optimization problems that can be formualated as SDP s. A typical eigenvalue optimization problem is to create a matrix S := B − k  wi Ai i=1 given symmetric data matrices B and Ai , i = 1, . . . , k, using weights w1 , . . . , wk , in such a way to minimize the difference between the largest and smallest eigenvalue of S. This problem can be written down as: EOP : minimize λmax (S) − λmin (S) w, S S=B− s.t. k
i=1 wi Ai , where λmin (S) and λmax (S) denote the smallest and the largest eigenvalue of S, respectively. We now show how to convert this problem into an SDP . Recall that S can be factored into S = QDQT where Q is an orthonormal matrix (i.e., Q−1 = QT ) and D is a diagonal matrix consisting of the eigenvalues of S. The conditions: λI S µI can be rewritten as: Q(λI)QT QDQT 26 Q(µI)QT . After premultiplying the above QT and postmultiplying Q, these conditions become: λI D µI which are equivalent to: λ ≤ λmin (S) and λmax (S) ≤ µ. Therefore EOP can be written as: EOP : minimize µ − λ w, S, µ, λ S=B− s.t. λI S k
i=1 wi Ai µI. This last problem is a semidefinite program. Using constructs such as those shown above, very many other types of eigenvalue optimization problems can be formulated as SDP s. 27 9 SDP in Control Theory A variety of control and system problems can be cast and solved as instances of SDP . However, this topic is beyond the scope of these notes. 10 Interior-point Methods for SDP At the heart of an interior-point method is a barrier function that exerts a repelling force from the boundary of the feasible region. For SDP , we need a barrier function whose values approach +∞ as points X approach n. the boundary of the semidefinite cone S+ n . Then X will have n eigenvalues, say λ (X), . . . , λ (X) Let X ∈ S+ 1 n (possibly counting multiplicities). We can characterize the interior of the semidefinite cone as follows: n intS+ = {X ∈ S n | λ1 (X) > 0, . . . , λn (X) > 0}. n A natural barrier function to use to repel X from the boundary of S+ then is − n  j=1 ln(λi (X)) = − ln( n  λi (X)) = − ln(det(X)). j=1 Consider the logarithmic barrier problem BSDP (θ) parameterized by the positive barrier parameter θ: 28 BSDP (θ) : minimize C • X − θ ln(det(X)) Ai • X = bi , i = 1, . . . , m, X  0. s.t. Let fθ (X) denote the objective function of BSDP (θ). Then it is not too difficult to derive: ∇fθ (X) = C − θX −1 , (1) and so the Karush-Kuhn-Tucker conditions for BSDP (θ) are: ⎧ ⎪ Ai • X = bi , i = 1, . . . , m, ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ X  0, ⎪ ⎪ ⎪ ⎪ m ⎪
⎪ ⎪ yi Ai . ⎩ C − θX −1 = (2) i=1 Because X is symmetric, we can factorize X into X = LLT . We then can define S = θX −1 = θL−T L−1 , which implies 1 T L SL = I, θ 29 and we can rewrite the Karush-Kuhn-Tucker conditions as: ⎧ Ai • X = bi , i = 1, . . . , m, ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ T ⎪ ⎪ X  0, X = LL ⎪ ⎨ (3) m
⎪ ⎪ ⎪ yi Ai + S = C ⎪ ⎪ ⎪ i=1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ I − 1θ LT SL = 0. From the equations of (3) it follows that if (X, y, S) is a solution of (3), then X is feasible for SDP , (y, S) is feasible for SDD, and the resulting duality gap is S•X = n  n  Sij Xij = i=1 j=1 n  (SX)jj = j=1 n  θ = nθ. j=1 This suggests that we try solving BSDP (θ) for a variety of values of θ as θ → 0. However, we cannot usually solve (3) exactly, because the fourth equation group is not linear in the variables. We will instead define a “β-approximate solution” of the Karush-Kuhn-Tucker conditions (3). Before doing so, we introduce the following norm on matrices, called the Frobenius norm: M  := √   n  n  Mij2 . M •M = i=1 j=1 For some important properties of the Frobenius norm, see the last subsection 30 of this section. A “β-approximate solution” of BSDP (θ) is defined as any solution (X, y, S) of ⎧ ⎪ Ai • X = bi , i = 1, . . . , m, ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ X  0, X = LLT ⎪ ⎪ ⎨ m
⎪ ⎪ ⎪ yi Ai + S = C ⎪ ⎪ ⎪ i=1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ I − 1θ LT SL ≤ β. (4) ¯ y, ¯ is a β-approximate solution of BSDP (θ) and Lemma 10.1 If (X, ¯ S) ¯ is feasible for SDP , (¯ β < 1, then X y, S̄) is feasible for SDD, and the duality gap satisfies: nθ(1 − β) ≤ C • X − m  ¯ • S¯ ≤ nθ(1 + β). yi bi = X (5) i=1 Proof: Primal feasibility is obvious. To prove dual feasibility, we need to show that S¯  0. To see this, define 1 ¯T ¯ ¯ SL R=I− L θ and note that R ≤ β < 1. Rearranging (6), we obtain (6) ¯ −1  0 ¯ −T (I − R)L S¯ = θL ¯ • S¯ = trace(XS̄) ¯ = because R < 1 implies that I − R  0. We also have X T T ¯ ¯ S̄) = trace(L ¯ SL) ¯ ¯ = θtrace(I − R) = θ(n − trace(R)). However, trace(LL √ |trace(R)| ≤ nR ≤ nβ, whereby we obtain ¯ • S¯ ≤ nθ(1 + β). nθ(1 − β) ≤ X 31 q.e.d. 10.1 The Algorithm Based on the analysis just presented, we are motivated to develop the following algorithm: Step 0 . Initialization. Data is (X 0 , y 0 , S 0 , θ0 ). k = 0. Assume that (X 0 , y 0 , S 0 ) is a β-approximate solution of BSDP (θ0 ) for some known value of β that satisfies β < 1. ¯ y, ¯ = (X k , y k , S k ), θ = θk . Step 1. Set Current values. (X, ¯ S)
Step 2. Shrink θ. Set θ = αθ for some α ∈ (0, 1). In fact, it will be appropriate to set √ β−β √ α=1− √ β+ n Step 3. Compute Newton Direction and Multipliers. Compute

¯ by factoring X ¯ =L ¯L ¯ T and the Newton step D for BSDP (θ ) at X = X solving the following system of equations in the variables (D, y): ⎧ m

⎪ ¯ −1 ¯ −1 ¯ −1 =
yi Ai ⎪ ⎨ C − θ X + θ X DX ⎪ ⎪ ⎩ i=1 Ai • D = 0, i = 1, . . . , m.

Denote the solution to this system by (D , y ). Step 4. Update All Values.
¯ + D
X =X 32 (7)
S =C− m 
yi Ai i=1


Step 5. Reset Counter and Continue. (X k+1 , y k+1 , S k+1 ) = (X , y , S ).
θk+1 = θ . k ← k + 1. Go to Step 1. Figure 1 shows a picture of the algorithm: θ=0 θ = 1/10 θ = 10 θ = 100 Figure 1: A conceptual picture of the interior-point algorithm. Some of the unresolved issues regarding this algorithm include: • how to set the fractional decrease parameter α
• the derivation of the Newton step D and the multipliers y 33
• whether or not successive iterative values (X k , y k , S k ) are β-approximate solutions to BSDP (θk ), and • how to get the method started in the first place. 10.2 The Newton Step ¯ is a feasible solution to BSDP (θ): Suppose that X BSDP (θ) : minimize C • X − θ ln(det(X)) s.t. Ai • X = bi , i = 1, . . . , m, X  0. Let us denote the objective function of BSDP (θ) by fθ (X), i.e., fθ (X) = C • X − θ ln(det(X)). Then we can derive: ¯ = C − θX ¯ −1 ∇fθ (X) ¯ can be derived and the quadratic approximation of BSDP (θ) at X = X as: ¯ + (C − X ¯ −1 ) • (X − X) ¯ + 1 θX ¯ −1 (X − X) ¯ •X ¯ −1 (X − X) ¯ minimize fθ (X) 2 X s.t. Ai • X = bi , i = 1, . . . , m. 34 ¯ this is equivalent to: Letting D = X − X, ¯ −1 ) • D + 1 θX ¯ −1 D • X ¯ −1 D minimize (C − θX 2 D s.t. Ai • D = 0, i = 1, . . . , m. The solution to this program will be the Newton direction. The KarushKuhn-Tucker conditions for this program are necessary and sufficient, and are: ⎧ m ⎪ ¯ −1 + θX ¯ −1 DX ¯ −1 =
yi Ai ⎪ C − θX ⎨ ⎪ ⎪ ⎩ i=1 Ai • D = 0, (8) i = 1, . . . , m.

These equations are called the Normal Equations. Let D and y denote the solution to the Normal Equations. Note in particular from the first equation


in (8) that D must be symmetric. Suppose that (D , y ) is the (unique) solution of the Normal Equations (8). We obtain the new value of the primal variable X by taking the Newton step, i.e.,
¯ + D
. X =X We can produce new values of the dual variables (y, S) by setting the new

value of y to be y and by setting S = C − have that m
i=1
y Ai . Using (8), then, we
¯ −1 − θX ¯ −1 D
X ¯ −1 . S = θX 35 (9) We have the following very powerful convergence theorem which demonstrates the quadratic convergence of Newton’s method for this problem, with an explicit guarantee of the range in which quadratic convergence takes place. Theorem 10.1 (Explicit Quadratic Convergence of Newton’s Method). ¯ y, ¯ is a β-approximate solution of BSDP (θ) and β < 1. Suppose that (X, ¯ S)

Let (D , y ) be the solution to the Normal Equations (8), and let
¯ + D
X =X and
¯ −1 − θX ¯ −1 D
X ¯ −1 . S = θX


Then (X , y , S ) is a β 2 -approximate solution of BSDP (θ). ¯ satisfies: Proof: Our current point X ¯T  0 ¯ = bi , i = 1, . . . , m, X ¯ =L ¯L Ai • X m  ȳi Ai + S¯ = C i=1 1 ¯T ¯ ¯ I − L SL ≤ β < 1. θ

Furthermore the Newton direction D and multipliers y satisfy:
Ai • D = 0, i = 1, . . . , m m 

yi Ai + S = C i=1
¯ −T )L ¯ + D
= L(I ¯ +L ¯ −1 D
L ¯T X =X
¯ −1 − θX ¯ −1 D
X ¯ −T (I − L ¯ −1 D
L ¯ −1 . ¯ −1 = θL ¯ −T )L S = θX ¯ −1 D
L ¯ −T  ≤ β. It turns out that this is the We will first show that L crucial fact from which everything will follow nicely. To prove this, note that m  i=1 ȳi Ai + S¯ = C = m  i=1

yi Ai + S = m  i=1 36
¯ −T (I − L ¯ −1 D
L ¯ −1 . ¯ −T )L yi Ai + θL
Taking the inner product with D yields:
¯ −1 • D
− θL ¯ −1 D
L ¯ −T L ¯ −1 • D
, ¯ −T L ¯ −T L S¯ • D = θL which we can rewrite as: ¯ −1 D
L ¯ −1 D
L ¯ −1 D
L ¯¯ • L ¯ −1 D
L ¯ −T = θI • L ¯ −T − θL ¯ −T • L ¯ −T , ¯ T SL L which we finally rewrite as:   ¯ −T 2 = I − 1 L ¯¯ • L ¯ −T . ¯ T SL ¯ −1 D
L ¯ −1 D
L L θ Invoking the Cauchy-Schwartz inequality we obtain: ¯ T SLL ¯ −1 D
L ¯ −T 2 ≤ I − 1 L ¯ ¯ ¯ −1 D
L ¯ −T  ≤ βL ¯ −T , ¯ −1 D
L L θ ¯ −1 D
L ¯ −T  ≤ β. from which we see that L It therefore follows that
¯ −T )L ¯ +L ¯ −1 D
L ¯T  0 X = L(I and
¯ −T )L ¯ −T (I − L ¯ −1 D
L ¯ −1  0, S = θL ¯ −T  ≤ β < 1, which guarantees that I ± L ¯ −T  0. ¯ −1 D
L ¯ −1 D
L since L Next, factorize ¯ −1 D
L ¯ −T = M 2 , I +L (where M = M T ) and note that
¯ ML ¯ T = L
(L
)T X = LM
¯ . Then note that where we define L = LM 1
1 ¯T
¯
I − (L )T S  L = I − M L S LM θ θ 1 ¯ T ¯ −T ¯ −1 D
L ¯ −1 )LM ¯ = I − M (I − L ¯ −1 D
L ¯ −T )L ¯ −T )M (θL (I − L = I − ML θ 37 ¯ −T )M = I−M M +M (M M −I)M = (I−M M )(I−M M ) ¯ −1 D
L = I−M M +M (L ¯ −T )(L ¯ −T ). ¯ −1 D
L ¯ −1 D
L = (L From this we next have: 1

¯ −1 D
L ¯ −1 D
L ¯ −1 D
L ¯ −T )(L ¯ −T ) ≤ (L ¯ −T )2 ≤ β 2 . I − (L )T S  L  = (L θ


This shows that (X , y , S ) is a β 2 -approximate solution of BSDP (θ). q.e.d. 10.3 Complexity Analysis of the Algorithm ¯ y, ¯ is a βTheorem 10.2 (Relaxation Theorem). Suppose that (X, ¯ S) approximate solution of BSDP (θ) and β < 1. Let √ β−β √ α=1− √ β+ n
¯ y, ¯ is a √β-approximate solution of BSDP (θ
). ¯ S) and let θ = αθ. Then (X, ¯ = bi , i = 1, . . . , m, X ¯  0, and ¯ y, ¯ satisfies Ai • X Proof: The triplet (X, ¯ S) m
ȳi Ai + S¯ = C, and so it remains to show that i=1     1 T ¯ SL ¯ ¯ − I  ≤ β, 
L θ  ¯ = LL ¯ ¯ T . We have where X           1 T      ¯ T SL ¯ T SL ¯ SL ¯ ¯ − I =  1 L ¯ ¯ − I =  1 1 L ¯ ¯ − I − 1 − 1 I 
L θ   αθ  α θ α        1 − α 1  1 ¯T ¯ ¯     I L SL − I  +  ≤ α θ α  √   ≤ β + α 1−α √ β+ n √ − n n = α α   √ √ = β + n − n = β. 38 q.e.d. Theorem 10.3 (Convergence Theorem). Suppose that (X 0 , y 0 , S 0 ) is a β-approximate solution of BSDP (θ0 ) and β < 1. Then for all k = 1, 2, 3, ..., (X k , y k , S k ) is a β-approximate solution of BSDP (θk ). Proof: By induction, suppose that the theorem is true for iterates 0, 1, 2, ..., k. Then (X k , y k , S k ) is a β-approximate solution of BSDP (θk ). From the Relaxation Theorem, (X k , y k , S k ) is a of BSDP (θk+1 ) where θk+1 = αθk . √ β-approximate solution From the Quadratic Convergence Theorem, (X k+1 , y k+1 , S k+1 ) is a β-approximate solution of BSDP (θk+1 ). Therefore, by induction, the theorem is true for all values of k. q.e.d. Figure 2 shows a better picture of the algorithm: Theorem 10.4 (Complexity Theorem). Suppose that (X 0 , y 0 , S 0 ) is a β = 14 -approximate solution of BSDP (θ0 ). In order to obtain primal and dual feasible solutions (X k , y k , S k ) with a duality gap of at most , one needs to run the algorithm for at most  1.25 X 0 • S 0 k = 6 n ln 0.75  √ iterations. 39  θ = 80 x^ ~ x θ = 90 _ x θ = 100 Figure 2: Another picture of the interior-point algorithm. 40 Proof: Let k be as defined above. Note that √ 1 −1 1 1 β−β √ ≤1− √ . √ = 1 −  2 √4  = 1 − α=1− √ 1 2+4 n 6 n β+ n + n 2 Therefore  1 θ ≤ 1− √ 6 n k k θ0 . This implies that C •X − k m   bi yik i=1 1 = X • S ≤ θ n(1 + β) ≤ 1 − √ 6 n k k  1 ≤ 1− √ 6 n k k (1.25n) X 0 • S0 0.75n k (1.25nθ0 ) , from (5). Taking logarithms, we obtain ln C • X − k m  i=1  bi yik 1 ≤ k ln 1 − √ 6 n   1.25 0 X • S0 + ln 0.75 −k 1.25 0 X • S0 ≤ √ + ln 6 n 0.75 1.25 X 0 • S 0 ≤ − ln 0.75  Therefore C • X k − q.e.d. 10.4 m
i=1      1.25 0 X • S 0 = ln(). + ln 0.75 bi yik ≤ . How to Start the Method from a Strictly Feasible Point The algorithm and its performance relies on having a starting point (X 0 , y 0 , S 0 ) that is a β-approximate solution of the problem BSDP (θ0 ). In this subsection, we show how to obtain such a starting point, given a positive definite feasible solution X 0 of SDP . 41 We suppose that we are given a target value θ0 of the barrier parameter, and we are given X = X 0 that is feasible for BSDP (θ0 ), that is, Ai • X 0 = bi , i = 1, . . . , m, and X 0  0. We will attempt to approximately solve BSDP (θ0 ) starting at X = X 0 , using the Newton direction at each iteration. The formal statement of the algorithm is as follows: Step 0 . Initialization. Data is (X 0 , θ0 ). k = 0. Assume that X 0 satisfies Ai • X 0 = bi , i = 1, . . . , m, X 0  0. ¯ =L ¯L ¯T . ¯ = X k . Factor X Step 1. Set Current values. X Step 2. Compute Newton Direction and Multipliers. Compute the
¯ by solving the following system Newton step D for BSDP (θ0 ) at X = X of equations in the variables (D, y): ⎧ m 0 ¯ −1 0 ¯ −1 ⎪ ¯ −1 =
yi Ai ⎪ ⎨ C − θ X + θ X DX ⎪ ⎪ ⎩ i=1 Ai • D = 0, (10) i = 1, . . . , m.

Denote the solution to this system by (D , y ). Set
S =C− m 
yi Ai . i=1 ¯ −T  ≤ 1 , stop. In this ¯ −1 D
L Step 3. Test the Current Point. If L 4 ¯ is a 1 -approximate solution of BSDP (θ0 ), along with the dual valcase, X 4

ues (y , S ). Step 4. Update Primal Point.
¯ + αD X =X 42
where α= 0.2 −1 ¯ ¯ −T  . L D
L
¯ + αD ). Alternatively, α can be computed by a line-search of fθ0 (X
Step 5. Reset Counter and Continue. X k+1 ← X , k ← k + 1. Go to Step 1. The following proposition validates Step 3 of the algorithm:

Proposition 10.1 Suppose that (D , y ) is the solution of the Normal equa¯ for the given value θ0 of the barrier parameter, tions (10) for the point X and that ¯ −T  ≤ 1 . ¯ −1 D
L L 4 1 ¯ is a -approximate solution of BSDP (θ0 ). Then X 4 Proof: We must exhibit values (y, S) that satisfy     ¯ ≤ 1 . ¯T SL I − 1 L   θ0 4

m
i=1 yi Ai + S = C and
Let (D , y ) solve the Normal equations (10), and let S = C − we have from (10) that I− m
i=1
yi Ai . Then 1 ¯ T 0 ¯ −T ¯ −1 ¯ −T ¯ −1
¯ −T ¯ −1 ¯ 1 ¯T
¯ ¯ −T , ¯ −1 D
L L S L = I− 0 L (θ (L L −L L D L L ))L = L 0 θ θ whereby I − 1 ¯T
¯ ¯ −1 D
L ¯ −T  ≤ 1 . L S L = L 0 θ 4 q.e.d. The next proposition shows that whenever the algorithm proceeds to Step 4, then the objective function fθ0 (X) decreases by at least 0.025θ0 : 43 ¯ satisfies Ai • X ¯ = bi , i = 1, . . . , m, and Proposition 10.2 Suppose that X

¯  0. Suppose that (D , y ) is the solution of the Normal equations (10) X ¯ for a given value θ0 of the barrier parameter, and that for the point X ¯ −1 D
L ¯ −T  > 1 . L 4 Then for all γ ∈ [0, 1),  γ
¯+ fθ0 X
¯ −T D −1 ¯ L D L   ¯ + θ0 −γL ¯ −1 D
L ¯ −T  + ≤ fθ0 (X) γ2 . 2(1 − γ) In particular,  0.2
¯+ fθ0 X
¯ −T D −1 ¯ L D L   ¯ − 0.025θ0 . ≤ fθ0 (X) (11) In order to prove this proposition, we will need two powerful facts about the logarithm function: Fact 1. Suppose that |x| ≤ δ < 1. Then ln(1 + x) ≥ x − x2 . 2(1 − δ) Proof: We have: ln(1 + x) = x − x2 2 + x3 3 − x4 4 + ... ≥ x− |x|2 2 − |x|3 3 − |x|4 4 − ... ≥ x− |x|2 2 − |x|3 2 − |x|4 2 − ... = x− x2 2 1 + |x| + |x|2 + |x|3 + . . . = x− x2 2(1−|x|) ≥ x− x2 2(1−δ) . ! 44 q.e.d. Fact 2. Suppose that R ∈ S n and that R ≤ γ < 1. Then ln(det(I + R)) ≥ I • R − γ2 . 2(1 − γ) Proof: Factor R = QDQT where Q is orthonormal and D"is a diagonal matrix of the eigenvalues of R. Then first note that R = then have: ln(det(I + R)) = ln(det(I + QDQT )) = ln(det(I + D)) = ≥ n
ln(1 + Djj ) j=1  n
j=1 Djj − = I •D− 2 Djj 2(1−γ) R2 2(1−γ) ≥ I • QT RQ − = I •R− Proof of Proposition 10.2. Let α= and notice that γ ¯ −1 D
L ¯ −T  , L ¯ −1 D
L ¯ −T  = γ. αL 45 γ2 2(1−γ) γ2 2(1−γ) q.e.d.  n
j=1 2 . We Djj Then  ¯+ fθ0 X
γ ¯ −1 D
L̄−T  D L   ¯ + αD
= fθ0 X  ¯ −T )L ¯ + αC • D
− θ0 ln(det(L(I ¯ + αL ¯ −1 D
L ¯ T )) = C •X ¯ + αC • D
− θ0 ln(det(I + αL ¯ −1 D
L ¯ −T )) ¯ − θ0 ln(det(X)) = C •X ¯ + αC • D
− θ0 αI • L ¯ −1 D
L ¯ −T + θ0 γ 2 ≤ fθ0 (X) 2(1−γ) ¯ −1 • D
+ θ0 γ 2 ¯ + αC • D
− θ0 αL ¯ −T L = fθ0 (X) 2(1−γ) !
2 ¯ −1 • D + θ0 γ . ¯ + α C − θ0 X = fθ0 (X) 2(1−γ)

Now, (D , y ) solve the Normal equations: ⎧ m 0 ¯ −1 0 ¯ −1
¯ −1 =
y
A ⎪ ⎪ i i ⎨ C −θ X +θ X D X ⎪ ⎪ ⎩ i=1
Ai • D = 0, (12) i = 1, . . . , m. Taking the inner product of both sides of the first equation above with D and rearranging yields: ¯ −1 D
= −(C − θ0 X ¯ −1 D
• X ¯ −1 ) • D
. θ0 X Substituting this in our inequality above yields: 46
 ¯+ fθ0 X
γ ¯ −1 D
L̄−T  D L  ¯ − αθ0 X ¯ −1 D
+ θ0 γ 2 ¯ −1 D
• X ≤ fθ0 (X) 2(1−γ) ¯ − αθ0 L ¯ −1 D
L ¯ −T • L ¯ −1 D
L ¯ −T + θ0 γ 2 = fθ0 (X) 2(1−γ) ¯ − αθ0 L ¯ −1 D
L ¯ −T 2 + θ0 γ 2 = fθ0 (X) 2(1−γ) ¯ − θ0 γL ¯ −1 D
L ¯ −T  + θ0 γ 2 = fθ0 (X) 2(1−γ)  ¯ + θ0 −γL ¯ −1 D
L ¯ −T  + = fθ0 (X) ¯ −1 D
L ¯ −T  > Subsituting γ = 0.2 and and L q.e.d. 1 4 γ2 2(1−γ)  . yields the final result. Last of all, we prove a bound on the number of iterations that the algorithm will need in order to find a 14 -approximate solution of BSDP (θ0 ): Proposition 10.3 Suppose that X 0 satisfies Ai • X 0 = bi , i = 1, . . . , m, and X 0  0. Let θ0 be given and let fθ∗0 be the optimal objective function value of BSDP (θ0 ). Then the algorithm initiated at X 0 will find a 14 -approximate solution of BSDP (θ0 ) in at most  fθ0 (X 0 ) − fθ∗0 k= 0.025θ0  iterations. Proof: This follows immediately from (11). Each iteration that is not a 14 approximate solution decreases the objective function fθ0 (X) of BSDP (θ0 ) by at least 0.025θ0 . Therefore, there cannot be more than  fθ0 (X 0 ) − fθ∗0 0.025θ0 47  iterations that are not 14 -approximate solutions of BSDP (θ0 ). q.e.d. 10.5 Some Properties of the Frobenius Norm Proposition 10.4 If M ∈ S n , then " 1. M  = n
(λj (M ))2 , where λ1 (M ), λ2 (M ), . . . , λn (M ) is an enu- j=1 meration of the n eigenvalues of M . 2. If λ is any eigenvalue of M , then |λ| ≤ M . √ 3. |trace(M )| ≤ nM . 4. If M  < 1, then I + M  0. Proof: We can factorize M = QDQT where Q is orthonormal and D is a diagonal matrix of the eigenvalues of M . Then √ M  = M • M = QDQT • QDQT = trace(QDQT QDQT ) =    n trace(DD) =  (λj (M ))2 . trace(QT QDQT QD) = j=1 This proves the first two assertions. To prove the third assertion, note that trace(M ) = trace(QDQT ) = trace(QT QD) = trace(D) = n    n √  √  λj (M ) ≤ n (λj (M ))2 = nM . j=1 j=1
To prove the fourth assertion, let λ be an eigenvalue of I + M . Then
λ = 1 + λ where λ is an eigenvalue of M . However, from the second asser
tion, λ = 1 + λ ≥ 1 − M  > 0, and so M  0. q.e.d. 48 Proposition 10.5 If A, B ∈ S n , then AB ≤ AB. Proof: We have   n  n    AB = n  i=1 j=1 k=1   n  n  ≤ n  i=1 j=1 k=1    n n   2 = Aik ⎛ ⎝ i=1 k=1 2 Aik Bkj n  A2ik 2 Bkj k=1 n  n  ⎞ 2 ⎠ = AB. Bkj j=1 k=1 q.e.d. 11 Issues in Solving SDP using the Ellipsoid Al­ gorithm To see how the ellipsoid algorithm is used to solve a semidefinite program, assume for convenience that the format of the problem is that of the dual problem SDD. Then the feasible region of the problem can be written as: # F = (y1 , . . . , ym ) ∈  m |C− m  $ yi Ai  0 , i=1 and the objective function is gion in m .
m i=1 yi bi . Note that F is just a convex re- Recall that at any iteration of the ellipsoid algorithm, the set of solutions of SDD is known to lie in the current ellipsoid, and the center of the current 49 ellipsoid is, say, y¯ = (¯ y1 , . . . , ȳm ). If y¯ ∈ F , then we perform an optimality

m cut of the form m y i=1 i bi ≥ i=1 ȳi bi , and use standard formulas to update the ellipsoid and its new center. If y¯ ∈ / F , then we perform a feasibility cut ¯ ∈ m such that h ¯T y > h ¯ T ȳ for all y ∈ F . by computing a vector h There are four issues that must be resolved in order to implement the above version of the ellipsoid algorithm to solve SDD: 1. Testing if ȳ ∈ F . This is done by computing the matrix S¯ = C −
m ¯ ¯ i=1 ȳi Ai . If S  0, then ȳ ∈ F . Testing if the matrix S  0 can be done by computing an exact Cholesky factorization of S̄, which takes O(n3 ) operations, assuming that computing square roots can be performed (exactly) in one operation. 2. Computing a feasibility cut. As above, testing if S¯  0 can be computed in O(n3 ) operations. If S¯  0, then again assuming exact arith¯v < 0 in O(n3 ) opermetic, we can find an n-vector v¯ such that v¯T S¯ ¯ is computed by the ations as well. Then the feasiblity cut vector h formula: ¯ i = v¯T Ai v, ¯ i = 1, . . . , m, h whose computation requires O(mn2 ) operations. Notice that for any y ∈ F , that ¯= yT h m  yi v̄ T Ai v̄ = v¯T i=1 ≤ v̄ T C v¯ < v̄ T m  m  yi Ai v̄ i=1 ȳi Ai v̄ = i=1 m  ¯ ȳi v̄ T Ai v̄ = y¯T h, i=1 ¯ indeed provides a feasibility cut for F at y = y. thereby showing h ¯ 3. Starting the ellipsoid algorithm. We need to determine an upper bound R on the distance of some optimal solution y ∗ from the origin. This 50 cannot be done by examining the input length of the data, as is the case in linear programming. One needs to know some special information about the specific problem at hand in order to determine R before solving the semidefinite program. 4. Stopping the ellipsoid algorithm. Suppose that we seek an -optimal solution of SDD. In order to prove a complexity bound on the number of iterations needed to find an -optimal solution, we need to know beforehand the radius r of a Euclidean ball that is contained in the set of -optimal solutions of SDD. The value of r also cannot be determined by examining the input length of the data, as is the case in linear programming. One needs to know some special information about the specific problem at hand in order to determine r before solving the semidefinite program. 12 Current Research in SDP There are many very active research areas in semidefinite programming in nonlinear (convex) programming, in combinatorial optimization, and in control theory. In the area of convex analysis, recent research topics include the geometry and the boundary structure of SDP feasible regions (including notions of degeneracy) and research related to the computational complexity of SDP such as decidability questions, certificates of infeasibility, and duality theory. In the area of combinatorial optimization, there has been much research on the practical and the theoretical use of SDP relaxations of hard combinatorial optimization problems. As regards interior point methods, there are a host of research issues, mostly involving the development of different interior point algorithms and their properties, including rates of convergence, performance guarantees, etc. 13 Computational State of the Art of SDP Because SDP has so many applications, and because interior point methods show so much promise, perhaps the most exciting area of research on SDP has to do with computation and implementation of interior point algorithms 51 for solving SDP . Much research has focused on the practical efficiency of interior point methods for SDP . However, in the research to date, computational issues have arisen that are much more complex than those for linear programming, and these computational issues are only beginning to be well-understood. They probably stem from a variety of factors, including the fact that SDP is not guaranteed to have strictly complementary optimal solutions (as is the case in linear programming). Finally, because SDP is such a new field, there is no representative suite of practical problems on which to test algorithms, i.e., there is no equivalent version of the netlib suite of industrial linear programming problems. A good website for semidefinite programming is: http://www.zib.de/helmberg/semidef.html. 14 Exercises 1. For a (square) matrix M ∈ IRn×n , define trace(M ) = n
j=1 Mjj , and for two matrices A, B ∈ IRk×l define A • B := k  l  Aij Bij . i=1 j=1 Prove that: (a) A • B = trace(AT B). (b) trace(M N ) = trace(N M ). k×k denote the cone of positive semi-definite symmetric matrices, 2. Let S+ k×k namely S+ = {X ∈ S k×k | v T Xv ≥ 0 for all v ∈ n }. Considering  k×k k×k as a cone, prove that S+ S+ is self-dual. 52 ∗ n×n k×k = S+ , thus showing that S+ 3. Consider the problem: (P) : minimized dT Qd dT M d = 1 . s.t. If M  0, show that (P) is equivalent to: (S) : minimizeX Q•X M •X =1 s.t. X0. What is the SDP dual of (S)? 4. Prove that X  xxT if and only if  X xT x 1  0. 5. Let K := {X ∈ S k×k | X  0} and define K ∗ := {S ∈ S k×k | S • X ≥ 0 for all X  0} . Prove that K ∗ = K. 6. Let λmin denote the smallest eigenvalue of the symmetric matrix Q. Show that the following three optimization problems each have optimal objective function value equal to λmin . (P1) : minimized dT Qd s.t. dT Id = 1 . (P2) : maximizeλ λ s.t. 53 Q  λI . Q•X (P3) : minimizeX I •X =1 X0. s.t. 7. Suppose that Q  0. Prove the following: xT Qx + q T x + c ≤ 0 if and only if there exists W for which c 1 2q 1 T 2q Q • 1 xT x W ≤ 0 and 1 xT x W 0. Hint: use the property shown in Exercise 4.  8. Consider the matrix M= A BT B C  , where A, B are symmetric matrices and A is nonsingular. Prove that M  0 if and only if C − B T A−1 B  0. 54