Two versions of the library

SDPSL comes in two versions, whose functionality is basically the same. One version uses double-precision floating-point numbers, the other uses arbitrary-precision floating-point numbers via the GMP library. To use the double-precision version, one includes the header sdpsl.h. All names defined by the library are then part of the namespace sdp. To use the GMP version, one includes the header mpfsdpsl.h, and then all names are part of the namespace mpf_sdp.

The main difference between both versions lies in the type Number, defined in the file number.h. In the double-precision version, Number is defined as double. In the GMP version, it is defined as mpf_class. Changing the type of Number changes the numbers used as coefficients in polynomials, as entries in the matrices that define constraints of the semidefinite programs, etc.

See Also: number.h

The basics of problem generation

Variables

Before you can generate a problem, you have to define its variable matrices. These are instances of the class SDPVariable. Each variable has a unique index, assigned to it automatically, which is later used to identify it. Each variable also has an associated size, which is the dimension of the variable matrix.

Variables can be created as follows:

SDPVariable X(1); // A 1-by-1 variable

SDPVariable X(5); // A 5-by-5 variable

Note that variables are merely placeholders, they do not in fact contain a matrix. Once a variable is created, a new index is assigned to it, and even if it is destroyed this index is not assigned to any other variable. Variables can later be used to generate linear and polynomial constraints and also to set the objective function of the problem.

Linear combinations

The left-hand side of a linear constraint of a semidefinite program, or the objective function of the program, is an expression of the kind

\[ \sum_{j=1}^r \langle A_j, X_j \rangle, \]

where the $A_j$ are numerical matrices and the $X_j$ are the variables of the problem.

In SDPSL, such expressions are represented by instances of SDPLinearCombination, which is basically a map that associates numerical matrices (instances of Matrix<Number>) to variables.

Such maps can be created using the multiplication and addition operators. Multiplying an instance of Matrix<Number> by an instance of SDPVariable gives an instance of SDPLinearCombination associating the matrix to the variable. Different instances of SDPLinearCombination can be added, meaning that the maps are joined together.

using namespace sdp;  // We use the double-precision version
SDPVariable X(5), Y(7);
Matrix<double> A(5), B(7);
SDPLinearCombination C = A * X + B * Y;

SDPLinearCombination keeps a reference to the matrices. This means that, if the matrices are changed, then the combination itself also changes. These are not standard C++ references, but rather shared references. This means that the data of the matrices will be kept as long as the instance of SDPLinearCombination exists. The data of a certain matrix will be destroyed upon destruction of the instance of SDPLinearCombination if no other shared references to it exist.

Addition of instances of SDPLinearCombination does not mean that matrices are added. If the combinations being added share a variable, then an exception is thrown.

Since we keep references for the matrices used, copying an instance of SDPLinearCombination does not mean that matrices are copied. However, the internal map associating matrices to variables is copied. So in cases when many variables are used this can be an expensive operation, but recall one of our design assumptions is that problems have few variables.

See Also: shared_ref, const_shared_ref, Matrix

Linear constraints

Linear constraints are represented by instances of the class SDPLinearConstraint. Such instances contain an instance of SDPLinearCombination, representing the left-hand side of the constraint, the relation of the constraint (one of $=$, $\leq$, or $\geq$), and the right-hand side.

Instances of SDPLinearConstraint can be created from instances of SDPLinearCombination and numbers using the operators ==, <=, and >=.

using namespace sdp;  // We use the double-precision version
SDPVariable X(5), Y(7);
Matrix<double> A(5), B(7);
SDPLinearConstraint C = A * X + B * Y <= 5.0;

Constraint containers and problems

Linear constraints, that is, instances of SDPLinearConstraint, can be added to problems via so-called constraint containers, instances of SDPConstraintContainer. These are objects that are used to transform instances of SDPLinearConstraint into raw data that is later used to generate solver input, for example. Such raw data is also stored by the container.

More than one container can be used at the same time to add different constraints; this allows for parallelism in constraint generation, since different containers can be used in different threads for example.

Semidefinite programs are represented by instances of SDPProblem. Problems contain an objective function, which can be set by providing an instance of SDPLinearCombination, and a list of constraint containers specifying the constraints it contains.

