6.7 Semidefinite Optimization

Semidefinite optimization is a generalization of conic optimization, allowing the use of matrix variables belonging to the convex cone of positive semidefinite matrices

\[\PSD^r = \left\lbrace X \in \Symm^r: z^T X z \geq 0, \quad \forall z \in \real^r \right\rbrace,\]

where \(\Symm^r\) is the set of \(r \times r\) real-valued symmetric matrices.

MOSEK can solve semidefinite optimization problems stated in the primal form,

(6.19)\[\begin{split}\begin{array}{lccccll} \mbox{minimize} & & & \sum_{j=0}^{p-1} \left\langle \barC_j, \barX_j \right\rangle + \sum_{j=0}^{n-1} c_j x_j + c^f & & &\\ \mbox{subject to} & l_i^c & \leq & \sum_{j=0}^{p-1} \left\langle \barA_{ij}, \barX_j \right\rangle + \sum_{j=0}^{n-1} a_{ij} x_j & \leq & u_i^c, & i = 0, \ldots, m-1,\\ & & & \sum_{j=0}^{p-1} \left\langle \barF_{ij}, \barX_j \right\rangle + \sum_{j=0}^{n-1} f_{ij} x_j + g_i & \in & \K_{i}, & i = 0, \ldots, q-1,\\ & l_j^x & \leq & x_j & \leq & u_j^x, & j = 0, \ldots, n-1,\\ & & & x \in \K, \barX_j \in \PSD^{r_j}, & & & j = 0, \ldots, p-1 \end{array}\end{split}\]

where the problem has \(p\) symmetric positive semidefinite variables \(\barX_j\in \PSD^{r_j}\) of dimension \(r_j\). The symmetric coefficient matrices \(\barC_j\in \Symm^{r_j}\) and \(\barA_{i,j}\in \Symm^{r_j}\) are used to specify PSD terms in the linear objective and the linear constraints, respectively. The symmetric coefficient matrices \(\barF_{i,j}\in \Symm^{r_j}\) are used to specify PSD terms in the affine conic constraints. Note that \(q\) ((6.19)) is the total dimension of all the cones, i.e. \(q=\text{dim}(\K_1 \times \ldots \times \K_k)\), given there are \(k\) ACCs. We use standard notation for the matrix inner product, i.e., for \(A,B\in \real^{m\times n}\) we have

\[\left\langle A,B \right\rangle := \sum_{i=0}^{m-1} \sum_{j=0}^{n-1} A_{ij} B_{ij}.\]

In addition to the primal form presented above, semidefinite problems can be expressed in their dual form. Constraints in this form are usually called linear matrix inequalities (LMIs). LMIs can be easily specified in MOSEK using the vectorized positive semidefinite cone which is defined as:

  • Vectorized semidefinite domain:

    \[\PSD^{d,\mathrm{vec}} = \left\{(x_1,\ldots,x_{d(d+1)/2})\in \real^n~:~ \mathrm{sMat}(x)\in\PSD^d\right\},\]

    where \(n=d(d+1)/2\) and,

    \[\begin{split}\mathrm{sMat}(x) = \left[\begin{array}{cccc}x_1 & x_2/\sqrt{2} & \cdots & x_{d}/\sqrt{2} \\ x_2/\sqrt{2} & x_{d+1} & \cdots & x_{2d-1}/\sqrt{2} \\ \cdots & \cdots & \cdots & \cdots \\ x_{d}/\sqrt{2} & x_{2d-1}/\sqrt{2} & \cdots & x_{d(d+1)/2}\end{array}\right],\end{split}\]

    or equivalently

    \[\PSD^{d,\mathrm{vec}} = \left\{\mathrm{sVec}(X)~:~X\in\PSD^d\right\},\]

    where

    \[\mathrm{sVec}(X) = (X_{11},\sqrt{2}X_{21},\ldots,\sqrt{2}X_{d1},X_{22},\sqrt{2}X_{32},\ldots,X_{dd}).\]

In other words, the domain consists of vectorizations of the lower-triangular part of a positive semidefinite matrix, with the non-diagonal elements additionally rescaled. LMIs can be expressed by restricting appropriate affine expressions to this cone type.

