# NAG Library Routine Document

## 1Purpose

e05usf is designed to find the global minimum of an arbitrary smooth sum of squares function subject to constraints (which may include simple bounds on the variables, linear constraints and smooth nonlinear constraints) by generating a number of different starting points and performing a local search from each using sequential quadratic programming.

## 2Specification

Fortran Interface
 Subroutine e05usf ( m, n, a, lda, bl, bu, y, npts, x, ldx, nb, objf, f, fjac, iter, c, ldc, cjac, opts, info,
 Integer, Intent (In) :: m, n, nclin, ncnln, lda, npts, ldx, nb, ldfjac, sdfjac, ldc, ldcjac, sdcjac, ldclda, listat Integer, Intent (Inout) :: istate(listat,*), iopts(740), iuser(*), ifail Integer, Intent (Out) :: iter(nb), info(nb) Real (Kind=nag_wp), Intent (In) :: a(lda,*), bl(n+nclin+ncnln), bu(n+nclin+ncnln), y(m) Real (Kind=nag_wp), Intent (Inout) :: x(ldx,*), f(m,*), fjac(ldfjac,sdfjac,*), c(ldc,*), cjac(ldcjac,sdcjac,*), clamda(ldclda,*), opts(485), ruser(*) Real (Kind=nag_wp), Intent (Out) :: objf(nb) Logical, Intent (In) :: repeat1 External :: confun, objfun, start
#include <nagmk26.h>
 void e05usf_ (const Integer *m, const Integer *n, const Integer *nclin, const Integer *ncnln, const double a[], const Integer *lda, const double bl[], const double bu[], const double y[], void (NAG_CALL *confun)(Integer *mode, const Integer *ncnln, const Integer *n, const Integer *ldcjsl, const Integer needc[], const double x[], double c[], double cjsl[], const Integer *nstate, Integer iuser[], double ruser[]),void (NAG_CALL *objfun)(Integer *mode, const Integer *m, const Integer *n, const Integer *ldfjsl, const Integer *needfi, const double x[], double f[], double fjsl[], const Integer *nstate, Integer iuser[], double ruser[]),const Integer *npts, double x[], const Integer *ldx, void (NAG_CALL *start)(const Integer *npts, double quas[], const Integer *n, const logical *repeat1, const double bl[], const double bu[], Integer iuser[], double ruser[], Integer *mode),const logical *repeat1, const Integer *nb, double objf[], double f[], double fjac[], const Integer *ldfjac, const Integer *sdfjac, Integer iter[], double c[], const Integer *ldc, double cjac[], const Integer *ldcjac, const Integer *sdcjac, double clamda[], const Integer *ldclda, Integer istate[], const Integer *listat, Integer iopts[], double opts[], Integer iuser[], double ruser[], Integer info[], Integer *ifail)
Before calling e05usf, the optional parameter arrays iopts and opts must be initialized for use with e05usf by calling e05zkf with optstr set to ‘Initialize = e05usf’. Optional parameters may subsequently be specified by calling e05zkf before the call to e05usf.

## 3Description

The local minimization method is e04usa. The problem is assumed to be stated in the following form:
 $minimize x∈Rn ​ ​ Fx = 12 ∑ i=1 m yi- fi x 2 subject to l≤ x ALx cx ≤u,$ (1)
where $F\left(x\right)$ (the objective function) is a nonlinear function which can be represented as the sum of squares of $m$ subfunctions $\left({y}_{1}-{f}_{1}\left(x\right)\right),\left({y}_{2}-{f}_{2}\left(x\right)\right),\dots ,\left({y}_{m}-{f}_{m}\left(x\right)\right)$, the ${y}_{i}$ are constant, ${A}_{L}$ is an ${n}_{L}$ by $n$ constant linear constraint matrix, and $c\left(x\right)$ is an ${n}_{N}$ element vector of nonlinear constraint functions. (The matrix ${A}_{L}$ and the vector $c\left(x\right)$ may be empty.) The objective function and the constraint functions are assumed to be smooth, i.e., at least twice-continuously differentiable. (This routine will usually solve (1) if any isolated discontinuities are away from the solution.)
e05usf solves a user-specified number of local optimization problems with different starting points. You may specify the starting points via the subroutine start. If a random number generator is used to generate the starting points then the argument repeat1 allows you to specify whether a repeatable set of points are generated or whether different starting points are generated on different calls. The resulting local minima are ordered and the best nb results returned in order of ascending values of the resulting objective function values at the minima. Thus the value returned in position $1$ will be the best result obtained. If a sufficiently high number of different points are chosen then this is likely to be the global minimum.
Gill P E, Murray W and Wright M H (1981) Practical Optimization Academic Press
Hock W and Schittkowski K (1981) Test Examples for Nonlinear Programming Codes. Lecture Notes in Economics and Mathematical Systems 187 Springer–Verlag

## 5Arguments

1:     $\mathbf{m}$ – IntegerInput
On entry: $m$, the number of subfunctions associated with $F\left(x\right)$.
Constraint: ${\mathbf{m}}>0$.
2:     $\mathbf{n}$ – IntegerInput
On entry: $n$, the number of variables.
Constraint: ${\mathbf{n}}>0$.
3:     $\mathbf{nclin}$ – IntegerInput
On entry: ${n}_{L}$, the number of general linear constraints.
Constraint: ${\mathbf{nclin}}\ge 0$.
4:     $\mathbf{ncnln}$ – IntegerInput
On entry: ${n}_{N}$, the number of nonlinear constraints.
Constraint: ${\mathbf{ncnln}}\ge 0$.
5:     $\mathbf{a}\left({\mathbf{lda}},*\right)$ – Real (Kind=nag_wp) arrayInput
Note: the second dimension of the array a must be at least ${\mathbf{n}}$ if ${\mathbf{nclin}}>0$, and at least $1$ otherwise.
On entry: the matrix ${A}_{L}$ of general linear constraints in (1). That is, the $\mathit{i}$th row contains the coefficients of the $\mathit{i}$th general linear constraint, for $\mathit{i}=1,2,\dots ,{\mathbf{nclin}}$.
If ${\mathbf{nclin}}=0$, the array a is not referenced.
6:     $\mathbf{lda}$ – IntegerInput
On entry: the first dimension of the array a as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{lda}}\ge {\mathbf{nclin}}$.
7:     $\mathbf{bl}\left({\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}\right)$ – Real (Kind=nag_wp) arrayInput
8:     $\mathbf{bu}\left({\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: bl must contain the lower bounds and bu the upper bounds for all the constraints in the following order. The first $n$ elements of each array must contain the bounds on the variables, the next ${n}_{L}$ elements the bounds for the general linear constraints (if any) and the next ${n}_{N}$ elements the bounds for the general nonlinear constraints (if any). To specify a nonexistent lower bound (i.e., ${l}_{j}=-\infty$), set ${\mathbf{bl}}\left(j\right)\le -\mathit{bigbnd}$, and to specify a nonexistent upper bound (i.e., ${u}_{j}=+\infty$), set ${\mathbf{bu}}\left(j\right)\ge \mathit{bigbnd}$; the default value of $\mathit{bigbnd}$ is ${10}^{20}$, but this may be changed by the optional parameter Infinite Bound Size. To specify the $j$th constraint as an equality, set ${\mathbf{bl}}\left(j\right)={\mathbf{bu}}\left(j\right)=\beta$, say, where $\left|\beta \right|<\mathit{bigbnd}$.
Constraints:
• ${\mathbf{bl}}\left(\mathit{j}\right)\le {\mathbf{bu}}\left(\mathit{j}\right)$, for $\mathit{j}=1,2,\dots ,{\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}$;
• if ${\mathbf{bl}}\left(j\right)={\mathbf{bu}}\left(j\right)=\beta$, $\left|\beta \right|<\mathit{bigbnd}$.
9:     $\mathbf{y}\left({\mathbf{m}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: the coefficients of the constant vector $y$ of the objective function.
10:   $\mathbf{confun}$ – Subroutine, supplied by the NAG Library or the user.External Procedure
confun must calculate the vector $c\left(x\right)$ of nonlinear constraint functions and (optionally) its Jacobian ($\text{}=\frac{\partial c}{\partial x}$) for a specified $n$-element vector $x$. If there are no nonlinear constraints (i.e., ${\mathbf{ncnln}}=0$), confun will never be called by e05usf and confun may be the dummy routine e05udm. (e05udm is included in the NAG Library.) If there are nonlinear constraints, the first call to confun will occur before the first call to objfun.
The specification of confun is:
Fortran Interface
 Subroutine confun ( mode, n, x, c, cjsl,
 Integer, Intent (In) :: ncnln, n, ldcjsl, needc(ncnln), nstate Integer, Intent (Inout) :: mode, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(n) Real (Kind=nag_wp), Intent (Inout) :: cjsl(ldcjsl,*), ruser(*) Real (Kind=nag_wp), Intent (Out) :: c(ncnln)
#include <nagmk26.h>
 void confun (Integer *mode, const Integer *ncnln, const Integer *n, const Integer *ldcjsl, const Integer needc[], const double x[], double c[], double cjsl[], const Integer *nstate, Integer iuser[], double ruser[])
1:     $\mathbf{mode}$ – IntegerInput/Output
On entry: indicates which values must be assigned during each call of confun. Only the following values need be assigned, for each value of $i$ such that ${\mathbf{needc}}\left(i\right)>0$:
${\mathbf{mode}}=0$
${\mathbf{c}}\left(i\right)$, the $i$th nonlinear constraint.
${\mathbf{mode}}=1$
All available elements in the $i$th row of cjsl.
${\mathbf{mode}}=2$
${\mathbf{c}}\left(i\right)$ and all available elements in the $i$th row of cjsl.
On exit: may be set to a negative value if you wish to abandon the solution to the current local minimization problem. In this case e05usf will move to the next local minimization problem.
2:     $\mathbf{ncnln}$ – IntegerInput
On entry: ${n}_{N}$, the number of nonlinear constraints.
3:     $\mathbf{n}$ – IntegerInput
On entry: $n$, the number of variables.
4:     $\mathbf{ldcjsl}$ – IntegerInput
On entry: ldcjsl is the first dimension of the array cjsl.
5:     $\mathbf{needc}\left({\mathbf{ncnln}}\right)$ – Integer arrayInput
On entry: the indices of the elements of c and/or cjsl that must be evaluated by confun. If ${\mathbf{needc}}\left(i\right)>0$, ${\mathbf{c}}\left(i\right)$ and/or the available elements of the $i$th row of cjsl (see argument mode) must be evaluated at $x$.
6:     $\mathbf{x}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: $x$, the vector of variables at which the constraint functions and/or the available elements of the constraint Jacobian are to be evaluated.
7:     $\mathbf{c}\left({\mathbf{ncnln}}\right)$ – Real (Kind=nag_wp) arrayOutput
On exit: if ${\mathbf{needc}}\left(i\right)>0$ and ${\mathbf{mode}}=0$ or $2$, ${\mathbf{c}}\left(i\right)$ must contain the value of ${c}_{i}\left(x\right)$. The remaining elements of c, corresponding to the non-positive elements of needc, need not be set.
8:     $\mathbf{cjsl}\left({\mathbf{ldcjsl}},*\right)$ – Real (Kind=nag_wp) arrayInput/Output
cjsl may be regarded as a two-dimensional ‘slice’ of the three-dimensional array cjac of e05usf.
On entry: unless ${\mathbf{Derivative Level}}=2$ or $3$, the elements of cjsl are set to special values which enable e05usf to detect whether they are changed by confun.
On exit: if ${\mathbf{needc}}\left(i\right)>0$ and ${\mathbf{mode}}=1$ or $2$, the $i$th row of cjsl must contain the available elements of the vector $\nabla {c}_{i}$ given by
 $∇ci= ∂ci ∂x1 , ∂ci ∂x2 ,…, ∂ci ∂xn T,$
where $\frac{\partial {c}_{i}}{\partial {x}_{j}}$ is the partial derivative of the $i$th constraint with respect to the $j$th variable, evaluated at the point $x$. See also the argument nstate. The remaining rows of cjsl, corresponding to non-positive elements of needc, need not be set.
If all elements of the constraint Jacobian are known (i.e., ${\mathbf{Derivative Level}}=2$ or $3$; note the default is ${\mathbf{Derivative Level}}=3$), any constant elements may be assigned to cjsl one time only at the start of each local optimization. An element of cjsl that is not subsequently assigned in confun will retain its initial value throughout the local optimization. Constant elements may be loaded into cjsl during the first call to confun for the local optimization (signalled by the value ${\mathbf{nstate}}=1$). The ability to preload constants is useful when many Jacobian elements are identically zero, in which case cjsl may be initialized to zero and nonzero elements may be reset by confun.
Note that constant nonzero elements do affect the values of the constraints. Thus, if ${\mathbf{cjsl}}\left(i,j\right)$ is set to a constant value, it need not be reset in subsequent calls to confun, but the value ${\mathbf{cjsl}}\left(i,j\right)×{\mathbf{x}}\left(j\right)$ must nonetheless be added to ${\mathbf{c}}\left(i\right)$. For example, if ${\mathbf{cjsl}}\left(1,1\right)=2$ and ${\mathbf{cjsl}}\left(1,2\right)=-5$ then the term $2×{\mathbf{x}}\left(1\right)-5×{\mathbf{x}}\left(2\right)$ must be included in the definition of ${\mathbf{c}}\left(1\right)$.
It must be emphasized that, if ${\mathbf{Derivative Level}}=0$ or $1$, unassigned elements of cjsl are not treated as constant; they are estimated by finite differences, at nontrivial expense. If you do not supply a value for the optional parameter Difference Interval, an interval for each element of $x$ is computed automatically at the start of each local optimization. The automatic procedure can usually identify constant elements of cjsl, which are then computed once only by finite differences.
9:     $\mathbf{nstate}$ – IntegerInput
On entry: if ${\mathbf{nstate}}=1$ then e05usf is calling confun for the first time on the current local optimization problem. This argument setting allows you to save computation time if certain data must be read or calculated only once.
10:   $\mathbf{iuser}\left(*\right)$ – Integer arrayUser Workspace
11:   $\mathbf{ruser}\left(*\right)$ – Real (Kind=nag_wp) arrayUser Workspace
confun is called with the arguments iuser and ruser as supplied to e05usf. You should use the arrays iuser and ruser to supply information to confun.
confun must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e05usf is called. Arguments denoted as Input must not be changed by this procedure.
Note: confun should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e05usf. If your code inadvertently does return any NaNs or infinities, e05usf is likely to produce unexpected results.
confun should be tested separately before being used in conjunction with e05usf. See also the description of the optional parameter Verify.
11:   $\mathbf{objfun}$ – Subroutine, supplied by the user.External Procedure
objfun must calculate either the $i$th element of the vector $f\left(x\right)={\left({f}_{1}\left(x\right),{f}_{2}\left(x\right),\dots ,{f}_{m}\left(x\right)\right)}^{\mathrm{T}}$ or all $m$ elements of $f\left(x\right)$ and (optionally) its Jacobian ($\text{}=\frac{\partial f}{\partial x}$) for a specified $n$-element vector $x$.
The specification of objfun is:
Fortran Interface
 Subroutine objfun ( mode, m, n, x, f, fjsl,
 Integer, Intent (In) :: m, n, ldfjsl, needfi, nstate Integer, Intent (Inout) :: mode, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(n) Real (Kind=nag_wp), Intent (Inout) :: fjsl(ldfjsl,*), ruser(*) Real (Kind=nag_wp), Intent (Out) :: f(m)
#include <nagmk26.h>
 void objfun (Integer *mode, const Integer *m, const Integer *n, const Integer *ldfjsl, const Integer *needfi, const double x[], double f[], double fjsl[], const Integer *nstate, Integer iuser[], double ruser[])
1:     $\mathbf{mode}$ – IntegerInput/Output
On entry: indicates which values must be assigned during each call of objfun. Only the following values need be assigned:
${\mathbf{mode}}=0$ and ${\mathbf{needfi}}=i$, where $i>0$
${\mathbf{f}}\left(i\right)$.
${\mathbf{mode}}=0$ and ${\mathbf{needfi}}<0$
f.
${\mathbf{mode}}=1$ and ${\mathbf{needfi}}<0$
All available elements of fjsl.
${\mathbf{mode}}=2$ and ${\mathbf{needfi}}<0$
f and all available elements of fjsl.
On exit: may be set to a negative value if you wish to abandon the solution to the current local minimization problem. In this case e05usf will move to the next local minimization problem.
2:     $\mathbf{m}$ – IntegerInput
On entry: $m$, the number of subfunctions.
3:     $\mathbf{n}$ – IntegerInput
On entry: $n$, the number of variables.
4:     $\mathbf{ldfjsl}$ – IntegerInput
On entry: ldfjsl is the first dimension of the array fjsl.
5:     $\mathbf{needfi}$ – IntegerInput
On entry: if ${\mathbf{needfi}}=i>0$, only the $i$th element of $f\left(x\right)$ needs to be evaluated at $x$; the remaining elements need not be set. This can result in significant computational savings when $m\gg n$.
6:     $\mathbf{x}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: $x$, the vector of variables at which the objective function and/or all available elements of its gradient are to be evaluated.
7:     $\mathbf{f}\left({\mathbf{m}}\right)$ – Real (Kind=nag_wp) arrayOutput
On exit: if ${\mathbf{mode}}=0$ and ${\mathbf{needfi}}=i>0$, ${\mathbf{f}}\left(i\right)$ must contain the value of ${f}_{i}$ at $x$.
If ${\mathbf{mode}}=0$ or $2$ and ${\mathbf{needfi}}<0$, ${\mathbf{f}}\left(\mathit{i}\right)$ must contain the value of ${f}_{\mathit{i}}$ at $x$, for $\mathit{i}=1,2,\dots ,m$.
8:     $\mathbf{fjsl}\left({\mathbf{ldfjsl}},*\right)$ – Real (Kind=nag_wp) arrayInput/Output
fjsl may be regarded as a two-dimensional ‘slice’ of the three-dimensional array fjac of e05usf.
On entry: is set to a special value.
On exit: if ${\mathbf{mode}}=1$ or $2$ and ${\mathbf{needfi}}<0$, the $i$th row of fjsl must contain the available elements of the vector $\nabla {f}_{i}$ given by
 $∇fi = ∂fi/ ∂x1, ∂fi/ ∂x2 ,…, ∂fi/ ∂xn T ,$
evaluated at the point $x$. See also the argument nstate.
9:     $\mathbf{nstate}$ – IntegerInput
On entry: if ${\mathbf{nstate}}=1$ then e05usf is calling objfun for the first time on the current local optimization problem. This argument setting allows you to save computation time if certain data must be read or calculated only once.
10:   $\mathbf{iuser}\left(*\right)$ – Integer arrayUser Workspace
11:   $\mathbf{ruser}\left(*\right)$ – Real (Kind=nag_wp) arrayUser Workspace
objfun is called with the arguments iuser and ruser as supplied to e05usf. You should use the arrays iuser and ruser to supply information to objfun.
objfun must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e05usf is called. Arguments denoted as Input must not be changed by this procedure.
Note: objfun should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e05usf. If your code inadvertently does return any NaNs or infinities, e05usf is likely to produce unexpected results.
objfun should be tested separately before being used in conjunction with e05usf. See also the description of the optional parameter Verify.
12:   $\mathbf{npts}$ – IntegerInput
On entry: the number of different starting points to be generated and used. The more points used, the more likely that the best returned solution will be a global minimum.
Constraint: $1\le {\mathbf{nb}}\le {\mathbf{npts}}$.
13:   $\mathbf{x}\left({\mathbf{ldx}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the second dimension of the array x must be at least ${\mathbf{nb}}$.
On exit: ${\mathbf{x}}\left(\mathit{j},i\right)$ contains the final estimate of the $i$th solution, for $\mathit{j}=1,2,\dots ,{\mathbf{n}}$.
14:   $\mathbf{ldx}$ – IntegerInput
On entry: the first dimension of the array x as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{ldx}}\ge {\mathbf{n}}$.
15:   $\mathbf{start}$ – Subroutine, supplied by the NAG Library or the user.External Procedure
start must calculate the npts starting points to be used by the local optimizer. If you do not wish to write a routine specific to your problem then e05ucz may be used as the actual argument. e05ucz is supplied in the NAG Library and uses the NAG quasi-random number generators to distribute starting points uniformly across the domain. It is affected by the value of repeat1.
The specification of start is:
Fortran Interface
 Subroutine start ( npts, quas, n, bl, bu, mode)
 Integer, Intent (In) :: npts, n Integer, Intent (Inout) :: iuser(*), mode Real (Kind=nag_wp), Intent (In) :: bl(n), bu(n) Real (Kind=nag_wp), Intent (Inout) :: quas(n,npts), ruser(*) Logical, Intent (In) :: repeat1
#include <nagmk26.h>
 void start (const Integer *npts, double quas[], const Integer *n, const logical *repeat1, const double bl[], const double bu[], Integer iuser[], double ruser[], Integer *mode)
1:     $\mathbf{npts}$ – IntegerInput
On entry: indicates the number of starting points.
2:     $\mathbf{quas}\left({\mathbf{n}},{\mathbf{npts}}\right)$ – Real (Kind=nag_wp) arrayInput/Output
On entry: all elements of quas will have been set to zero, so only nonzero values need be set subsequently.
On exit: must contain the starting points for the npts local minimizations, i.e., ${\mathbf{quas}}\left(j,i\right)$ must contain the $j$th component of the $i$th starting point.
3:     $\mathbf{n}$ – IntegerInput
On entry: the number of variables.
4:     $\mathbf{repeat1}$ – LogicalInput
On entry: specifies whether a repeatable or non-repeatable sequence of points are to be generated.
5:     $\mathbf{bl}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: the lower bounds on the variables. These may be used to ensure that the starting points generated in some sense ‘cover’ the region, but there is no requirement that a starting point be feasible.
6:     $\mathbf{bu}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) arrayInput
On entry: the upper bounds on the variables. (See bl.)
7:     $\mathbf{iuser}\left(*\right)$ – Integer arrayUser Workspace
8:     $\mathbf{ruser}\left(*\right)$ – Real (Kind=nag_wp) arrayUser Workspace
start is called with the arguments iuser and ruser as supplied to e05usf. You should use the arrays iuser and ruser to supply information to start.
9:     $\mathbf{mode}$ – IntegerInput/Output
On entry: mode will contain $0$.
On exit: if you set mode to a negative value then e05usf will terminate immediately with ${\mathbf{ifail}}={\mathbf{9}}$.
start must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e05usf is called. Arguments denoted as Input must not be changed by this procedure.
Note: start should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e05usf. If your code inadvertently does return any NaNs or infinities, e05usf is likely to produce unexpected results.
16:   $\mathbf{repeat1}$ – LogicalInput
On entry: is passed as an argument to start and may be used to initialize a random number generator to a repeatable, or non-repeatable, sequence. See Section 9 for more detail.
17:   $\mathbf{nb}$ – IntegerInput
On entry: the number of solutions to be returned. The routine saves up to nb local minima ordered by increasing value of the final objective function. If the defining criterion for ‘best solution’ is only that the value of the objective function is as small as possible then nb should be set to $1$. However, if you want to look at other solutions that may have desirable properties then setting ${\mathbf{nb}}>1$ will produce nb local minima, ordered by increasing value of their objective functions at the minima.
Constraint: $1\le {\mathbf{nb}}\le {\mathbf{npts}}$.
18:   $\mathbf{objf}\left({\mathbf{nb}}\right)$ – Real (Kind=nag_wp) arrayOutput
On exit: ${\mathbf{objf}}\left(i\right)$ contains the value of the objective function at the final iterate for the $i$th solution.
19:   $\mathbf{f}\left({\mathbf{m}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the second dimension of the array f must be at least ${\mathbf{nb}}$.
On exit: ${\mathbf{f}}\left(\mathit{j},i\right)$ contains the value of the $\mathit{j}$th function ${f}_{j}$ at the final iterate, for $\mathit{j}=1,2,\dots ,{\mathbf{m}}$, for the $\mathit{i}$th solution, for $\mathit{i}=1,2,\dots ,{\mathbf{nb}}$.
20:   $\mathbf{fjac}\left({\mathbf{ldfjac}},{\mathbf{sdfjac}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the last dimension of the array fjac must be at least ${\mathbf{nb}}$.
On exit: for the $i$th returned solution, the Jacobian matrix of the functions ${f}_{1},{f}_{2},\dots ,{f}_{m}$ at the final iterate, i.e., ${\mathbf{fjac}}\left(\mathit{k},\mathit{j},\mathit{i}\right)$ contains the partial derivative of the $\mathit{k}$th function with respect to the $\mathit{j}$th variable, for $\mathit{k}=1,2,\dots ,{\mathbf{m}}$, $\mathit{j}=1,2,\dots ,{\mathbf{n}}$ and $\mathit{i}=1,2,\dots ,{\mathbf{nb}}$. (See also the discussion of argument fjsl under objfun.)
21:   $\mathbf{ldfjac}$ – IntegerInput
On entry: the first dimension of the array fjac as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{ldfjac}}\ge {\mathbf{m}}$.
22:   $\mathbf{sdfjac}$ – IntegerInput
On entry: the second dimension of the array fjac as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{sdfjac}}\ge {\mathbf{n}}$.
23:   $\mathbf{iter}\left({\mathbf{nb}}\right)$ – Integer arrayOutput
On exit: ${\mathbf{iter}}\left(i\right)$ contains the number of major iterations performed to obtain the $i$th solution. If less than nb solutions are returned then ${\mathbf{iter}}\left({\mathbf{nb}}\right)$ contains the number of starting points that have resulted in a converged solution. If this is close to npts then this might be indicative that fewer than nb local minima exist.
24:   $\mathbf{c}\left({\mathbf{ldc}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the second dimension of the array c must be at least ${\mathbf{nb}}$.
On exit: if ${\mathbf{ncnln}}>0$, ${\mathbf{c}}\left(\mathit{j},\mathit{i}\right)$ contains the value of the $\mathit{j}$th nonlinear constraint function ${c}_{\mathit{j}}$ at the final iterate, for the $\mathit{i}$th solution, for $\mathit{j}=1,2,\dots ,{\mathbf{ncnln}}$.
If ${\mathbf{ncnln}}=0$, the array c is not referenced.
25:   $\mathbf{ldc}$ – IntegerInput
On entry: the first dimension of the array c as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{ldc}}\ge {\mathbf{ncnln}}$.
26:   $\mathbf{cjac}\left({\mathbf{ldcjac}},{\mathbf{sdcjac}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the last dimension of the array cjac must be at least ${\mathbf{nb}}$.
On exit: if ${\mathbf{ncnln}}>0$, cjac contains the Jacobian matrices of the nonlinear constraint functions at the final iterate for each of the returned solutions, i.e., ${\mathbf{cjac}}\left(\mathit{k},\mathit{j},i\right)$ contains the partial derivative of the $\mathit{k}$th constraint function with respect to the $\mathit{j}$th variable, for $\mathit{k}=1,2,\dots ,{\mathbf{ncnln}}$ and $\mathit{j}=1,2,\dots ,{\mathbf{n}}$, for the $i$th solution. (See the discussion of argument cjsl under confun.)
If ${\mathbf{ncnln}}=0$, the array cjac is not referenced.
27:   $\mathbf{ldcjac}$ – IntegerInput
On entry: the first dimension of the array cjac as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{ldcjac}}\ge {\mathbf{ncnln}}$.
28:   $\mathbf{sdcjac}$ – IntegerInput
On entry: the second dimension of the array cjac as declared in the (sub)program from which e05usf is called.
Constraint: if ${\mathbf{ncnln}}>0$, ${\mathbf{sdcjac}}\ge {\mathbf{n}}$.
29:   $\mathbf{clamda}\left({\mathbf{ldclda}},*\right)$ – Real (Kind=nag_wp) arrayOutput
Note: the second dimension of the array clamda must be at least ${\mathbf{nb}}$.
On exit: the values of the QP multipliers from the last QP subproblem solved for the $i$th solution. ${\mathbf{clamda}}\left(j,i\right)$ should be non-negative if ${\mathbf{istate}}\left(j,i\right)=1$ and non-positive if ${\mathbf{istate}}\left(j,i\right)=2$.
30:   $\mathbf{ldclda}$ – IntegerInput
On entry: the first dimension of the array clamda as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{ldclda}}\ge {\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}$.
31:   $\mathbf{istate}\left({\mathbf{listat}},*\right)$ – Integer arrayOutput
Note: the second dimension of the array istate must be at least ${\mathbf{nb}}$.
On exit: ${\mathbf{istate}}\left(j,i\right)$ contains the status of the constraints in the QP working set for the $i$th solution. The significance of each possible value of ${\mathbf{istate}}\left(j,i\right)$ is as follows:
 ${\mathbf{istate}}\left(j,i\right)$ Meaning $\phantom{-}0$ The constraint is satisfied to within the feasibility tolerance, but is not in the QP working set. $\phantom{-}1$ This inequality constraint is included in the QP working set at its lower bound. $\phantom{-}2$ This inequality constraint is included in the QP working set at its upper bound. $\phantom{-}3$ This constraint is included in the QP working set as an equality. This value of istate can occur only when ${\mathbf{bl}}\left(j\right)={\mathbf{bu}}\left(j\right)$.
32:   $\mathbf{listat}$ – IntegerInput
On entry: the first dimension of the array istate as declared in the (sub)program from which e05usf is called.
Constraint: ${\mathbf{listat}}\ge {\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}$.
33:   $\mathbf{iopts}\left(740\right)$ – Integer arrayCommunication Array
34:   $\mathbf{opts}\left(485\right)$ – Real (Kind=nag_wp) arrayCommunication Array
The arrays iopts and opts must not be altered between calls to any of the routines e05usf and e05zkf.
35:   $\mathbf{iuser}\left(*\right)$ – Integer arrayUser Workspace
36:   $\mathbf{ruser}\left(*\right)$ – Real (Kind=nag_wp) arrayUser Workspace
iuser and ruser are not used by e05usf, but are passed directly to confun, objfun and start and may be used to pass information to these routines.
With SMP-enabled versions of e05usf the arrays iuser and ruser provided are classified as OpenMP shared memory. Use of iuser and ruser has to take account of this in order to preserve thread safety whenever information is written back to either of these arrays.
37:   $\mathbf{info}\left({\mathbf{nb}}\right)$ – Integer arrayOutput
On exit: if ${\mathbf{ifail}}={\mathbf{0}}$, ${\mathbf{info}}\left(i\right)$ does not contain an error value returned by e04usf/e04usa .
If ${\mathbf{ifail}}={\mathbf{8}}$ on exit, then not all nb solutions have been found, and ${\mathbf{info}}\left({\mathbf{nb}}\right)$ contains the number of solutions actually found.
38:   $\mathbf{ifail}$ – IntegerInput/Output
On entry: ifail must be set to $0$, . If you are unfamiliar with this argument you should refer to Section 3.4 in How to Use the NAG Library and its Documentation for details.
For environments where it might be inappropriate to halt program execution when an error is detected, the value  is recommended. If the output of error messages is undesirable, then the value $1$ is recommended. Otherwise, because for this routine the values of the output arguments may be useful even if ${\mathbf{ifail}}\ne {\mathbf{0}}$ on exit, the recommended value is $-1$. When the value  is used it is essential to test the value of ifail on exit.
On exit: ${\mathbf{ifail}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see Section 6).

## 6Error Indicators and Warnings

If on entry ${\mathbf{ifail}}=0$ or $-1$, explanatory error messages are output on the current error message unit (as defined by x04aaf).
Note: e05usf may return useful information for one or more of the following detected errors or warnings.
Errors or warnings detected by the routine:
${\mathbf{ifail}}=1$
An input value is incorrect. One or more of the following requirements are violated:
• On entry, ${\mathbf{bl}}\left(i\right)>{\mathbf{bu}}\left(i\right)$: $i=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{bl}}\left(i\right)\le {\mathbf{bu}}\left(i\right)$, for all $i$.
• On entry, ${\mathbf{lda}}=〈\mathit{\text{value}}〉$ and ${\mathbf{nclin}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{lda}}\ge {\mathbf{nclin}}$.
• On entry, ${\mathbf{ldc}}=〈\mathit{\text{value}}〉$ and ${\mathbf{ncnln}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ldc}}\ge {\mathbf{ncnln}}$.
• On entry, ${\mathbf{ldcjac}}=〈\mathit{\text{value}}〉$ and ${\mathbf{ncnln}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ldcjac}}\ge {\mathbf{ncnln}}$.
• On entry, ${\mathbf{ldclda}}=〈\mathit{\text{value}}〉$, ${\mathbf{n}}=〈\mathit{\text{value}}〉$, ${\mathbf{nclin}}=〈\mathit{\text{value}}〉$ and ${\mathbf{ncnln}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ldclda}}\ge {\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}$.
• On entry, ${\mathbf{ldfjac}}=〈\mathit{\text{value}}〉$ and ${\mathbf{m}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ldfjac}}\ge {\mathbf{m}}$.
• On entry, ${\mathbf{ldx}}=〈\mathit{\text{value}}〉$ and ${\mathbf{n}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ldx}}\ge {\mathbf{n}}$.
• On entry, ${\mathbf{listat}}=〈\mathit{\text{value}}〉$, ${\mathbf{n}}=〈\mathit{\text{value}}〉$, ${\mathbf{nclin}}=〈\mathit{\text{value}}〉$ and ${\mathbf{ncnln}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{listat}}\ge {\mathbf{n}}+{\mathbf{nclin}}+{\mathbf{ncnln}}$.
• On entry, ${\mathbf{m}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{m}}>0$.
• On entry, ${\mathbf{n}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{n}}>0$.
• On entry, ${\mathbf{nb}}=〈\mathit{\text{value}}〉$ and ${\mathbf{npts}}=〈\mathit{\text{value}}〉$.
Constraint: $1\le {\mathbf{nb}}\le {\mathbf{npts}}$.
• On entry, ${\mathbf{nclin}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{nclin}}\ge 0$.
• On entry, ${\mathbf{ncnln}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ncnln}}\ge 0$.
• On entry, ${\mathbf{ncnln}}>0$, ${\mathbf{sdcjac}}=〈\mathit{\text{value}}〉$ and ${\mathbf{n}}=〈\mathit{\text{value}}〉$.
Constraint: if ${\mathbf{ncnln}}>0$, ${\mathbf{sdcjac}}\ge {\mathbf{n}}$.
• On entry, ${\mathbf{sdfjac}}=〈\mathit{\text{value}}〉$ and ${\mathbf{n}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{sdfjac}}\ge {\mathbf{n}}$.
${\mathbf{ifail}}=2$
e05usf has terminated without finding any solutions. The majority of calls to the local optimizer have failed to find a feasible point for the linear constraints and bounds, which means that either no feasible point exists for the given value of the optional parameter Linear Feasibility Tolerance (default value $\sqrt{\mathit{macheps}}$, where $\mathit{macheps}$ is the machine precision), or no feasible point could be found in the number of iterations specified by the optional parameter Minor Iteration Limit. You should check that there are no constraint redundancies. If the data for the constraints are accurate only to an absolute precision $\sigma$, you should ensure that the value of the optional parameter Linear Feasibility Tolerance is greater than $\sigma$. For example, if all elements of ${A}_{L}$ are of order unity and are accurate to only three decimal places, Linear Feasibility Tolerance should be at least ${10}^{-3}$.
No solution obtained. Linear constraints may be infeasible.
${\mathbf{ifail}}=3$
e05usf has failed to find any solutions. The majority of local optimizations could not find a feasible point for the nonlinear constraints. The problem may have no feasible solution. This behaviour will occur if there is no feasible point for the nonlinear constraints. (However, there is no general test that can determine whether a feasible point exists for a set of nonlinear constraints.)
No solution obtained. Nonlinear constraints may be infeasible.
${\mathbf{ifail}}=4$
e05usf has failed to find any solutions. The majority of local optimizations have failed because the limiting number of iterations have been reached.
No solution obtained. Many potential solutions reach iteration limit.
${\mathbf{ifail}}=7$
The user-supplied derivatives of the objective function and/or nonlinear constraints appear to be incorrect.
Large errors were found in the derivatives of the objective function and/or nonlinear constraints. This value of ifail will occur if the verification process indicated that at least one gradient or Jacobian element had no correct figures. You should refer to or enable the printed output to determine which elements are suspected to be in error.
As a first-step, you should check that the code for the objective and constraint values is correct – for example, by computing the function at a point where the correct value is known. However, care should be taken that the chosen point fully tests the evaluation of the function. It is remarkable how often the values $x=0$ or $x=1$ are used to test function evaluation procedures, and how often the special properties of these numbers make the test meaningless.
Gradient checking will be ineffective if the objective function uses information computed by the constraints, since they are not necessarily computed before each function evaluation.
Errors in programming the function may be quite subtle in that the function value is ‘almost’ correct. For example, the function may not be accurate to full precision because of the inaccurate calculation of a subsidiary quantity, or the limited accuracy of data upon which the function depends. A common error on machines where numerical calculations are usually performed in double precision is to include even one single precision constant in the calculation of the function; since some compilers do not convert such constants to double precision, half the correct figures may be lost by such a seemingly trivial error.
${\mathbf{ifail}}=8$
Only $〈\mathit{\text{value}}〉$ solutions obtained.
Not all nb solutions have been found. ${\mathbf{info}}\left({\mathbf{nb}}\right)$ contains the number actually found.
${\mathbf{ifail}}=9$
User terminated computation from start procedure: ${\mathbf{mode}}=〈\mathit{\text{value}}〉$.
If e05ucz has been used as an actual argument for start then the message displayed, when ${\mathbf{ifail}}={\mathbf{0}}$ or $-{\mathbf{1}}$ on entry to e05usf, will have the following meaning:
 $998$ failure to allocate space, a smaller value of NPTS should be tried. $997$ an internal error has occurred. Please contact NAG for assistance.
${\mathbf{ifail}}=10$
Failed to initialize optional parameter arrays.
${\mathbf{ifail}}=-99$
See Section 3.9 in How to Use the NAG Library and its Documentation for further information.
${\mathbf{ifail}}=-399$
Your licence key may have expired or may not have been installed correctly.
See Section 3.8 in How to Use the NAG Library and its Documentation for further information.
${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.
See Section 3.7 in How to Use the NAG Library and its Documentation for further information.

## 7Accuracy

If ${\mathbf{ifail}}={\mathbf{0}}$ on exit and the value of ${\mathbf{info}}\left(i\right)=0$, then the vector returned in the array x for solution $i$ is an estimate of the solution to an accuracy of approximately Optimality Tolerance.

## 8Parallelism and Performance

e05usf is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library. In these implementations, this routine may make calls to the user-supplied functions from within an OpenMP parallel region. Thus OpenMP directives within the user functions can only be used if you are compiling the user-supplied function and linking the executable in accordance with the instructions in the Users' Note for your implementation. The user workspace arrays iuser and ruser are classified as OpenMP shared memory and use of iuser and ruser has to take account of this in order to preserve thread safety whenever information is written back to either of these arrays. If at all possible, it is recommended that these arrays are only used to supply read-only data to the user functions when a multithreaded implementation is being used.
e05usf makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the Users' Note for your implementation for any additional implementation-specific information.

You should be wary of requesting much intermediate output from the local optimizer, since large volumes may be produced if npts is large.
The auxiliary routine e05ucz smakes use of the NAG quasi-random Sobol generator (g05ylf and g05ymf). If e05ucz is used as the actual argument for start and ${\mathbf{repeat1}}=\mathrm{.FALSE.}$ then a randomly chosen value for iskip is used, otherwise iskip is set to $100$. If repeat1 is set to .FALSE. and the program is executed several times, each time producing the same best answer, then there is increased probability that this answer is a global minimum. However, if it is important that identical results be obtained on successive runs, then repeat1 should be set to .TRUE..

### 9.1Description of the Printed Output

See Section 9.1 in e04usf/e04usa.

## 10Example

This example is based on Problem 57 in Hock and Schittkowski (1981) and involves the minimization of the sum of squares function
 $Fx = 12 ∑ i=1 44 yi - fi x 2 ,$
where
 $fi x = x1 + 0.49-x1 e -x2 ai- 8$
and
 $i yi ai i yi ai 1 0.49 8 23 0.41 22 2 0.49 8 24 0.40 22 3 0.48 10 25 0.42 24 4 0.47 10 26 0.40 24 5 0.48 10 27 0.40 24 6 0.47 10 28 0.41 26 7 0.46 12 29 0.40 26 8 0.46 12 30 0.41 26 9 0.45 12 31 0.41 28 10 0.43 12 32 0.40 28 11 0.45 14 33 0.40 30 12 0.43 14 34 0.40 30 13 0.43 14 35 0.38 30 14 0.44 16 36 0.41 32 15 0.43 16 37 0.40 32 16 0.43 16 38 0.40 34 17 0.46 18 39 0.41 36 18 0.45 18 40 0.38 36 19 0.42 20 41 0.40 38 20 0.42 20 42 0.40 38 21 0.43 20 43 0.39 40 22 0.41 22 44 0.39 42$
subject to the bounds
 $x1≥-0.4 x2≥-4.0$
to the general linear constraint
 $x1+x2≥1.0$
and to the nonlinear constraint
 $0.49x2-x1 x2≥ 0.09 .$
The optimal solution (to five figures) is
 $x*=0.41995,1.28484T,$
and $F\left({x}^{*}\right)=0.01423$. The nonlinear constraint is active at the solution.
The document for e04uqf/e04uqa includes an example program to solve the same problem using some of the optional parameters described in Section 12.

### 10.1Program Text

Program Text (e05usfe.f90)

### 10.2Program Data

Program Data (e05usfe.d)

### 10.3Program Results

Program Results (e05usfe.r)

## 11Algorithmic Details

See Section 11 in e04usf/e04usa.

## 12Optional Parameters

Several optional parameters in e05usf define choices in the problem specification or the algorithm logic. In order to reduce the number of formal arguments of e05usf these optional parameters have associated default values that are appropriate for most problems. Therefore you need only specify those optional parameters whose values are to be different from their default values.
Optional parameters may be specified by calling e05zkf before a call to e05usf. Before calling e05usf, the optional parameter arrays iopts and opts must be initialized for use with e05usf by calling e05zkf with optstr set to ‘Initialize = e05usf’.
All optional parameters not specified are set to their default values. Optional parameters specified are unaltered by e05usf (unless they define invalid values) and so remain in effect for subsequent calls to e05usf.
See Section 12 in e04usf/e04usa for full details.
The Warm Start option of e04usf/e04usa is not a valid option for use with e05usf.
e05usf supports two options that are distinct from those of e04usf/e04usa:
 Punch Unit $i$ Default $\text{}=6$
This option allows you to send information arising from an appropriate setting of Out_Level to be sent to the Fortran unit number defined by Punch Unit. If you wish this file to be different to the standard output unit ($6$) where other output is displayed then this file should be attached by calling x04acf prior to calling e05usf.
 Out_Level $i$ Default $\text{}=0$
This option defines the amount of extra information to be sent to the Fortran unit number defined by Punch Unit. The possible choices for $i$ are the following:
 $i$ Meaning 0 No extra output. 1 Updated solutions only. This is useful during long runs to observe progress. 2 Successful start points only. This is useful to save the starting points that gave rise to the final solution. 3 Both updated solutions and successful start points.

## 13Description of Monitoring Information

See Section 13 in e04usf/e04usa.