e04vh solves sparse linear and nonlinear programming problems.

Syntax

C#
public static void e04vh( int start, int nf, int n, int nxname, int nfname, double objadd, int objrow, string prob, E04..::..E04VH_USRFUN usrfun, int[] iafun, int[] javar, int nea, double[] a, int[] igfun, int[] jgvar, int neg, double[] xlow, double[] xupp, string[] xnames, double[] flow, double[] fupp, string[] fnames, double[] x, int[] xstate, double[] xmul, double[] f, int[] fstate, double[] fmul, ref int ns, out int ninf, out double sinf, E04..::..e04vhOptions options, out int ifail )

public static void e04vh(
	int start,
	int nf,
	int n,
	int nxname,
	int nfname,
	double objadd,
	int objrow,
	string prob,
	E04..::..E04VH_USRFUN usrfun,
	int[] iafun,
	int[] javar,
	int nea,
	double[] a,
	int[] igfun,
	int[] jgvar,
	int neg,
	double[] xlow,
	double[] xupp,
	string[] xnames,
	double[] flow,
	double[] fupp,
	string[] fnames,
	double[] x,
	int[] xstate,
	double[] xmul,
	double[] f,
	int[] fstate,
	double[] fmul,
	ref int ns,
	out int ninf,
	out double sinf,
	E04..::..e04vhOptions options,
	out int ifail
)

Visual Basic
Public Shared Sub e04vh ( _ start As Integer, _ nf As Integer, _ n As Integer, _ nxname As Integer, _ nfname As Integer, _ objadd As Double, _ objrow As Integer, _ prob As String, _ usrfun As E04..::..E04VH_USRFUN, _ iafun As Integer(), _ javar As Integer(), _ nea As Integer, _ a As Double(), _ igfun As Integer(), _ jgvar As Integer(), _ neg As Integer, _ xlow As Double(), _ xupp As Double(), _ xnames As String(), _ flow As Double(), _ fupp As Double(), _ fnames As String(), _ x As Double(), _ xstate As Integer(), _ xmul As Double(), _ f As Double(), _ fstate As Integer(), _ fmul As Double(), _ ByRef ns As Integer, _ <OutAttribute> ByRef ninf As Integer, _ <OutAttribute> ByRef sinf As Double, _ options As E04..::..e04vhOptions, _ <OutAttribute> ByRef ifail As Integer _ )

Visual Basic

Public Shared Sub e04vh ( _
	start As Integer, _
	nf As Integer, _
	n As Integer, _
	nxname As Integer, _
	nfname As Integer, _
	objadd As Double, _
	objrow As Integer, _
	prob As String, _
	usrfun As E04..::..E04VH_USRFUN, _
	iafun As Integer(), _
	javar As Integer(), _
	nea As Integer, _
	a As Double(), _
	igfun As Integer(), _
	jgvar As Integer(), _
	neg As Integer, _
	xlow As Double(), _
	xupp As Double(), _
	xnames As String(), _
	flow As Double(), _
	fupp As Double(), _
	fnames As String(), _
	x As Double(), _
	xstate As Integer(), _
	xmul As Double(), _
	f As Double(), _
	fstate As Integer(), _
	fmul As Double(), _
	ByRef ns As Integer, _
	<OutAttribute> ByRef ninf As Integer, _
	<OutAttribute> ByRef sinf As Double, _
	options As E04..::..e04vhOptions, _
	<OutAttribute> ByRef ifail As Integer _
)

Visual C++
public: static void e04vh( int start, int nf, int n, int nxname, int nfname, double objadd, int objrow, String^ prob, E04..::..E04VH_USRFUN^ usrfun, array<int>^ iafun, array<int>^ javar, int nea, array<double>^ a, array<int>^ igfun, array<int>^ jgvar, int neg, array<double>^ xlow, array<double>^ xupp, array<String^>^ xnames, array<double>^ flow, array<double>^ fupp, array<String^>^ fnames, array<double>^ x, array<int>^ xstate, array<double>^ xmul, array<double>^ f, array<int>^ fstate, array<double>^ fmul, int% ns, [OutAttribute] int% ninf, [OutAttribute] double% sinf, E04..::..e04vhOptions^ options, [OutAttribute] int% ifail )

Visual C++

public:
static void e04vh(
	int start, 
	int nf, 
	int n, 
	int nxname, 
	int nfname, 
	double objadd, 
	int objrow, 
	String^ prob, 
	E04..::..E04VH_USRFUN^ usrfun, 
	array<int>^ iafun, 
	array<int>^ javar, 
	int nea, 
	array<double>^ a, 
	array<int>^ igfun, 
	array<int>^ jgvar, 
	int neg, 
	array<double>^ xlow, 
	array<double>^ xupp, 
	array<String^>^ xnames, 
	array<double>^ flow, 
	array<double>^ fupp, 
	array<String^>^ fnames, 
	array<double>^ x, 
	array<int>^ xstate, 
	array<double>^ xmul, 
	array<double>^ f, 
	array<int>^ fstate, 
	array<double>^ fmul, 
	int% ns, 
	[OutAttribute] int% ninf, 
	[OutAttribute] double% sinf, 
	E04..::..e04vhOptions^ options, 
	[OutAttribute] int% ifail
)

F#
static member e04vh : start : int * nf : int * n : int * nxname : int * nfname : int * objadd : float * objrow : int * prob : string * usrfun : E04..::..E04VH_USRFUN * iafun : int[] * javar : int[] * nea : int * a : float[] * igfun : int[] * jgvar : int[] * neg : int * xlow : float[] * xupp : float[] * xnames : string[] * flow : float[] * fupp : float[] * fnames : string[] * x : float[] * xstate : int[] * xmul : float[] * f : float[] * fstate : int[] * fmul : float[] * ns : int byref * ninf : int byref * sinf : float byref * options : E04..::..e04vhOptions * ifail : int byref -> unit

static member e04vh : 
        start : int * 
        nf : int * 
        n : int * 
        nxname : int * 
        nfname : int * 
        objadd : float * 
        objrow : int * 
        prob : string * 
        usrfun : E04..::..E04VH_USRFUN * 
        iafun : int[] * 
        javar : int[] * 
        nea : int * 
        a : float[] * 
        igfun : int[] * 
        jgvar : int[] * 
        neg : int * 
        xlow : float[] * 
        xupp : float[] * 
        xnames : string[] * 
        flow : float[] * 
        fupp : float[] * 
        fnames : string[] * 
        x : float[] * 
        xstate : int[] * 
        xmul : float[] * 
        f : float[] * 
        fstate : int[] * 
        fmul : float[] * 
        ns : int byref * 
        ninf : int byref * 
        sinf : float byref * 
        options : E04..::..e04vhOptions * 
        ifail : int byref -> unit

Parameters

start

Type: System..::..Int32

On entry: indicates how a starting point is to be obtained.

$start = 0$: Requests that the Crash procedure be used, unless a Basis file is provided via optional parameters Old Basis File, Insert File or Load File.
$start = 1$: Is the same as $start = 0$ but is more meaningful when a Basis file is given.
$start = 2$: Means that xstate and fstate define a valid starting point (probably from an earlier call, though not necessarily).

Constraint:

start = 0

1

2

nf: Type: System..::..Int32
On entry: $nf$ , the number of problem functions in $F (x)$ , including the objective function (if any) and the linear and nonlinear constraints. Upper and lower bounds on $x$ can be defined using the parameters xlow and xupp and should not be included in $F$ .

Constraint: $nf > 0$ .

n: Type: System..::..Int32
On entry: $n$ , the number of variables.

Constraint: $n > 0$ .

nxname

Type: System..::..Int32

On entry: the number of names provided in the array xnames.

$nxname = 1$: There are no names provided and generic names will be used in the output.
$nxname = n$: Names for all variables must be provided and will be used in the output.

Constraint:

nxname = 1 ​ or ​ n

nfname

Type: System..::..Int32

On entry: the number of names provided in the array fnames.

$nfname = 1$: There are no names provided and generic names will be used in the output.
$nfname = nf$: Names for all functions must be provided and will be used in the output.

Constraint:

nfname = 1 ​ or ​ nf

Note: if

nxname = 1

then nfname must also be

1

(and vice versa). Similarly, if

nxname = n

then nfname must be nf (and vice versa).

objadd: Type: System..::..Double
On entry: is a constant that will be added to the objective row $F_{obj}$ for printing purposes. Typically, $objadd = 0.0E+0$ .

objrow: Type: System..::..Int32
On entry: says which row of $F (x)$ is to act as the objective function. If there is no such row, set $objrow = 0$ . Then e04vh will seek a feasible point such that $l_{F} \leq F (x) \leq u_{F}$ and $l_{x} \leq x \leq u_{x}$ .

Constraint: $1 \leq objrow \leq nf or objrow = 0$ (or a feasible point problem).

prob: Type: System..::..String
On entry: is an $8$ -character name for the problem. prob is used in the printed solution and in some methods that output Basis files. A blank name may be used.

usrfun

Type: NagLibrary..::..E04..::..E04VH_USRFUN

usrfun must define the nonlinear portion

f (x)

of the problem functions

F (x) = f (x) + A x

, along with its gradient elements

G_{i j} (x) = \frac{\partial f_{i} (x)}{\partial x_{j}}

. (A dummy method is needed even if

f \equiv 0

and all functions are linear.)

In general, usrfun should return all function and gradient values on every entry except perhaps the last. This provides maximum reliability and corresponds to the default option setting,

Derivative Option = 1

The elements of

G (x)

are stored in the array

g (1 : usrfun_leng)

in the order specified by the input arrays igfun and jgvar.

In practice it is often convenient not to code gradients. e04vh is able to estimate them by finite differences, using a call to usrfun for each variable

x_{j}

for which some

\frac{\partial f_{i} (x)}{\partial x_{j}}

needs to be estimated. However, this reduces the reliability of the optimization algorithm, and it can be very expensive if there are many such variables

x_{j}

As a compromise, e04vh allows you to code as many gradients as you like. This option is implemented as follows. Just before usrfun is called, each element of the derivative array g is initialized to a specific value. On exit, any element retaining that value must be estimated by finite differences.

Some rules of thumb follow:

(i)	for maximum reliability, compute all gradients;
(ii)	if the gradients are expensive to compute, specify optional parameter Nonderivative Linesearch and use the value of the input parameter needg to avoid computing them on certain entries. (There is no need to compute gradients if $needg = 0$ on entry to usrfun.);
(iii)	if not all gradients are known, you must specify $Derivative Option = 0$ . You should still compute as many gradients as you can. (It often happens that some of them are constant or zero.);
(iv)	again, if the known gradients are expensive, don't compute them if $needg = 0$ on entry to usrfun;
(v)	use the input parameter status to test for special actions on the first or last entries;
(vi)	while usrfun is being developed, use the optional parameter Verify Level to check the computation of gradients that are supposedly known;
(vii)	usrfun is not called until the linear constraints and bounds on $x$ are satisfied. This helps confine $x$ to regions where the functions $f_{i} (x)$ are likely to be defined. However, be aware of the optional parameter Minor Feasibility Tolerance if the functions have singularities on the constraint boundaries;
(viii)	set $status = - 1$ if some of the functions are undefined. The linesearch will shorten the step and try again;
(ix)	set $status \leq - 2$ if you want e04vh to stop.

A delegate of type E04VH_USRFUN.

iafun: Type: array<System..::..Int32>[]()[][]
An array of size [lena]
On entry: define the coordinates $(i, j)$ and values $A_{i j}$ of the nonzero elements of the linear part $A$ of the function $F (x) = f (x) + A x$ .
In particular, nea triples $(iafun [k - 1], javar [k - 1], a [k - 1])$ define the row and column indices $i = iafun [k - 1]$ and $j = javar [k - 1]$ of the element $A_{i j} = a [k - 1]$ .

The coordinates may define the elements of $A$ in any order.

javar: Type: array<System..::..Int32>[]()[][]
An array of size [lena]
On entry: define the coordinates $(i, j)$ and values $A_{i j}$ of the nonzero elements of the linear part $A$ of the function $F (x) = f (x) + A x$ .
In particular, nea triples $(iafun [k - 1], javar [k - 1], a [k - 1])$ define the row and column indices $i = iafun [k - 1]$ and $j = javar [k - 1]$ of the element $A_{i j} = a [k - 1]$ .

The coordinates may define the elements of $A$ in any order.

nea: Type: System..::..Int32
On entry: is the number of nonzero entries in $A$ such that $F (x) = f (x) + A x$ .

Constraint: $0 \leq nea \leq lena$ .

a: Type: array<System..::..Double>[]()[][]
An array of size [lena]
On entry: define the coordinates $(i, j)$ and values $A_{i j}$ of the nonzero elements of the linear part $A$ of the function $F (x) = f (x) + A x$ .
In particular, nea triples $(iafun [k - 1], javar [k - 1], a [k - 1])$ define the row and column indices $i = iafun [k - 1]$ and $j = javar [k - 1]$ of the element $A_{i j} = a [k - 1]$ .

The coordinates may define the elements of $A$ in any order.

igfun: Type: array<System..::..Int32>[]()[][]
An array of size [leng]
On entry: define the coordinates $(i, j)$ of the nonzero elements of $G$ , the nonlinear part of the derivative $J (x) = G (x) + A$ of the function $F (x) = f (x) + A x$ . e04vj may be used to define these two arrays.
The coordinates can define the elements of $G$ in any order. However, usrfun must define the actual elements of g in exactly the same order as defined by the coordinates $(igfun, jgvar)$ .

jgvar: Type: array<System..::..Int32>[]()[][]
An array of size [leng]
On entry: define the coordinates $(i, j)$ of the nonzero elements of $G$ , the nonlinear part of the derivative $J (x) = G (x) + A$ of the function $F (x) = f (x) + A x$ . e04vj may be used to define these two arrays.
The coordinates can define the elements of $G$ in any order. However, usrfun must define the actual elements of g in exactly the same order as defined by the coordinates $(igfun, jgvar)$ .

neg: Type: System..::..Int32
On entry: the number of nonzero entries in $G$ .

Constraint: $0 \leq neg \leq leng$ .

xlow: Type: array<System..::..Double>[]()[][]
An array of size [n]
On entry: contain the lower and upper bounds $l_{x}$ and $u_{x}$ on the variables $x$ .
To specify a nonexistent lower bound ${[l_{x}]}_{j} = - \infty$ , set $xlow [j - 1] \leq - bigbnd$ , where $bigbnd$ is the optional parameter Infinite Bound Size. To specify a nonexistent upper bound ${[u_{x}]}_{j} = \infty$ , set $xupp [j - 1] \geq bigbnd$ .

To fix the $j$ th variable at $x_{j} = β$ , where $|β| < bigbnd$ , set $xlow [j - 1] = xupp [j - 1] = β$ .

Constraint: $xlow [i - 1] \leq xupp [i - 1]$ , for $i = 1, 2, \dots, n$ .

xupp: Type: array<System..::..Double>[]()[][]
An array of size [n]
On entry: contain the lower and upper bounds $l_{x}$ and $u_{x}$ on the variables $x$ .
To specify a nonexistent lower bound ${[l_{x}]}_{j} = - \infty$ , set $xlow [j - 1] \leq - bigbnd$ , where $bigbnd$ is the optional parameter Infinite Bound Size. To specify a nonexistent upper bound ${[u_{x}]}_{j} = \infty$ , set $xupp [j - 1] \geq bigbnd$ .

To fix the $j$ th variable at $x_{j} = β$ , where $|β| < bigbnd$ , set $xlow [j - 1] = xupp [j - 1] = β$ .

Constraint: $xlow [i - 1] \leq xupp [i - 1]$ , for $i = 1, 2, \dots, n$ .

xnames: Type: array<System..::..String>[]()[][]
An array of size [nxname]
On entry: the optional names for the variables.
If $nxname = 1$ , xnames is not referenced and default names will be used for output.

If $nxname = n$ , $xnames [j - 1]$ should contain the $8$ -character name of the $j$ th variable.

flow: Type: array<System..::..Double>[]()[][]
An array of size [nf]
On entry: contain the lower and upper bounds $l_{F}$ and $u_{F}$ on $F (x)$ .
To specify a nonexistent lower bound ${[l_{F}]}_{i} = - \infty$ , set $flow [i - 1] \leq - bigbnd$ . For a nonexistent upper bound ${[u_{F}]}_{i} = \infty$ , set $fupp [i - 1] \geq bigbnd$ .

To make the $i$ th constraint an equality at $F_{i} = β$ , where $|β| < bigbnd$ , set $flow [i - 1] = fupp [i - 1] = β$ .

Constraint: $flow [i - 1] \leq fupp [i - 1]$ , for $i = 1, 2, \dots, n$ .

fupp: Type: array<System..::..Double>[]()[][]
An array of size [nf]
On entry: contain the lower and upper bounds $l_{F}$ and $u_{F}$ on $F (x)$ .
To specify a nonexistent lower bound ${[l_{F}]}_{i} = - \infty$ , set $flow [i - 1] \leq - bigbnd$ . For a nonexistent upper bound ${[u_{F}]}_{i} = \infty$ , set $fupp [i - 1] \geq bigbnd$ .

To make the $i$ th constraint an equality at $F_{i} = β$ , where $|β| < bigbnd$ , set $flow [i - 1] = fupp [i - 1] = β$ .

Constraint: $flow [i - 1] \leq fupp [i - 1]$ , for $i = 1, 2, \dots, n$ .

fnames: Type: array<System..::..String>[]()[][]
An array of size [nfname]
On entry: the optional names for the problem functions.
If $nfname = 1$ , fnames is not referenced and default names will be used for output.

If $nfname = nf$ , $fnames [i - 1]$ should contain the $8$ -character name of the $i$ th row of $F$ .

x: Type: array<System..::..Double>[]()[][]
An array of size [n]
On entry: an initial estimate of the variables $x$ . See the following description of xstate.
On exit: the final values of the variable $x$ .

xstate

Type: array<System..::..Int32>[]()[][]

An array of size [n]

On entry: the initial state for each variable

x

start = 0

1

and no basis information is provided (the optional parameters Old Basis File, Insert File and Load File are all set to

0

; the default) x and xstate must be defined.

If nothing special is known about the problem, or if there is no wish to provide special information, you may set

x [j - 1] = 0.0

xstate [j - 1] = 0

, for all

j = 1, 2, \dots, n

. If you set

x [j - 1] = xlow [j - 1]

set

xstate [j - 1] = 4

; if you set

x [j - 1] = xupp [j - 1]

then set

xstate [j - 1] = 5

. In this case a Crash procedure is used to select an initial basis.

start = 0

1

and basis information is provided (at least one of the optional parameters Old Basis File, Insert File and Load File is nonzero) x and xstate need not be set.

start = 2

(Warm Start), x and xstate must be set (probably from a previous call). In this case

xstate [j - 1]

must be

0

1

2 ​ or ​ 3

, for

j = 1, 2, \dots, n

On exit: the final state of the variables.

$xstate [j - 1]$	State of variable $j$	Usual value of $x [j - 1]$
0	nonbasic	$xlow [j - 1]$
1	nonbasic	$xupp [j - 1]$
2	superbasic	Between $xlow [j - 1]$ and $xupp [j - 1]$
3	basic	Between $xlow [j - 1]$ and $xupp [j - 1]$

Basic and superbasic variables may be outside their bounds by as much as the optional parameter Minor Feasibility Tolerance. Note that if scaling is specified, the feasibility tolerance applies to the variables of the scaled problem. In this case, the variables of the original problem may be as much as

0.1

outside their bounds, but this is unlikely unless the problem is very badly scaled. Check the value of Primal infeasibility output to the unit number associated with the optional parameter Print File.

Very occasionally some nonbasic variables may be outside their bounds by as much as the optional parameter Minor Feasibility Tolerance, and there may be some nonbasics for which

x [j - 1]

lies strictly between its bounds.

ninf > 0

, some basic and superbasic variables may be outside their bounds by an arbitrary amount (bounded by sinf if scaling was not used).

Constraint:

0 \leq xstate [j - 1] \leq 5, for ​ j = 1, 2, \dots, n

xmul: Type: array<System..::..Double>[]()[][]
An array of size [n]
On exit: the vector of the dual variables (Lagrange multipliers) for the simple bounds $l_{x} \leq x \leq u_{x}$ .

f: Type: array<System..::..Double>[]()[][]
An array of size [nf]
On entry: an initial value for the problem functions $F$ . See the following description of fstate.
On exit: the final values for the problem functions $F$ (the values $F$ at the final point x).

fstate

Type: array<System..::..Int32>[]()[][]

An array of size [nf]

On entry: the initial state for the problem functions

F

start = 0

1

