Main Content

A second-order cone programming problem has the form

$$\underset{x}{\mathrm{min}}{f}^{T}x$$

subject to the constraints

$$\begin{array}{c}\Vert {A}_{\text{sc}}(i)\cdot x-{b}_{\text{sc}}(i)\Vert \le {d}_{\text{sc}}^{T}(i)\cdot x-\gamma (i)\\ A\cdot x\le b\\ \text{Aeq}\cdot x=\text{beq}\\ \text{lb}\le x\le \text{ub}.\end{array}$$

*f*, *x*, *b*,
*beq*, *lb*, and *ub* are
vectors, and *A* and *Aeq* are matrices. For each
*i*, the matrix
*A*_{sc}(*i*), the vectors
*b*_{sc}(*i*) and
*d*_{sc}(*i*), and the
scalar *γ*(*i*) are in a second-order cone
constraint that you create using `secondordercone`

.

In other words, the problem has a linear objective function and linear constraints, as well as a set of second-order cone constraints of the form $$\Vert {A}_{\text{sc}}(i)\cdot x-{b}_{\text{sc}}(i)\Vert \le {d}_{\text{sc}}^{T}(i)\cdot x-\gamma (i)$$.

`coneprog`

AlgorithmThe `coneprog`

solver uses the algorithm described in Andersen, Roos, and
Terlaky [1]. This method is
an interior-point algorithm similar to the Interior-Point linprog Algorithm.

The algorithm starts by placing the problem in *standard
form*. The algorithm adds nonnegative slack variables so that the
problem has the form

$$\underset{x}{\mathrm{min}}{f}^{T}x$$

subject to the constraints

$$\begin{array}{l}A\cdot x=b\\ x\in K.\end{array}$$

The solver expands the sizes of the linear coefficient vector
*f* and linear constraint matrix *A* to
account for the slack variables.

The region *K* is the cross product of *Lorentz
cones*
Equation 1 and the
nonnegative orthant. To convert each convex cone

$$\Vert {A}_{\text{sc}}(i)\cdot x-{b}_{\text{sc}}(i)\Vert \le {d}_{\text{sc}}^{T}(i)\cdot x-\gamma (i)$$

to a Lorentz cone Equation 1, create a
column vector of variables *t*_{1},
*t*_{2}, …,
*t*_{n+1}:

$$\begin{array}{c}{t}_{1}={d}^{T}x-\gamma \\ {t}_{2:(n+1)}={A}_{\text{sc}}x-{b}_{\text{sc}}.\end{array}$$

Here, the number of variables *n* for each cone
*i* is the number of rows in
*A*_{sc}(*i*). By its
definition, the variable vector *t* satisfies the
inequality

$$\Vert {t}_{2:(n+1)}\Vert \le {t}_{1}.$$ | (1) |

Equation 1 is the
definition of a Lorentz cone in (*n*+1) variables. The
variables *t* appear in the problem in place of the variables
*x* in the convex region *K*.

Internally, the algorithm also uses a *rotated Lorentz
cone* in the reformulation of cone constraints, but this topic
does not address that case. For details, see Andersen, Roos, and Terlaky [1].

When adding slack variables, the algorithm negates variables, as needed, and adds appropriate constants so that:

Variables with only one bound have a lower bound of zero.

Variables with two bounds have a lower bound of zero and, using a slack variable, have no upper bound.

Variables without bounds are placed in a Lorentz cone with a slack variable as the constrained variable. This slack variable is not part of any other expression, objective or constraint.

The dual cone is

$${K}_{*}=\left\{s\text{\hspace{0.17em}}:\text{\hspace{0.17em}}{s}^{T}x\ge 0\text{\hspace{0.17em}}\forall x\in K\right\}.$$

The dual problem is

$$\underset{y}{\mathrm{max}}{b}^{T}y$$

such that

$${A}^{T}y+s=f$$

for some

$$s\in {K}_{*}.$$

A dual optimal solution is a point (*y*,*s*)
that satisfies the dual constraints and maximizes the dual objective.

To handle potentially infeasible or unbounded problems, the algorithm adds two
more variables *τ* and *κ* and formulates the
problem as homogeneous (equal to zero) and self-dual.

$$\begin{array}{c}Ax-b\tau =0\\ {A}^{T}y+s-f\tau =0\\ -{f}^{T}x+{b}^{T}y-\kappa =0\end{array}$$ | (2) |

along with the constraints

$$(x;\tau )\in \overline{K},\text{\hspace{1em}}(s;\kappa )\in {\overline{K}}_{*}\text{\hspace{0.17em}}.$$ | (3) |

Here, $$\overline{K}$$ is the cone *K* adjoined with the nonnegative
real line, which is the space for (*x*;*τ*).
Similarly $${\overline{K}}_{*}$$ is the cone $${K}_{*}$$ adjoined with the nonnegative real line, which is the space
for (*s*;*κ*). In this formulation, the
following lemma shows that *τ* is the scaling for feasible
solutions, and *κ* is the indicator of an infeasible
problem.

