The VMCON Algorithm

VMCON aims to solve a nonlinearly constrained problem (general nonlinear programming problem), specifically, minimising an objective function subject to several nonlinear constraints.

The Problem

Also called the general non-linear programming problem, we want to minimise an objective function (figure of merit), $f(\overrightarrow{x})$, subject to some non-linear constraints:

$$\begin{array}{r}\begin{array}{c}{c}_i(\overrightarrow{x})=0,i=1,...,k\\ c_i(\overrightarrow{x})\ge 0,i=k+1,...,m\end{array}\end{array}$$

where the objective and constraint functions are all $n$-dimensional.

Several facts are worth noting about problem formulation:

1. To maximise some function $g(\overrightarrow{x})$ we would minimise the objective $f(\overrightarrow{x})=-g(\overrightarrow{x})$.

2. To constrain the solution such that $h(\overrightarrow{x})=a$ we would apply the constraint $h(\overrightarrow{x})-a=0$.

3. To constrain the solution such that $h(\overrightarrow{x})\ge a$ we would apply the constraint $h(\overrightarrow{x})-a\ge 0$.

4. To constrain the solution such that $h(\overrightarrow{x})\le 0$ we would apply the constraint $-h(\overrightarrow{x})\ge 0$.

5. To constrain the solution such that $h(\overrightarrow{x})\le a$ we would apply the constraint $a-h(\overrightarrow{x})\ge 0$.

The Lagrangian

VMCON is an augmented Lagrangian solver which means it uses the Lagrange multipliers, $\overrightarrow{\lambda }$, of some solution to characterise the quality of the solution.

Of specific note is the Lagrangian function:

$$L(\overrightarrow{x},\overrightarrow{\lambda })=f(\overrightarrow{x})-\sum _{i=1}^m{\overrightarrow{\lambda }}_i{c}_i(\overrightarrow{x})$$

and the derivative of the Lagrangian function, with respect to $\overrightarrow{x}$:

(1)$${\mathrm{\nabla }}_{X}L(\overrightarrow{x},\overrightarrow{\lambda })=\mathrm{\nabla }f(\overrightarrow{x})-\sum _{i=1}^m{\overrightarrow{\lambda }}_i\mathrm{\nabla }{c}_i(\overrightarrow{x})$$

Initialisation of VMCON

VMCON is initialised with:

• The objective function to minimise, $f(\overrightarrow{x})$, as described above.

• The constraints $c_i(\overrightarrow{x}),i=1,...,m$, as described above.

• An initial sample point ${\overrightarrow{x}}_0$.

• $\mathbf{B}$: the initial Hessian approximation matrix, usually the identity matrix.

• $ϵ$: the “user-supplied error tolerance”.

It should be noted that $\mathbf{B}$ will need to be of dimension $d\times d$ where $d=\max (n,m)$.

We also set the iteration number to 1, $j=1$.

The Quadratic Programming Problem

The Quadratic Programming Probelm (QPP) will also be known as the Quadratic Sub-Problem (QSP) because it forms only a part of the VMCON algorithm–with the other half being the Augmented Lagrangian.

The QPP provides the search direction ${\delta }_j$ which is a vector upon which ${\overrightarrow{x}}_j$ will lay. Solving the QPP also provides the Lagrange multipliers, ${\lambda }_j$.

The quadratic program to be minimised on iteration $j$ is:

$$Q(\delta )=f({\overrightarrow{x}}_{j-1})+{\delta }^T\mathrm{\nabla }f({\overrightarrow{x}}_{j-1})+\frac{1}{2}{\delta }^{T}B\delta$$

subject to

$$\begin{array}{r}\begin{array}{c}\mathrm{\nabla }{c}_i({\overrightarrow{x}}_{j-1}{)}^T\delta +c_i({\overrightarrow{x}}_{j-1})=0,i=1,...,k\\ \mathrm{\nabla }{c}_i({\overrightarrow{x}}_{j-1}{)}^T\delta +c_i({\overrightarrow{x}}_{j-1})\ge 0,i=k+1,...,m\end{array}\end{array}$$

The Convergence Test

The convergence test is performed on the $j$’th iteration after the QSP. The convergence test is the sum of two terms:

• The predicted change in magnitude of the objective function.

• The complimentary error; where the complimentary error being 0 would mean that a specific constraint is at equality or the Lagrange multipliers are 0.

This is encapsulated in the equation:

$$|\mathrm{\nabla }f({\overrightarrow{x}}_{j-1}{)}^T\cdot {\delta }_j|+\sum _{i=1}^m|{\lambda }_{j,i}{c}_i({\overrightarrow{x}}_{j-1})|<ϵ$$

