2015-05-27

Guide to CVXOPT's quadprog() for row-major and/or MATLAB-speaking minds

I am using CVXOPT now. My mind is very row-major, because I mostly program in C and Python. And my mind is also MATLAB-style, because my friend send me a code in MATLAB. Now I need to convert his code to CVXOPT in Python. Here are some thing you wanna pay attention.

1. Matrixes in CVXOPT are column-major


By the code:
G=cvxopt.matrix([[1.,1],[-1,2],[2,1],[-1,0],[0,-1]])
I thought I created a 5-by-2 matrix (5 horizontal rows and 2 vertical columns). However, I got a 2-by-5 matrix.

In [72]: G
Out[72]: <2x5 matrix, tc='d'>

In [73]: print G
[ 1.00e+00 -1.00e+00  2.00e+00 -1.00e+00  0.00e+00]
[ 1.00e+00  2.00e+00  1.00e+00  0.00e+00 -1.00e+00]

2. The matrixes must be of type float.


Make sure your matrixes are of type float, even though your coefficients are integers.

Earlier I created a matrix:
In [53]: q=cvxopt.matrix([[-2],[-6]])

Later I got the error like
TypeError: 'q' must be a 'd' matrix with one column

This is because the matrix is not of type float:
In [54]: q
Out[54]: <1x2 matrix, tc='i'>

To fix, simply make at least one element float, e.g., "1." Then the type code will become d.
In [89]: q
Out[89]: <1x2 matrix, tc='d'>


3. Generating the matrixes needed for your optimization problem.


First of all, the cvxopt formulates quadratic programming so much simpler than MATLAB. Here is how cvxopt formulates the problem:

$$ \begin{array}{ll}
\mbox{minimize} & (1/2) x^TPx + q^T x \\
\mbox{subject to} & G x \preceq h, \\
& Ax = b.
\end{array}
$$

The boundaries for all variables can be defined as part of the matrix G - just use eyes. E.g., the lower boundary constraints

$$x_1 > 1, x_2>2, x_3> 3
$$
is just
$$
\begin{bmatrix}
-1 & 0 & 0 \\
0 & -1 & 0 \\
0 & 0 & -1
\end{bmatrix}
\begin{bmatrix}
x_1\\
x_2 \\
x_ 3
\end{bmatrix}
< \begin{bmatrix} -1\\ -2\\ -3 \end{bmatrix} $$ Note that in either MATLAB's or cvxopt's formulation, the linear terms are always on the side of "smaller than or equal to". Hence, the minus one above. MATLAB formulations have two additional column arrays for upper and lower boundaries: $$ \begin{array}{ll} \mbox{minimize} & (1/2) x^THx + f^T x \\ \mbox{subject to} & A x \preceq b, \\ & A_{eq}x = b_{eq}, \\ & lb \preceq x \preceq ub. \end{array} $$ Apparently, there are lazy people who wants to have a straightforward conversion from MATLAB's formulation to cvxopt's. However, I am not sure whether that translation is out-dated because cvxopt changed their APIs. Here is an easier one.

Let the MATLAB version be
x = quadprog(H,f,A,b,Aeq,beq,lb,ub, ...)
This is how you do it in cvxopt:
import numpy
import cvxopt
import cvxopt.solvers
n = H.shape[1]   # n is the number of variables 

P = H
q = f
G = numpy.vstack([A, -numpy.eye(n), numpy.eye(n)])
h = numpy.vstack([b, -lb, ub])
A = Aeq
b = beq

sol = cvxopt.solvers.qp(cvxopt.matrix(P), cvxopt.matrix(q), cvxopt.matrix(G), cvxopt.matrix(h), cvxopt.matrix(A), cvxopt.matrix(b))
x = sol['x']

If you don't have lb nor/and ub, no need to have it/them in G and h.

Put things together

So, let's see an example. The example in MATLAB Optimization Toolbox can be solved in cvxopt as follows (in iPython shell):

In [81]: G=cvxopt.matrix([[1.,1],[-1,2],[2,1],[-1,0],[0,-1]])

In [82]: P=cvxopt.matrix([[1,-1],[-1,2.]])

In [83]: q=cvxopt.matrix([[-2,-6.]])

In [84]: h=cvxopt.matrix([2,2,3,0,0.])

In [85]: Solv=cvxopt.solvers.qp(P.T, q, G.T, h)
     pcost       dcost       gap    pres   dres
 0: -1.0389e+01 -8.2778e+00  2e+01  9e-01  1e+00
 1: -7.2856e+00 -9.9286e+00  3e+00  1e-16  4e-16
 2: -8.1161e+00 -8.6188e+00  5e-01  8e-17  2e-16
 3: -8.2068e+00 -8.2359e+00  3e-02  6e-17  3e-15
 4: -8.2220e+00 -8.2224e+00  3e-04  5e-17  1e-15
 5: -8.2222e+00 -8.2222e+00  3e-06  8e-17  6e-16
Optimal solution found.

In [86]: print Solv['x']
[ 6.67e-01]
[ 1.33e+00]

1 comment:

Anonymous said...

Thank you for this guide! I have tried to debug my code by myself but in vain. This topic is really helpful.