6.9 Quadratic Optimization

MOSEK can solve quadratic and quadratically constrained problems, as long as they are convex. This class of problems can be formulated as follows:

(6.28)$$\begin{array}{r}\begin{array}{lrcccll}\text{minimize}& & & \frac{1}{2}{x}^T{Q}^{o}x+c^{T}x+c^f& & & \\ \text{subject to}& l_k^c& \le & \frac{1}{2}{x}^T{Q}^{k}x+\sum _{j=0}^{n-1}{a}_{k,j}{x}_j& \le & u_k^c,& k=0,\dots ,m-1,\\ & l_j^x& \le & x_j& \le & u_j^x,& j=0,\dots ,n-1.\end{array}\end{array}$$

Without loss of generality it is assumed that $Q^o$ and $Q^k$ are all symmetric because

$$x^{T}Qx=\frac{1}{2}{x}^T(Q+Q^T)x.$$

This implies that a non-symmetric $Q$ can be replaced by the symmetric matrix $\frac{1}{2}(Q+Q^T)$.

The problem is required to be convex. More precisely, the matrix $Q^o$ must be positive semi-definite and the $k$th constraint must be of the form

(6.29)$$l_k^c\le \frac{1}{2}{x}^T{Q}^{k}x+\sum _{j=0}^{n-1}{a}_{k,j}{x}_j$$

with a negative semi-definite $Q^k$ or of the form

$$\frac{1}{2}{x}^T{Q}^{k}x+\sum _{j=0}^{n-1}{a}_{k,j}{x}_j\le u_k^c.$$

with a positive semi-definite $Q^k$. This implies that quadratic equalities are not allowed. Specifying a non-convex problem will result in an error when the optimizer is called.

A matrix is positive semidefinite if all the eigenvalues of $Q$ are nonnegative. An alternative statement of the positive semidefinite requirement is

$$x^{T}Qx\ge 0,\mathrm{\forall }x.$$

If the convexity (i.e. semidefiniteness) conditions are not met MOSEK will not produce reliable results or work at all.

6.9.1 Example: Quadratic Objective

We look at a small problem with linear constraints and quadratic objective:

(6.30)$$\begin{array}{r}\begin{array}{lll}\text{minimize}& & x_1^2+0.1x_2^2+x_3^2-x_1{x}_3-x_2\\ \text{subject to}& 1\le & x_1+x_2+x_3\\ & 0\le & x.\end{array}\end{array}$$

The matrix formulation of (6.30) has:

$$\begin{array}{r}{Q}^o=\left[\begin{array}{ccc}2& 0& -1\\ 0& 0.2& 0\\ -1& 0& 2\end{array}\right],c=\left[\begin{array}{c}0\\ -1\\ 0\end{array}\right],A=\left[\begin{array}{ccc}1& 1& 1\end{array}\right],\end{array}$$

with the bounds:

$$\begin{array}{r}{l}^c=1,u^c=\mathrm{\infty },l^x=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]\text{ and }{u}^x=\left[\begin{array}{c}\mathrm{\infty }\\ \mathrm{\infty }\\ \mathrm{\infty }\end{array}\right]\end{array}$$

Please note the explicit $\frac{1}{2}$ in the objective function of (6.28) which implies that diagonal elements must be doubled in $Q$, i.e. $Q_{11}=2$ even though $1$ is the coefficient in front of $x_1^2$ in (6.30).

Using mosekopt

In Listing 6.17 we show how to use mosekopt to solve problem (6.30). This is the preferred way.

Listing 6.17 How to solve problem (6.30) using mosekopt. Click here to download.

function qo2()

clear prob;

% c vector.
prob.c = [0 -1 0]';

% Define the data.

% First the lower triangular part of q in the objective 
% is specified in a sparse format. The format is:
%
%   Q(prob.qosubi(t),prob.qosubj(t)) = prob.qoval(t), t=1,...,4

prob.qosubi = [ 1  3 2   3]';
prob.qosubj = [ 1  1 2   3]';
prob.qoval  = [ 2 -1 0.2 2]';

% a, the constraint matrix
subi  = ones(3,1);
subj  = 1:3;
valij = ones(3,1);

prob.a = sparse(subi,subj,valij);

% Lower bounds of constraints.
prob.blc  = [1.0]';

% Upper bounds of constraints.
prob.buc  = [inf]';

% Lower bounds of variables.
prob.blx  = sparse(3,1);

% Upper bounds of variables.
prob.bux = [];   % There are no bounds.

[r,res] = mosekopt('minimize',prob);

