The OPTLP Procedure

The Interior Point Algorithm

The interior point solver in PROC OPTLP implements an infeasible primal-dual predictor-corrector interior point algorithm. To illustrate the algorithm and the concepts of duality and dual infeasibility, consider the following LP formulation (the primal):

$\begin{array}{rl} \mbox{min} & \mathbf{c}^\mr {T} \mathbf{x} \\ \mbox{subject to} & \mathbf{A} \mathbf{x} \ge \mathbf{b} \\ & \mathbf{x} \ge \mathbf{0} \end{array}$

The corresponding dual formulation is as follows:

$\begin{array}{rccccl} \mbox{max} & \mathbf{b}^\mr {T} \mathbf{y} & & & & \\ \mbox{subject to} & \mathbf{A}^\mr {T} \mathbf{y} & + & \mathbf{w} & = & \mathbf{c} \\ & & & \mathbf{y} & \ge & \mathbf{0} \\ & & & \mathbf{w} & \ge & \mathbf{0} \end{array}$

where $\mathbf{y} \in \mathbb {R}^ m$ refers to the vector of dual variables and $\mathbf{w} \in \mathbb {R}^ n$ refers to the vector of dual slack variables.

The dual formulation makes an important contribution to the certificate of optimality for the primal formulation. The primal and dual constraints combined with complementarity conditions define the first-order optimality conditions, also known as KKT (Karush-Kuhn-Tucker) conditions, which can be stated as follows:

$\begin{array}{rclr} \mathbf{A}\mathbf{x} - \mathbf{s} & = & \mathbf{b} & \mr {(primal~ feasibility)} \\ \mathbf{A}^\mr {T} \mathbf{y} + \mathbf{w} & = & \mathbf{c} & \mr {(dual~ feasibility)}\\ \mathbf{W} \mathbf{X} \mathbf{e} & = & \mathbf{0} & \mr {(complementarity)}\\ \mathbf{S} \mathbf{Y} \mathbf{e} & = & \mathbf{0} & \mr {(complementarity)}\\ \mathbf{x},\, \mathbf{y},\, \mathbf{w}, \, \mathbf{s} & \ge & \mathbf{0} & \\ \end{array}$

where $\mathbf{e} \equiv (1, \ldots , 1)^\mr {T}$ of appropriate dimension and $\mathbf{s} \in \mathbb {R}^ m$ is the vector of primal slack variables.

Note: Slack variables (the $s$ vector) are automatically introduced by the solver when necessary; it is therefore recommended that you not introduce any slack variables explicitly. This enables the solver to handle slack variables much more efficiently.

The letters $\mathbf{X},$ $\mathbf{Y},$ $\mathbf{W},$ and $\mathbf{S}$ denote matrices with corresponding $x,$ $y,$ $w,$ and $s$ on the main diagonal and zero elsewhere, as in the following example:

$\mathbf{X} \equiv \left[ \begin{array}{cccc} x_1 & 0 & \cdots & 0 \\ 0 & x_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & x_ n \end{array} \right]$

If $(\mathbf{x}^\ast , \mathbf{y}^\ast , \mathbf{w}^\ast , \mathbf{s}^\ast )$ is a solution of the previously defined system of equations that represent the KKT conditions, then $\mathbf{x}^\ast$ is also an optimal solution to the original LP model.

At each iteration the interior point algorithm solves a large, sparse system of linear equations,

$\left[ \begin{array}{cc} \mathbf{Y}^{-1}\mathbf{S} & \mathbf{A} \\ \mathbf{A}^\mr {T} & - \mathbf{X}^{-1}\mathbf{W} \end{array}\right] \left[ \begin{array}{c} \Delta \mathbf{y} \\ \Delta \mathbf{x} \end{array} \right] = \left[ \begin{array}{c} \Xi \\ \Theta \\ \end{array} \right]$

where $\Delta \mathbf{x}$ and $\Delta \mathbf{y}$ denote the vector of search directions in the primal and dual spaces, respectively, and $\Theta$ and $\Xi$ constitute the vector of the right-hand sides.

The preceding system is known as the reduced KKT system. PROC OPTLP uses a preconditioned quasi-minimum residual algorithm to solve this system of equations efficiently.

An important feature of the interior point solver is that it takes full advantage of the sparsity in the constraint matrix, thereby enabling it to efficiently solve large-scale linear programs.

The interior point algorithm works simultaneously in the primal and dual spaces. It attains optimality when both primal and dual feasibility are achieved and when complementarity conditions hold. Therefore, it is of interest to observe the following four measures where $\| v\| _2$ is the Euclidean norm of the vector $v$ :

relative primal infeasibility measure $\alpha$ :

$\alpha = \frac{\| \mathbf{Ax} - \mathbf{b} - \mathbf{s}\| _2}{\| \mathbf{b}\| _2 + 1}$
relative dual infeasibility measure $\beta$ :

$\beta = \frac{\| \mathbf{c} - \mathbf{A}^\textrm {T} \mathbf{y} - \mathbf{w} \| _2}{\| \mathbf{c}\| _2 + 1}$
relative duality gap $\delta$ :

$\delta = \frac{|\mathbf{c}^\textrm {T}\mathbf{x} - \mathbf{b}^\textrm {T}\mathbf{y}|}{|\mathbf{c}^\textrm {T}\mathbf{x}| + 1}$
absolute complementarity $\gamma$ :

$\gamma =\sum _{i=1}^ n x_ iw_ i + \sum _{i=1}^ m y_ i s_ i$

These measures are displayed in the iteration log.