For other types of cones supported by MOSEK, see Sec. 13.6 (Supported domains) and the other tutorials in this chapter. Different cone types can appear together in one optimization problem.

We demonstrate the setup of semidefinite variables and their coefficient matrices in the following examples:

6.7.1 Example SDO1

We consider the simple optimization problem with semidefinite and conic quadratic constraints:

(6.20)\[\begin{split}\begin{array} {llcc} \mbox{minimize} & \left\langle \left[ \begin{array} {ccc} 2 & 1 & 0 \\ 1 & 2 & 1 \\ 0 & 1 & 2 \end{array} \right], \barX \right\rangle + x_0 & & \\ \mbox{subject to} & \left\langle \left[ \begin{array} {ccc} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{array} \right], \barX \right\rangle + x_0 & = & 1, \\ & \left\langle \left[ \begin{array}{ccc} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{array} \right], \barX \right\rangle + x_1 + x_2 & = & 1/2, \\ & x_0 \geq \sqrt{{x_1}^2 + {x_2}^2}, & \barX \succeq 0, & \end{array}\end{split}\]

The problem description contains a 3-dimensional symmetric semidefinite variable which can be written explicitly as:

\[\begin{split}\barX = \left[ \begin{array} {ccc} \barX_{00} & \barX_{10} & \barX_{20} \\ \barX_{10} & \barX_{11} & \barX_{21} \\ \barX_{20} & \barX_{21} & \barX_{22} \end{array} \right] \in \PSD^3,\end{split}\]

and an affine conic constraint (ACC) \((x_0, x_1, x_2) \in \Q^3\). The objective is to minimize

\[2(\barX_{00} + \barX_{10} + \barX_{11} + \barX_{21} + \barX_{22}) + x_0,\]

subject to the two linear constraints

\[\begin{split}\begin{array}{ccc} \barX_{00} + \barX_{11} + \barX_{22} + x_0 & = & 1, \\ \barX_{00} + \barX_{11} + \barX_{22} + 2(\barX_{10} + \barX_{20} + \barX_{21}) + x_1 + x_2 & = & 1/2. \end{array}\end{split}\]

Setting up the linear and conic part

The linear and conic parts (constraints, variables, objective, ACC) are set up using the methods described in the relevant tutorials; Sec. 6.1 (Linear Optimization), Sec. 6.2 (From Linear to Conic Optimization). Here we only discuss the aspects directly involving semidefinite variables.

Appending semidefinite variables

The dimensions of semidefinite variables are passed in prob$bardim.

Coefficients of semidefinite terms.

Every term of the form \((\barA_{i,j})_{k,l}(\barX_j)_{k,l}\) is determined by four indices \((i,j,k,l)\) and a coefficient value \(v=(\barA_{i,j})_{k,l}\). Here \(i\) is the number of the constraint in which the term appears, \(j\) is the index of the semidefinite variable it involves and \((k,l)\) is the position in that variable. This data is passed in the structure prob$barA. Note that only the lower triangular part should be specified explicitly, that is one always has \(k\geq l\).

Semidefinite terms \((\barC_j)_{k,l}(\barX_j)_{k,l}\) of the objective are specified in the same way in prob$barc but only include \((j,k,l)\) and \(v\).

Source code

Listing 6.8 R implementation of model (6.20). Click here to download.
library("Rmosek")

getbarvarMatrix <- function(barvar, bardim)
{
    N <- as.integer(bardim)
    new("dspMatrix", x=barvar, uplo="L", Dim=c(N,N))
}