Problems can be used, for instance, to generate input for the solver. One way to do this is through the class SDPAWriter, which can be used to generate solver input in sparse SDPA format.

#include <sdp.h>
using namespace sdp;  // We use the double-precision version
int main()
{
  SDPVariable X(5), Y(7);
  Matrix<double> A(5), B(7);
  Matrix<double> J(5, 5, 1);  // An all-ones matrix, 5-by-5
  /* (We fill in matrices A and B somehow...) */
  /* Making a problem */
  SDPConstraintContainer cc;
  cc.add_constraint(A * X + B * Y <= 7);  // Adds a linear constraint
  SDPProblem problem;
  problem.set_objective(J * X);  // Sets the objective function
  problem << cc;  // Adds constraint container to the problem
  cout << SDPAWriter(problem);  // Outputs the problem in sparse SDPA format
  return 0;
}

In the above example, notice that not all variables have to be always present in a linear combination, such as when we specify the objective function.

Constraint containers copy the data from the instance of SDPLinearConstraint that is provided to add_constraint. So after this operation the constraint that was provided can be changed without altering the data in the container. This allows one to create a single instance of SDPLinearConstraint and, by changing the matrices used before adding the constraint each time to a container, create different constraints.

Finally, SDPProblem keeps shared references (see Linear combinations above) to constraint containers, so even if a container is destroyed before the problem to which it was added, its data will not be lost.

Using masks for efficiency

Sometimes, matrices provided as parts of constraints are sparse, and this can be exploited so as to make the inclusion of the constraint faster. This can be done via masks, objects that tell which positions of a matrix are nonzero.

An instance of MatrixMask contains a matrix bit-mask. Bits can be turned on and off to symbolize that the entries of an associated matrix are to be considered nonzeros or not. Masks can be associated to matrices when linear combinations are created, as shown below.

using namespace sdp;  // We use the double-precision version
SDPVariable X(5);
Matrix<double> J(5, 5, 1);  // A 5-by-5 all-ones matrix
MatrixMask mask(5);
mask.turn_on(1, 3);
SDPLinearConstraint C = SDPMaskedMatrix(J, mask) * X == 3;

Here, the all-ones matrix J is used in the constraint C, but only entry (1, 3) of this matrix is considered nonzero, because of the mask we specified. The next section contains another example of how masks can be used to speed-up constraint generation.

The Lovász theta number of the pentagon

Say we wish to compute the Lovász theta number of the pentagon[1], or 5-cycle. This means that we wish to solve the following semidefinite programming problem:

\[ \begin{array}{rl} {\rm maximize}&\langle J, X\rangle\\ &\langle I, X\rangle = 1,\\ &X_{12} = X_{23} = X_{34} = X_{45} = X_{51} = 0,\\ &\text{$X \in \mathbb{R}^{5 \times 5}$ is positive semidefinite}. \end{array} \]

Here, $J$ is the all-one matrix. The complete code using SDPSL to generate this problem in sparse SDPA format is below.

 /*
   theta.cpp
 */
 
 #include <iostream>
 #include <sdpsl.h>
 
 using namespace std;
 using namespace sdp;
 
 int main()
 {
   SDPVariable X(5);
   SDPConstraintContainer cc;
 
   /* We first generate one constraint for each edge of the 5-cycle,
      using masks for efficiency. 
   */
   MatrixMask mask(5);
   Matrix<double> J(5, 5, 1);  // A 5 x 5 all-ones matrix
   SDPLinearConstraint edge_constraint = SDPMaskedMatrix(J, mask) * X == 0;
 
   for (int i = 0; i < 5; i++) {
     int j = (i + 1) % 5;
 
     mask.turn_on(i, j);
     cc.add_constraint(edge_constraint);
     mask.turn_off(i, j);
   }
 
   /* The trace == 1 constraint */
   cc.add_constraint(Matrix<double>::identity(5) * X == 1);
 
   /* We add the constraints to the problem */
   SDPProblem problem;
 
   problem << cc;
 
   /* We set the objective function */
   problem.set_objective(J * X);
 
   /* And we output the SDPA file */
   cout << SDPAWriter(problem);
 
   return 0;
 }

