PMI class examples

In the following, several examples of PMI classes with their definitions in Matlab are given.

Example 1

This PMI class has

problem variables

and

problem parameters

. There are

residual parameters. There are

constraints with

. Since all of the constraints have dimension one, this is a POP class.

Matlab code snippet with the definition of the above PMI problem class:

syms x1 x2 a1 a2 a3 real;

problem.vars = [x1, x2];
problem.ppars = [a1, a2, a3];
problem.obj = -(x1 - 1)^2 - (x1 - x2)^2 - (x2 - a1)^2;
problem.cons = { 2 - a2^3 - (x1 - 1)^2 >= 0, ...
               1 - (x1 - x2)^2 >= 0, 1 - (x2 - a3)^2 >= 0 };

Example 2

This PMI class has

problem variables

and

problem parameters

. We have no residual parameters, hence

. Finally, we have

constraints with

being the size of the single polynomial matrix constraint.

Matlab code snippet with the definition of the above PMI problem class:

syms x1 x2 a1 a2 real;

problem.vars = [x1, x2];
problem.ppars = [a1, a2];
problem.obj = -x1^2 - x2^2;
problem.cons = {[1 - a1*x1*x2, x1; x1, a2 - x1^2 - x2^2] >= 0};

Example 3

In this example, we have

problem variables

and

problem parameters

. We have

constraints with

and

. Note that since the problem variables

appear as degree one monomials only, this problem is in fact an LMI class that does not have to be solved using an additional LMI relaxation. GpoSolver will recognize LMI classes and will solve them directly.

Matlab code snippet with the definition of the above PMI problem class:

syms x1 x2 x3 a1 a2 a3 a4 a5 real;

problem.vars = [x1, x2];
problem.ppars = [a1, a2, a3, a4, a5];
problem.obj = x2;
problem.cons = { a1 <= x1, a2 >= x1, a3 <= x3, a4 >= x3, ...
  [1, x1, a5^2*x2; x1, 1, x3; a5^2*x2, x3, 1] >= 0 };

Example 4

Here, we have

problem variables

problem parameters

, and, for the first time, also

residual parameters

. We have

constraint.

Matlab code snippet with the definition of the above PMI problem class:

syms x1 x2 a1 a2 a3 a4 b1 b2 b3 b4 real;

problem.vars = [x1 x2];
problem.ppars = [a1 a2 a3 a4];
problem.rpars = [b1 b2 b3 b4];
problem.obj = (a1^2*x1^2 - b1*x2 - b2)^2 + (a2^2*x1^2 - b3*x2 - b4)^2;
problem.cons = {a3^2 * x1^2 + a4^2 * x2^2 == 1};

C++ code generation

After a PMI class is defined via problem data structure, the C++ code with the LMI relaxation definition can be readily generated using function gpogenerator. The function will generate a header file and a source file with the same names to a directory of the user’s choosing. In the next, we will denote the directory root of the GpoSolver distribution as GPO_ROOT.

Let us execute function gpogenerator to produce the second order LMI relaxation of the PMI class from Example 4↑:

>> cd GPO_ROOT/matlab;
>> params.relax_order = 2;
>> params.classname = ’PopProblem’;
>> params.filename = ’pop_problem’;
>> gpogenerator(problem, params);
Optimization variables ... [x1, x2] 
Computing problem polynomial degree ... 4
Computing minimum relaxation degree ... 2
Computing moment matrix ... done
Computing localizing matrices ... done
Identifying monomials ... 15
Decomposing LMI constraints ... done
Decomposing objective function ... done
Generating pop_problem.h... done
Generating pop_problem.cpp... done
>>

We can see that the second order relaxation leads to a

moment matrix. Two files were generated into the current working directory GPO_ROOT/matlab: pop_problem.h and pop_problem.cpp.

Solving PMI instances in C++

The problem solving part of GpoSolver consists of a C++ template library based in GPO_ROOT/include directory. GpoSolver can work with the following three SDP solvers: CSDP, SDPA, and MOSEK. Both CSDP and SDPA are available under open source licenses. MOSEK is a commercial product, however, an academic license can be obtained free of charge. All of these solvers can be used and linked at the same time and a user can choose between them by simply instantiating objects of different names.

To use the library, one must first include at least one of the following header files:

gposolver/gposolver_csdp.h
gposolver/gposolver_sdpa.h
gposolver/gposolver_mosek.h