sdo1 <- function()
{
    # Specify the non-matrix variable part of the problem.
    prob <- list(sense="min")
    prob$c     <- c(1, 0, 0)
    prob$A <- sparseMatrix(i=c(1, 2, 2),
                           j=c(1, 2, 3),
                           x=c(1, 1, 1), dims=c(2, 3))
    prob$bc    <- rbind(blc=c(1, 0.5), 
                        buc=c(1, 0.5))
    prob$bx    <- rbind(blx=rep(-Inf,3), 
                        bux=rep( Inf,3))

    # NOTE: The F matrix is internally stored in the sparse
    #       triplet form. Use 'giveCsparse' or 'repr' option 
    #       in the sparseMatrix() call to construct the F 
    #       matrix directly in the sparse triplet form. 
    prob$F     <- sparseMatrix(i=c(1,2,3), 
                               j=c(1,2,3), 
                               x=c(1,1,1), 
                               dims = c(3,3))
    prob$g     <- c(1:3)*0
    prob$cones <- cbind(list("QUAD", 3, NULL))

    # Specify semidefinite matrix variables (one 3x3 block)
    prob$bardim <- c(3)

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barc':
    prob$barc$j <- c(1, 1, 1, 1, 1)
    prob$barc$k <- c(1, 2, 3, 2, 3)
    prob$barc$l <- c(1, 2, 3, 1, 2)
    prob$barc$v <- c(2, 2, 2, 1, 1)

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barA':
    prob$barA$i <- c(1, 1, 1, 2, 2, 2, 2, 2, 2)
    prob$barA$j <- c(1, 1, 1, 1, 1, 1, 1, 1, 1)
    prob$barA$k <- c(1, 2, 3, 1, 2, 3, 2, 3, 3)
    prob$barA$l <- c(1, 2, 3, 1, 2, 3, 1, 1, 2)
    prob$barA$v <- c(1, 1, 1, 1, 1, 1, 1, 1, 1)

    # Solve the problem
    r <- mosek(prob)

    # Print matrix variable and return the solution
    stopifnot(identical(r$response$code, 0))
    print( list(barx=getbarvarMatrix(r$sol$itr$barx[[1]], prob$bardim[1])) )
    r$sol
}

sdo1()

The numerical values of \(\barX_j\) are returned in the list r$sol$itr$barx; the \(j\)-th element of the list is the lower triangular part of each \(\barX_j\) stacked column-by-column into a numeric vector. Similarly, the dual semidefinite variables \(\barS_j\) are recovered through r$sol$itr$bars.

6.7.2 Example SDO2

We now demonstrate how to define more than one semidefinite variable using the following problem with two matrix variables and two types of constraints:

(6.21)\[\begin{split}\begin{array}{lrll} \mbox{minimize} & \langle C_1,\barX_1\rangle + \langle C_2,\barX_2\rangle & & \\ \mbox{subject to} & \langle A_1,\barX_1\rangle + \langle A_2,\barX_2\rangle & = & b, \\ & (\barX_2)_{01} & \leq & k, \\ & \barX_1, \barX_2 & \succeq & 0. \end{array}\end{split}\]

In our example \(\dim(\barX_1)=3\), \(\dim(\barX_2)=4\), \(b=23\), \(k=-3\) and

\[\begin{split}C_1= \left[\begin{array}{ccc} 1 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 6 \end{array}\right], A_1= \left[\begin{array}{ccc} 1 & 0 & 1 \\ 0 & 0 & 0 \\ 1 & 0 & 2 \end{array}\right],\end{split}\]
\[\begin{split}C_2= \left[\begin{array}{cccc} 1 & -3 & 0 & 0\\ -3 & 2 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 \\ \end{array}\right], A_2= \left[\begin{array}{cccc} 0 & 1 & 0 & 0\\ 1 & -1 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -3 \\ \end{array}\right],\end{split}\]

are constant symmetric matrices.

Note that this problem does not contain any scalar variables, but they could be added in the same fashion as in Sec. 6.7.1 (Example SDO1).

For explanations of data structures used in the example see Sec. 6.7.1 (Example SDO1). Note that the field bardim is used to specify that we have two semidefinite variables of dimensions 3 and 4.

The code representing the above problem is shown below.

Listing 6.9 Implementation of model (6.21). Click here to download.
library("Rmosek")

getbarvarMatrix <- function(barvar, bardim)
{
    N <- as.integer(bardim)
    new("dspMatrix", x=barvar, uplo="L", Dim=c(N,N))
}