% Display return code.
fprintf('Return code: %d\n',r);

% Display primal solution for the constraints.
res.sol.itr.xc'

% Display primal solution for the variables.
res.sol.itr.xx'

This sequence of commands looks much like the one that was used to solve the linear optimization example using mosekopt except that the definition of the $Q$ matrix in prob. mosekopt requires that $Q$ is specified in a sparse format. Indeed the vectors qosubi, qosubj, and qoval are used to specify the coefficients of $Q$ in the objective using the principle

$$Q_{\mathtt{qosubi}\mathtt{(}\mathtt{t}\mathtt{)}\mathtt{,}\mathtt{qosubj}\mathtt{(}\mathtt{t}\mathtt{)}}=\mathtt{qoval}\mathtt{(}\mathtt{t}\mathtt{)},\text{for}t=1,\dots ,\text{length}(\mathtt{qosubi}).$$

An important observation is that due to $Q$ being symmetric, only the lower triangular part of $Q$ should be specified.

Using mskqpopt

In Listing 6.18 we show how to use mskqpopt to solve problem (6.30).

Listing 6.18 Function solving problem (6.30) using mskqpopt. Click here to download.

function qo1()

% Set up Q.
q     = [[2 0 -1];[0 0.2 0];[-1 0 2]];

% Set up the linear part of the problem.
c     = [0 -1 0]';
a     = ones(1,3);
blc   = [1.0];
buc   = [inf];
blx   = sparse(3,1);
bux   = [];

% Optimize the problem.
[res] = mskqpopt(q,c,a,blc,buc,blx,bux);

% Show the primal solution.
res.sol.itr.xx

It should be clear that the format for calling mskqpopt is very similar to calling msklpopt except that the $Q$ matrix is included as the first argument of the call. Similarly, the solution can be inspected by viewing the res.sol field.

6.9.2 Example: Quadratic constraints

In this section we show how to solve a problem with quadratic constraints. Please note that quadratic constraints are subject to the convexity requirement (6.29).

Consider the problem:

$$\begin{array}{r}\begin{array}{lcccl}\text{minimize}& & & x_1^2+0.1x_2^2+x_3^2-x_1{x}_3-x_2& \\ \text{subject to}& 1& \le & x_1+x_2+x_3-x_1^2-x_2^2-0.1x_3^2+0.2x_1{x}_3,& \\ & & & x\ge 0.& \end{array}\end{array}$$

This is equivalent to

(6.31)$$\begin{array}{r}\begin{array}{lccl}\text{minimize}& \frac{1}{2}{x}^T{Q}^{o}x+c^{T}x& & \\ \text{subject to}& \frac{1}{2}{x}^T{Q}^{0}x+Ax& \ge & b,\\ & x\ge 0,\end{array}\end{array}$$

where

$$\begin{array}{r}{Q}^o=\left[\begin{array}{ccc}2& 0& -1\\ 0& 0.2& 0\\ -1& 0& 2\end{array}\right],c={\left[\begin{array}{ccc}0& -1& 0\end{array}\right]}^T,A=\left[\begin{array}{ccc}1& 1& 1\end{array}\right],b=1.\end{array}$$

$$\begin{array}{r}{Q}^0=\left[\begin{array}{ccc}-2& 0& 0.2\\ 0& -2& 0\\ 0.2& 0& -0.2\end{array}\right].\end{array}$$

The linear parts and quadratic objective are set up the way described in the previous tutorial.

Setting up quadratic constraints

Listing 6.19 Script implementing problem (6.31). Click here to download.

function qcqo1()
clear prob;

% Specify the linear objective terms.
prob.c      = [0, -1, 0];

% Specify the quadratic terms of the constraints.
prob.qcsubk = [1     1    1   1  ]';
prob.qcsubi = [1     2    3   3  ]';
prob.qcsubj = [1     2    3   1  ]';
prob.qcval  = [-2.0 -2.0 -0.2 0.2]';

% Specify the quadratic terms of the objective.
prob.qosubi = [1     2    3    3  ]';
prob.qosubj = [1     2    3    1  ]';
prob.qoval  = [2.0   0.2  2.0 -1.0]';

% Specify the linear constraint matrix
prob.a      = [1 1 1];

% Specify the lower bounds
prob.blc    = [1];
prob.blx    = zeros(3,1);

[r,res]     = mosekopt('minimize',prob);

% Display the solution.
fprintf('\nx:');
fprintf(' %-.4e',res.sol.itr.xx');
fprintf('\n||x||: %-.4e',norm(res.sol.itr.xx));