depending on the SDP solvers available. Also, the header file generated by gpogenerator must be included. Note that every class of the GpoSolver library lives in the namespace GpoSolver. Next, an object of GpoSolverCsdp, GpoSolverSdpa, or GpoSolverMOSEK class is instantiated. All of these classes are templated and take gpogenerator generated class as the only template parameter. Finally, to solve a PMI instance, one just needs to call a virtual method solve declared in the common ancestor class GpoSolverBase:

GpoSolverStatus GpoSolverBase::solve(
  const int num_res, 
  const double *rvals, 
  const double *pvals, 
  GpoSolverSolutions &sols
);

Let’s continue with the results of the previous section. There, two files pop_problem.cpp and pop_problem.h were generated by gpogenerator into the GPO_ROOT/matlab directory; the C++ class name was chosen as ‘PopProblem’. The following shows an example of a PMI instance solving code based on the CSDP solver. A full working version of this solver can be found in GPO_ROOT/examples/pop_solver.cpp.

First, let’s make the necessary includes.

#include <gposolver/gposolver_csdp.h>
#include "pop_problem.h"

Next, we will be using the GpoSolver namespace:

using namespace GpoSolver;

Let’s say we want to solve a PMI instance with problem parameters

and

residual parameter blocks

and

. To do that, we can define the following two arrays:

double pvals[] = {1, 2, 3, 4};
double rvals[][4] = {{1, 2, 2, 3}, {2, 1, 1, 2}};

Now, we define the main function, construct the solver object using the PopProblem class, and call solve method:

int main(int argc, const char* argv[]) { 
  GpoSolverCsdp<PopProblem> gposolver;
  GpoSolverStatus status;
  GpoSolverSolutions sols;
  status = gposolver.solve(2, (double*) rvals, (double*) pvals, 
                           sols);

To interpret the solution, we first check the status variable:

  if (sols.size() == 0)
    return -1;
  else if (status == SUCCESS_GLOPT)
    cout << "Global optimality certified numerically" << endl;
  else if (status == SUCCESS)
    cout << "Alas, global optimality could not be certified" 
         << endl;

Now, we just print out the solutions returned in the sols list:

  int c = 1;
  for (GpoSolverSolutions::iterator it = sols.begin(); 
         it != sols.end(); it++)
    {
      std::vector<double> &sol = *(it);
        cout << "Solution " << c++ << ": ";
        for (int i = 0; i < sol.size(); i++)
          cout << sol[i] << " ";
        cout << endl;
    }

Finally, we close the main function:

  return 0;
}

After compiling and linking, the resulting executable should give the following output:

GpoSolver: 2 solution(s) extracted
Global optimality certified numerically
Solution 1: -0.271539 -0.145 
Solution 2: 0.271539 -0.145

Solving PMI instances in Matlab

Even though it may seem to go against the main idea behind GpoSolver, it is possible to solve PMI instances directly from Matlab using function gposolve. This feature becomes useful when testing correctness of a PMI class definition. However, this approach cannot be recommended as a GloptiPoly substitution, since gposolve is much slower. This is because gposolve manipulates parametrized problems using Symbolic Math Toolbox whereas GloptiPoly can work with concrete problem instances from the beginning. The gposolve function requires SeDuMi as its underlying SDP solver. The gposolve function again resides in the directory GPO_ROOT/matlab.

Before solving a PMI instance, we need to construct an LMI relaxation

. This is done using function gporelax:

 >> relax = gporelax(problem, order);

Given a PMI class definition problem and a relaxation order , gporelax will return data structure relax - the LMI relaxation

. Relaxation relax is then an input parameter of function gposolve:

 >> [status, solutions] = gposolve(relax, rvals, pvals, params);

Let us once again use the PMI class from Example 4↑ and let us assume that the variable problem already contains the PMI class definition. As in the previous section, we would like to solve the problem instance with problem parameters

and

residual parameter blocks

and

using the second order relaxation:

>> relax = gporelax(problem, 2);
Optimization variables ... [x1, x2] 
Computing problem polynomial degree ... 4
Computing minimum relaxation degree ... 2
Computing moment matrix ... done
Computing localizing matrices ... done
Identifying monomials ... 15
Decomposing LMI constraints ... done
Decomposing objective function ... done
>> rvals = [1, 2, 2, 3; 2, 1, 1, 2];
>> pvals = [1, 2, 3, 4];
>> [status, solutions] = gposolve(relax, rvals, pvals);
>> disp(solutions)
   -0.2715   -0.1450
    0.2715   -0.1450
>>

Notice that we did not use gposolve parameter params. As expected, we obtained two solutions

and

, as in the previous section.

PMI class examples

C++ code generation

Solving PMI instances in C++

Solving PMI instances in Matlab

Navigate