The workings of the code should by now be clear; only a couple of observations are in place. The constraints relating to the edges of the pentagon have all the same structure. So notice how we create only one instance of SDPLinearConstraint, edge_constraint, for all of them, and use it repeatedly to add constraints to the container. Since the mask is being changed at each iteration of the loop of line 23, at each iteration a different constraint is added.

Finally, in line 32 we add the trace constraint. Here, Matrix<double>::identity(5) returns a 5-by-5 identity matrix. The matrix returned by the call is a temporary, but this is no matter, because of the shared references mechanism. In fact, one could even have the following code:

SDPVariable X(5);
SDPLinearConstraint C = Matrix<Number>::identity(5) * X == 1;
SDPConstraintContainer cc;
cc.add_constraint(C);              

Here, even though the temporary returned by Matrix<Number>::identity no longer exists when the constraint C is added to the container, there is no problem, because the data of the identity matrix is saved in the constraint via the shared references mechanism.

Fan's theorem

Fan's theorem provides a good opportunity to investigate a semidefinite programming problem in which two variable matrices are used. The theorem states that the optimal value of the following semidefinite programming problem equals the sum of the $k$ largest eigenvalues of the symmetric matrix $A$:

\[ \begin{array}{rl} {\rm maximize}&\langle A, X\rangle\\ &\langle I, X\rangle = k,\\ &Y = I - X,\\ &\text{$X$ and $Y$ are positive semidefinite.} \end{array} \]

This is a semidefinite programming problem with two variables, $X$ and $Y$. The following program uses SDPSL and Fan's theorem to compute the sum of the $k$ largest eigenvalues of a given symmetric matrix.

 /*
   fan.cpp
 */
 
 #include <iostream>
 #include <sdpsl.h>
 
 using namespace std;
 using namespace sdp;
 
 //
 // Prototypes
 //
 static void output_fan_sdp(const Matrix<double>& A, int k);
 
 /*
   main
   
   We apply function output_fan_sdp to a sample matrix.
 */
 int main()
 {
   Matrix<double> A({ { -2,  4,  3, -5, -1, -6 }, 
                      {  4,  8,  0,  0,  6,  3 }, 
                      {  3,  0,  0,  0, -5,  3 }, 
                      { -5,  0,  0, -2, -2, -4 }, 
                      { -1,  6, -5, -2, -2,  1 }, 
                      { -6,  3,  3, -4,  1, -2 } });
   
   output_fan_sdp(A, 3);  // sum of 3 largest eigenvalues
 
   return 0;
 }
 
 /*
   output_fan_sdp
 
   Outputs an SDP whose optimal value is the sum of the k largest
   eigenvalues of symmetric matrix A.
 */
 void output_fan_sdp(const Matrix<double>& A, int k)
 {
   int n = A.num_rows();
   SDPVariable X(n), Y(n);
   SDPConstraintContainer cc;
 
   /* First the trace constraint for X */
   cc.add_constraint(Matrix<double>::identity(n) * X == k);
 
   /* Then the correspondence between X and Y */
   Matrix<double> E(n, n, 1);
   MatrixMask E_mask(n);
 
   for (int i = 0; i < n; i++)
     for (int j = i; j < n; j++) {
       double rhs = i == j ? 1 : 0;
 
       E_mask.turn_on(i, j);
       cc.add_constraint(SDPMaskedMatrix(E, E_mask) * X
                         + SDPMaskedMatrix(E, E_mask) * Y == rhs);
       E_mask.turn_off(i, j);      
     }
 
   /* Then we create the problem and output it */
   SDPProblem problem;
 
   problem.set_objective(A * X);
   problem << cc;
 
   cout << SDPAWriter(problem);
 }

Environments

Finally, it should be mentioned that some parameters regulate the behavior of some of the classes of SDPSL. Such parameters are part of the environment, and their values can be changed by the user.

One such parameter is output_eps. It is used, for instance, by SDPProblem and SDPConstraintContainer, to decide which numbers are nonzero: only numbers whose absolute value is greater than output_eps are considered nonzero.

Environment parameters can be specified through the class SDPEnvironment. They can be globally specified, or also only locally. Classes such as SDPProblem and SDPConstraintContainer take instances of SDPEnvironment when constructed, so that custom environments can be supplied. If no instance is supplied, then the default environment is used.

See Also: SDPEnvironment

Next chapter: Polynomials and polynomial constraints

Table of Contents