**Lemma** ([1] Lemma
2.1)

Let (*x*, *τ*, *y*,
*s*, *κ*) be a feasible solution of Equation 2 along with
the constraints in Equation 3.

*x*+^{T}s*τκ*= 0.If

*τ*> 0, then (*x*,*y*,*s*)/*τ*is a primal-dual optimal solution of the standard form second-order cone problem.If

*κ*> 0, then at least one of these strict inequalities holds:*b*> 0^{T}y*f*< 0.^{T}xIf the first inequality holds, then the standard form, primal second-order cone problem is infeasible. If the second inequality holds, then the standard form, dual second-order cone problem is infeasible.

In summary, for feasible problems, the variable *τ* scales
the solution between the original standard form problem and the homogeneous
self-dual problem. For infeasible problems, the final iterate
(*x*, *y*, *s*,
*τ*, *κ*) provides a certificate of
infeasibility for the original standard form problem.

The start point for the iterations is the feasible point:

*x*= 1 for each nonnegative variable, 1 for the first variable in each Lorentz cone, and 0 otherwise.*y*= 0.*s*= (1,0,…,0) for each cone, 1 for each nonnegative variable.*τ*= 1.*κ*= 1.

The algorithm attempts to follow the *central path*,
which is the parameterized solution to the following equations for
*γ* decreasing from 1 toward 0.

$$\begin{array}{c}Ax-b\tau =\gamma \left(A{x}_{0}-b{\tau}_{0}\right)\\ {A}^{T}y+s-c\tau =\gamma \left({A}^{T}{y}_{0}+{s}_{0}-f{\tau}_{0}\right)\\ -{f}^{T}x+{b}^{T}y-\kappa =\gamma \left(-{f}^{T}{x}_{0}+{b}^{T}{y}_{0}-{\kappa}_{0}\right)\\ XSe=\gamma {\mu}_{0}e\\ \tau \kappa =\gamma {\mu}_{0}.\end{array}$$ | (4) |

Each variable with a 0 subscript indicates the start point of the variable.

The variables

*X*and*S*are*arrow head*matrices formed from the*x*and*s*vectors, respectively. For a vector*x*= [*x*_{1},*x*_{2},…,*x*], the arrow head matrix_{n}*X*has the definition$$X=\text{mat}(x)=\left[\begin{array}{cc}{x}_{1}& {x}_{2:n}^{T}\\ {x}_{2:n}& {x}_{1}I\end{array}\right].$$

By its definition,

*X*is symmetric.The variable

*e*is the vector with a 1 in each cone coordinate corresponding to the*x*_{1}Lorentz cone coordinate.The variable

*μ*_{0}has the definition$${\mu}_{0}=\frac{{x}_{0}^{T}{s}_{0}+{\tau}_{0}{\kappa}_{0}}{k+1},$$

where

*k*is the number of nonzero elements in*x*_{0}.

The central path begins at the start point and ends at an optimal solution to the homogeneous self-dual problem.

Andersen, Roos, and Terlaky [1] show in
Lemma 3.1 that the complementarity condition
*x ^{T}s* = 0, where

$${X}_{i}{S}_{i}{e}_{i}={S}_{i}{X}_{i}{e}_{i}=0$$

for every cone *i*. Here
*X _{i}* =
mat(

To obtain points near the central path as the parameter *γ*
decreases from 1 toward 0, the algorithm uses Newton's method. The variables to
find are labeled (*x*, *τ*,
*y*, *s*, *κ*). Let
*d _{x}* represent the search
direction for the

$$\begin{array}{c}A{d}_{x}-b{d}_{\tau}=(\gamma -1)\left(A{x}_{0}-b{\tau}_{0}\right)\\ {A}^{T}{d}_{y}+{d}_{s}-f{d}_{\tau}=(\gamma -1)\left({A}^{T}{y}_{0}+{s}_{0}-f{\tau}_{0}\right)\\ -{f}^{T}{d}_{x}+{b}^{T}{d}_{y}-{d}_{\kappa}=(\gamma -1)\left(-{f}^{T}{x}_{0}+{b}^{T}{y}_{0}-\kappa \right)\\ {X}_{0}{d}_{s}+{S}_{0}{d}_{x}=-{X}_{0}{S}_{0}e+\gamma {\mu}_{0}e\\ {\tau}_{0}{d}_{\kappa}+{\kappa}_{0}{d}_{t}au=-{\tau}_{0}{\kappa}_{0}+\gamma {\mu}_{0}.\end{array}$$

The algorithm obtains its next point by taking a step in the
*d* direction.

$$\left[\begin{array}{c}{x}_{1}\\ {\tau}_{1}\\ {y}_{1}\\ \begin{array}{l}{s}_{1}\\ {\kappa}_{1}\end{array}\end{array}\right]=\left[\begin{array}{c}{x}_{0}\\ {\tau}_{0}\\ {y}_{0}\\ \begin{array}{l}{s}_{0}\\ {\kappa}_{0}\end{array}\end{array}\right]+\alpha \left[\begin{array}{c}{d}_{x}\\ {d}_{\tau}\\ {d}_{y}\\ \begin{array}{l}{d}_{s}\\ {d}_{\kappa}\end{array}\end{array}\right]$$

for some step $$\alpha \in [0,1]$$.

For both numerical stability and accelerated convergence, the algorithm scales the step according to a suggestion in Nesterov and Todd [5]. Also, the algorithm corrects the step according to a variant of Mehrotra's predictor-corrector [4]. (For further details, see Andersen, Roos, and Terlaky [1].)

At each iteration *k*, the algorithm computes three relative
convergence measures:

Primal infeasibility

$${\text{Infeas}}_{\text{Primal}}^{k}=\frac{\Vert A{x}_{k}-b{\tau}_{k}\Vert}{\mathrm{max}\left(1,\Vert A{x}_{0}-b{\tau}_{0}\Vert \right)}.$$

Dual infeasibility

$${\text{Infeas}}_{\text{Dual}}^{k}=\frac{\Vert {A}^{T}{y}_{k}+{s}_{k}-f{\tau}_{k}\Vert}{\mathrm{max}\left(1,\Vert {A}^{T}{y}_{0}+{s}_{0}-f{\tau}_{0}\Vert \right)}.$$

Gap infeasibility

$${\text{Infeas}}_{\text{Gap}}^{k}=\frac{\left|-{f}^{T}{x}_{k}+{b}^{T}{y}_{k}-{\kappa}_{k}\right|}{\mathrm{max}\left(1,\left|-{f}^{T}{x}_{0}+{b}^{T}{y}_{0}-{\kappa}_{0}\right|\right)}.$$

You can view these three statistics at the command line by specifying iterative display.

options = optimoptions('coneprog','Display','iter');

All three should approach zero when the problem is feasible and the solver
converges. For a feasible problem, the variable
*κ _{k}* approaches zero, and the
variable

One stopping condition is somewhat related to the gap infeasibility. The stopping condition is when the following optimality measure decreases below the optimality tolerance.

$${\text{Optimality}}^{k}=\frac{\left|{f}^{T}{x}_{k}-{b}^{T}{y}_{k}\right|}{{\tau}_{k}+\left|{b}^{T}{y}_{k}\right|}=\frac{\left|{f}^{T}{x}_{k}/{\tau}_{k}-{b}^{T}{y}_{k}/{\tau}_{k}\right|}{1+\left|{b}^{T}{y}_{k}/{\tau}_{k}\right|}.$$

This statistic measures the precision of the objective value.

The solver also stops and declares the problem to be infeasible under the
following conditions. The three relative infeasibility measures are less than
*c* = `ConstraintTolerance`

, and

$${\tau}_{k}\le c\text{\hspace{0.17em}}\mathrm{max}\left(1,{\kappa}_{k}\right).$$

If *b ^{T}y_{k}*
> 0, then the solver declares that the primal problem is
infeasible. If

The algorithm also stops when

$${\mu}_{k}\le c{\mu}_{0}$$

and

$${\tau}_{k}\le c\text{\hspace{0.17em}}\mathrm{max}\left(1,{\kappa}_{k}\right).$$

In this case, `coneprog`

reports that the problem is
numerically unstable (exit flag `-10`

).

The remaining stopping condition occurs when at least one infeasibility
measure is greater than `ConstraintTolerance`

and the computed
step size is too small. In this case, `coneprog`

reports that
the search direction became too small and no further progress could be made
(exit flag `-7`

).

[1] Andersen, E. D., C. Roos,
and T. Terlaky. *On implementing a primal-dual interior-point method for
conic quadratic optimization.* Math. Program., Ser. B **95**, pp. 249–277 (2003). https://doi.org/10.1007/s10107-002-0349-3

[2] Ben-Tal, Aharon, and
Arkadi Nemirovski. *Convex Optimization in Engineering: Modeling,
Analysis, Algorithms.* (1998). Available at https://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.455.2733&rep=rep1&type=pdf.

[3] Luo, Zhi-Quan, Jos F.
Sturm, and Shuzhong Zhang. *Duality and Self-Duality for Conic Convex
Programming.* (1996). Available at https://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.48.6432

[4] Mehrotra, Sanjay. “On the
Implementation of a Primal-Dual Interior Point Method.” *SIAM Journal on
Optimization 2*, no. 4 (November 1992): 575–601. https://doi.org/10.1137/0802028.

[5] Nesterov, Yu. E., and M.
J. Todd. “Self-Scaled Barriers and Interior-Point Methods for Convex Programming.”
*Mathematics of Operations Research 22*, no. 1 (February
1997): 1–42. https://doi.org/10.1287/moor.22.1.1.