sdo2 <- function()
{
    # Sample data in sparse lower-triangular triplet form
    C1_k <- c(1, 3);
    C1_l <- c(1, 3);
    C1_v <- c(1, 6);
    A1_k <- c(1, 3, 3);
    A1_l <- c(1, 1, 3);
    A1_v <- c(1, 1, 2);
    C2_k <- c(1, 2, 2, 3);
    C2_l <- c(1, 1, 2, 3);
    C2_v <- c(1, -3, 2, 1);
    A2_k <- c(2, 2, 4);
    A2_l <- c(1, 2, 4);
    A2_v <- c(1, -1, -3);
    b <- 23.0;
    k <- -3.0;
 
    # Specify the dimensions of the problem
    prob <- list(sense="min");
    # Two constraints
    prob$A <- Matrix(nrow=2, ncol=0);                   # 2 constraints
    prob$c <- numeric(0);
    prob$bx <- rbind(blc=numeric(0), 
                     buc=numeric(0));

    # Dimensions of semidefinite matrix variables 
    prob$bardim <- c(3, 4);
    # Constraint bounds    
    prob$bc <- rbind(blc=c(b, -Inf), 
                     buc=c(b, k));

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barc':
    prob$barc$j <- c(rep(1, length(C1_v)), 
                     rep(2, length(C2_v)));         # Which PSD variable (j)
    prob$barc$k <- c(C1_k, C2_k);                   # Entries: (k,l)->v
    prob$barc$l <- c(C1_l, C2_l);
    prob$barc$v <- c(C1_v, C2_v);

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barA':
    prob$barA$i <- c(rep(1, length(A1_v)+length(A2_v)), 2);   # Which constraint (i)
    prob$barA$j <- c(rep(1, length(A1_v)),
                     rep(2, length(A2_v)),
                     2);                                      # Which PSD variable (j)
    prob$barA$k <- c(A1_k, A2_k, 2);                          # Entries: (k,l)->v
    prob$barA$l <- c(A1_l, A2_l, 1);
    prob$barA$v <- c(A1_v, A2_v, 0.5);

    # Solve the problem
    r <- mosek(prob);

    # Print matrix variable and return the solution
    stopifnot(identical(r$response$code, 0));
    print( list(X1=getbarvarMatrix(r$sol$itr$barx[[1]], prob$bardim[1])) );
    print( list(X2=getbarvarMatrix(r$sol$itr$barx[[2]], prob$bardim[2])) );

}

sdo2();

The numerical values of \(\barX_j\) are returned in the list r$sol$itr$barx; the \(j\)-th element of the list is the lower triangular part of each \(\barX_j\) stacked column-by-column into a numeric vector. Similarly, the dual semidefinite variables \(\barS_j\) are recovered through r$sol$itr$bars.

6.7.3 Example SDO_LMI: Linear matrix inequalities and the vectorized semidefinite domain

The standard form of a semidefinite problem is usually either based on semidefinite variables (primal form) or on linear matrix inequalities (dual form). However, MOSEK allows mixing of these two forms, as shown in (6.22)

(6.22)\[\begin{split}\begin{array} {llcc} \mbox{minimize} & \left\langle \left[ \begin{array} {cc} 1 & 0 \\ 0 & 1 \end{array} \right], \barX \right\rangle + x_0 + x_1 + 1 & & \\ \mbox{subject to} & \left\langle \left[ \begin{array} {cc} 0 & 1 \\ 1 & 0 \end{array} \right], \barX \right\rangle - x_0 - x_1 & \in & \real_{\geq 0}^1, \\ & x_0 \left[ \begin{array}{cc} 0 & 1 \\ 1 & 3 \end{array} \right] + x_1 \left[ \begin{array}{cc} 3 & 1 \\ 1 & 0 \end{array} \right] - \left[ \begin{array}{cc} 1 & 0 \\ 0 & 1 \end{array} \right] & \succeq & 0, \\ & \barX \succeq 0. & & \end{array}\end{split}\]

The first affine expression is restricted to a linear domain and could also be modelled as a linear constraint (instead of an ACC). The lower triangular part of the linear matrix inequality (second constraint) can be vectorized and restricted to the "MSK_DOMAIN_SVEC_PSD_CONE". This allows us to express the constraints in (6.22) as the affine conic constraints shown in (6.23).

(6.23)\[\begin{split}\begin{array}{ccccccc} \left\langle\left[\begin{array}{cc}0&1\\1&0\end{array}\right],\barX \right\rangle & + & \left[\begin{array}{cc}-1&-1\end{array}\right] x & + & \left[\begin{array}{c}0\end{array}\right] & \in & \real_{\geq 0}^1, \\ & & \left[\begin{array}{cc}0&3\\ \sqrt{2}&\sqrt{2}\\3&0\end{array}\right] x & + & \left[\begin{array}{c}-1\\0\\-1\end{array}\right] & \in & \PSD^{3,\mathrm{vec}} \end{array}\end{split}\]

