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 ) 
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 ) 
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 
Parameters
 start
 Type: System..::..Int32On entry: indicates how a starting point is to be obtained.
 ${\mathbf{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.
 ${\mathbf{start}}=1$
 Is the same as ${\mathbf{start}}=0$ but is more meaningful when a Basis file is given.
 ${\mathbf{start}}=2$
 Means that xstate and fstate define a valid starting point (probably from an earlier call, though not necessarily).
Constraint: ${\mathbf{start}}=0$, $1$ or $2$.
 nf
 Type: System..::..Int32On entry: $\mathit{nf}$, the number of problem functions in $F\left(x\right)$, 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: ${\mathbf{nf}}>0$.
 n
 Type: System..::..Int32On entry: $n$, the number of variables.Constraint: ${\mathbf{n}}>0$.
 nxname
 Type: System..::..Int32On entry: the number of names provided in the array xnames.
 ${\mathbf{nxname}}=1$
 There are no names provided and generic names will be used in the output.
 ${\mathbf{nxname}}={\mathbf{n}}$
 Names for all variables must be provided and will be used in the output.
Constraint: ${\mathbf{nxname}}=1\text{ or}{\mathbf{n}}$.
 nfname
 Type: System..::..Int32On entry: the number of names provided in the array fnames.
 ${\mathbf{nfname}}=1$
 There are no names provided and generic names will be used in the output.
 ${\mathbf{nfname}}={\mathbf{nf}}$
 Names for all functions must be provided and will be used in the output.
 objadd
 Type: System..::..DoubleOn entry: is a constant that will be added to the objective row ${F}_{\mathrm{obj}}$ for printing purposes. Typically, ${\mathbf{objadd}}=\text{0.0E+0}$.
 objrow
 Type: System..::..Int32On entry: says which row of $F\left(x\right)$ is to act as the objective function. If there is no such row, set ${\mathbf{objrow}}=0$. Then e04vh will seek a feasible point such that ${l}_{F}\le F\left(x\right)\le {u}_{F}$ and ${l}_{x}\le x\le {u}_{x}$.Constraint: $1\le {\mathbf{objrow}}\le {\mathbf{nf}}\text{ or}{\mathbf{objrow}}=0$ (or a feasible point problem).
 prob
 Type: System..::..StringOn 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_USRFUNusrfun must define the nonlinear portion $f\left(x\right)$ of the problem functions $F\left(x\right)=f\left(x\right)+Ax$, along with its gradient elements ${G}_{ij}\left(x\right)=\frac{\partial {f}_{i}\left(x\right)}{\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, ${\mathbf{Derivative\; Option}}=1$.The elements of $G\left(x\right)$ are stored in the array ${\mathbf{g}}\left(1:\mathbf{usrfun\_leng}\right)$ 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}\left(x\right)}{\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 ${\mathbf{needg}}=0$ on entry to usrfun.); (iii) if not all gradients are known, you must specify ${\mathbf{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 ${\mathbf{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}\left(x\right)$ 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 ${\mathbf{status}}=1$ if some of the functions are undefined. The linesearch will shorten the step and try again; (ix) set ${\mathbf{status}}\le 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 $\left(i,j\right)$ and values ${A}_{ij}$ of the nonzero elements of the linear part $A$ of the function $F\left(x\right)=f\left(x\right)+Ax$.In particular, nea triples $\left({\mathbf{iafun}}\left[k1\right],{\mathbf{javar}}\left[k1\right],{\mathbf{a}}\left[k1\right]\right)$ define the row and column indices $i={\mathbf{iafun}}\left[k1\right]$ and $j={\mathbf{javar}}\left[k1\right]$ of the element ${A}_{ij}={\mathbf{a}}\left[k1\right]$.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 $\left(i,j\right)$ and values ${A}_{ij}$ of the nonzero elements of the linear part $A$ of the function $F\left(x\right)=f\left(x\right)+Ax$.In particular, nea triples $\left({\mathbf{iafun}}\left[k1\right],{\mathbf{javar}}\left[k1\right],{\mathbf{a}}\left[k1\right]\right)$ define the row and column indices $i={\mathbf{iafun}}\left[k1\right]$ and $j={\mathbf{javar}}\left[k1\right]$ of the element ${A}_{ij}={\mathbf{a}}\left[k1\right]$.The coordinates may define the elements of $A$ in any order.
 nea
 Type: System..::..Int32On entry: is the number of nonzero entries in $A$ such that $F\left(x\right)=f\left(x\right)+Ax$.Constraint: $0\le {\mathbf{nea}}\le {\mathbf{lena}}$.
 a
 Type: array<System..::..Double>[]()[][]An array of size [lena]On entry: define the coordinates $\left(i,j\right)$ and values ${A}_{ij}$ of the nonzero elements of the linear part $A$ of the function $F\left(x\right)=f\left(x\right)+Ax$.In particular, nea triples $\left({\mathbf{iafun}}\left[k1\right],{\mathbf{javar}}\left[k1\right],{\mathbf{a}}\left[k1\right]\right)$ define the row and column indices $i={\mathbf{iafun}}\left[k1\right]$ and $j={\mathbf{javar}}\left[k1\right]$ of the element ${A}_{ij}={\mathbf{a}}\left[k1\right]$.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 $\left(i,j\right)$ of the nonzero elements of $G$, the nonlinear part of the derivative $J\left(x\right)=G\left(x\right)+A$ of the function $F\left(x\right)=f\left(x\right)+Ax$. e04vj may be used to define these two arrays.
 jgvar
 Type: array<System..::..Int32>[]()[][]An array of size [leng]On entry: define the coordinates $\left(i,j\right)$ of the nonzero elements of $G$, the nonlinear part of the derivative $J\left(x\right)=G\left(x\right)+A$ of the function $F\left(x\right)=f\left(x\right)+Ax$. e04vj may be used to define these two arrays.
 neg
 Type: System..::..Int32On entry: the number of nonzero entries in $G$.Constraint: $0\le {\mathbf{neg}}\le {\mathbf{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 ${\left[{l}_{x}\right]}_{j}=\infty $, set ${\mathbf{xlow}}\left[j1\right]\le \mathit{bigbnd}$, where $\mathit{bigbnd}$ is the optional parameter Infinite Bound Size. To specify a nonexistent upper bound ${\left[{u}_{x}\right]}_{j}=\infty $, set ${\mathbf{xupp}}\left[j1\right]\ge \mathit{bigbnd}$.To fix the $j$th variable at ${x}_{j}=\beta $, where $\left\beta \right<\mathit{bigbnd}$, set ${\mathbf{xlow}}\left[j1\right]={\mathbf{xupp}}\left[j1\right]=\beta $.Constraint: ${\mathbf{xlow}}\left[\mathit{i}1\right]\le {\mathbf{xupp}}\left[\mathit{i}1\right]$, for $\mathit{i}=1,2,\dots ,{\mathbf{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 ${\left[{l}_{x}\right]}_{j}=\infty $, set ${\mathbf{xlow}}\left[j1\right]\le \mathit{bigbnd}$, where $\mathit{bigbnd}$ is the optional parameter Infinite Bound Size. To specify a nonexistent upper bound ${\left[{u}_{x}\right]}_{j}=\infty $, set ${\mathbf{xupp}}\left[j1\right]\ge \mathit{bigbnd}$.To fix the $j$th variable at ${x}_{j}=\beta $, where $\left\beta \right<\mathit{bigbnd}$, set ${\mathbf{xlow}}\left[j1\right]={\mathbf{xupp}}\left[j1\right]=\beta $.Constraint: ${\mathbf{xlow}}\left[\mathit{i}1\right]\le {\mathbf{xupp}}\left[\mathit{i}1\right]$, for $\mathit{i}=1,2,\dots ,{\mathbf{n}}$.
 xnames
 Type: array<System..::..String>[]()[][]An array of size [nxname]On entry: the optional names for the variables.If ${\mathbf{nxname}}=1$, xnames is not referenced and default names will be used for output.If ${\mathbf{nxname}}={\mathbf{n}}$, ${\mathbf{xnames}}\left[j1\right]$ 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\left(x\right)$.To specify a nonexistent lower bound ${\left[{l}_{F}\right]}_{i}=\infty $, set ${\mathbf{flow}}\left[i1\right]\le \mathit{bigbnd}$. For a nonexistent upper bound ${\left[{u}_{F}\right]}_{i}=\infty $, set ${\mathbf{fupp}}\left[i1\right]\ge \mathit{bigbnd}$.To make the $i$th constraint an equality at ${F}_{i}=\beta $, where $\left\beta \right<\mathit{bigbnd}$, set ${\mathbf{flow}}\left[i1\right]={\mathbf{fupp}}\left[i1\right]=\beta $.Constraint: ${\mathbf{flow}}\left[\mathit{i}1\right]\le {\mathbf{fupp}}\left[\mathit{i}1\right]$, for $\mathit{i}=1,2,\dots ,{\mathbf{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\left(x\right)$.To specify a nonexistent lower bound ${\left[{l}_{F}\right]}_{i}=\infty $, set ${\mathbf{flow}}\left[i1\right]\le \mathit{bigbnd}$. For a nonexistent upper bound ${\left[{u}_{F}\right]}_{i}=\infty $, set ${\mathbf{fupp}}\left[i1\right]\ge \mathit{bigbnd}$.To make the $i$th constraint an equality at ${F}_{i}=\beta $, where $\left\beta \right<\mathit{bigbnd}$, set ${\mathbf{flow}}\left[i1\right]={\mathbf{fupp}}\left[i1\right]=\beta $.Constraint: ${\mathbf{flow}}\left[\mathit{i}1\right]\le {\mathbf{fupp}}\left[\mathit{i}1\right]$, for $\mathit{i}=1,2,\dots ,{\mathbf{n}}$.
 fnames
 Type: array<System..::..String>[]()[][]An array of size [nfname]On entry: the optional names for the problem functions.If ${\mathbf{nfname}}=1$, fnames is not referenced and default names will be used for output.If ${\mathbf{nfname}}={\mathbf{nf}}$, ${\mathbf{fnames}}\left[i1\right]$ 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$.If ${\mathbf{start}}=0$ or $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 ${\mathbf{x}}\left[j1\right]=0.0$, ${\mathbf{xstate}}\left[j1\right]=0$, for all $j=1,2,\dots ,{\mathbf{n}}$. If you set ${\mathbf{x}}\left[j1\right]={\mathbf{xlow}}\left[j1\right]$ set ${\mathbf{xstate}}\left[j1\right]=4$; if you set ${\mathbf{x}}\left[j1\right]={\mathbf{xupp}}\left[j1\right]$ then set ${\mathbf{xstate}}\left[j1\right]=5$. In this case a Crash procedure is used to select an initial basis.If ${\mathbf{start}}=0$ or $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.On exit: the final state of the variables.
${\mathbf{xstate}}\left[j1\right]$ State of variable $j$ Usual value of ${\mathbf{x}}\left[j1\right]$ 0 nonbasic ${\mathbf{xlow}}\left[j1\right]$ 1 nonbasic ${\mathbf{xupp}}\left[j1\right]$ 2 superbasic Between ${\mathbf{xlow}}\left[j1\right]$ and ${\mathbf{xupp}}\left[j1\right]$ 3 basic Between ${\mathbf{xlow}}\left[j1\right]$ and ${\mathbf{xupp}}\left[j1\right]$ 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 ${\mathbf{x}}\left[j1\right]$ lies strictly between its bounds.If ${\mathbf{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\le {\mathbf{xstate}}\left[j1\right]\le 5\text{, for}j=1,2,\dots ,{\mathbf{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}\le x\le {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.
 fstate
 Type: array<System..::..Int32>[]()[][]An array of size [nf]On entry: the initial state for the problem functions $F$.If ${\mathbf{start}}=0$ or $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 ${\mathbf{f}}\left[i1\right]=0.0$, ${\mathbf{fstate}}\left[i1\right]=0$, for all $i=1,2,\dots ,{\mathbf{nf}}$. Less trivially, to say that the optimal value of function ${\mathbf{f}}\left[i1\right]$ will probably be equal to one of its bounds, set ${\mathbf{f}}\left[i1\right]={\mathbf{flow}}\left[i1\right]$ and ${\mathbf{fstate}}\left[i1\right]=4$ or ${\mathbf{f}}\left[i1\right]={\mathbf{fupp}}\left[i1\right]$ and ${\mathbf{fstate}}\left[i1\right]=5$ as appropriate. In this case a Crash procedure is used to select an initial basis.If ${\mathbf{start}}=0$ or $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.On exit: the final state of the variables. The elements of fstate have the following meaning:
${\mathbf{fstate}}\left[i1\right]$ State of the corresponding
slack variableUsual value of ${\mathbf{f}}\left[i1\right]$ 0 nonbasic ${\mathbf{flow}}\left[i1\right]$ 1 nonbasic ${\mathbf{fupp}}\left[i1\right]$ 2 superbasic Between ${\mathbf{flow}}\left[i1\right]$ and ${\mathbf{fupp}}\left[i1\right]$ 3 basic Between ${\mathbf{flow}}\left[i1\right]$ and ${\mathbf{fupp}}\left[i1\right]$ 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 ${\mathbf{f}}\left[i1\right]$ lies strictly between its bounds.If ${\mathbf{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\le {\mathbf{fstate}}\left[i1\right]\le 5\text{, for}i=1,2,\dots ,{\mathbf{nf}}$.
 fmul
 Type: array<System..::..Double>[]()[][]An array of size [nf]On entry: an estimate of $\gamma $, the vector of Lagrange multipliers (shadow prices) for the constraints ${l}_{F}\le F\left(x\right)\le {u}_{F}$. All nf components must be defined. If nothing is known about $\gamma $, set ${\mathbf{fmul}}\left[\mathit{i}1\right]=0.0$, for $\mathit{i}=1,2,\dots ,{\mathbf{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}\le F\left(x\right)\le {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\left(x\right)$ 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\left(x\right)$ 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..::..e04vhOptionsAn Object of type E04.e04vhOptions. Used to configure optional parameters to this method.
 ifail
 Type: System..::..Int32%On exit: ${\mathbf{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 largescale linear and quadratic programming and for linearly constrained optimization, as well as for general nonlinear programs of the form
where $x$ is an $n$vector of variables, $l$ and $u$ are constant lower and upper bounds, ${f}_{0}\left(x\right)$ is a smooth scalar objective function, ${A}_{L}$ is a sparse matrix, and $f\left(x\right)$ is a vector of smooth nonlinear constraint functions $\left\{{f}_{i}\left(x\right)\right\}$. The optional parameter Maximize specifies that ${f}_{0}\left(x\right)$ should be maximized instead of minimized.
$$\underset{x}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}\text{}{f}_{0}\left(x\right)\text{\hspace{1em} subject to}l\le \left(\begin{array}{c}x\\ f\left(x\right)\\ {A}_{L}x\end{array}\right)\le u\text{,}$$  (1) 
Ideally, the first derivatives (gradients) of ${f}_{0}\left(x\right)$ and ${f}_{i}\left(x\right)$ should be known and coded by you. If only some of the gradients are known, e04vh estimates the missing ones by finite differences.
If ${f}_{0}\left(x\right)$ is linear and $f\left(x\right)$ 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 limitedmemory quasiNewton approximations to the Hessian of the Lagrangian. The merit function for steplength 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
where the upper and lower bounds are constant, $F\left(x\right)$ is a vector of smooth linear and nonlinear constraint functions $\left\{{F}_{i}\left(x\right)\right\}$, and ${F}_{\mathrm{obj}}\left(x\right)$ 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).
$$\underset{x}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}{F}_{\mathrm{obj}}\left(x\right)\text{\hspace{1em} subject to}{l}_{x}\le x\le {u}_{x}\text{, \hspace{1em}}{l}_{F}\le F\left(x\right)\le {u}_{F}\text{,}$$  (2) 
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$ or $u$ should be set to special values that are treated as $\infty $ or $+\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\left(x\right)$ 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\left(x\right)$ is of the form
where ${f}_{i}\left(x\right)$ is a nonlinear function (possibly zero) and the elements ${A}_{ij}$ are constant. The $nf$ by $n$ Jacobian of $F\left(x\right)$ is the sum of two sparse matrices of the same size: ${F}^{\prime}\left(x\right)=G\left(x\right)+A$, where $G\left(x\right)={f}^{\prime}\left(x\right)$ and $A$ is the matrix with elements $\left\{{A}_{ij}\right\}$. The two matrices must be nonoverlapping in the sense that each element of the Jacobian ${F}^{\prime}\left(x\right)=G\left(x\right)+A$ comes from $G\left(x\right)$ or $A$, but not both. The element cannot be split between $G\left(x\right)$ and $A$.
$${F}_{i}\left(x\right)={f}_{i}\left(x\right)+\sum _{j=1}^{n}{A}_{ij}{x}_{j}\text{,}$$ 
For example, the function
can be written as
in which case
can be written as ${F}^{\prime}\left(x\right)={f}^{\prime}\left(x\right)+A=G\left(x\right)+A$, where
$$F\left(x\right)=\left(\begin{array}{c}3{x}_{1}+{e}^{{x}_{2}}{x}_{4}+{x}_{2}^{2}+4{x}_{4}{x}_{3}+{x}_{5}\\ {x}_{2}+{x}_{3}^{2}+\mathrm{sin}\u200a{x}_{4}3{x}_{5}\\ {x}_{1}{x}_{3}\end{array}\right)$$ 
$$F\left(x\right)=f\left(x\right)+Ax=\left(\begin{array}{c}{e}^{{x}_{2}}{x}_{4}+{x}_{2}^{2}+4{x}_{4}\\ {x}_{3}^{2}+\mathrm{sin}\u200a{x}_{4}\\ 0\end{array}\right)+\left(\begin{array}{c}3{x}_{1}{x}_{3}+{x}_{5}\\ {x}_{2}3{x}_{5}\\ {x}_{1}{x}_{3}\end{array}\right)\text{,}$$ 
$${F}^{\prime}\left(x\right)=\left(\begin{array}{ccccc}3& {e}^{{x}_{2}}{x}_{4}+2{x}_{2}& 1& {e}^{{x}_{2}}+4& \phantom{}1\\ 0& 1& 2{x}_{3}& \mathrm{cos}\u200a{x}_{4}& 3\\ 1& 0& 1& 0& \phantom{}0\end{array}\right)$$ 
$$G\left(x\right)=\left(\begin{array}{ccccc}0& {e}^{{x}_{2}}{x}_{4}+2{x}_{2}& 0& {e}^{{x}_{2}}+4& 0\\ 0& 0& 2{x}_{3}& \mathrm{cos}\u200a{x}_{4}& 0\\ 0& 0& 0& 0& 0\end{array}\right)\text{, \hspace{1em}}A=\left(\begin{array}{ccccc}3& 0& 1& 0& \phantom{}1\\ 0& 1& \phantom{}0& 0& 3\\ 1& 0& 1& 0& \phantom{}0\end{array}\right)\text{.}$$ 
Note: the element ${e}^{{x}_{2}}+4$ of ${F}^{\prime}\left(x\right)$ appears in $G\left(x\right)$ and is not split between $G\left(x\right)$ 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 $\left(i,j,{A}_{ij}\right)$ in the arrays iafun, javar and a. The sparsity pattern $G$ is entered as pairs $\left(i,j\right)$ in the arrays igfun and jgvar. The corresponding entries ${G}_{ij}$ (any that are known) are assigned to appropriate array elements ${\mathbf{g}}\left(k\right)$ 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 ${\mathbf{Derivative\; Option}}=0$ is specified and usrfun does not provide any gradients.
Throughout this document the symbol $\epsilon $ 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) Largescale 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 Largescale 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 largescale linear and quadratic programming Report NA 051 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 largescale linear nonlinear programming Report NA 052 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 862 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) Largescale 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 8320R 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:
 ${\mathbf{ifail}}=1$
 ${\mathbf{ifail}}=2$
 An input parameter is invalid. The output message provides more details of the invalid argument.
 ${\mathbf{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.
 ${\mathbf{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 $Axs=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 socalled ‘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.
 ${\mathbf{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.
 ${\mathbf{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.
 ${\mathbf{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 ${\mathbf{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}$
4. An $LU$ 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 $Axs=0$ sufficiently well. This probably means that the current basis is very illconditioned. If there are some linear constraints and variables, try ${\mathbf{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.
 ${\mathbf{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 ${\mathbf{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 forwarddifference estimate is therefore wrong), you can specify ${\mathbf{Verify\; Level}}=0$ to prevent individual elements from being checked. However, the optimization procedure may have difficulty.
 ${\mathbf{ifail}}=9$
 Undefined usersupplied function.You have indicated that the problem functions are undefined by assigning the value ${\mathbf{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.
 ${\mathbf{ifail}}=10$
 User requested termination.
 ${\mathbf{ifail}}=11$
 Internal memory allocation failed when attempting to obtain the required workspace. Please contact NAG.
 ${\mathbf{ifail}}=12$
 Internal memory allocation was insufficient. Please contact NAG.
 ${\mathbf{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.
 ${\mathbf{ifail}}=14$
 An unexpected error has occurred. Set the optional parameter Print File and examine the output carefully for further information.
 ${\mathbf{ifail}}=9000$
 An error occured, see message report.
 ${\mathbf{ifail}}=1000$
 The array lengths are not the same for arrays
 ${\mathbf{ifail}}=1000$
 The array lengths are not the same for arrays
 ${\mathbf{ifail}}=8000$
 Negative dimension for array $\u2329\mathit{\text{value}}\u232a$
 ${\mathbf{ifail}}=6000$
 Invalid Parameters $\u2329\mathit{\text{value}}\u232a$
Accuracy
If the value of the optional parameter Major Optimality Tolerance is set to ${10}^{d}$ ($\text{default value}=\sqrt{\epsilon}$) and ${\mathbf{ifail}}={0}$ on exit, then the final value of $f\left(x\right)$ 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
If ${\mathbf{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\le {A}_{L}x\le u$. The $i$th constraint is therefore of the form
where ${\nu}_{i}$ is the $i$th row of ${A}_{L}$.
$$\alpha \le {\nu}_{i}x\le \beta \text{,}$$ 
Internally, the constraints take the form ${A}_{L}xs=0$, where $s$ is the set of slack variables (which satisfy the bounds $l\le s\le 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 $\alpha \le {f}_{i}\left(x\right)+{\nu}_{i}x\le \beta $ are treated similarly, except that the row activity and degree of infeasibility are computed directly from ${f}_{i}\left(x\right)+{\nu}_{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 $\alpha $ and $\beta $. The various states possible are as follows:
A key is sometimes printed before State.
Note that unless the optional parameter ${\mathbf{Scale\; Option}}=0$ is specified, the tests for assigning a key are applied to the variables of the scaled problem.


Activity  is the value of ${\nu}_{i}x$ (or ${f}_{i}\left(x\right)+{\nu}_{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  is $\alpha $, the lower bound on the row.  
Upper Limit  is $\beta $, the upper bound on the row.  
Dual Activity  is the value of the dual variable ${\pi}_{i}$ (the Lagrange multiplier for the $i$th constraint). The full vector $\pi $ always satisfies ${B}^{\mathrm{T}}\pi ={g}_{B}$, where $B$ is the current basis matrix and ${g}_{B}$ contains the associated gradients for the current objective function. For FP problems, ${\pi}_{i}$ is set to zero.  
i  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 $\alpha \le {x}_{j}\le \beta $. 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 $\alpha $ and $\beta $. The various states possible are as follows:
A key is sometimes printed before State.
Note that unless the optional parameter ${\mathbf{Scale\; Option}}=0$ is specified, the tests for assigning a key are applied to the variables of the scaled problem.


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 ${\mathbf{xlow}}\left[j1\right]\le \mathit{infbnd}$.  
Upper Limit  is the upper bound specified for the variable. None indicates that ${\mathbf{xupp}}\left[j1\right]\ge \mathit{infbnd}$.  
Reduced Gradnt  is the value of the reduced gradient ${d}_{j}={g}_{j}{\pi}^{\mathrm{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
subject to the bounds
to the nonlinear constraints
and to the linear constraints
The initial point, which is infeasible, is
and $f\left({x}_{0}\right)=0$.
$$f\left(x\right)={10}^{6}{x}_{3}^{3}+\frac{2}{3}\times {10}^{6}{x}_{4}^{3}+3{x}_{3}+2{x}_{4}$$ 
$$\begin{array}{c}0.55\le {x}_{1}\le 0.55\text{,}\\ 0.55\le {x}_{2}\le 0.55\text{,}\\ \phantom{0.0}0\le {x}_{3}\le 1200\text{,}\\ \phantom{0.0}0\le {x}_{4}\le 1200\text{,}\end{array}$$ 
$$\begin{array}{ccc}\hfill 1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{2}0.25\right){x}_{3}& \hfill =\hfill & \hfill 894.8\text{,}\\ \hfill 1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{1}{x}_{2}0.25\right){x}_{4}& \hfill =\hfill & \hfill 894.8\text{,}\\ \hfill 1000\mathrm{sin}\left({x}_{2}0.25\right)+1000\mathrm{sin}\left({x}_{2}{x}_{1}0.25\right)& \hfill =\hfill & \hfill 1294.8\text{,}\end{array}$$ 
$$\begin{array}{c}{x}_{1}+{x}_{2}\ge 0.55\text{,}\\ \phantom{}{x}_{1}{x}_{2}\ge 0.55\text{.}\end{array}$$ 
$${x}_{0}={\left(\begin{array}{cccc}0\text{,}& 0\text{,}& 0\text{,}& 0\end{array}\right)}^{\mathrm{T}}\text{,}$$ 
The optimal solution (to five figures) is
and $f\left({x}^{*}\right)=5126.4$. All the nonlinear constraints are active at the solution.
$${x}^{*}={\left(0.11887,0.39623,679.94,1026.0\right)}^{\mathrm{T}}\text{,}$$ 
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
where
and
$$F=\left(\begin{array}{c}1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{2}0.25\right){x}_{3}\\ 1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{1}{x}_{2}0.25\right){x}_{4}\\ 1000\mathrm{sin}\left({x}_{2}0.25\right)+1000\mathrm{sin}\left({x}_{2}{x}_{1}0.25\right)\\ {x}_{1}+{x}_{2}\\ \phantom{}{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{array}\right)=f+{A}_{L}x$$ 
$$f=\left(\begin{array}{c}1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{2}0.25\right)\\ 1000\mathrm{sin}\left({x}_{1}0.25\right)+1000\mathrm{sin}\left({x}_{1}{x}_{2}0.25\right)\\ 1000\mathrm{sin}\left({x}_{2}0.25\right)+1000\mathrm{sin}\left({x}_{2}{x}_{1}0.25\right)\\ 0\\ 0\\ {10}^{6}{x}_{3}^{3}+\frac{2}{3}\times {10}^{6}{x}_{4}^{3}\end{array}\right)$$ 
$${A}_{L}=\left(\begin{array}{cccc}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{array}\right)\text{.}$$ 
The nonzero elements of ${A}_{L}$ need to be stored in the triples $\left({\mathbf{iafun}}\left[k1\right],{\mathbf{javar}}\left[k1\right],{\mathbf{a}}\left[k1\right]\right)$ in any order. For example
$k$  $1$  $2$  $3$  $4$  $5$  $6$  $7$  $8$ 
${\mathbf{iafun}}\left[k1\right]$  $1$  $2$  $4$  $4$  $5$  $5$  $6$  $6$ 
${\mathbf{javar}}\left[k1\right]$  $3$  $4$  $1$  $2$  $1$  $2$  $3$  $4$ 
${\mathbf{a}}\left[k1\right]$  $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}$ or ${f}_{5}$ as there is no nonlinear part in ${F}_{4}$ or ${F}_{5}$.
The nonzero entries of the Jacobian of $f$ are
$$\begin{array}{cc}\frac{\partial {f}_{1}}{\partial {x}_{1}}& =1000\mathrm{cos}\left({x}_{1}0.25\right)\\ \frac{\partial {f}_{1}}{\partial {x}_{2}}& =1000\mathrm{cos}\left({x}_{2}0.25\right)\\ \frac{\partial {f}_{2}}{\partial {x}_{1}}& =1000\mathrm{cos}\left({x}_{1}0.25\right)+1000\mathrm{cos}\left({x}_{1}{x}_{2}0.25\right)\\ \frac{\partial {f}_{2}}{\partial {x}_{2}}& =1000\mathrm{cos}\left({x}_{1}{x}_{2}0.25\right)\\ \frac{\partial {f}_{3}}{\partial {x}_{1}}& =1000\mathrm{cos}\left({x}_{2}{x}_{1}0.25\right)\\ \frac{\partial {f}_{3}}{\partial {x}_{2}}& =1000\mathrm{cos}\left({x}_{2}0.25\right)+1000\mathrm{cos}\left({x}_{2}{x}_{1}0.25\right)\\ \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{array}$$ 
So the arrays igfun and jgvar must contain:
and usrfun must return in ${\mathbf{g}}\left[k1\right]$ the value of $\frac{\partial {f}_{i}}{\partial {x}_{j}}$, where $i={\mathbf{igfun}}\left[k1\right]$ and $j={\mathbf{jgvar}}\left[k1\right]$.
$k$  $1$  $2$  $3$  $4$  $5$  $6$  $7$  $8$ 
${\mathbf{igfun}}\left[k1\right]$  $1$  $1$  $2$  $2$  $3$  $3$  $6$  $6$ 
${\mathbf{jgvar}}\left[k1\right]$  $1$  $2$  $1$  $2$  $1$  $2$  $3$  $4$ 
Example program (C#): e04vhe.cs
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\left(x\right)$ 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={\left({s}_{1},{s}_{2},\dots ,{s}_{m}\right)}^{\mathrm{T}}$. For example, the linear constraint $5\le 2{x}_{1}+3{x}_{2}\le \infty $ is replaced by $2{x}_{1}+3{x}_{2}{s}_{1}=0$ together with the bounded slack $5\le {s}_{1}\le \infty $. The minimization problem (1) can therefore be written in the equivalent form
$$\underset{x,s}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}{f}_{0}\left(x\right)\text{\hspace{1em}subject to}\left(\begin{array}{c}f\left(x\right)\\ {A}_{L}x\end{array}\right)s=0,\hspace{1em}l\le \left(\begin{array}{c}x\\ s\end{array}\right)\le u\text{.}$$  (3) 
The general constraints become the equalities $f\left(x\right){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 $\left\{{x}_{k}\right\}$ that satisfy the linear constraints and converge to a point that satisfies the nonlinear constraints and the firstorder 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
where ${f}^{\prime}\left({x}_{k}\right)$ denotes the Jacobian matrix, whose elements are the first derivatives of $f\left(x\right)$ evaluated at ${x}_{k}$. The QP constraints therefore comprise the $m$ linear constraints
where $x$ and $s$ are bounded above and below by $u$ and $l$ as before. If the $m$ by $n$ matrix $A$ and $m$vector $b$ are defined as
then the QP subproblem can be written as
where $q\left(x,{x}_{k}\right)$ is a quadratic approximation to a modified Lagrangian function (see Gill et al. (2002)). The matrix ${H}_{k}$ is a quasiNewton 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.
$$f\left({x}_{k}\right)+{f}^{\prime}\left({x}_{k}\right)\left(x{x}_{k}\right){s}_{N}=0\text{,}$$  (4) 
$$\begin{array}{ccc}\hfill {f}^{\prime}\left({x}_{k}\right)x{s}_{N}& \hfill =\hfill & f\left({x}_{k}\right)+{f}^{\prime}\left({x}_{k}\right){x}_{k}\text{,}\hfill \\ \hfill {A}_{L}x\phantom{{s}_{N}}{s}_{L}& \hfill =\hfill & 0\text{,}\hfill \end{array}$$  (5) 
$$A=\left(\begin{array}{c}{f}^{\prime}\left({x}_{k}\right)\\ {A}_{L}\end{array}\right)\text{\hspace{1em}and\hspace{1em}}b=\left(\begin{array}{c}f\left({x}_{k}\right)+{f}^{\prime}\left({x}_{k}\right){x}_{k}\\ 0\end{array}\right)\text{,}$$  (6) 
$$\underset{x,s}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}q\left(x,{x}_{k}\right)={g}_{k}^{\mathrm{T}}\left(x{x}_{k}\right)+\frac{1}{2}{\left(x{x}_{k}\right)}^{\mathrm{T}}{H}_{k}\left(x{x}_{k}\right)\text{\hspace{1em}subject to}Axs=b\text{,\hspace{1em}}l\le \left(\begin{array}{c}x\\ s\end{array}\right)\le u\text{,}$$  (7) 
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 reducedHessian activeset method implemented as a reducedgradient method. At each minor iteration, the constraints $Axs=b$ are partitioned into the form
where the basis matrix $B$ is square and nonsingular, and the matrices $S$ and $N$ are the remaining columns of $\left(\begin{array}{cc}A& I\end{array}\right)$. 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 $\left(x,s\right)$ continues to satisfy $Axs=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.
$$B{x}_{B}+S{x}_{S}+N{x}_{N}=b\text{,}$$  (8) 
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 $Axs=b$ are the dual variables $\pi $. Similarly, each variable in $\left(x,s\right)$ has an associated reduced gradient ${d}_{j}$. The reduced gradients for the variables $x$ are the quantities $g{A}^{\mathrm{T}}\pi $, where $g$ is the gradient of the QP objective, and the reduced gradients for the slacks are the dual variables $\pi $. The QP subproblem is optimal if ${d}_{j}\ge 0$ for all nonbasic variables at their lower bounds, ${d}_{j}\le 0$ for all nonbasic variables at their upper bounds, and ${d}_{j}=0$ for other variables, including superbasics. In practice, an approximate QP solution $\left({\hat{x}}_{k},{\hat{s}}_{k},{\hat{\pi}}_{k}\right)$ 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
where $D$ is a diagonal matrix of penalty parameters $\left({D}_{ii}\ge 0\right)$, and $\pi $ now refers to dual variables for the nonlinear constraints in (1). If $\left({x}_{k},{s}_{k},{\pi}_{k}\right)$ denotes the current solution estimate and $\left({\hat{x}}_{k},{\hat{s}}_{k},{\hat{\pi}}_{k}\right)$ denotes the QP solution, the linesearch determines a step ${\alpha}_{k}$ $\left(0<{\alpha}_{k}\le 1\right)$ such that the new point
gives a sufficient decrease in the merit function $\mathcal{M}$. When necessary, the penalties in $D$ are increased by the minimumnorm perturbation that ensures descent for $\mathcal{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)).
$$\mathcal{M}\left(x,s,\pi \right)={f}_{0}\left(x\right){\pi}^{\mathrm{T}}\left(f\left(x\right){s}_{N}\right)+\frac{1}{2}{\left(f\left(x\right){s}_{N}\right)}^{\mathrm{T}}D\left(f\left(x\right){s}_{N}\right)\text{,}$$  (9) 
$$\left(\begin{array}{c}{x}_{k+1}\\ {s}_{k+1}\\ {\pi}_{k+1}\end{array}\right)=\left(\begin{array}{c}{x}_{k}\\ {s}_{k}\\ {\pi}_{k}\end{array}\right)+{\alpha}_{k}\left(\begin{array}{c}{\hat{x}}_{k}{x}_{k}\\ {\hat{s}}_{k}{s}_{k}\\ {\hat{\pi}}_{k}{\pi}_{k}\end{array}\right)$$  (10) 
Treatment of Constraint Infeasibilities
e04vh makes explicit allowance for infeasible constraints. First, infeasible linear constraints are detected by solving the linear program
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 ${\ell}_{1}$norm of the linear constraint violations. In the linear programming literature, the approach is called elastic programming.)
$$\underset{x,v,w}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}{e}^{\mathrm{T}}\left(v+w\right)\text{\hspace{1em}subject to}l\le \left(\begin{array}{c}x\\ {A}_{L}xv+w\end{array}\right)\le u\text{,\hspace{1em}}v\ge 0\text{,\hspace{1em}}w\ge 0\text{,}$$  (11) 
The linear constraints are infeasible if the optimal solution of (11) has $v\ne 0$ or $w\ne 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 $\pi $ for the nonlinear constraints become large), e04vh enters ‘elastic’ mode and thereafter solves the problem
where $\gamma $ is a nonnegative parameter (the elastic weight), and ${f}_{0}\left(x\right)+\gamma {e}^{\mathrm{T}}\left(v+w\right)$ is called a composite objective (the ${\ell}_{1}$ penalty function for the nonlinear constraints).
$$\underset{x,v,w}{\mathrm{minimize}}\phantom{\rule{0.25em}{0ex}}{f}_{0}\left(x\right)+\gamma {e}^{\mathrm{T}}\left(v+w\right)\text{\hspace{1em}subject to}l\le \left(\begin{array}{c}x\\ f\left(x\right)v+w\\ {A}_{L}x\end{array}\right)\le u\text{,\hspace{1em}}v\ge 0\text{,\hspace{1em}}w\ge 0\text{,}$$  (12) 
The value of $\gamma $ may increase automatically by multiples of $10$ if the optimal $v$ and $w$ continue to be nonzero. If $\gamma $ 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 $\gamma $ 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 ${\mathbf{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 $\alpha $ taken along the current search direction $p$. The variables $x$ have just been changed to $x+\alpha p$. On reasonably wellbehaved 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}_{\mathrm{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}_{\mathrm{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=LU$ is computed at the start of the first minor iteration. At this stage, $\mathtt{L+U}=\mathtt{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 uppertriangular 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 $LU$ 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 $LU$ factorization of the rectangular matrix ${B}_{S}={\left(B\text{}S\right)}^{\mathrm{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}^{\mathrm{T}}R$, itself an estimate of ${Z}^{\mathrm{T}}HZ$, 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}^{\mathrm{T}}R$. condHz gives a rough indication of whether or not the optimization procedure is having difficulty. If $\epsilon $ is the relative machine precision being used, the SQP algorithm will make slow progress if condHz becomes as large as ${\epsilon}^{1/2}\approx {10}^{8}$, and will probably fail to find a better solution if condHz reaches ${\epsilon}^{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 resolve 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 ${\mathbf{status}}=1$ on exit from usrfun, indicating that the linesearch needed to be done with a smaller value of the step length $\alpha $. 
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 quasiNewton 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 illconditioned 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 selfscaled 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
If ${\mathbf{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}_{\mathtt{pp}}$ or ${I}_{\mathtt{pp}}$, the $\mathtt{pp}$th segments of the constraint matrix $\left(\begin{array}{cc}A& I\end{array}\right)$.
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}{\pi}^{\mathrm{T}}{a}_{j}$ for $j=\mathtt{jq}$, where ${g}_{j}$ is the gradient of the current objective function, $\pi $ is the vector of dual variables for the QP subproblem, and ${a}_{j}$ is the $j$th column of $\left(\begin{array}{cc}A& I\end{array}\right)$.
Note that the reduced cost is the $1$norm of the reducedgradient vector at the start of the iteration, just after the pricing procedure. 
LPstep or QPstep  is the step length $\alpha $ taken along the current search direction $p$. The variables $x$ have just been changed to $x+\alpha p$. Write Step to stand for LPStep or QPStep, depending on the problem. If a variable is made superbasic during the current iteration ($\mathtt{+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 $\mathtt{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 reducedgradient 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 $By={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=LU$, 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 uppertriangular 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
If ${\mathbf{Major\; Print\; Level}}\ge 10$ and System Information Yes has been specified, the following items are output to the Print file when ${\mathbf{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., ${\mathbf{xstate}}\left[j1\right]=3$ for some $j\le 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
If ${\mathbf{Major\; Print\; Level}}\ge 10$, the first seven items listed below are output to the Print file whenever the basis $B$ or the rectangular matrix ${B}_{S}={\left(B\text{}S\right)}^{\mathrm{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 $LU$ factorization of $B$ or ${B}_{S}$, where $PL{P}^{\mathrm{T}}$ and $PUQ$ are lower and upper triangular matrices, for some permutation matrices $P$ and $Q$. Stability is ensured as described under optional parameter LU Factor Tolerance.
If ${\mathbf{Minor\; Print\; Level}}\ge 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:


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 $LU$ factorization.


m  is the number of rows in $B$ or ${B}_{S}$.  
n  is the number of columns in $B$ or ${B}_{S}$. Preceded by ‘=’ or ‘>’ respectively.  
Elems  is the number of nonzero elements in $B$ or ${B}_{S}$.  
Amax  is the largest nonzero in $B$ or ${B}_{S}$.  
Density  is the percentage nonzero density of $B$ or ${B}_{S}$.  
Merit/MerRP/MerCP  Merit is the average Markowitz merit count for the elements chosen to be the diagonals of $PUQ$. Each merit count is defined to be $\left(c1\right)\left(r1\right)$ 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=LU$, 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 uppertriangular 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$ or $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$ or ${B}_{S}$.  
Utri  is the number of triangular rows of $B$ or ${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 $\mathtt{Umax}/\mathtt{Amax}$, which ideally should not be substantially larger than $10.0$ or $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$ or $2.0$, say (but bigger than $1.0$).
As long as Lmax is not large (say $5.0$ or less), $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(\mathtt{Amax},\mathtt{Umax}\right)/\mathtt{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$ or ${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 $LU$ factorization. (Values much larger than Amax indicate instability.) Akmax is not printed if LU Partial Pivoting is selected.  
Agrwth  is the ratio $\mathtt{Akmax}/\mathtt{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$ or ${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 $PUQ$.  
DUmin  is the smallest diagonal of $PUQ$.  
condU  the ratio $\mathtt{DUmax}/\mathtt{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 (carriagecontrol) 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 selfcontained 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
If ${\mathbf{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. 