The Line Search

The line search helps to mitigate poor initial conditions. It does this by searching in a line along the ‘search direction’ $\delta$ such that:

$${\overrightarrow{x}}_j={\overrightarrow{x}}_{j-1}+{\alpha }_j{\overrightarrow{\delta }}_j$$

$\alpha$ is found via the minimisation of:

$$\varphi (\alpha )=f({\overrightarrow{x}}_j)+\sum _{i=1}^k{\overrightarrow{\mu }}_{j,i}|c_i({\overrightarrow{x}}_j)|+\sum _{i=k+1}^m{\overrightarrow{\mu }}_{j,i}|min(c_i(0,{\overrightarrow{x}}_j))|$$

On the $j$ th iteration, ${\overrightarrow{\mu }}_{j,i}$ is a 1D vector which contains $i=1,...,m$ elements. On the first iteration:

$${\overrightarrow{\mu }}_1=|{\overrightarrow{\lambda }}_1|$$

On subsequent iterations:

$${\overrightarrow{\mu }}_j=max[|{\overrightarrow{\lambda }}_0|,\frac{1}{2}({\overrightarrow{\mu }}_{j-1}+|{\overrightarrow{\lambda }}_j|)]$$

The line search iterates for a maximum of 10 steps and exits if the chosen value of $\alpha$ satisfies either the Armijo condition:

$$\varphi (\alpha )\le \varphi (0)+0.1\alpha (\varphi (1)-\varphi (0))$$

or the so-called Kovari condition, which was an ad-hoc break condition in the PROCESS implementation of VMCON, therefore does not appear in the paper:

$$\varphi (\alpha )>\varphi (0)$$

Once the line search exits, we have found our optimal value and ${\alpha }_j=\alpha$.

On each iteration of the line search, we revise $\alpha$ using a quadratic approximation:

$$\alpha =min(0.1\alpha ,\frac{-{\alpha }^2}{\varphi (\alpha )-\varphi (0)-\alpha (\varphi (1)-\varphi (0))})$$

The Broyden-Fletcher-Goldfarb-Shanno (BFGS) Quasi-Newton Update

The final stage of an iteration of the VMCON optimiser is to update the Hessian approximation via a BFGS update.

For an unconstrained problem, we use the following differences to update $\mathbf{B}$:

$$\overrightarrow{\xi }={\overrightarrow{x}}_j-{\overrightarrow{x}}_{j-1}$$

$$\overrightarrow{\gamma }={\mathrm{\nabla }}_{X}L({\overrightarrow{x}}_j,{\overrightarrow{\lambda }}_j)-{\mathrm{\nabla }}_{X}L({\overrightarrow{x}}_{j-1},{\overrightarrow{\lambda }}_j)$$

which is calculated using (1).

Since we have a constrained problem, we define a further quantity:

$$\overrightarrow{\eta }=\theta \overrightarrow{\gamma }+(1-\theta )\mathbf{B}\overrightarrow{\xi }$$

where

$$\begin{array}{r}\theta =\{\begin{array}{ll}1,& \text{if }{\overrightarrow{\xi }}^T\overrightarrow{\gamma }\ge 0.2{\overrightarrow{\xi }}^T\mathbf{B}\overrightarrow{\xi }\\ \frac{0.8{\overrightarrow{\xi }}^T\mathbf{B}\overrightarrow{\xi }}{{\overrightarrow{\xi }}^T\mathbf{B}\overrightarrow{\xi }-{\overrightarrow{\xi }}^T\overrightarrow{\gamma }},& \text{otherwise}\end{array}\end{array}$$

The definition of $\overrightarrow{\eta }$ ensures $\mathbf{B}$ remains positive semi-definite, which is a prereqesite to solving the QSP.

We can then perform the BFGS update:

$${\mathbf{B}}_{\mathbf{NEW}}=\mathbf{B}-\frac{\mathbf{B}\overrightarrow{\xi }{\overrightarrow{\xi }}^T\mathbf{B}}{{\overrightarrow{\xi }}^T\mathbf{B}\overrightarrow{\xi }}+\frac{\overrightarrow{\eta }{\overrightarrow{\eta }}^T}{{\overrightarrow{\xi }}^T\overrightarrow{\eta }}$$

The VMCON Algorithm

This page covers the mathematics and theory behind the VMCON algorithm. For completeness, the following flow diagram demonstrates how the algorithm is implemented at a high level.