Vectorization of the LMI is performed as explained in Sec. 13.6 (Supported domains).

Setting up the linear part

The linear parts (objective, constraints, variables) and the semidefinite terms in the linear expressions are defined exactly as shown in the previous examples.

Setting up the affine conic constraints with semidefinite terms

To define affine conic constraints, we set prob$F and prob$g to the values that are shown in (6.23). The coefficient for the semidefinite variable is defined by setting prob$barF equal to the symmetric matrix shown in (6.23).

    prob$barF$i <- c(1, 1)
    prob$barF$j <- c(1, 1)
    prob$barF$k <- c(1, 2)
    prob$barF$l <- c(1, 1)
    prob$barF$v <- c(0, 1)

The domains are specified as columns in prob$cones where the first row selects the "type" of cone, such as "MSK_DOMAIN_RPLUS" or "MSK_DOMAIN_SVEC_PSD_CONE" (note that RPLUS and SVEC_PSD_CONE are valid aliases for each domain, respectively). The second row specifies the dimension ("dim"), here set to 1 and 3, respectively. The third row sets the parameters ("conepar") for parametric domains but here it is set to NULL.

    prob$F     <- sparseMatrix(i=c(1, 1, 2, 3, 3, 4),
                               j=c(1, 2, 2, 1, 2, 1),
                               x=c(-1, -1, 3, sqrt(2), sqrt(2), 3), dims=c(4, 2))
    prob$g     <- c(0, -1, 0, -1)
    prob$cones <- matrix(list(), nrow=3, ncol=2)
    prob$cones[,1] <- list("RPLUS", 1, NULL)
    prob$cones[,2] <- list("SVEC_PSD_CONE", 3, NULL)

Source code

Listing 6.10 Source code solving problem (6.22). Click here to download.
library("Rmosek")

getbarvarMatrix <- function(barvar, bardim)
{
    N <- as.integer(bardim)
    new("dspMatrix", x=barvar, uplo="L", Dim=c(N,N))
}

sdo_lmi <- function()
{
    # Specify the non-matrix variable part of the problem.
    prob       <- list(sense="min")
    prob$c     <- c(1, 1)
    prob$c0    <- 1

    # Specify variable bounds
    prob$bx    <- rbind(blx=rep(-Inf,2), bux=rep( Inf,2))
    # The following two entries must always be defined, even if set to zero.
    prob$A  <- Matrix(c(0, 0), nrow=1, sparse=TRUE)
    prob$bc <- rbind(blc=rep(-Inf,1),buc=rep(Inf, 1))

    prob$F     <- sparseMatrix(i=c(1, 1, 2, 3, 3, 4),
                               j=c(1, 2, 2, 1, 2, 1),
                               x=c(-1, -1, 3, sqrt(2), sqrt(2), 3), dims=c(4, 2))
    prob$g     <- c(0, -1, 0, -1)
    prob$cones <- matrix(list(), nrow=3, ncol=2)
    prob$cones[,1] <- list("RPLUS", 1, NULL)
    prob$cones[,2] <- list("SVEC_PSD_CONE", 3, NULL)

    # Specify semidefinite matrix variables (one 2x2 block)
    prob$bardim <- c(2)

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barc':
    prob$barc$j <- c(1, 1, 1)
    prob$barc$k <- c(1, 2, 2)
    prob$barc$l <- c(1, 1, 2)
    prob$barc$v <- c(1, 0, 1)

    # Block triplet format specifying the lower triangular part 
    # of the symmetric coefficient matrix 'barF' for the ACC:
    prob$barF$i <- c(1, 1)
    prob$barF$j <- c(1, 1)
    prob$barF$k <- c(1, 2)
    prob$barF$l <- c(1, 1)
    prob$barF$v <- c(0, 1)

    # Solve the problem
    r <- mosek(prob)

    # Print matrix variable and return the solution
    stopifnot(identical(r$response$code, 0))
    print( list(barx=getbarvarMatrix(r$sol$itr$barx[[1]], prob$bardim[1]), xx=r$sol$itr$xx) )
    r$sol
}

sdo_lmi()