and no basis information is provided (the optional parameters Old Basis File, Insert File and Load File are all set to

0

; the default, f and fstate must be defined.

If nothing special is known about the problem, or if there is no wish to provide special information, you may set

f [i - 1] = 0.0

fstate [i - 1] = 0

, for all

i = 1, 2, \dots, nf

. Less trivially, to say that the optimal value of function

f [i - 1]

will probably be equal to one of its bounds, set

f [i - 1] = flow [i - 1]

and

fstate [i - 1] = 4

f [i - 1] = fupp [i - 1]

and

fstate [i - 1] = 5

as appropriate. In this case a Crash procedure is used to select an initial basis.

start = 0

1

and basis information is provided (at least one of the optional parameters Old Basis File, Insert File and Load File is nonzero), f and fstate need not be set.

start = 2

(Warm Start), f and fstate must be set (probably from a previous call). In this case

fstate [i - 1]

must be

0

1

2

3

, for

i = 1, 2, \dots, nf

On exit: the final state of the variables. The elements of fstate have the following meaning:

$fstate [i - 1]$	State of the corresponding slack variable	Usual value of $f [i - 1]$
0	nonbasic	$flow [i - 1]$
1	nonbasic	$fupp [i - 1]$
2	superbasic	Between $flow [i - 1]$ and $fupp [i - 1]$
3	basic	Between $flow [i - 1]$ and $fupp [i - 1]$

Basic and superbasic slack variables may lead to the corresponding functions being outside their bounds by as much as the optional parameter Minor Feasibility Tolerance.

Very occasionally some functions may be outside their bounds by as much as the optional parameter Minor Feasibility Tolerance, and there may be some nonbasics for which

f [i - 1]

lies strictly between its bounds.

ninf > 0

, some basic and superbasic variables may be outside their bounds by an arbitrary amount (bounded by sinf if scaling was not used).

Constraint:

0 \leq fstate [i - 1] \leq 5, for ​ i = 1, 2, \dots, nf

fmul: Type: array<System..::..Double>[]()[][]
An array of size [nf]
On entry: an estimate of $γ$ , the vector of Lagrange multipliers (shadow prices) for the constraints $l_{F} \leq F (x) \leq u_{F}$ . All nf components must be defined. If nothing is known about $γ$ , set $fmul [i - 1] = 0.0$ , for $i = 1, 2, \dots, nf$ . For warm start use the values from a previous call.
On exit: the vector of the dual variables (Lagrange multipliers) for the general constraints $l_{F} \leq F (x) \leq u_{F}$

ns: Type: System..::..Int32%
On entry: the number of superbasic variables. ns need not be specified for cold starts, but should retain its value from a previous call when warm start is used.
On exit: the final number of superbasic variables.

ninf: Type: System..::..Int32%
On exit: are the number and the sum of the infeasibilities of constraints that lie outside one of their bounds by more than the optional parameter Minor Feasibility Tolerance before the solution is unscaled.
If any linear constraints are infeasible, $x$ minimizes the sum of the infeasibilities of the linear constraints subject to the upper and lower bounds being satisfied. In this case ninf gives the number of variables and linear constraints lying outside their upper or lower bounds. The nonlinear constraints are not evaluated.

Otherwise, $x$ minimizes the sum of infeasibilities of the nonlinear constraints subject to the linear constraints and upper and lower bounds being satisfied. In this case ninf gives the number of components of $F (x)$ lying outside their bounds by more than the optional parameter Minor Feasibility Tolerance. Again this is before the solution is unscaled.

sinf: Type: System..::..Double%
On exit: are the number and the sum of the infeasibilities of constraints that lie outside one of their bounds by more than the optional parameter Minor Feasibility Tolerance before the solution is unscaled.
If any linear constraints are infeasible, $x$ minimizes the sum of the infeasibilities of the linear constraints subject to the upper and lower bounds being satisfied. In this case ninf gives the number of variables and linear constraints lying outside their upper or lower bounds. The nonlinear constraints are not evaluated.

Otherwise, $x$ minimizes the sum of infeasibilities of the nonlinear constraints subject to the linear constraints and upper and lower bounds being satisfied. In this case ninf gives the number of components of $F (x)$ lying outside their bounds by more than the optional parameter Minor Feasibility Tolerance. Again this is before the solution is unscaled.

options: Type: NagLibrary..::..E04..::..e04vhOptions
An Object of type E04.e04vhOptions. Used to configure optional parameters to this method.

ifail: Type: System..::..Int32%
On exit: $ifail = 0$ unless the method detects an error or a warning has been flagged (see [Error Indicators and Warnings]).

Description

e04vh is designed to minimize a linear or nonlinear function subject to bounds on the variables and sparse linear or nonlinear constraints. It is suitable for large-scale linear and quadratic programming and for linearly constrained optimization, as well as for general nonlinear programs of the form

\underset{x}{minimize} f_{0} (x) subject to ​ l \leq (\begin{matrix} x \\ f (x) \\ A_{L} x \end{matrix}) \leq u,

(1)

where

x

is an

n

-vector of variables,

l

and

u

are constant lower and upper bounds,

f_{0} (x)

is a smooth scalar objective function,

A_{L}

is a sparse matrix, and

f (x)

is a vector of smooth nonlinear constraint functions

\{f_{i} (x)\}

. The optional parameter Maximize specifies that

f_{0} (x)

should be maximized instead of minimized.

Ideally, the first derivatives (gradients) of

f_{0} (x)

and

f_{i} (x)

should be known and coded by you. If only some of the gradients are known, e04vh estimates the missing ones by finite differences.

f_{0} (x)

is linear and

f (x)

is absent, (1) is a linear program (LP) and e04vh applies the primal simplex method (see Dantzig (1963)). Sparse basis factors are maintained by LUSOL (see Gill et al. (1987)) as in MINOS (see Murtagh and Saunders (1995)).

If only the objective is nonlinear, the problem is linearly constrained (LC) and tends to solve more easily than the general case with nonlinear constraints (NC). For both nonlinear cases, e04vh applies a sparse sequential quadratic programming (SQP) method (see Gill et al. (2002)), using limited-memory quasi-Newton approximations to the Hessian of the Lagrangian. The merit function for step-length control is an augmented Lagrangian, as in the dense SQP solver e04wd (see Gill et al. (1986) and Gill et al. (1992)).

e04vh is suitable for nonlinear problems with thousands of constraints and variables, and is most efficient if only some of the variables enter nonlinearly, or there are relatively few degrees of freedom at a solution (i.e., many constraints are active). However, there is no limit on the number of degrees of freedom.

e04vh allows linear and nonlinear constraints and variables to be entered in an arbitrary order, and uses one method to define all the nonlinear functions.

The optimization problem is assumed to be in the form

\underset{x}{minimize} F_{obj} (x) subject to ​ l_{x} \leq x \leq u_{x}, l_{F} \leq F (x) \leq u_{F},

(2)

where the upper and lower bounds are constant,

F (x)

is a vector of smooth linear and nonlinear constraint functions

\{F_{i} (x)\}

, and

F_{obj} (x)

is one of the components of

F

to be minimized, as specified by the input parameter objrow. e04vh reorders the variables and constraints so that the problem is in the form (1).

Upper and lower bounds are specified for all variables and functions. The

j

th constraint may be defined as an equality by setting

l_{j} = u_{j}

. If certain bounds are not present, the associated elements of

l

u

should be set to special values that are treated as

- \infty

+ \infty

. Free variables and free constraints (‘free rows’) have both bounds infinite.

In general, the components of

F

are structured in the sense that they are formed from sums of linear and nonlinear functions of just some of the variables. This structure can be exploited by e04vh.

In many cases, the vector

F (x)

is a sum of linear and nonlinear functions. e04vh allows these terms to be specified separately, so that the linear part is defined just once by the input arguments iafun, javar and a. Only the nonlinear part is recomputed at each

x

Suppose that each component of

F (x)

is of the form

F_{i} (x) = f_{i} (x) + \sum_{j = 1}^{n} A_{i j} x_{j},

where

f_{i} (x)

is a nonlinear function (possibly zero) and the elements

A_{i j}

are constant. The

n f

n

Jacobian of

F (x)

is the sum of two sparse matrices of the same size:

F^{'} (x) = G (x) + A

, where

G (x) = f^{'} (x)

and

A

is the matrix with elements

\{A_{i j}\}

. The two matrices must be non-overlapping in the sense that each element of the Jacobian

F^{'} (x) = G (x) + A

comes from

G (x)

A

, but not both. The element cannot be split between

G (x)

and

A

For example, the function

F (x) = (\begin{matrix} 3 x_{1} + e^{x_{2}} x_{4} + x_{2}^{2} + 4 x_{4} - x_{3} + x_{5} \\ x_{2} + x_{3}^{2} + \sin x_{4} - 3 x_{5} \\ x_{1} - x_{3} \end{matrix})

can be written as

F (x) = f (x) + A x = (\begin{matrix} e^{x_{2}} x_{4} + x_{2}^{2} + 4 x_{4} \\ x_{3}^{2} + \sin x_{4} \\ 0 \end{matrix}) + (\begin{matrix} 3 x_{1} - x_{3} + x_{5} \\ x_{2} - 3 x_{5} \\ x_{1} - x_{3} \end{matrix}),

in which case

F^{'} (x) = (\begin{matrix} 3 & e^{x_{2}} x_{4} + 2 x_{2} & - 1 & e^{x_{2}} + 4 & 1 \\ 0 & 1 & 2 x_{3} & \cos x_{4} & - 3 \\ 1 & 0 & - 1 & 0 & 0 \end{matrix})

can be written as

F^{'} (x) = f^{'} (x) + A = G (x) + A

, where

G (x) = (\begin{matrix} 0 & e^{x_{2}} x_{4} + 2 x_{2} & 0 & e^{x_{2}} + 4 & 0 \\ 0 & 0 & 2 x_{3} & \cos x_{4} & 0 \\ 0 & 0 & 0 & 0 & 0 \end{matrix}), A = (\begin{matrix} 3 & 0 & - 1 & 0 & 1 \\ 0 & 1 & 0 & 0 & - 3 \\ 1 & 0 & - 1 & 0 & 0 \end{matrix}) .

Note: the element

e^{x_{2}} + 4

F^{'} (x)

appears in

G (x)

and is not split between

G (x)

and

A

although it contains a linear term.

The nonzero elements of

A

and

G

are provided to e04vh in coordinate form. The elements of

A

are entered as triples

(i, j, A_{i j})

in the arrays iafun, javar and a. The sparsity pattern

G

is entered as pairs

(i, j)

in the arrays igfun and jgvar. The corresponding entries

G_{i j}

(any that are known) are assigned to appropriate array elements

g (k)

in usrfun.

The elements of

A

and

G

may be stored in any order. Duplicate entries are ignored. igfun and jgvar may be defined automatically by method e04vj when

Derivative Option = 0

is specified and usrfun does not provide any gradients.

Throughout this document the symbol

ε

is used to represent the machine precision (see x02aj).

e04vh is based on SNOPTA, which is part of the SNOPT package described in Gill et al. (2005b).

References

Dantzig G B (1963) Linear Programming and Extensions Princeton University Press

Eldersveld S K (1991) Large-scale sequential quadratic programming algorithms PhD Thesis Department of Operations Research, Stanford University, Stanford

Fourer R (1982) Solving staircase linear programs by the simplex method Math. Programming 23 274–313

Gill P E, Murray W and Saunders M A (2002) SNOPT: An SQP Algorithm for Large-scale Constrained Optimization 12 979–1006 SIAM J. Optim.

Gill P E, Murray W and Saunders M A (2005a) Users' guide for SQOPT 7: a Fortran package for large-scale linear and quadratic programming Report NA 05-1 Department of Mathematics, University of California, San Diego http://www.ccom.ucsd.edu/~peg/papers/sqdoc7.pdf

Gill P E, Murray W and Saunders M A (2005b) Users' guide for SNOPT 7.1: a Fortran package for large-scale linear nonlinear programming Report NA 05-2 Department of Mathematics, University of California, San Diego http://www.ccom.ucsd.edu/~peg/papers/sndoc7.pdf

Gill P E, Murray W, Saunders M A and Wright M H (1986) Users' guide for NPSOL (Version 4.0): a Fortran package for nonlinear programming Report SOL 86-2 Department of Operations Research, Stanford University

Gill P E, Murray W, Saunders M A and Wright M H (1987) Maintaining LU factors of a general sparse matrix Linear Algebra and its Applics. 88/89 239–270

Gill P E, Murray W, Saunders M A and Wright M H (1992) Some theoretical properties of an augmented Lagrangian merit function Advances in Optimization and Parallel Computing (ed P M Pardalos) 101–128 North Holland

Hock W and Schittkowski K (1981) Test Examples for Nonlinear Programming Codes. Lecture Notes in Economics and Mathematical Systems 187 Springer–Verlag

Murtagh B A and Saunders M A (1978) Large-scale linearly constrained optimization 14 41–72 Math. Programming

Murtagh B A and Saunders M A (1982) A projected Lagrangian algorithm and its implementation for sparse nonlinear constraints Math. Program. Stud. 16 84–118

Murtagh B A and Saunders M A (1995) MINOS 5.4 users' guide Report SOL 83-20R Department of Operations Research, Stanford University

Error Indicators and Warnings

Note: e04vh may return useful information for one or more of the following detected errors or warnings.

Errors or warnings detected by the method:

$ifail = 1$: The initialization method (E04VGF not in this release) has not been called or at least one of lencw, leniw and lenrw is less than $600$ .

$ifail = 2$: An input parameter is invalid. The output message provides more details of the invalid argument.

$ifail = 3$: Requested accuracy could not be achieved.
A feasible solution has been found, but the requested accuracy in the dual infeasibilities could not be achieved. An abnormal termination has occurred, but e04vh is within $10^{- 2}$ of satisfying the Major Optimality Tolerance. Check that the Major Optimality Tolerance is not too small.

$ifail = 4$: The problem appears to be infeasible.
When the constraints are linear, this message is based on a relatively reliable indicator of infeasibility. Feasibility is measured with respect to the upper and lower bounds on the variables and slacks. Among all the points satisfying the general constraints $A x - s = 0$ (see (6) and (7) in [Major Iterations]), there is apparently no point that satisfies the bounds on $x$ and $s$ . Violations as small as the Minor Feasibility Tolerance are ignored, but at least one component of $x$ or $s$ violates a bound by more than the tolerance.
When nonlinear constraints are present, infeasibility is much harder to recognize correctly. Even if a feasible solution exists, the current linearization of the constraints may not contain a feasible point. In an attempt to deal with this situation, when solving each QP subproblem, e04vh is prepared to relax the bounds on the slacks associated with nonlinear rows.
If a QP subproblem proves to be infeasible or unbounded (or if the Lagrange multiplier estimates for the nonlinear constraints become large), e04vh enters so-called ‘nonlinear elastic’ mode. The subproblem includes the original QP objective and the sum of the infeasibilities – suitably weighted using the optional parameter Elastic Weight. In elastic mode, some of the bounds on the nonlinear rows are ‘elastic’ – i.e., they are allowed to violate their specific bounds. Variables subject to elastic bounds are known as elastic variables. An elastic variable is free to violate one or both of its original upper or lower bounds. If the original problem has a feasible solution and the elastic weight is sufficiently large, a feasible point eventually will be obtained for the perturbed constraints, and optimization can continue on the subproblem. If the nonlinear problem has no feasible solution, e04vh will tend to determine a ‘good’ infeasible point if the elastic weight is sufficiently large. (If the elastic weight were infinite, e04vh would locally minimize the nonlinear constraint violations subject to the linear constraints and bounds.)
Unfortunately, even though e04vh locally minimizes the nonlinear constraint violations, there may still exist other regions in which the nonlinear constraints are satisfied. Wherever possible, nonlinear constraints should be defined in such a way that feasible points are known to exist when the constraints are linearized.

$ifail = 5$: The problem appears to be unbounded (or badly scaled).
For linear problems, unboundedness is detected by the simplex method when a nonbasic variable can be increased or decreased by an arbitrary amount without causing a basic variable to violate a bound. Consider adding an upper or lower bound to the variable. Also, examine the constraints that have nonzeros in the associated column, to see if they have been formulated as intended.
Very rarely, the scaling of the problem could be so poor that numerical error will give an erroneous indication of unboundedness. Consider using the optional parameter Scale Option.
For nonlinear problems, e04vh monitors both the size of the current objective function and the size of the change in the variables at each step. If either of these is very large (as judged by the unbounded parameters (see [Major Iteration Log])), the problem is terminated and declared unbounded. To avoid large function values, it may be necessary to impose bounds on some of the variables in order to keep them away from singularities in the nonlinear functions.
The message may indicate an abnormal termination while enforcing the limit on the constraint violations. This exit implies that the objective is not bounded below in the feasible region defined by expanding the bounds by the value of the Violation Limit.

$ifail = 6$: Iteration limit reached.
Either the Iterations Limit or the Major Iterations Limit was exceeded before the required solution could be found. Check the iteration log to be sure that progress was being made. If so, and if you caused a basis file to be saved by using the optional parameter New Basis File, consider restarting the run using the optional parameter Old Basis File to see whether further progress can be made. If you have no basis file available, you might rerun the problem after increasing the optional parameters Minor Iterations Limit and/or Major Iterations Limit.
If none of the above limits have been reached, this error may mean that the problem appears to be more nonlinear than anticipated. The current set of basic and superbasic variables have been optimized as much as possible and a pricing operation (where a nonbasic variable is selected to become superbasic) is necessary to continue, but it can't continue as the number of superbasic variables has already reached the limit specified by the optional parameter Superbasics Limit. In general, raise the Superbasics Limit $s$ by a reasonable amount, bearing in mind the storage needed for the reduced Hessian.

$ifail = 7$

Numerical difficulties have been encountered and no further progress can be made.

Several circumstances could lead to this exit.

1.	usrfun could be returning accurate function values but inaccurate gradients (or vice versa). This is the most likely cause. Study the comments given for $ifail = 8$ , and do your best to ensure that the coding is correct.
2.	The function and gradient values could be consistent, but their precision could be too low. For example, accidental use of a low precision data type when a higher precision was intended would lead to a relative function precision of about $10^{- 6}$ instead of something like $10^{- 15}$ . The default Major Optimality Tolerance of $2 \times 10^{- 6}$ would need to be raised to about $10^{- 3}$ for optimality to be declared (at a rather suboptimal point). Of course, it is better to revise the function coding to obtain as much precision as economically possible.
3.	If function values are obtained from an expensive iterative process, they may be accurate to rather few significant figures, and gradients will probably not be available. One should specify Function Precision $t$ Major Optimality Tolerance $\sqrt{t}$ but even then, if $t$ is as large as $10^{- 5}$ or $10^{- 6}$ (only $5$ or $6$ significant figures), the same exit condition may occur. At present the only remedy is to increase the accuracy of the function calculation.
4.	An $L U$ factorization of the basis has just been obtained and used to recompute the basic variables $x_{B}$ , given the present values of the superbasic and nonbasic variables. A step of ‘iterative refinement’ has also been applied to increase the accuracy of $x_{B}$ . However, a row check has revealed that the resulting solution does not satisfy the current constraints $A x - s = 0$ sufficiently well. This probably means that the current basis is very ill-conditioned. If there are some linear constraints and variables, try $Scale Option = 1$ if scaling has not yet been used. For certain highly structured basis matrices (notably those with band structure), a systematic growth may occur in the factor $U$ . Consult the description of Umax and Growth in [Basis Factorization Statistics] and set the LU Factor Tolerance to $2.0$ (or possibly even smaller, but not less than $1.0$ ).
5.	The first factorization attempt will have found the basis to be structurally or numerically singular. (Some diagonals of the triangular matrix $U$ were respectively zero or smaller than a certain tolerance.) The associated variables are replaced by slacks and the modified basis is refactorized, but singularity persists. This must mean that the problem is badly scaled, or the LU Factor Tolerance is too much larger than $1.0$ . This is highly unlikely to occur.

$ifail = 8$: Derivative appears to be incorrect.
A check has been made on some elements of the Jacobian as returned in the parameter g of usrfun. At least one value disagrees remarkably with its associated forward difference estimate (the relative difference between the computed and estimated values is $1.0$ or more). This exit is a safeguard, since e04vh will usually fail to make progress when the computed gradients are seriously inaccurate. In the process it may expend considerable effort before terminating with $ifail = 7$ .
Check the function and Jacobian computation very carefully in usrfun. A simple omission could explain everything. If a component is very large, then give serious thought to scaling the function or the nonlinear variables.
If you feel certain that the computed Jacobian is correct (and that the forward-difference estimate is therefore wrong), you can specify $Verify Level = 0$ to prevent individual elements from being checked. However, the optimization procedure may have difficulty.

$ifail = 9$

Undefined user-supplied function.

You have indicated that the problem functions are undefined by assigning the value

status = - 1

on exit from usrfun. e04vh attempts to evaluate the problem functions closer to a point at which the functions are already known to be defined. This exit occurs if e04vh is unable to find a point at which the functions are defined. This will occur in the case of:

–	undefined functions with no recovery possible;
–	undefined functions at the first point;
–	undefined functions at the first feasible point; or
–	undefined functions when checking derivatives.

$ifail = 10$: User requested termination.
You have indicated the wish to terminate solution of the current problem by setting status to a value $< - 1$ on exit from usrfun.

$ifail = 11$: Internal memory allocation failed when attempting to obtain the required workspace. Please contact NAG.

$ifail = 12$: Internal memory allocation was insufficient. Please contact NAG.

$ifail = 13$: An error has occurred in the basis package, perhaps indicating incorrect setup of arrays. Set the optional parameter Print File and examine the output carefully for further information.

$ifail = 14$: An unexpected error has occurred. Set the optional parameter Print File and examine the output carefully for further information.

$ifail = -9000$: An error occured, see message report.
$ifail = -1000$: The array lengths are not the same for arrays
$ifail = -1000$: The array lengths are not the same for arrays
$ifail = -8000$: Negative dimension for array $〈value〉$
$ifail = -6000$: Invalid Parameters $〈value〉$

Accuracy

If the value of the optional parameter Major Optimality Tolerance is set to

10^{- d}

(

default value = \sqrt{ε}

) and

ifail = 0

on exit, then the final value of

f (x)

should have approximately

d

correct significant digits.

Parallelism and Performance

None.

Further Comments

This section describes the final output produced by e04vh. Intermediate and other output are given in [Description of Monitoring Information].

The Final Output

Print File > 0

, the final output, including a listing of the status of every variable and constraint will be sent to the channel numbers associated with Print File. The following describes the output for each constraint (row) and variable (column).

The ROWS section

General linear constraints take the form

l \leq A_{L} x \leq u

. The

i

th constraint is therefore of the form

α \leq ν_{i} x \leq β,

where

ν_{i}

is the

i

th row of

A_{L}

Internally, the constraints take the form

A_{L} x - s = 0

, where

s

is the set of slack variables (which satisfy the bounds

l \leq s \leq u

). For the

i

th row it is the slack variable

s_{i}

that is directly available and it is sometimes convenient to refer to its state. Nonlinear constraints

α \leq f_{i} (x) + ν_{i} x \leq β

are treated similarly, except that the row activity and degree of infeasibility are computed directly from

f_{i} (x) + ν_{i} x

, rather than

s_{i}

A full stop (.) is printed for any numerical value that is exactly zero.

Label

Description

Number

is the value of

n + i

. (This is used internally to refer to

s_{i}

in the intermediate output.)

Row

gives the name of the

i

th row.

State

the state of the

i

th row relative to the bounds

α

and

β

. The various states possible are as follows:

LL	the row is at its lower limit, $α$ .
UL	the row is at its upper limit, $β$ .
EQ	the limits are the same ( $α = β$ ).
FR	$s_{i}$ is nonbasic and currently zero, even though it is free to take any value between its bounds $α$ and $β$ .
BS	$s_{i}$ is basic.
SBS	$s_{i}$ is superbasic.

A key is sometimes printed before State. Note that unless the optional parameter

Scale Option = 0

is specified, the tests for assigning a key are applied to the variables of the scaled problem.

A	Alternative optimum possible. The variable is nonbasic, but its reduced gradient is essentially zero. This means that if the variable were allowed to start moving away from its bound, there would be no change in the value of the objective function. The values of the other free variables might change, giving a genuine alternative solution. However, if there are any degenerate variables (labelled D), the actual change might prove to be zero, since one of them could encounter a bound immediately. In either case, the values of the Lagrange multipliers might also change.
D	Degenerate. The variable is basic or superbasic, but it is equal (or very close) to one of its bounds.
I	Infeasible. The variable is basic or superbasic and is currently violating one of its bounds by more than the value of the Feasibility Tolerance.
N	Not precisely optimal. The variable is nonbasic or superbasic. If the value of the reduced gradient for the variable exceeds the value of the optional parameter Major Optimality Tolerance, the solution would not be declared optimal because the reduced gradient for the variable would not be considered negligible.

Activity

is the value of

ν_{i} x

(or

f_{i} (x) + ν_{i} x

for nonlinear rows) at the final iterate.

Slack Activity

is the value by which the row differs from its nearest bound. (For the free row (if any), it is set to Activity.)

Lower Limit

α

, the lower bound on the row.

Upper Limit

β

, the upper bound on the row.

Dual Activity

is the value of the dual variable

π_{i}

(the Lagrange multiplier for the

i

th constraint). The full vector

π

always satisfies

B^{T} π = g_{B}

, where

B

is the current basis matrix and

g_{B}

contains the associated gradients for the current objective function. For FP problems,

π_{i}

is set to zero.

gives the index

i

of the

i

th row.

The COLUMNS section

Let the

j

th component of

x

be the variable

x_{j}

and assume that it satisfies the bounds

α \leq x_{j} \leq β

. A fullstop (.) is printed for any numerical value that is exactly zero.

Label

Description

Number

is the column number

j

. (This is used internally to refer to

x_{j}

in the intermediate output.)

Column

gives the name of

x_{j}

State

the state of

x_{j}

relative to the bounds

α

and

β

. The various states possible are as follows:

LL	$x_{j}$ is nonbasic at its lower limit, $α$ .
UL	$x_{j}$ is nonbasic at its upper limit, $β$ .
EQ	$x_{j}$ is nonbasic and fixed at the value $α = β$ .
FR	$x_{j}$ is nonbasic at some value strictly between its bounds: $α < x_{j} < β$ .
BS	$x_{j}$ is basic. Usually $α < x_{j} < β$ .
SBS	$x_{j}$ is superbasic. Usually $α < x_{j} < β$ .

A key is sometimes printed before State. Note that unless the optional parameter

Scale Option = 0

is specified, the tests for assigning a key are applied to the variables of the scaled problem.

A	Alternative optimum possible. The variable is nonbasic, but its reduced gradient is essentially zero. This means that if the variable were allowed to start moving away from its bound, there would be no change in the value of the objective function. The values of the other free variables might change, giving a genuine alternative solution. However, if there are any degenerate variables (labelled D), the actual change might prove to be zero, since one of them could encounter a bound immediately. In either case, the values of the Lagrange multipliers might also change.
D	Degenerate. The variable is basic or superbasic, but it is equal (or very close) to one of its bounds.
I	Infeasible. The variable is basic or superbasic and is currently violating one of its bounds by more than the value of the Feasibility Tolerance.
N	Not precisely optimal. The variable is nonbasic or superbasic. If the value of the reduced gradient for the variable exceeds the value of the optional parameter Major Optimality Tolerance, the solution would not be declared optimal because the reduced gradient for the variable would not be considered negligible.

Activity

is the value of

x_{j}

at the final iterate.

Obj Gradient

is the value of

g_{j}

at the final iterate. For FP problems,

g_{j}

is set to zero.

Lower Limit

is the lower bound specified for the variable. None indicates that

xlow [j - 1] \leq - infbnd

Upper Limit

is the upper bound specified for the variable. None indicates that

xupp [j - 1] \geq infbnd

Reduced Gradnt

is the value of the reduced gradient

d_{j} = g_{j} - π^{T} a_{j}

where

a_{j}

is the

j

th column of the constraint matrix. For FP problems,

d_{j}

is set to zero.

m + j

is the value of

m + j

Note that movement off a constraint (as opposed to a variable moving away from its bound) can be interpreted as allowing the entry in the Slack Activity column to become positive.

Numerical values are output with a fixed number of digits; they are not guaranteed to be accurate to this precision.

Example

This example is a reformulation of Problem 74 from Hock and Schittkowski (1981) and involves the minimization of the nonlinear function

f (x) = 10^{- 6} x_{3}^{3} + \frac{2}{3} \times 10^{- 6} x_{4}^{3} + 3 x_{3} + 2 x_{4}

subject to the bounds

\begin{matrix} - 0.55 \leq x_{1} \leq 0.55, \\ - 0.55 \leq x_{2} \leq 0.55, \\ 0 \leq x_{3} \leq 1200, \\ 0 \leq x_{4} \leq 1200, \end{matrix}

to the nonlinear constraints

\begin{matrix} 1000 \sin (- x_{1} - 0.25) + 1000 \sin (- x_{2} - 0.25) - x_{3} & = & - 894.8, \\ 1000 \sin (x_{1} - 0.25) + 1000 \sin (x_{1} - x_{2} - 0.25) - x_{4} & = & - 894.8, \\ 1000 \sin (x_{2} - 0.25) + 1000 \sin (x_{2} - x_{1} - 0.25) & = & - 1294.8, \end{matrix}

and to the linear constraints

\begin{matrix} - x_{1} + x_{2} \geq - 0.55, \\ x_{1} - x_{2} \geq - 0.55 . \end{matrix}

The initial point, which is infeasible, is

x_{0} = {(\begin{matrix} 0, & 0, & 0, & 0 \end{matrix})}^{T},

and

f (x_{0}) = 0

The optimal solution (to five figures) is

x^{*} = {(0.11887, - 0.39623, 679.94, 1026.0)}^{T},

and

f (x^{*}) = 5126.4

. All the nonlinear constraints are active at the solution.

The formulation of the problem combines the constraints and the objective into a single vector (

F

) which is split into linear part (

A_{L} x

) and a nonlinear part (

f

). For example we could write

F = (\begin{matrix} 1000 \sin (- x_{1} - 0.25) + 1000 \sin (- x_{2} - 0.25) - x_{3} \\ 1000 \sin (x_{1} - 0.25) + 1000 \sin (x_{1} - x_{2} - 0.25) - x_{4} \\ 1000 \sin (x_{2} - 0.25) + 1000 \sin (x_{2} - x_{1} - 0.25) \\ - x_{1} + x_{2} \\ x_{1} - x_{2} \\ 10^{- 6} x_{3}^{3} + \frac{2}{3} \times 10^{- 6} x_{4}^{3} + 3 x_{3} + 2 x_{4} \end{matrix}) = f + A_{L} x

where

f = (\begin{matrix} 1000 \sin (- x_{1} - 0.25) + 1000 \sin (- x_{2} - 0.25) \\ 1000 \sin (x_{1} - 0.25) + 1000 \sin (x_{1} - x_{2} - 0.25) \\ 1000 \sin (x_{2} - 0.25) + 1000 \sin (x_{2} - x_{1} - 0.25) \\ 0 \\ 0 \\ 10^{- 6} x_{3}^{3} + \frac{2}{3} \times 10^{- 6} x_{4}^{3} \end{matrix})

and

A_{L} = (\begin{matrix} 0 & 0 & - 1 & 0 \\ 0 & 0 & 0 & - 1 \\ 0 & 0 & 0 & 0 \\ - 1 & 1 & 0 & 0 \\ 1 & - 1 & 0 & 0 \\ 0 & 0 & 3 & 2 \end{matrix}) .

The nonzero elements of

A_{L}

need to be stored in the triples

(iafun [k - 1], javar [k - 1], a [k - 1])

in any order. For example

$k$	$1$	$2$	$3$	$4$	$5$	$6$	$7$	$8$
$iafun [k - 1]$	$1$	$2$	$4$	$4$	$5$	$5$	$6$	$6$
$javar [k - 1]$	$3$	$4$	$1$	$2$	$1$	$2$	$3$	$4$
$a [k - 1]$	$- 1$	$- 1$	$- 1$	$1$	$1$	$- 1$	$3$	$2$

The nonlinear functions

f

and the Jacobian need to be supplied in usrfun. Please note that there is no need to assign any value to

f_{4}

f_{5}

as there is no nonlinear part in

F_{4}

F_{5}

The nonzero entries of the Jacobian of

f

are

\begin{matrix} \frac{\partial f_{1}}{\partial x_{1}} & = - 1000 \cos (- x_{1} - 0.25) \\ \frac{\partial f_{1}}{\partial x_{2}} & = - 1000 \cos (- x_{2} - 0.25) \\ \frac{\partial f_{2}}{\partial x_{1}} & = 1000 \cos (x_{1} - 0.25) + 1000 \cos (x_{1} - x_{2} - 0.25) \\ \frac{\partial f_{2}}{\partial x_{2}} & = - 1000 \cos (x_{1} - x_{2} - 0.25) \\ \frac{\partial f_{3}}{\partial x_{1}} & = - 1000 \cos (x_{2} - x_{1} - 0.25) \\ \frac{\partial f_{3}}{\partial x_{2}} & = 1000 \cos (x_{2} - 0.25) + 1000 \cos (x_{2} - x_{1} - 0.25) \\ \frac{\partial f_{6}}{\partial x_{3}} & = {3 \times 10}^{- 6} x_{3}^{2} \\ \frac{\partial f_{6}}{\partial x_{4}} & = {2 \times 10}^{- 6} x_{4}^{2} \end{matrix}

So the arrays igfun and jgvar must contain:

$k$	$1$	$2$	$3$	$4$	$5$	$6$	$7$	$8$
$igfun [k - 1]$	$1$	$1$	$2$	$2$	$3$	$3$	$6$	$6$
$jgvar [k - 1]$	$1$	$2$	$1$	$2$	$1$	$2$	$3$	$4$

and usrfun must return in

g [k - 1]

the value of

\frac{\partial f_{i}}{\partial x_{j}}

, where

i = igfun [k - 1]

and

j = jgvar [k - 1]

Example program (C#): e04vhe.cs

Example program data: e04vhe.d

Example program results: e04vhe.r

Algorithmic Details

Here we summarise the main features of the SQP algorithm used in e04vh and introduce some terminology used in the description of the method and its arguments. The SQP algorithm is fully described in Gill et al. (2002).

Constraints and Slack Variables

Problem (1) contains

n

variables in

x

. Let

m

be the number of components of

f (x)

and

A_{L} x

combined. The upper and lower bounds on those terms define the general constraints of the problem. e04vh converts the general constraints to equalities by introducing a set of slack variables

s = {(s_{1}, s_{2}, \dots, s_{m})}^{T}

. For example, the linear constraint

5 \leq 2 x_{1} + 3 x_{2} \leq \infty

is replaced by

2 x_{1} + 3 x_{2} - s_{1} = 0

together with the bounded slack

5 \leq s_{1} \leq \infty

. The minimization problem (1) can therefore be written in the equivalent form

\underset{x, s}{minimize} f_{0} (x) subject to ​ (\begin{matrix} f (x) \\ A_{L} x \end{matrix}) - s = 0, l \leq (\begin{matrix} x \\ s \end{matrix}) \leq u .

(3)

The general constraints become the equalities

f (x) - s_{N} = 0

and

A_{L} x - s_{L} = 0

, where

s_{L}

and

s_{N}

are the linear and nonlinear slacks.

Major Iterations

The basic structure of the SQP algorithm involves major and minor iterations. The major iterations generate a sequence of iterates

\{x_{k}\}

that satisfy the linear constraints and converge to a point that satisfies the nonlinear constraints and the first-order conditions for optimality. At each iterate

x_{k}

a QP subproblem is used to generate a search direction towards the next iterate

x_{k + 1}

. The constraints of the subproblem are formed from the linear constraints

A_{L} x - s_{L} = 0

and the linearized constraint

f (x_{k}) + f^{'} (x_{k}) (x - x_{k}) - s_{N} = 0,

(4)

where

f^{'} (x_{k})

denotes the Jacobian matrix, whose elements are the first derivatives of

f (x)

evaluated at

x_{k}

. The QP constraints therefore comprise the

m

linear constraints

\begin{matrix} f^{'} (x_{k}) x - s_{N} & = & - f (x_{k}) + f^{'} (x_{k}) x_{k}, \\ A_{L} x - s_{L} & = & 0, \end{matrix}

(5)

where

x

and

s

are bounded above and below by

u

and

l

as before. If the

m

n

matrix

A

and

m

-vector

b

are defined as

A = (\begin{matrix} f^{'} (x_{k}) \\ A_{L} \end{matrix}) and b = (\begin{matrix} - f (x_{k}) + f^{'} (x_{k}) x_{k} \\ 0 \end{matrix}),

(6)

then the QP subproblem can be written as

\underset{x, s}{minimize} q (x, x_{k}) = g_{k}^{T} (x - x_{k}) + \frac{1}{2} {(x - x_{k})}^{T} H_{k} (x - x_{k}) subject to ​ A x - s = b, l \leq (\begin{matrix} x \\ s \end{matrix}) \leq u,

(7)

where

q (x, x_{k})

is a quadratic approximation to a modified Lagrangian function (see Gill et al. (2002)). The matrix

H_{k}

is a quasi-Newton approximation to the Hessian of the Lagrangian. A BFGS update is applied after each major iteration. If some of the variables enter the Lagrangian linearly the Hessian will have some zero rows and columns. If the nonlinear variables appear first, then only the

n_{1}

rows and columns of the Hessian need to be approximated, where

n_{1}

is the number of nonlinear variables. This quantity is determined by the implicit values of the number of nonlinear objective and Jacobian variables determined after the constraints and variables are reordered.

Minor Iterations

Solving the QP subproblem is itself an iterative procedure. Here, the iterations of the QP solver e04nq form the minor iterations of the SQP method. e04nq uses a reduced-Hessian active-set method implemented as a reduced-gradient method. At each minor iteration, the constraints

A x - s = b

are partitioned into the form

B x_{B} + S x_{S} + N x_{N} = b,

(8)

where the basis matrix

B

is square and nonsingular, and the matrices

S

and

N

are the remaining columns of

(\begin{matrix} A & - I \end{matrix})

. The vectors

x_{B}

x_{S}

and

x_{N}

are the associated basic, superbasic and nonbasic variables respectively; they are a permutation of the elements of

x

and

s

. At a QP subproblem, the basic and superbasic variables will lie somewhere between their bounds, while the nonbasic variables will normally be equal to one of their bounds. At each iteration,

x_{S}

is regarded as a set of independent variables that are free to move in any desired direction, namely one that will improve the value of the QP objective (or the sum of infeasibilities). The basic variables are then adjusted in order to ensure that

(x, s)

continues to satisfy

A x - s = b

. The number of superbasic variables (

n_{S}

, say) therefore indicates the number of degrees of freedom remaining after the constraints have been satisfied. In broad terms,

n_{S}

is a measure of how nonlinear the problem is. In particular,

n_{S}

will always be zero for LP problems.

If it appears that no improvement can be made with the current definition of

B

S

and

N

, a nonbasic variable is selected to be added to

S

, and the process is repeated with the value of

n_{S}

increased by one. At all stages, if a basic or superbasic variable encounters one of its bounds, the variable is made nonbasic and the value of

n_{S}

is decreased by one.

Associated with each of the

m

equality constraints

A x - s = b

are the dual variables

π

. Similarly, each variable in

(x, s)

has an associated reduced gradient

d_{j}

. The reduced gradients for the variables

x

are the quantities

g - A^{T} π

, where

g

is the gradient of the QP objective, and the reduced gradients for the slacks are the dual variables

π

. The QP subproblem is optimal if

d_{j} \geq 0

for all nonbasic variables at their lower bounds,

d_{j} \leq 0

for all nonbasic variables at their upper bounds, and

d_{j} = 0

for other variables, including superbasics. In practice, an approximate QP solution

({\hat{x}}_{k}, {\hat{s}}_{k}, {\hat{π}}_{k})

is found by relaxing these conditions.

The Merit Function

After a QP subproblem has been solved, new estimates of the solution are computed using a linesearch on the augmented Lagrangian merit function

M (x, s, π) = f_{0} (x) - π^{T} (f (x) - s_{N}) + \frac{1}{2} {(f (x) - s_{N})}^{T} D (f (x) - s_{N}),

(9)

where

D

is a diagonal matrix of penalty parameters

(D_{i i} \geq 0)

, and

π

now refers to dual variables for the nonlinear constraints in (1). If

(x_{k}, s_{k}, π_{k})

denotes the current solution estimate and

({\hat{x}}_{k}, {\hat{s}}_{k}, {\hat{π}}_{k})

denotes the QP solution, the linesearch determines a step

α_{k}

(0 < α_{k} \leq 1)

such that the new point

(\begin{matrix} x_{k + 1} \\ s_{k + 1} \\ π_{k + 1} \end{matrix}) = (\begin{matrix} x_{k} \\ s_{k} \\ π_{k} \end{matrix}) + α_{k} (\begin{matrix} {\hat{x}}_{k} - x_{k} \\ {\hat{s}}_{k} - s_{k} \\ {\hat{π}}_{k} - π_{k} \end{matrix})

(10)

gives a sufficient decrease in the merit function

M

. When necessary, the penalties in

D

are increased by the minimum-norm perturbation that ensures descent for

M

(see Gill et al. (1992)). The value of

s_{N}

is adjusted to minimize the merit function as a function of

s

before the solution of the QP subproblem (see Gill et al. (1986) and Eldersveld (1991)).

Treatment of Constraint Infeasibilities

e04vh makes explicit allowance for infeasible constraints. First, infeasible linear constraints are detected by solving the linear program

\underset{x, v, w}{minimize} e^{T} (v + w) subject to ​ l \leq (\begin{matrix} x \\ A_{L} x - v + w \end{matrix}) \leq u, v \geq 0, w \geq 0,

(11)

where

e

is a vector of ones, and the nonlinear constraint bounds are temporarily excluded from

l

and

u

. This is equivalent to minimizing the sum of the general linear constraint violations subject to the bounds on

x

. (The sum is the

ℓ_{1}

-norm of the linear constraint violations. In the linear programming literature, the approach is called elastic programming.)

The linear constraints are infeasible if the optimal solution of (11) has

v \neq 0

w \neq 0

. e04vh then terminates without computing the nonlinear functions.

Otherwise, all subsequent iterates satisfy the linear constraints. (Such a strategy allows linear constraints to be used to define a region in which the functions can be safely evaluated.) e04vh proceeds to solve nonlinear problems as given, using search directions obtained from the sequence of QP subproblems (see (7)).

If a QP subproblem proves to be infeasible or unbounded (or if the dual variables

π

for the nonlinear constraints become large), e04vh enters ‘elastic’ mode and thereafter solves the problem

\underset{x, v, w}{minimize} f_{0} (x) + γ e^{T} (v + w) subject to ​ l \leq (\begin{matrix} x \\ f (x) - v + w \\ A_{L} x \end{matrix}) \leq u, v \geq 0, w \geq 0,

(12)

where

γ

is a non-negative parameter (the elastic weight), and

f_{0} (x) + γ e^{T} (v + w)

is called a composite objective (the

ℓ_{1}

penalty function for the nonlinear constraints).

The value of

γ

may increase automatically by multiples of

10

if the optimal

v

and

w

continue to be nonzero. If

γ

is sufficiently large, this is equivalent to minimizing the sum of the nonlinear constraint violations subject to the linear constraints and bounds.

The initial value of

γ

is controlled by the optional parameter Elastic Weight.

Description of Monitoring Information

e04vh produces monitoring information, statistical information and information about the solution. [The Final Output] contains details of the final output information sent to the unit specified by the optional parameter Print File. This section contains other details of output information.

Major Iteration Log

This section describes the output to unit Print File if

Major Print Level > 0

. One line of information is output every

k

th major iteration, where

k

is Print Frequency.

Label	Description
Itns	is the cumulative number of minor iterations.
Major	is the current major iteration number.
Minors	is the number of iterations required by both the feasibility and optimality phases of the QP subproblem. Generally, Minors will be $1$ in the later iterations, since theoretical analysis predicts that the correct active set will be identified near the solution (see [Algorithmic Details]).
Step	is the step length $α$ taken along the current search direction $p$ . The variables $x$ have just been changed to $x + α p$ . On reasonably well-behaved problems, the unit step will be taken as the solution is approached.
nCon	the number of times usrfun has been called to evaluate the nonlinear problem functions. Evaluations needed for the estimation of the derivatives by finite differences are not included. nCon is printed as a guide to the amount of work required for the linesearch.
Feasible	is the value of $v_{\max}$ (see (13)), the maximum component of the scaled nonlinear constraint residual (see optional parameter Major Feasibility Tolerance). The solution is regarded as acceptably feasible if Feasible is less than the Major Feasibility Tolerance. In this case, the entry is contained in parentheses. If the constraints are linear, all iterates are feasible and this entry is not printed.
Optimal	is the value of $c_{\max}$ (see (14)), the maximum complementary gap (see optional parameter Major Optimality Tolerance). It is an estimate of the degree of nonoptimality of the reduced costs. Both Feasible and Optimal are small in the neighbourhood of a solution.
MeritFunction	is the value of the augmented Lagrangian merit function (see (8)). This function will decrease at each iteration unless it was necessary to increase the penalty parameters (see [The Merit Function]). As the solution is approached, MeritFunction will converge to the value of the objective at the solution. In elastic mode, the merit function is a composite function involving the constraint violations weighted by the elastic weight. If the constraints are linear, this item is labelled Objective, the value of the objective function. It will decrease monotonically to its optimal value.
L+U	is the number of nonzeros representing the basis factors $L$ and $U$ on completion of the QP subproblem. If nonlinear constraints are present, the basis factorization $B = L U$ is computed at the start of the first minor iteration. At this stage, $L+U = lenL+lenU$ , where lenL (see [Basis Factorization Statistics]) is the number of subdiagonal elements in the columns of a lower triangular matrix and lenU (see [Basis Factorization Statistics]) is the number of diagonal and superdiagonal elements in the rows of an upper-triangular matrix. As columns of $B$ are replaced during the minor iterations, L+U may fluctuate up or down but, in general, will tend to increase. As the solution is approached and the minor iterations decrease towards zero, L+U will reflect the number of nonzeros in the $L U$ factors at the start of the QP subproblem. If the constraints are linear, refactorization is subject only to the Factorization Frequency, and L+U will tend to increase between factorizations.
BSwap	is the number of columns of the basis matrix $B$ that were swapped with columns of $S$ to improve the condition of $B$ . The swaps are determined by an $L U$ factorization of the rectangular matrix $B_{S} = {(B S)}^{T}$ with stability being favoured more than sparsity.
nS	is the current number of superbasic variables.
condHz	is an estimate of the condition number of $R^{T} R$ , itself an estimate of $Z^{T} H Z$ , the reduced Hessian of the Lagrangian. The condition number is the square of the ratio of the largest and smallest diagonals of the upper triangular matrix $R$ , this being a lower bound on the condition number of $R^{T} R$ . condHz gives a rough indication of whether or not the optimization procedure is having difficulty. If $ε$ is the relative machine precision being used, the SQP algorithm will make slow progress if condHz becomes as large as $ε^{- 1 / 2} \approx 10^{8}$ , and will probably fail to find a better solution if condHz reaches $ε^{- 3 / 4} \approx 10^{12}$ . To guard against high values of condHz, attention should be given to the scaling of the variables and the constraints. In some cases it may be necessary to add upper or lower bounds to certain variables to keep them a reasonable distance from singularities in the nonlinear functions or their derivatives.
Penalty	is the Euclidean norm of the vector of penalty parameters used in the augmented Lagrangian merit function (not printed if there are no nonlinear constraints).

The summary line may include additional code characters that indicate what happened during the course of the major iteration. These will follow the separator ‘_’ in the output

Label	Description
c	central differences have been used to compute the unknown components of the objective and constraint gradients. A switch to central differences is made if either the linesearch gives a small step, or $x$ is close to being optimal. In some cases, it may be necessary to re-solve the QP subproblem with the central difference gradient and Jacobian.
d	during the linesearch it was necessary to decrease the step in order to obtain a maximum constraint violation conforming to the value of the optional parameter Violation Limit.
D	you set $status = - 1$ on exit from usrfun, indicating that the linesearch needed to be done with a smaller value of the step length $α$ .
l	the norm wise change in the variables was limited by the value of the optional parameter Major Step Limit. If this output occurs repeatedly during later iterations, it may be worthwhile increasing the value of the optional parameter Major Step Limit.
i	if e04vh is not in elastic mode, an i signifies that the QP subproblem is infeasible. This event triggers the start of nonlinear elastic mode, which remains in effect for all subsequent iterations. Once in elastic mode, the QP subproblems are associated with the elastic problem (12) (see [Treatment of Constraint Infeasibilities]). If e04vh is already in elastic mode, an i indicates that the minimizer of the elastic subproblem does not satisfy the linearized constraints. (In this case, a feasible point for the usual QP subproblem may or may not exist.)
M	an extra evaluation of the problem functions was needed to define an acceptable positive definite quasi-Newton update to the Lagrangian Hessian. This modification is only done when there are nonlinear constraints.
m	this is the same as M except that it was also necessary to modify the update to include an augmented Lagrangian term.
n	no positive definite BFGS update could be found. The approximate Hessian is unchanged from the previous iteration.
R	the approximate Hessian has been reset by discarding all but the diagonal elements. This reset will be forced periodically by the Hessian Frequency and Hessian Updates keywords. However, it may also be necessary to reset an ill-conditioned Hessian from time to time.
r	the approximate Hessian was reset after ten consecutive major iterations in which no BFGS update could be made. The diagonals of the approximate Hessian are retained if at least one update has been done since the last reset. Otherwise, the approximate Hessian is reset to the identity matrix.
s	a self-scaled BFGS update was performed. This update is used when the Hessian approximation is diagonal, and hence always follows a Hessian reset.
t	the minor iterations were terminated because of the Minor Iterations Limit.
T	the minor iterations were terminated because of the New Superbasics Limit.
u	the QP subproblem was unbounded.
w	a weak solution of the QP subproblem was found.
z	the Superbasics Limit was reached.

Minor Iteration Log

Minor Print Level > 0

, one line of information is output to the Print file every

k

th minor iteration, where

k

is the specified Print Frequency. A heading is printed before the first such line following a basis factorization. The heading contains the items described below. In this description, a pricing operation is the process by which a nonbasic variable is selected to become superbasic (in addition to those already in the superbasic set). The selected variable is denoted by jq. Variable jq often becomes basic immediately. Otherwise it remains superbasic, unless it reaches its opposite bound and returns to the nonbasic set.

If Partial Price is in effect, variable jq is selected from

A_{pp}

I_{pp}

, the

pp

th segments of the constraint matrix

(\begin{matrix} A & - I \end{matrix})

Label	Description
Itn	the current iteration number.
LPmult or QPmult	is the reduced cost (or reduced gradient) of the variable jq selected by the pricing procedure at the start of the present iteration. Algebraically, the reduced gradient is $d_{j} = g_{j} - π^{T} a_{j}$ for $j = jq$ , where $g_{j}$ is the gradient of the current objective function, $π$ is the vector of dual variables for the QP subproblem, and $a_{j}$ is the $j$ th column of $(\begin{matrix} A & - I \end{matrix})$ . Note that the reduced cost is the $1$ -norm of the reduced-gradient vector at the start of the iteration, just after the pricing procedure.
LPstep or QPstep	is the step length $α$ taken along the current search direction $p$ . The variables $x$ have just been changed to $x + α p$ . Write Step to stand for LPStep or QPStep, depending on the problem. If a variable is made superbasic during the current iteration ( $+SBS > 0$ ), Step will be the step to the nearest bound. During Phase 2, the step can be greater than one only if the reduced Hessian is not positive definite.
nInf	is the number of infeasibilities after the present iteration. This number will not increase unless the iterations are in elastic mode.
SumInf	is the sum of infeasibilities after the present iteration, if $nInf > 0$ . The value usually decreases at each nonzero Step, but if it decreases by $2$ or more, SumInf may occasionally increase.
rgNorm	is the norm of the reduced-gradient vector at the start of the iteration. (It is the norm of the vector with elements $d_{j}$ for variables $j$ in the superbasic set.) During Phase 2 this norm will be approximately zero after a unit step. (The heading is not printed if the problem is linear.)
LPobjective or QPobjective	the QP objective function after the present iteration. In elastic mode, the heading is changed to Elastic QPobj. In either case, the value printed decreases monotonically.
+SBS	is the variable jq selected by the pricing operation to be added to the superbasic set.
-SBS	is the superbasic variable chosen to become nonbasic.
-BS	is the basis variable removed (if any) to become nonbasic.
Pivot	if column $a_{q}$ replaces the $r$ th column of the basis $B$ , Pivot is the $r$ th element of a vector $y$ satisfying $B y = a_{q}$ . Wherever possible, Step is chosen to avoid extremely small values of Pivot (since they cause the basis to be nearly singular). In rare cases, it may be necessary to increase the Pivot Tolerance to exclude very small elements of $y$ from consideration during the computation of Step.
L+U	is the number of nonzeros representing the basis factors $L$ and $U$ . Immediately after a basis factorization $B = L U$ , L+U is lenL+lenU, the number of subdiagonal elements in the columns of a lower triangular matrix and the number of diagonal and superdiagonal elements in the rows of an upper-triangular matrix. Further nonzeros are added to L when various columns of $B$ are later replaced. As columns of $B$ are replaced, the matrix $U$ is maintained explicitly (in sparse form). The value of L will steadily increase, whereas the value of U may fluctuate up or down. Thus the value of L+U may fluctuate up or down (in general, it will tend to increase).
ncp	is the number of compressions required to recover storage in the data structure for $U$ . This includes the number of compressions needed during the previous basis factorization.
nS	is the current number of superbasic variables. (The heading is not printed if the problem is linear.)
condHz	see [Major Iteration Log]. (The heading is not printed if the problem is linear.)

Crash Statistics

Major Print Level \geq 10

and System Information Yes has been specified, the following items are output to the Print file when

start = 0

and no Basis file is loaded. They refer to the number of columns that the Crash procedure selects during selected passes through

A

while searching for a triangular basis matrix.

Label	Description
Slacks	is the number of slacks selected initially.
Free cols	is the number of free columns in the basis.
Preferred	is the number of ‘preferred’ columns in the basis (i.e., $xstate [j - 1] = 3$ for some $j \leq n$ ).
Unit	is the number of unit columns in the basis.
Double	is the number of columns in the basis containing $2$ nonzeros.
Triangle	is the number of triangular columns in the basis.
Pad	is the number of slacks used to pad the basis (to make it a nonsingular triangle).

Basis Factorization Statistics

Major Print Level \geq 10

, the first seven items listed below are output to the Print file whenever the basis

B

or the rectangular matrix

B_{S} = {(B S)}^{T}

is factorized before solution of the next QP subproblem (see [Description of the Optional Parameters]).

Note that

B_{S}

may be factorized at the start of just some of the major iterations. It is immediately followed by a factorization of

B

itself.

Gaussian elimination is used to compute a sparse

L U

factorization of

B

B_{S}

, where

P L P^{T}

and

P U Q

are lower and upper triangular matrices, for some permutation matrices

P

and

Q

. Stability is ensured as described under optional parameter LU Factor Tolerance.

Minor Print Level \geq 10

, the same items are printed during the QP solution whenever the current

B

is factorized. In addition, if System Information Yes has been specified, the entries from Elems onwards are also printed.

Label

Description

Factor

the number of factorizations since the start of the run.

Demand

a code giving the reason for the present factorization, as follows:

Code	Meaning
0	First $L U$ factorization.
1	The number of updates reached the Factorization Frequency.
2	The nonzeros in the updated factors have increased significantly.
7	Not enough storage to update factors.
10	Row residuals are too large (see the description of the optional parameter Check Frequency).
11	Ill-conditioning has caused inconsistent results.

Itn

is the current minor iteration number.

Nonlin

is the number of nonlinear variables in the current basis

B

Linear

is the number of linear variables in

B

Slacks

is the number of slack variables in

B

B, BR, BS or BT factorize

is the type of

L U

factorization.

B	periodic factorization of the basis $B$ .
BR	more careful rank-revealing factorization of $B$ using threshold rook pivoting. This occurs mainly at the start, if the first basis factors seem singular or ill-conditioned. Followed by a normal B factorize.
BS	$B_{S}$ is factorized to choose a well-conditioned $B$ from the current $(B S)$ . Followed by a normal B factorize.
BT	same as BS except the current $B$ is tried first and accepted if it appears to be not much more ill-conditioned than after the previous BS factorize.

is the number of rows in

B

B_{S}

is the number of columns in

B

B_{S}

. Preceded by ‘=’ or ‘>’ respectively.

Elems

is the number of nonzero elements in

B

B_{S}

Amax

is the largest nonzero in

B

B_{S}

Density

is the percentage nonzero density of

B

B_{S}

Merit/MerRP/MerCP

Merit is the average Markowitz merit count for the elements chosen to be the diagonals of

P U Q

. Each merit count is defined to be

(c - 1) (r - 1)

where

c

and

r

are the number of nonzeros in the column and row containing the element at the time it is selected to be the next diagonal. Merit is the average of n such quantities. It gives an indication of how much work was required to preserve sparsity during the factorization. If LU Complete Pivoting or LU Rook Pivoting has been selected, this heading is changed to MerCP, respectively MerRP.

lenL

is the number of nonzeros in

L

L+U

is the number of nonzeros representing the basis factors

L

and

U

. Immediately after a basis factorization

B = L U

, this is lenL+lenU, the number of subdiagonal elements in the columns of a lower triangular matrix and the number of diagonal and superdiagonal elements in the rows of an upper-triangular matrix. Further nonzeros are added to L when various columns of

B

are later replaced. As columns of

B

are replaced, the matrix

U

is maintained explicitly (in sparse form). The value of L will steadily increase, whereas the value of U may fluctuate up or down. Thus the value of L+U may fluctuate up or down (in general, it will tend to increase).

Cmpressns

is the number of times the data structure holding the partially factored matrix needed to be compressed to recover unused storage. Ideally this number should be zero. If it is more than

3

4

, the amount of workspace available to e04vh should be increased for efficiency.

Incres

is the percentage increase in the number of nonzeros in

L

and

U

relative to the number of nonzeros in

B

B_{S}

Utri

is the number of triangular rows of

B

B_{S}

at the top of

U

lenU

the number of nonzeros in

U

, including its diagonals.

Ltol

is the largest subdiagonal element allowed in

L

. This is the specified LU Factor Tolerance or a smaller value that is currently being used for greater stability.

Umax

the maximum nonzero element in

U

Ugrwth

is the ratio

Umax / Amax

, which ideally should not be substantially larger than

10.0

100.0

. If it is orders of magnitude larger, it may be advisable to reduce the LU Factor Tolerance to

5.0

4.0

3.0

2.0

, say (but bigger than

1.0

As long as Lmax is not large (say

5.0

or less),

\max (Amax, Umax) / DUmin

gives an estimate of the condition number

B

. If this is extremely large, the basis is nearly singular. Slacks are used to replace suspect columns of

B

and the modified basis is refactored.

Ltri

is the number of triangular columns of

B

B_{S}

at the left of

L

dense1

is the number of columns remaining when the density of the basis matrix being factorized reached

0.3

Lmax

is the actual maximum subdiagonal element in

L

(bounded by Ltol).

Akmax

is the largest nonzero generated at any stage of the

L U

factorization. (Values much larger than Amax indicate instability.) Akmax is not printed if LU Partial Pivoting is selected.

Agrwth

is the ratio

Akmax / Amax

. Values much larger than

100

(say) indicate instability. Agrwth is not printed if LU Partial Pivoting is selected.

bump

is the size of the block to be factorized nontrivially after the triangular rows and columns of

B

B_{S}

have been removed.

dense2

is the number of columns remaining when the density of the basis matrix being factorized reached

0.6

. (The Markowitz pivot strategy searches fewer columns at that stage.)

DUmax

is the largest diagonal of

P U Q

DUmin

is the smallest diagonal of

P U Q

condU

the ratio

DUmax / DUmin

, which estimates the condition number of

U

(and of

B

if Ltol is less than

5.0

, say).

The Solution File

At the end of a run, the final solution may be output as a Solution file, according to Solution File. Some header information appears first to identify the problem and the final state of the optimization procedure. A ROWS section and a COLUMNS section then follow, giving one line of information for each row and column. The format used is similar to certain commercial systems, though there is no industry standard.

In general, numerical values are output with format f16.5. The maximum record length is

111

characters, including the first (carriage-control) character.

To reduce clutter, a full stop (.) is printed for any numerical value that is exactly zero. The values

\pm 1

are also printed specially as

1.0

and

- 1.0

. Infinite bounds (

\pm 10^{20}

or larger) are printed as None.

A Solution file is intended to be read from disk by a self-contained program that extracts and saves certain values as required for possible further computation. Typically, the first

14

records would be ignored. Each subsequent record may be read using adapted to suit the occasion. The end of the ROWS section is marked by a record that starts with a

1

and is otherwise blank. If this and the next

4

records are skipped, the COLUMNS section can then be read under the same format. (There should be no need for backspace statements.)

A full description of the ROWS section and the COLUMNS section is given in [The ROWS section] and [The COLUMNS section].

The Summary File

Summary File > 0

, the following information is output to the unit number associated with Summary File. (It is a brief summary of the output directed to unit Print File):

–	the optional parameters supplied via the option setting methods, if any;
–	the Basis file loaded, if any;
–	a brief major iteration log (see [Major Iteration Log]);
–	a brief minor iteration log (see [Minor Iteration Log]);
–	the exit condition, ifail;
–	a summary of the final iterate.

Syntax

Parameters

Description

References

Error Indicators and Warnings

Accuracy

Parallelism and Performance

Further Comments

The Final Output

The ROWS section

The COLUMNS section

Example

Algorithmic Details

Constraints and Slack Variables

Major Iterations

Minor Iterations

The Merit Function

Treatment of Constraint Infeasibilities

Description of Monitoring Information

Major Iteration Log

Minor Iteration Log

Crash Statistics

Basis Factorization Statistics

The Solution File

The Summary File

See Also