# NAG FL Interfacee04srf (handle_​solve_​ssqp)

Note: this routine uses optional parameters to define choices in the problem specification and in the details of the algorithm. If you wish to use default settings for all of the optional parameters, you need only read Sections 1 to 10 of this document. If, however, you wish to reset some or all of the settings please refer to Section 11 for a detailed description of the algorithm and to Section 12 for a detailed description of the specification of the optional parameters.
Note: please be advised that this routine is classed as ‘experimental’ and its interface may be developed further in the future. Please see Section 4 in How to Use the NAG Library for further information.

## ▸▿ Contents

Settings help

FL Name Style:

FL Specification Language:

## 1Purpose

e04srf is a solver from the NAG optimization modelling suite for large-scale Nonlinear Programming (NLP) problems based on an active-set Sequential Quadratic Programming (SQP) method, using limited-memory quasi-Newton approximations to the Hessian of the Lagrangian.

## 2Specification

Fortran Interface
 Subroutine e04srf ( hess, nvar, x, nnzu, u,
 Integer, Intent (In) :: nvar, nnzu Integer, Intent (Inout) :: iuser(*), ifail Real (Kind=nag_wp), Intent (Inout) :: x(nvar), u(nnzu), ruser(*) Real (Kind=nag_wp), Intent (Out) :: rinfo(100), stats(100) Type (c_ptr), Intent (In) :: handle, cpuser External :: objfun, objgrd, confun, congrd, hess, monit
#include <nag.h>
 void e04srf_ (void **handle, void (NAG_CALL *objfun)(const Integer *nvar, const double x[], double *fx, Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *objgrd)(const Integer *nvar, const double x[], const Integer *nnzfd, double fdx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *confun)(const Integer *nvar, const double x[], const Integer *ncnln, double gx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *congrd)(const Integer *nvar, const double x[], const Integer *nnzgd, double gdx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *hess)(const Integer *nvar, const double x[], const Integer *ncnln, const Integer *idf, const double *sigma, const double lambda[], const Integer *nnzh, double hx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *monit)(const Integer *nvar, const double x[], const Integer *nnzu, const double u[], Integer *inform, const double rinfo[], const double stats[], Integer iuser[], double ruser[], void **cpuser),const Integer *nvar, double x[], const Integer *nnzu, double u[], double rinfo[], double stats[], Integer iuser[], double ruser[], void **cpuser, Integer *ifail)
The routine may be called by the names e04srf or nagf_opt_handle_solve_ssqp.

## 3Description

e04srf is designed to minimize (maximize) a linear or nonlinear function possibly subject to a set of constraints including bounds on the variables and sparse linear and nonlinear constraints. It is suitable for large-scale general nonlinear programming (NLP) problems, as well as for linearly constrained optimization. It is particularly 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.
e04srf solves a nonlinear programming problem in the following form
 $minimize x∈ℝn f(x) subject to lg≤g(x)≤ug, 12 xTQix + riTx + si≤0 , 1≤i≤mQ, lB≤Bx≤uB, lx≤x≤ux,$ (1)
where
• $n$ is the number of the decision variables,
• ${m}_{g}$ is the number of the nonlinear constraints and $g\left(x\right)$, ${l}_{g}$ and ${u}_{g}$ are ${m}_{g}$-dimensional vectors,
• ${m}_{Q}$ is the number of the quadratic constraints,
• ${m}_{B}$ is the number of the linear constraints and $B$ is an ${m}_{B}×n$ matrix, ${l}_{B}$ and ${u}_{B}$ are ${m}_{B}$-dimensional vectors,
• there are $n$ box constraints and ${l}_{x}$ and ${u}_{x}$ are $n$-dimensional vectors.
e04srf serves as a solver for problems stored as a handle. The handle points to an internal data structure which defines the problem and serves as a means of communication for routines in the NAG optimization modelling suite. First, the problem handle is initialized, typically by calling e04raf. Then some of the components of the model need to be defined by calling one or more of the following routines based on the component type:
• the objective function $f\left(x\right)$: linear (e04ref or e04rff), quadratic (e04rff, e04rsf or e04rtf) or nonlinear (e04rgf),
• bound constraints on the variables ${l}_{x},{u}_{x}$ (e04rhf),
• one or more blocks of linear constraints ${l}_{B},B,{u}_{B}$ (e04rjf),
• quadratic constraints (e04rsf and e04rtf for quadratics in a factorized form),
• general nonlinear constraints ${l}_{g},g\left(x\right),{u}_{g}$ (e04rkf).
Once the problem is fully described, the handle may be passed to the solver e04srf. When the handle is no longer needed, e04rzf should be called to destroy it and deallocate the memory held within. There are other ways that the problem model can be defined and/or updated, see Section 3.1 in the E04 Chapter Introduction for more details about the NAG optimization modelling suite.
The algorithm behaviour can be modified by various optional parameters (see Section 12) which can be set by e04zmf and e04zpf anytime between the initialization of the handle and a call to the solver. Once the solver has finished, options or the model may be modified for the next solve. The solver may be called repeatedly with various starting points and/or optional parameters. Option getter e04znf can be called to retrieve the current value of any option.
The optional parameter Task may be used to switch the problem to maximization.
In certain situations the solver can be sped up by warm starting, see Section 9.4 and optional parameter SSQP Start Type for details.
Several options may have a significant impact on the performance of the solver. Even if the defaults were chosen to suit the majority of problems, it is recommended that you experiment in order to find the most suitable set of options for a particular problem, see Sections 11 and 12 for further details.

### 3.1Note regarding nonlinear functions

In the presence of a nonlinear objective function or constraint function, the structure (sparsity pattern) of these functions are defined by e04rgf and e04rkf respectively; the function values and the first derivatives (gradients) are computed at requested points by the corresponding supplied functions objfun, objgrd, confun and congrd. Ideally, the first derivatives should be known and coded but e04srf estimates any missing ones by finite differences, see SSQP Estimate Derivatives.
The solver can strongly benefit from the correct identification of linear variables; these are the variables which appear only in linear expressions in the whole optimization model. The solver cannot autodetect if any variable present in nonlinear constraints or the nonlinear objective is linear thus there is a mechanism to declare the variable as linear by e04rcf and it is highly recommended to do so whenever possible. See Section 3.1 in e04rcf for a detailed explanation and examples. The main advantage is that linear variables have constant first partial derivatives and must be zero in second derivatives, thus the internal approximation of the Hessian doesn't need to take them into account which speeds up the solver.

### 3.2Alternative solvers

e04srf is a generic solver which can solve a great variety of problems, however, a dedicated solver might work better in the following cases. See Section 4 in the E04 Chapter Introduction.
• If $f\left(x\right)$ is linear, ${m}_{Q}=0$, and $g\left(x\right)$ is absent, (1) is a linear program (LP) and e04mtf is recommended.
• If the problem is a convex quadratic or quadratically constrained problem, try e04ptf.
• If the problem is unconstrained or only box-constrained, see e04kff, or e04jdf if the problem is small and derivatives cannot be provided.
• A general nonlinear problem lacking the sparse structure (which is typically the case for problems with a relatively small number of variables and constraints) might be better handled by e04ucf/​e04uca or e04wdf.
• An alternative solver for the same class of problems but based on a different algorithm (interior point method) is e04stf, see Section 2.5.4 in the E04 Chapter Introduction for a discussion of the advantages of each approach.

## 4References

Dantzig G B (1963) Linear Programming and Extensions Princeton University Press
Eldersveld S K (1991) Large-scale sequential quadratic programming algorithms PhD Thesis Department of Operations Research, Stanford University, Stanford
Fourer R (1982) Solving staircase linear programs by the simplex method Math. Programming 23 274–313
Gill P E, Murray W and Saunders M A (2002) SNOPT: An SQP Algorithm for Large-scale Constrained Optimization 12 979–1006 SIAM J. Optim.
Gill P E, Murray W and Saunders M A (2015) Users' guide for SNOPT 7.5: Software for large-scale linear nonlinear programming Report SNOPT Department of Mathematics, University of California, San Diego https://www.ccom.ucsd.edu/~peg/papers/sndoc7.pdf
Gill P E, Murray W and Saunders M A (2016) Users' guide for SQOPT 7.5: Software for large-scale linear and quadratic programming Report SQOPT Department of Mathematics, University of California, San Diego https://www.ccom.ucsd.edu/~peg/papers/sqdoc7.pdf
Gill P E, Murray W, Saunders M A and Wright M H (1986) Users' guide for NPSOL (Version 4.0): a Fortran package for nonlinear programming Report SOL 86-2 Department of Operations Research, Stanford University
Gill P E, Murray W, Saunders M A and Wright M H (1987) Maintaining $LU$ factors of a general sparse matrix Linear Algebra and its Applics. 88/89 239–270
Gill P E, Murray W, Saunders M A and Wright M H (1992) Some theoretical properties of an augmented Lagrangian merit function Advances in Optimization and Parallel Computing (ed P M Pardalos) 101–128 North Holland
Hock W and Schittkowski K (1981) Test Examples for Nonlinear Programming Codes. Lecture Notes in Economics and Mathematical Systems 187 Springer–Verlag
Murtagh B A and Saunders M A (1978) Large-scale linearly constrained optimization 14 41–72 Math. Programming
Murtagh B A and Saunders M A (1982) A projected Lagrangian algorithm and its implementation for sparse nonlinear constraints Math. Program. Stud. 16 84–118
Murtagh B A and Saunders M A (1995) MINOS 5.4 users' guide Report SOL 83-20R Department of Operations Research, Stanford University

## 5Arguments

1: $\mathbf{handle}$Type (c_ptr) Input
On entry: the handle to the problem. It needs to be initialized (e.g., by e04raf) and to hold a problem formulation compatible with e04srf. It must not be changed between calls to the NAG optimization modelling suite.
2: $\mathbf{objfun}$Subroutine, supplied by the NAG Library or the user. External Procedure
objfun must calculate the value of the nonlinear objective function $f\left(x\right)$ at a specified point $x$. If there is no nonlinear objective (i.e., e04rgf was not called), objfun will never be called by e04srf and objfun may be the dummy routine e04srv included in the NAG Library.
The specification of objfun is:
Fortran Interface
 Subroutine objfun ( nvar, x, fx,
 Integer, Intent (In) :: nvar Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar) Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Real (Kind=nag_wp), Intent (Out) :: fx Type (c_ptr), Intent (In) :: cpuser
 void objfun (const Integer *nvar, const double x[], double *fx, Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the objective function is to be evaluated.
3: $\mathbf{fx}$Real (Kind=nag_wp) Output
On exit: the value of the objective function at $x$.
4: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
On exit: may be used to indicate that the function cannot be evaluated at the requested point $x$ by setting ${\mathbf{inform}}<0$. The algorithm will try to recover if the objective cannot be evaluated; if recovery is not possible it will stop with ${\mathbf{ifail}}={\mathbf{21}}$ or ${\mathbf{ifail}}={\mathbf{25}}$.
5: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
6: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
7: $\mathbf{cpuser}$Type (c_ptr) User Workspace
objfun is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to objfun.
objfun must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
3: $\mathbf{objgrd}$Subroutine, supplied by the NAG Library or the user. External Procedure
objgrd must calculate the values of the nonlinear objective function gradient $\frac{\partial f}{\partial x}$ at a specified point $x$. If there is no nonlinear objective (i.e., e04rgf was not called), objgrd will never be called by e04srf and objgrd may be the dummy subroutine e04srw (included in the NAG Library.)
If the optional parameter ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{YES}$ or $\mathrm{ONLYOBJ}$, then after returning from objgrd the gradient vector is checked for missing entries which you have not supplied. Missing entries will be estimated using a finite difference method, see optional parameter SSQP Estimate Derivatives description for more details.
Note: objgrd should not return floating-point NaN (Not a Number) or infinity values, if your code inadvertently returns any NaNs or infinities, e04srf is likely to terminate with ${\mathbf{ifail}}={\mathbf{21}}$, ${\mathbf{24}}$ or ${\mathbf{25}}$.
The specification of objgrd is:
Fortran Interface
 Subroutine objgrd ( nvar, x, fdx,
 Integer, Intent (In) :: nvar, nnzfd Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar) Real (Kind=nag_wp), Intent (Inout) :: fdx(nnzfd), ruser(*) Type (c_ptr), Intent (In) :: cpuser
 void objgrd (const Integer *nvar, const double x[], const Integer *nnzfd, double fdx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the objective function gradient is to be evaluated.
3: $\mathbf{nnzfd}$Integer Input
On entry: the number of nonzero elements in the sparse gradient vector of the objective function, as was set in a previous call to e04rgf.
4: $\mathbf{fdx}\left({\mathbf{nnzfd}}\right)$Real (Kind=nag_wp) array Input/Output
On entry: the elements should only be assigned and not referenced.
On exit: the values of the nonzero elements in the sparse gradient vector of the objective function, in the order specified by idxfd in a previous call to e04rgf. ${\mathbf{fdx}}\left(\mathit{i}\right)$ will store the gradient element $\frac{\partial f}{\partial {x}_{j}}$, where $j={\mathbf{idxfd}}\left(\mathit{i}\right)$.
5: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
If ${\mathbf{inform}}=100$, then this is a special presolve call. This call is used to fill-in the constant elements of the objective gradient, e.g., when there are linear variables in the model. Note that this call may occur with an iterate x that may not be feasible with respect to the linear constraints but it will always be feasible with respect to the variable bounds. Furthermore, this call only requires to return the constant elements of the gradient, $\frac{\partial f}{\partial {x}_{\mathit{i}}}$, where $f\left(x\right)$ is linear in ${x}_{\mathit{i}}$. If you do not provide these elements, subsequent calls to objfun with inform$=0$ may occur in the vicinity of x with the purpose of estimating the missing derivatives using finite differences.
On exit: may be used to inform that the gradient cannot be evaluated at the requested point $x$ by setting ${\mathbf{inform}}<0$. The algorithm will try to recover if the gradient cannot be evaluated; if recovery is not possible it will stop with ${\mathbf{ifail}}={\mathbf{21}}$ or ${\mathbf{ifail}}={\mathbf{25}}$.
6: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
7: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
8: $\mathbf{cpuser}$Type (c_ptr) User Workspace
objgrd is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to objgrd.
objgrd must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
4: $\mathbf{confun}$Subroutine, supplied by the NAG Library or the user. External Procedure
confun must calculate the values of the nonlinear constraints $g\left(x\right)$ at a specified point $x$. If there are no nonlinear constraints then confun will never be called by e04srf and it may be the dummy subroutine e04srx included in the NAG Library.
Note: confun should not return floating-point NaN (Not a Number) or infinity values, if your code inadvertently returns any NaNs or infinities, e04srf is likely to terminate with ${\mathbf{ifail}}={\mathbf{21}}$, ${\mathbf{24}}$ or ${\mathbf{25}}$.
The specification of confun is:
Fortran Interface
 Subroutine confun ( nvar, x, gx,
 Integer, Intent (In) :: nvar, ncnln Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar) Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Real (Kind=nag_wp), Intent (Out) :: gx(ncnln) Type (c_ptr), Intent (In) :: cpuser
 void confun (const Integer *nvar, const double x[], const Integer *ncnln, double gx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the constraint functions are to be evaluated.
3: $\mathbf{ncnln}$Integer Input
On entry: ${m}_{g}$, the number of nonlinear constraints, as specified in an earlier call to e04rkf.
4: $\mathbf{gx}\left({\mathbf{ncnln}}\right)$Real (Kind=nag_wp) array Output
On exit: the ${m}_{g}$ values of the nonlinear constraint functions at $x$.
5: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
If ${\mathbf{inform}}=100$, then this is a special presolve call. This call is used to quickly identify if the problem is infeasible, e.g., due to fixed variables. Note that this call can occur with an iterate x that may not be feasible with respect to the linear constraints but it will always be feasible with respect to the variable bounds.
On exit: may be used to inform that the constraint values cannot be evaluated at the requested point $x$ by setting ${\mathbf{inform}}<0$. The algorithm will try to recover if the gradient cannot be evaluated; if recovery is not possible it will stop with ${\mathbf{ifail}}={\mathbf{21}}$ or ${\mathbf{ifail}}={\mathbf{25}}$.
6: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
7: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
8: $\mathbf{cpuser}$Type (c_ptr) User Workspace
confun is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to confun.
confun must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
5: $\mathbf{congrd}$Subroutine, supplied by the NAG Library or the user. External Procedure
congrd must calculate the Jacobian of the nonlinear constraints $\frac{\partial g}{\partial x}$ at a specified point $x$. If there are no nonlinear constraints, congrd will never be called by e04srf and congrd may be the dummy subroutine e04sry included in the NAG Library.
If the optional parameter ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{YES}$ or $\mathrm{ONLYCON}$, then after returning from congrd the Jacobian is checked for missing entries which you have not supplied. Missing entries will be estimated using a finite difference method, see optional parameter SSQP Estimate Derivatives description for more details.
Note: congrd should not return floating-point NaN (Not a Number) or infinity values, if your code inadvertently returns any NaNs or infinities, e04srf is likely to terminate with ${\mathbf{ifail}}={\mathbf{21}}$, ${\mathbf{24}}$ or ${\mathbf{25}}$.
The specification of congrd is:
Fortran Interface
 Subroutine congrd ( nvar, x, gdx,
 Integer, Intent (In) :: nvar, nnzgd Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar) Real (Kind=nag_wp), Intent (Inout) :: gdx(nnzgd), ruser(*) Type (c_ptr), Intent (In) :: cpuser
 void congrd (const Integer *nvar, const double x[], const Integer *nnzgd, double gdx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the Jacobian of the constraint functions is to be evaluated.
3: $\mathbf{nnzgd}$Integer Input
On entry: is the number of nonzero elements in the sparse Jacobian of the constraint functions, as was set in a previous call to e04rkf.
4: $\mathbf{gdx}\left({\mathbf{nnzgd}}\right)$Real (Kind=nag_wp) array Input/Output
On entry: the elements should only be assigned and not referenced.
On exit: the nonzero values of the Jacobian of the nonlinear constraints, in the order specified by irowgd and icolgd in an earlier call to e04rkf. ${\mathbf{gdx}}\left(\mathit{i}\right)$ will be the gradient $\frac{\partial {g}_{j}}{\partial {x}_{k}}$, where $j={\mathbf{irowgd}}\left(\mathit{i}\right)$ and $k={\mathbf{icolgd}}\left(\mathit{i}\right)$.
5: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
If ${\mathbf{inform}}=100$, then this is a special presolve call. This call is used to fill-in the constant elements of the Jacobian, e.g., when there are linear variables in the model. Note that this call may occur with an iterate x that may not be feasible with respect to the linear constraints but it will always be feasible with respect to the variable bounds. Furthermore, this call only requires to return the constant elements of the Jacobian, $\frac{\partial {g}_{\mathit{j}}}{\partial {x}_{\mathit{i}}}$, where ${g}_{\mathit{j}}\left(x\right)$ is linear in ${x}_{\mathit{i}}$. If the user does not provide these elements, subsequent calls to confun with inform$=0$ may occur in the vicinity of x with the purpose of estimating the missing derivatives using finite differences.
On exit: may be used to inform that the Jacobian of the constraints cannot be evaluated at the requested point $x$ by setting ${\mathbf{inform}}<0$. The algorithm will try to recover if the gradient cannot be evaluated; if recovery is not possible it will stop with ${\mathbf{ifail}}={\mathbf{21}}$ or ${\mathbf{ifail}}={\mathbf{25}}$.
6: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
7: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
8: $\mathbf{cpuser}$Type (c_ptr) User Workspace
congrd is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to congrd.
congrd must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
6: $\mathbf{hess}$Subroutine, supplied by the NAG Library or the user. External Procedure
hess is reserved for future releases of the NAG Library which will allow calculation of the second derivatives. It will never be called in the current implementation and hess may be the dummy subroutine e04srz included in the NAG Library.
The specification of hess is:
Fortran Interface
 Subroutine hess ( nvar, x, idf, nnzh, hx,
 Integer, Intent (In) :: nvar, ncnln, idf, nnzh Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar), sigma, lambda(ncnln) Real (Kind=nag_wp), Intent (Inout) :: hx(nnzh), ruser(*) Type (c_ptr), Intent (In) :: cpuser
 void hess (const Integer *nvar, const double x[], const Integer *ncnln, const Integer *idf, const double *sigma, const double lambda[], const Integer *nnzh, double hx[], Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the Hessian functions are to be evaluated.
3: $\mathbf{ncnln}$Integer Input
On entry: ${m}_{g}$, the number of nonlinear constraints, as specified in an earlier call to e04rkf.
4: $\mathbf{idf}$Integer Input
On entry: specifies the quantities to be computed in hx.
${\mathbf{idf}}=-1$
The values of the Hessian of the Lagrangian will be computed in hx. This will be the case if e04rlf has been called with idf of the same value.
${\mathbf{idf}}=0$
The values of the Hessian of the objective function will be computed in hx. This will be the case if e04rlf has been called with idf of the same value.
${\mathbf{idf}}>0$
The values of the Hessian of the constraint function with index idf will be computed in hx. This will be the case if e04rlf has been called with idf of the same value.
5: $\mathbf{sigma}$Real (Kind=nag_wp) Input
On entry: if ${\mathbf{idf}}=-1$, the value of the $\sigma$ quantity in the definition of the Hessian of the Lagrangian. Otherwise, sigma should not be referenced.
6: $\mathbf{lambda}\left({\mathbf{ncnln}}\right)$Real (Kind=nag_wp) array Input
On entry: if ${\mathbf{idf}}=-1$, the values of the ${\lambda }_{i}$ quantities in the definition of the Hessian of the Lagrangian. Otherwise, lambda should not be referenced.
7: $\mathbf{nnzh}$Integer Input
On entry: the number of nonzero elements in the Hessian to be computed.
8: $\mathbf{hx}\left({\mathbf{nnzh}}\right)$Real (Kind=nag_wp) array Input/Output
On entry: the elements should only be assigned and not referenced.
On exit: the nonzero values of the requested Hessian evaluated at $x$. For each value of idf, the ordering of nonzeros must follow the sparsity structure registered in the handle by earlier calls to e04rlf through the arguments irowh and icolh.
9: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
On exit: must be set to a value describing the action to be taken by the solver on return from hess. Specifically, if the value is negative the solution of the current problem will terminate immediately with ${\mathbf{ifail}}={\mathbf{21}}$ or ${\mathbf{ifail}}={\mathbf{25}}$; otherwise, computations will continue.
10: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
11: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
12: $\mathbf{cpuser}$Type (c_ptr) User Workspace
hess is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to hess.
hess must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
7: $\mathbf{monit}$Subroutine, supplied by the NAG Library or the user. External Procedure
monit is provided to enable you to monitor the progress of the optimization and optionally to terminate the solver early if necessary. It is invoked at the end of every $i$th major iteration where $i$ is given by the SSQP Monitor Frequency (the default is $0$, monit is not called).
monit may be the dummy subroutine e04sru included in the NAG Library.
The specification of monit is:
Fortran Interface
 Subroutine monit ( nvar, x, nnzu, u,
 Integer, Intent (In) :: nvar, nnzu Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar), u(nnzu), rinfo(100), stats(100) Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Type (c_ptr), Intent (In) :: cpuser
 void monit (const Integer *nvar, const double x[], const Integer *nnzu, const double u[], Integer *inform, const double rinfo[], const double stats[], Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: ${x}^{i}$, the values of the decision variables $x$ at the current iteration.
3: $\mathbf{nnzu}$Integer Input
On entry: the dimension of array u.
4: $\mathbf{u}\left({\mathbf{nnzu}}\right)$Real (Kind=nag_wp) array Input
On entry: if ${\mathbf{nnzu}}>0$, u holds the values of Lagrange multipliers (dual variables, see Section 11.3) for the constraints at the current iteration. See Section 3.1 in e04stf for layout information.
5: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
On exit: may be used to request the solver to stop immediately. Specifically, if ${\mathbf{inform}}<0$ the solver will terminate immediately with ${\mathbf{ifail}}={\mathbf{20}}$; otherwise, the solver will proceed normally.
6: $\mathbf{rinfo}\left(100\right)$Real (Kind=nag_wp) array Input
On entry: error measures and various indicators at the end of the current iteration as described in rinfo.
7: $\mathbf{stats}\left(100\right)$Real (Kind=nag_wp) array Input
On entry: solver statistics at the end of the current iteration as described in stats.
8: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
9: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
10: $\mathbf{cpuser}$Type (c_ptr) User Workspace
monit is called with the arguments iuser, ruser and cpuser as supplied to e04srf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to monit.
monit must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e04srf is called. Arguments denoted as Input must not be changed by this procedure.
8: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
9: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input/Output
On entry: ${x}_{0}$, the initial estimates of the variables, $x$.
On exit: the final values of the variables, $x$.
10: $\mathbf{nnzu}$Integer Input
On entry: the number of Lagrange multipliers in array u.
If ${\mathbf{nnzu}}=0$, u will not be referenced; otherwise, it needs to match the dimension $q$ as explained in Section 3.1 in e04stf.
Constraints:
• ${\mathbf{nnzu}}\ge 0$;
• if ${\mathbf{nnzu}}>0$, ${\mathbf{nnzu}}=q$.
11: $\mathbf{u}\left({\mathbf{nnzu}}\right)$Real (Kind=nag_wp) array Input/Output
Note: if ${\mathbf{nnzu}}>0$, u holds Lagrange multipliers (dual variables) for the constraints. See Section 3.1 in e04stf for layout information. If ${\mathbf{nnzu}}=0$, u will not be referenced.
On entry: optionally provides the initial estimates of Lagrange multipliers, $\left({z}_{0},{\pi }_{0}\right)$. These values are only referenced if a warm start is requested (Section 9.4). If there are no initial estimates available, then set to zero.
On exit: the final value of Lagrange multipliers $\left({z}^{*},{\pi }^{*}\right)$, see Section 11.3.
12: $\mathbf{rinfo}\left(100\right)$Real (Kind=nag_wp) array Output
On exit: error measures and various indicators at the end of the final iteration as given in the table below:
$1$ Objective function value $f\left(x\right)$.
$2$ Constraint violation (primal infeasibility), see Section 7.
$3$ Dual infeasibility.
$4$$7$ Reserved for future use.
$8$ Step size for the primal variables, see Iteration log in Section 9.1.
$9$ Estimate of the condition number of ${R}^{\mathrm{T}}R$, itself an estimate of ${Z}^{\mathrm{T}}HZ$, the reduced Hessian of the Lagrangian. See Section 9.1.
$10$ Largest nonlinear constraint violation.
$11$ Largest nonlinear constraint violation, normalized using $1+‖{x}_{k}‖$, with ${x}_{k}$ the latest iterate. See (12).
$12$ Sum of the infeasibilities of constraints that lie outside one of their bounds by more than the optional parameter SSQP Minor Feasibility Tol before the solution is unscaled.
$13$ Flag indicating whether the current iterate satisfies the first-order Kuhn–Tucker conditions for optimality:
$\mathbit{j}$ Description
0 none of the conditions are satisfied;
1 condition on primal infeasibility satisfied, rinfo$\left(2\right)$ is within SSQP Major Feasibility Tol tolerance;
2 condition on dual infeasibility satisfied, rinfo$\left(3\right)$ is within SSQP Minor Feasibility Tol tolerance; and
3 both conditions are satisfied.
$14$ The size of the penalty parameter $D$ of (8). See Section 11.4.
$15$ Flag indicating whether the current iterate is feasible with respect to the linear constraints, if any. A nonzero value indicates feasibility.
$16$$100$ Reserved for future use.
13: $\mathbf{stats}\left(100\right)$Real (Kind=nag_wp) array Output
On exit: solver statistics at the end of the final iteration as given in the table below:
 $1$ Number of major iterations. $2$ Total number of minor iterations. $3$ Number of minor iterations performed on the last major iteration. $4$ Number of Hessian evaluations. $5$ Number of gradient evaluations. $6$ Number of infeasibilities of constraints that lie outside one of their bounds by more than the optional parameter SSQP Minor Feasibility Tol before the solution is unscaled. If any linear constraints are infeasible, it contains the number of variables and linear constraints lying outside their upper or lower bounds. The nonlinear constraints are not evaluated. Otherwise, it contains the number of components of $g\left(x\right)$ lying outside their bounds by more than the optional parameter SSQP Minor Feasibility Tol. Again this is before the solution is unscaled. $7$ Total time spent in presolver (including user-supplied function calls). $8$ Total time spent in solver (including user-supplied function calls). $9$ Total time spent evaluating constraint function. $10$ Total time spent evaluating constraint Jacobian. $11$ Total time spent evaluating objective function. $12$ Total time spent evaluating objective gradient. $13$ Total time spent evaluating Hessian. $14$–$18$ Reserved for future use. $19$ Number of function evaluations. $20$ Number of constraint function evaluations. $21$ Number of constraint Jacobian evaluations. $22$–$24$ Reserved for future use. $25$ Number of superbasic variables. $26$ Number of column swaps between $B$ and $S$. $27$ Number of elements in the $LU$ factorization of $B$. $28$–$29$ Reserved for future use. $30$ Number of slack variables. See Sections 11.1 and 9.2. $31$–$100$ Reserved for future use.
14: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
15: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
16: $\mathbf{cpuser}$Type (c_ptr) User Workspace
iuser, ruser and cpuser are not used by e04srf, but are passed directly to objfun, objgrd, confun, congrd, hess and monit and may be used to pass information to these routines. If you do not need to reference cpuser, it should be initialized to c_null_ptr.
17: $\mathbf{ifail}$Integer Input/Output
On entry: ifail must be set to $0$, $-1$ or $1$ to set behaviour on detection of an error; these values have no effect when no error is detected.
A value of $0$ causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of $-1$ means that an error message is printed while a value of $1$ means that it is not.
If halting is not appropriate, the value $-1$ or $1$ is recommended. If message printing is undesirable, then the value $1$ is recommended. Otherwise, the value $-1$ is recommended since useful values can be provided in some output arguments even when ${\mathbf{ifail}}\ne {\mathbf{0}}$ on exit. When the value $-\mathbf{1}$ or $\mathbf{1}$ 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).
Errors or warnings detected by the routine:
Note: in some cases e04srf may return useful information.
${\mathbf{ifail}}=1$
The supplied handle does not define a valid handle to the data structure for the NAG optimization modelling suite. It has not been properly initialized or it has been corrupted.
${\mathbf{ifail}}=2$
The problem is already being solved.
This solver does not support the model defined in the handle.
${\mathbf{ifail}}=4$
On entry, ${\mathbf{nvar}}=⟨\mathit{\text{value}}⟩$, expected $\mathrm{value}=⟨\mathit{\text{value}}⟩$.
Constraint: nvar must match the current number of variables of the model in the handle.
${\mathbf{ifail}}=5$
On entry, ${\mathbf{nnzu}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{nnzu}}=⟨\mathit{\text{value}}⟩$ or $0$.
On entry, ${\mathbf{nnzu}}=⟨\mathit{\text{value}}⟩$.
Constraint: no constraints present, so ${\mathbf{nnzu}}$ must be $0$.
${\mathbf{ifail}}=7$
The dummy confun routine was called but the problem requires these values. Please provide a proper confun routine.
The dummy congrd routine was called but the problem requires these derivatives. Please provide a proper congrd routine.
The dummy hess routine was called but the problem requires these derivatives. Please provide a proper hess routine.
The dummy objfun routine was called but the problem requires these values. Please provide a proper objfun routine.
The dummy objgrd routine was called but the problem requires these derivatives. Please provide a proper objgrd routine.
${\mathbf{ifail}}=20$
User requested termination during a monitoring step. inform was set to a negative value in monit.
${\mathbf{ifail}}=21$
The current starting point is unusable.
Either inform was set to a negative value within the user-supplied functions objfun, objgrd, confun, congrd or hess, or an Infinity or NaN was detected in values returned from them. Note that the user-provided starting point is made bound-constrain feasible (if applicable) before calling any of the user-supplied functions.
See also ${\mathbf{ifail}}={\mathbf{25}}$.
${\mathbf{ifail}}=22$
Maximum number of iterations reached.
Either the SSQP Iteration Limit, SSQP Major Iteration Limit or SSQP Minor Iteration Limit limit was exceeded before a solution with the requested accuracy could be found. Check the iteration log to be sure that progress was being made. If so, rerun the problem after increasing these optional parameters.
${\mathbf{ifail}}=23$
The solver terminated after the maximum time allowed was exhausted.
Maximum number of seconds exceeded. Use optional parameter Time Limit to change the limit.
${\mathbf{ifail}}=24$
The solver was terminated because no further progress could be achieved.
This can indicate that numerical difficulties have been encountered and little or no progress was made. It could also happen if the problem has been solved to the best numerical accuracy possible given the current scaling. Additionally, several other circumstances could lead to this error:
1. 1.objfun or confun could be returning accurate function values but objgrd or congrd inaccurate gradients (or vice versa). This is the most likely cause. Study the comments given for ${\mathbf{ifail}}={\mathbf{26}}$ as well as those given for Verify Derivatives and SSQP Estimate Derivatives. Ensure that the coding is correct.
2. 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 SSQP Major Optimality Tol of $2×{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. 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
but even then, if $t$ is as large as ${10}^{-5}$ or ${10}^{-6}$ (only $5$ or $6$ significant figures), the same exit condition may occur. At present, the only remedy is to increase the accuracy of the function calculation.
4. 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 $Ax-s=0$ sufficiently well. This probably means that the current basis is very ill-conditioned. If there are some linear constraints and variables, try ${\mathbf{SSQP Scale Option}}=\mathrm{VARS AND LC}$ if scaling has not yet been used.
5. 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 should mean that the problem is badly scaled.
${\mathbf{ifail}}=25$
Invalid value detected in user function and recovery failed.
Either inform was set to a negative value within the user-supplied functions objfun, objgrd, confun, congrd or hess, or an Infinity or NaN was detected in values returned from them. This will occur in the case when one or more of the problem functions are
• undefined in a region around the current iterate and no recovery is possible;
• undefined at the first feasible point; or
• undefined when checking derivatives.
In these cases, it may be convenient to add a linear constraint block to exclude the undefined region. See also ${\mathbf{ifail}}={\mathbf{21}}$.
• The solver is detecting irregular or badly scaled problem functions. See optional parameter SSQP Function Precision.
• The solver iterates are in a region where the QP Hessian is numerically indefinite and the optimization process was stopped.
${\mathbf{ifail}}=26$
User-provided derivatives are likely to be incorrect.
This error indicates that the user-supplied gradient vector fdx or Jacobian matrix gdx have at least one value that 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 e04srf 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}}={\mathbf{24}}$.
Check the function and Jacobian computation very carefully in objgrd and congrd. A simple omission could explain everything. If a component is very large, then give serious thought to scaling the function or the nonlinear variables. If you feel certain that the computed Jacobian is correct (and that the forward-difference estimate is, therefore, wrong), you can specify ${\mathbf{Verify Derivatives}}=\mathrm{NO}$ to prevent individual elements from being checked. However, the optimization procedure may have difficulties or even fail.
${\mathbf{ifail}}=27$
Maximum number of superbasic variables reached.
Increase value of the optional parameter SSQP Superbasics Limit.
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 cannot continue as the number of superbasic variables has already reached the limit specified by the optional parameter SSQP Superbasics Limit. In general, raise the limit by a reasonable amount, bearing in mind the storage needed for the reduced Hessian.
${\mathbf{ifail}}=50$
Problem was solved to an acceptable level, full accuracy was not achieved.
This indicates that the algorithm found a feasible solution, but the requested accuracy in the dual infeasibilities could not be achieved. The solution found is within ${10}^{-2}$ of satisfying the SSQP Major Feasibility Tol. This may happen if the SSQP Major Feasibility Tol tolerance is too small for the current problem, or if the problem is badly scaled.
${\mathbf{ifail}}=51$
The problem was found to be infeasible during preprocessing.
One or more of the constraints (or its part after preprocessing) violates the constraints.
${\mathbf{ifail}}=53$
The problem seems to be infeasible, the algorithm was stopped.
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 $Ax-s=0$ (see (5) and (6) in Section 11.2), there is apparently no point that satisfies the bounds on $x$ and $s$. Violations as small as the SSQP Minor Feasibility Tol 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, e04srf 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), e04srf enters a so-called ‘nonlinear elastic’ mode. The subproblem includes the original QP objective and the sum of the infeasibilities — suitably weighted. 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 will eventually be obtained for the perturbed constraints, and optimization can continue on the subproblem. If the nonlinear problem has no feasible solution, e04srf will tend to determine a ‘good’ infeasible point if the elastic weight is sufficiently large. If the elastic weight were infinite, e04srf would locally minimize the nonlinear constraint violations subject to the linear constraints and bounds.
Unfortunately, even though e04srf 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}}=54$
The problem seems to be unbounded and the algorithm was stopped.
This indicates that the solver found a direction where the objective value decreases arbitrarily. In the case of maximizing, the direction found increases arbitrarily the objective value. Cases where this happens are when a linear objective is provided with no constraints, or when the constraints defined do not bound the search direction.
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 SSQP Scale Option.
For nonlinear problems, e04srf 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, 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.
${\mathbf{ifail}}=90$
Memory limit is too small.
e04srf has terminated because the LU factorization exhausted the allocated memory limit. This indicates that the memory allocation heuristic underestimated the requirements. Please contact NAG.
${\mathbf{ifail}}=-99$
See Section 7 in the Introduction to the NAG Library FL Interface for further information.
${\mathbf{ifail}}=-399$
Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library FL Interface for further information.
${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.
See Section 9 in the Introduction to the NAG Library FL Interface for further information.

## 7Accuracy

If the value of the optional parameter SSQP Major Optimality Tol is set to ${10}^{-d}$ and ${\mathbf{ifail}}={\mathbf{0}}$ on exit, then the final value of $f\left(x\right)$ should have approximately $d$ correct significant digits.
e04srf returns with ${\mathbf{ifail}}={\mathbf{0}}$ if the iterates have converged to a point $x$ that satisfies the first-order Kuhn–Tucker (see Section 9.1.6) conditions to the accuracy requested by the optional parameter SSQP Major Optimality Tol, i.e., the projected gradient and active constraint residuals are negligible at $x$.
You should check whether the following four conditions are satisfied:
1. (i)the final value of rgNorm (see Section 9.1.6) is significantly less than that at the starting point;
2. (ii)during the final major iterations, the values of Step and Minors (see Section 9.1.5) are both one;
3. (iii)the last few values of both rgNorm and SumInf (see Section 9.1.6) become small at a fast linear rate; and
4. (iv)condHz (see Section 9.1.5) is small.
If all these conditions hold, $x$ is almost certainly a local minimum of (1).
A cautionary note about ‘Optimal solutions’. Some of the variables or slacks may lie outside their bounds more than desired, especially if scaling was requested. Max Primal infeas (printed at the end of the solve when Print Level and Monitoring Level$\ge 4$ or Print Level and Monitoring Level$\ge 4$) refers to the largest bound infeasibility and which variable is involved. If it is too large, consider restarting with a smaller SSQP Minor Feasibility Tol (say $10$ times smaller) and perhaps ${\mathbf{SSQP Scale Option}}=\mathrm{NONE}$.
Similarly, Max Dual infeas (also printed when Print Level or Monitoring Level$\ge 4$) indicates which variable is most likely to be at a nonoptimal value. Broadly speaking, if
 $Max Dual infeas / Max pi = 10-d ,$
where Max pi (printed at the end of a solve when Print Level and Monitoring Level$\ge 1$ or Print Level and Monitoring Level$\ge 1$) indicates the largest dual variable in magnitude, then the objective function would probably change in the $d$th significant digit if optimization could be continued. If $d$ seems too large, consider restarting with a smaller SSQP Major Optimality Tol.
Finally, Nonlinear constraint violn (printed at the end of a solve when Print Level and Monitoring Level$\ge 1$ or Print Level and Monitoring Level$\ge 1$) reports the maximum infeasibility for nonlinear rows. If it seems too large, consider restarting with a smaller SSQP Major Feasibility Tol.

## 8Parallelism and Performance

e04srf is not threaded in any implementation.

### 9.1Description of the Printed Output

The solver can print information to give an overview of the problem and of the progress of the computation. The output may be sent to two independent streams (files) which are set by optional parameters Print File and Monitoring File. Optional parameters Print Level, Monitoring Level, Print Solution and Print Options determine the exposed level of detail. This allows, for example, a detailed log file to be generated while the condensed information is displayed on the screen.
By default (${\mathbf{Print File}}=6$, ${\mathbf{Print Level}}=2$), the following sections are printed to the standard output:
• Optional parameters list (if ${\mathbf{Print Options}}=\mathrm{YES}$)
• Problem statistics
• Verification of derivatives (if ${\mathbf{Verify Derivatives}}=\mathrm{YES}$ or $\mathrm{CHEAP}$)
• Major iteration log
• Summary
• Solution (if ${\mathbf{Print Solution}}=\mathrm{YES}$)

The header is a message indicating the start of the solver. It should look like:
```-------------------------------------------------------------------
E04SR, SQP method for large-scale nonlinear optimization problems
-------------------------------------------------------------------```

#### 9.1.2Optional parameters list

The list shows all options of the solver, each displayed on one line. The output contains the option name, its current value and an indicator for how it was set. The options unchanged from the default setting are noted by ‘d’, options you set are noted by ‘U’, and options reset by the solver are noted by ‘S’. Note that the output format is compatible with the file format expected by e04zpf. The output might look as follows:
```Begin of Options
Print File                    =                   6     * d
Print Level                   =                   2     * U
Print Options                 =                 Yes     * d
Print Solution                =                 Yes     * U
Monitoring File               =                   9     * U
Monitoring Level              =                   3     * U

Infinite Bound Size           =         1.00000E+20     * d
Stats Time                    =                  No     * d
Time Limit                    =         1.00000E+06     * d
Verify Derivatives            =                 Yes     * U

...
Ssqp Estimate Derivatives     =                  No     * d
Ssqp Function Precision       =         3.00000E-13     * d
Ssqp Finite Diff Ctrl Interval=         6.69433E-05     * d
Ssqp Finite Diff Interval     =         5.47723E-07     * d
Ssqp Hessian                  =                Auto     * d
Ssqp Hessian Updates          =                  10     * d
Ssqp Iteration Limit          =               10000     * d
Ssqp Major Feasibility Tol    =         1.00001E-06     * d
Ssqp Major Iteration Limit    =                1000     * d
...
End of Options```

#### 9.1.3Problem statistics

Problem statistics
If ${\mathbf{Print Level}}\ge 2$, statistics on the problem are printed. It should look as follows:
```Problem Statistics
No of variables                  6
linear                         4
nonlinear                      2
free (unconstrained)           2
bounded                        2
No of lin. constraints           0
nonzeroes                      0
No of nln. constraints          12
nonzeroes                     15
Note: the list is showing that out of the six variables, four are detected as present in linear expressions (of the objective and constraint functions), while the remaining two are present nonlinearly.
It is followed by information regarding warm start (see Section 9.4) and statistics on the Jacobian matrix, the output may look similar to:
``` Warm start information loaded successfully from handle.
Handle [WARM START BASIS] data origin: solver
Handle [SLACK] data origin           : solver

Matrix statistics
Total      Normal        Free       Fixed     Bounded
Rows              4           1           1           2           0
Columns           4           2           2           0           0

The user has defined       5   out of       5   constraint gradients.
The user has defined       4   out of       4   objective  gradients.```

#### 9.1.4Verification of derivatives

If ${\mathbf{Verify Derivatives}}=\mathrm{CHEAP}$, then only a ‘cheap’ test will be performed on the user-supplied nonlinear objective gradient and constraint Jacobian. The output may be similar to:
``` Cheap test of user-supplied problem derivatives...
The constraint gradients seem to be OK.
-->  The largest discrepancy was    2.52E-07  in constraint     3 (NLC:      2)
The objective  gradients seem to be OK.
Gradient projected in one direction   6.70906812736E-01
Difference approximation              6.70906898365E-01```
If ${\mathbf{Verify Derivatives}}=\mathrm{YES}$, then the user-supplied nonlinear objective gradient and constraint Jacobian are individually verified. The output might look as follows:
``` Verification of user-supplied problem derivatives.
The constraint gradients seem to be OK.
-->  The largest discrepancy was    2.52E-07  in constraint     3 (NLC:      2)
The objective  gradients seem to be OK.
Gradient projected in one direction   6.70906812736E-01
Difference approximation              6.70906898365E-01
Column       x(j)        dx(j)    Element no.    Row        Derivative    Difference approxn
1 -7.06220960E-02  5.86E-07         1         2   -1.41244192E-01   -1.41243605E-01  ok
1 -7.06220960E-02  5.86E-07           Objective    2.68365409E+00    2.68365467E+00  ok
2  1.41244914E+00  1.32E-06         2         2    2.82489828E+00    2.82489960E+00  ok
3         3    1.12714152E+01    1.12714310E+01  ok
2  1.41244914E+00  1.32E-06           Objective    2.68365409E+00    2.68365441E+00  ok
...
-->  The largest relative error was    1.29E-06   in row     3 (NLC:      2),  column     2
4  objective gradients out of     1  through     4  seem to be OK.
-->  The largest relative error was    3.59E-07   in column     2```

#### 9.1.5Major Iteration Log

This section describes the output to unit Print File or Monitoring File if ${\mathbf{Print Level}}\text{​ or ​}{\mathbf{Monitoring Level}}>1$. One line of information is output every major iteration. A heading is printed before the first such line containing the items described below.
 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 Section 11. 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 well-behaved problems, the unit step will be taken as the solution is approached. nCon the number of times objfun or confun have 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}}$ from (12), the maximum component of the scaled nonlinear constraint residual (see optional parameter SSQP Major Feasibility Tol). The solution is regarded as acceptably feasible if Feasible is less than the SSQP Major Feasibility Tol. 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}}$ from (13), the maximum complementary gap (see optional parameter SSQP Major Optimality Tol). 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 from (7). This function will decrease at each iteration unless it was necessary to increase the penalty parameters, see Section 11.4. 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 is the number of subdiagonal elements in the columns of a lower triangular matrix and lenU is the number of diagonal and superdiagonal elements in the rows of an upper-triangular matrix, see Section 9.1.9. 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. 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 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 be under the ...Flags column.
 Flag Description + indicates that the log is from a major iteration. c central differences have been used to compute the unknown components of the objective and constraint gradients. A switch to central differences is made if either the linesearch gives a small step, or $x$ is close to being optimal. In some cases, it may be necessary to re-solve the QP subproblem with the central difference gradient and Jacobian. d during the linesearch it was necessary to decrease the step in order to obtain a maximum constraint violation. D ${\mathbf{inform}}<0$ on return from objfun, indicating that the linesearch needed to be done with a smaller value of the step length ${\alpha }_{k}$ of equation (9). l the norm wise change in the variables was limited by a safeguard upper bound. i if e04srf 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 (11) described in Section 11.5. If e04srf is already in elastic mode, an i indicates that the minimizer of the elastic subproblem does not satisfy the linearized constraints. (In this case, a feasible point for the usual QP subproblem may or may not exist.) M an extra evaluation of the problem functions was needed to define an acceptable positive definite quasi-Newton update to the Lagrangian Hessian. This modification is only done when there are nonlinear constraints. m this is the same as M except that it was also necessary to modify the update to include an augmented Lagrangian term. n no positive definite BFGS update could be found. The approximate Hessian is unchanged from the previous iteration. R the approximate Hessian has been reset by discarding all but the diagonal elements. This reset will be forced periodically by the SSQP Hessian Updates keyword. However, it may also be necessary to reset an ill-conditioned Hessian from time to time. r the approximate Hessian was reset after ten consecutive major iterations in which no BFGS update could be made. The diagonals of the approximate Hessian are retained if at least one update has been done since the last reset. Otherwise, the approximate Hessian is reset to the identity matrix. s a self-scaled BFGS update was performed. This update is used when the Hessian approximation is diagonal, and hence always follows a Hessian reset. t the minor iterations were terminated because of the SSQP Minor Iteration Limit. T the minor iterations were terminated because the limit on new superbasics was reached. u the QP subproblem was unbounded. w a weak solution of the QP subproblem was found. z the SSQP Superbasics Limit was reached.

#### 9.1.6Minor Iteration Log

Is only printed if ${\mathbf{Print Level}}\text{​ or ​}{\mathbf{Monitoring Level}}\ge 3$, one line of information is output every $k$th minor iteration, where $k$ is the specified SSQP Print Frequency, by default $k=1$. 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.
 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 reduced-gradient 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 $1$ 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 reduced-gradient vector at the start of the iteration. (It is the norm of the vector with elements ${d}_{j}$ for variables $j$ in the superbasic set.) During Phase 2 this norm will be approximately zero after a unit step. (The heading is not printed if the problem is linear.) LPobjective or QPobjective the QP objective function after the present iteration. In elastic mode, the heading is changed to Elastic QPobj. In either case, the value printed decreases monotonically. +SBS is the variable jq selected by the pricing operation to be added to the superbasic set. -SBS is the superbasic variable chosen to become nonbasic. -BS is the basis variable removed (if any) to become nonbasic. Pivot if column ${a}_{q}$ replaces the $r$th column of the basis $B$, Pivot is the $r$th element of a vector $y$ satisfying $By={a}_{q}$. Wherever possible, Step is chosen to avoid extremely small values of Pivot (since they cause the basis to be nearly singular). 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 upper-triangular matrix. Further nonzeros are added to L when various columns of $B$ are later replaced. As columns of $B$ are replaced, the matrix $U$ is maintained explicitly (in sparse form). The value of L will steadily increase, whereas the value of U may fluctuate up or down. Thus the value of L+U may fluctuate up or down (in general, it will tend to increase). ncp is the number of compressions required to recover storage in the data structure for $U$. This includes the number of compressions needed during the previous basis factorization. nS is the current number of superbasic variables. (The heading is not printed if the problem is linear.) condHz see Section 9.1.5. (The heading is not printed if the problem is linear.)

#### 9.1.7Summary

Once the solver finishes, a detailed summary is produced:
``` Problem                      Minimize
Max x                       2 1.4E+00   Max pi                      2 1.9E+01
Nonlinear constraint violn    5.2E-08

-------------------------------------------------------------------------------
Status: converged, an optimal solution was found
-------------------------------------------------------------------------------

Value of the objective               1.90012E+00
Primal infeasibility                 5.24140E-08
Dual infeasibility                   5.19087E-06
Penalty parameter                    9.04109E+00
Major iterations                              11
Minor iterations                              11
Objective function calls                      18
Constraint function calls                     16
Constraint Jacobian calls                     10
Hessian evaluations                            0
-------------------------------------------------------------------------------```
It starts with a status line of the overall result, followed by the final objective value as well as the primal and dual measures of infeasibility and the final value of the penalty parameter. If ${\mathbf{Print Level}}\text{​ and ​}{\mathbf{Monitoring Level}}\ge 2$, it will additionally report major and minor iteration count, objective or constraint function calls. Optionally, if ${\mathbf{Stats Time}}=\mathrm{YES}$ or $\mathrm{CPU}$, some statistics on timings is displayed. It might look as follows:
``` Timing                                        Seconds         Average
Total                                         18.32                  100.0%
Total time spent in pre-solve                  1.09                    5.9%
Total time spent in solver                    17.23                   94.1%
Time in objective function                     0.31            0.02    1.7%
Time in objective gradient                     1.24            0.12    6.8%
Time in constraint function                    0.45            0.03    2.5%
Time in constraint Jacobian                    2.68            0.27   14.6%
Time in Hessian callback                       0.00            0.00    0.0%```

#### 9.1.8Solution

If ${\mathbf{Print Solution}}=\mathrm{YES}$, the values of the primal variables x are printed, furthermore if the problem is constrained, also the dual variables (Lagrangian multipliers) u and the constraint bounds are reported. It might look as follows:
``` Primal variables:
idx   Lower bound       Value       Upper bound
1       -inf      -7.06221E-02         inf
2       -inf       1.41245E+00         inf
...

Box bounds dual variables:
idx   Lower bound       Value       Upper bound       Value
1       -inf       0.00000E+00         inf       2.92732E-06
2       -inf       0.00000E+00         inf       5.31212E-08
...

idx   Lower bound       Value       Upper bound       Value
1       -inf       0.00000E+00    0.00000E+00   -0.00000E+00
...

Nonlinear constraints dual variables:
idx   Lower bound       Value       Upper bound       Value
1   2.00000E+00    0.00000E+00    2.00000E+00    1.90001E+01
2   4.00000E+00    5.00000E+00    4.00000E+00    0.00000E+00
...```

#### 9.1.9Basis Factorization Statistics

If ${\mathbf{Print Level}}\text{​ or ​}{\mathbf{Monitoring Level}}\ge 4$, the first seven items listed below are printed whenever the basis $B$ or the rectangular matrix is factorized before finding the solution of the next QP subproblem. 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$.
If ${\mathbf{Print Level}}\text{​ or ​}{\mathbf{Monitoring Level}}\ge 4$, the same items are printed during the QP solution whenever the current $B$ is factorized.
A heading is printed containing the items described below.
Label Description
Factor the number of factorizations since the start of the run.
Demand a code giving the reason for the present factorization, as follows:
Code Meaning
0 First $LU$ factorization.
1 The number of updates reached the allowed limit.
2 The nonzeros in the updated factors have increased significantly.
7 Not enough storage to update factors.
10 Row residuals are too large.
11 Ill-conditioning has caused inconsistent results.
Itn is the current minor iteration number.
Nonlin is the number of nonlinear variables in the current basis $B$.
Linear is the number of linear variables in $B$.
Slacks is the number of slack variables in $B$.
B, BR, BS or BT factorize is the type of $LU$ factorization.
 B periodic factorization of the basis $B$. BR more careful rank-revealing factorization of $B$ using threshold rook pivoting. This occurs mainly at the start, if the first basis factors seem singular or ill-conditioned. Followed by a normal B factorize. BS ${B}_{S}$ is factorized to choose a well-conditioned $B$ from the current . Followed by a normal B factorize. BT same as BS except the current $B$ is tried first and accepted if it appears to be not much more ill-conditioned than after the previous BS factorize.
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(c-1\right)\left(r-1\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. Depending on the $LU$ pivoting strategy, this heading may change to MerCP or 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 upper-triangular matrix. Further nonzeros are added to L when various columns of $B$ are later replaced. As columns of $B$ are replaced, the matrix $U$ is maintained explicitly (in sparse form). The value of L will steadily increase, whereas the value of U may fluctuate up or down. Thus the value of L+U may fluctuate up or down (in general, it will tend to increase).
Cmpressns is the number of times the data structure holding the partially factored matrix needed to be compressed to recover unused storage. Ideally this number should be zero. If it is more than $3$ or $4$, the amount of workspace available to e04srf 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$.
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$.
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.)
Agrwth is the ratio $\mathtt{Akmax}/\mathtt{Amax}$. Values much larger than $100$ (say) indicate instability.
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 element of $PUQ$.
DUmin is the smallest diagonal element 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).

### 9.2Accessing Slack Variables

e04srf reformulates any inequality constraint in the model into an equality constraint. This is done implicitly by adding artificial slack variables to the problem, see Section 11.1.
e04srf stores in the handle under the label ‘SLACK’ the final value of the slack variables. It also retrieves this information from the handle when a warm start is requested using the optional parameter ${\mathbf{SSQP Start Type}}=\mathrm{WARM}$, see Section 9.4.
Slack variables can be stored or retrieved from the handle with e04rxf using ${\mathbf{cmdstr}}=\text{'SLACK'}$ and ${\mathbf{lrarr}}={\mathbf{stats}}\left(30\right)$, where stats$\left(30\right)$ is the correct number of slack variables and can be retrieved from a previous call to e04srf.
Note: any entry of the slack variables array with an index greater than the amount of defined variables and constraints will refer to the artificial slack variable associated with the objective function, this may be the case when linear variables are identified in the objective function.

### 9.3Retrieving and Storing a Basis

A basis refers to a partitioning of the primal and slack variables. This partitioning plays a fundamental role in the underlying LP and QP algorithms of e04srf, Section 11.3 provides a brief description.
e04srf stores in the handle under the label ‘BASIS’ (or ‘WARM START BASIS’) the final state of the primal and slack variables. It also retrieves this information from the handle when a warm start is requested using optional parameter ${\mathbf{SSQP Start Type}}=\mathrm{WARM}$, see Section 9.4.
The stored integer array is of length ${\mathbf{nvar}}+{\mathbf{stats}}\left(30\right)$ and the values describe the state of the primal variables x and the slacks as follows:
BASIS($\mathbit{i}$) State of variable $\mathbit{i}$ Usual value Note
-1 Ignored None of the enabled variables can have this value. It indicates that the basis element is ignored by the solver. Ignored elements occur when a primal variable or constraint is marked by you as disabled (see e04tcf) or presolve dropped a fixed constraint or primal variable.
0 Nonbasic lower bound
1 Nonbasic upper bound
2 Superbasic between lower and
upper bounds
3 Basic between lower and
upper bounds
4 Nonbasic lower bound Variable is ignored during the crash procedure and its associated slack variable needs to be set to the lower bound.
5 Nonbasic upper bound Variable is ignored during the crash procedure and its associated slack variable needs to be set to the upper bound.
The basis can be stored or retrieved from the handle with e04rwf using ${\mathbf{cmdstr}}=\text{'BASIS'}$ or $\text{'WARM START BASIS'}$ and ${\mathbf{liarr}}={\mathbf{nvar}}+{\mathbf{stats}}\left(30\right)$, where ${\mathbf{stats}}\left(30\right)$ is the correct number of slack variables and can be retrieved from a previous call to e04srf.
Note: basic and superbasic variables may be outside their bounds by as much as the SSQP Minor Feasibility Tol tolerance. If scaling is specified (using optional parameter SSQP Scale Option), 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 ‘Primal infeasibility’ line printed at the end of the solve, see Section 9.1.7.

### 9.4Warm Starting

Warm starting a problem refers to providing a starting point x and additional information used by the solver to start the optimization process, for example, by providing information on which variables are suberbasic or active and thus hinting on the possible final active-set or providing a good initial guess for the final values of the Lagrange multipliers.
In order to warm start e04srf, it is necessary to
1. (i)provide on the call to e04srf the initial guess ${x}_{0}$;
2. (ii)provide on the call to e04srf the initial guess for the Lagrange multipliers u. If nnzu$>0$ then the solver will access array u and so it must be provided. If nothing is known about the multipliers then u should be set to zero in the call to e04srf;
3. (iii)store in the handle (under the label ‘SLACK’) the slack variables array. See Section 9.2;
4. (iv)store in the handle (under the label ‘BASIS’ or ‘WARM START BASIS’) a valid basis vector of length nvar $+$ stats$\left(30\right)$. See Section 9.3;
5. (v)request the solver to attempt a warm start by setting optional parameter ${\mathbf{SSQP Start Type}}=\mathrm{WARM}$.
If optional parameter ${\mathbf{SSQP Start Type}}=\mathrm{WARM}$ but e04srf does not find the required information or it is inconsistent, then it will revert to a cold start.
Note: e04srf at exit (if the information is available) stores the slack variables and basis arrays into the handle under the labels ‘SLACK’ and ‘BASIS’. A next call to e04srf with the same handle along with ${\mathbf{SSQP Start Type}}=\mathrm{WARM}$, and the latest x and u, should trigger a warm start successfully. It will also notify the source of the warm starting information with a message similar to:
``` Warm start information loaded successfully from handle.
Handle [WARM START BASIS] data origin: solver
Handle [SLACK] data origin           : user```
Which indicates that the warm start information was successfully loaded. It also informs that you provided the slack variables array (using e04rxf with cmdstr$=$‘SLACK’) and that the basis information was provided by the solver itself, say, from a previous call to e04srf.

### 9.5Internal Changes

Internal changes have been made to this routine as follows:
• At Mark 28.4:
Updates that in general provide overall performance enhancements.
For details of all known issues which have been reported for the NAG Library please refer to the Known Issues.

## 10Example

This example involves the minimization of the nonlinear function
 $f(x) = (x1+x2+x3)2 + 3x3 + 5x4 + cos (0.01x1) - 1 ,$
subject to the bounds
 $0 ≤ x3, 0 ≤ x4,$
to the nonlinear constraints
 $x12 + x22 + x3 = 2, x24 + x4 = 4,$
and the linear constraint
 $2x1 + 4x2 ≥ 0 .$
The initial point is
 $x0 = (1,2,3,4)T ,$
and $f\left({x}_{0}\right)=65$. The optimal solution (to five significant figures) is
 $x* = (−0.070639,1.4124,0.0,0.019934) T , with ​ f(x*) = 1.9001 .$
Note: the variable ${x}_{4}$ appears only linearly in the example. The linear property of the variable can be hinted to the solver using the routine e04rcf, see Section 3.1 and Section 3.1 in e04rcf for further details.

### 10.1Program Text

Program Text (e04srfe.f90)

### 10.2Program Results

Program Results (e04srfe.r)

## 11Algorithmic Details

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

### 11.1Constraints and Slack Variables

Let $Q\left(x\right)$ be the ${m}_{Q}$-dimensional vector made up of the ${m}_{Q}$ quadratic constraints. Furthermore, let $m$ be the total number of components of the nonlinear constraints, $g\left(x\right)$, the quadratic constraints, $Q\left(x\right)$, and the linear constraints, $Bx$ combined, that is, $m={m}_{g}+{m}_{Q}+{m}_{B}$. The upper and lower bounds on those terms define the general constraints of the problem. e04srf 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
 $minimize x,s f(x) subject to ​ ( g(x) Q(x) Bx ) - s = 0 , l ≤ ( x s ) ≤ u .$ (2)
The general constraints become the equalities ${\left(g\left(x\right),Q\left(x\right)\right)}^{\mathrm{T}}-{s}_{N}=0$ and $Bx-{s}_{L}=0$, where ${s}_{L}$ and ${s}_{N}$ are the linear and nonlinear slacks.
Note: see Section 9.2 on how to access the slack variables stored in the handle.

### 11.2Major 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 first-order conditions for optimality. At each iterate ${x}_{k}$ a QP subproblem is used to generate a search direction towards the next iterate ${x}_{k+1}$.
Without loss of generality and for ease of exposition, suppose that the model to optimize consists of linear and nonlinear constraints only, that is, ${m}_{Q}=0$, then the constraints of the subproblem are formed from the linear constraints $Bx-{s}_{L}=0$ and the linearized constraint
 $g(xk) + g′ (xk) (x-xk) - sN = 0 ,$ (3)
where ${g}^{\prime }\left({x}_{k}\right)$ denotes the Jacobian matrix, whose elements are the first derivatives of $g\left(x\right)$ evaluated at ${x}_{k}$. The QP constraints, therefore, comprise the $m$ linear constraints
 $g′ (xk) x - sN = -g(xk) + g′ (xk) xk , Bx -sN - sL = 0 ,$ (4)
where $x$ and $s$ are bounded above and below by $u$ and $l$ as before. If the $m×n$ matrix $A$ and $m$-vector $b$ are defined as
 $A = ( g′ (xk) B ) and b = ( -g(xk) + g′ (xk) xk 0 ) ,$ (5)
then the QP subproblem can be written as
 $minimize x,s Q (x,xk) = gkT (x-xk) + 12 (x-xk)T Hk (x-xk) subject to ​ Ax - s = b , l ≤ ( x s ) ≤ u ,$ (6)
where $\mathcal{Q}\left(x,{x}_{k}\right)$ is a quadratic approximation to a modified Lagrangian function (see Gill et al. (2002)), $\mathcal{g}$ is the gradient of the Lagrangian, and the matrix ${H}_{k}$ is a quasi-Newton approximation to the Hessian of the Lagrangian. A BFGS update is applied after each major iteration. If some of the variables enter the Lagrangian linearly the Hessian will have some zero rows and columns. If the nonlinear variables appear first, then only the ${n}_{\text{NLN}}$ rows and columns of the Hessian need to be approximated, where ${n}_{\text{NLN}}$ is the number of nonlinear variables, see Section 3.1 and Section 3.1 in e04rcf.

### 11.3Minor Iterations

Solving the QP subproblem is itself an iterative procedure. Here, the iterations of the QP solver e04nqf form the minor iterations of the SQP method. e04nqf uses a reduced-Hessian active-set method implemented as a reduced-gradient method. At each minor iteration, the constraints $Ax-s=b$ are partitioned into the form
 $BxB + SxS + NxN = b ,$ (7)
where the basis matrix $\mathcal{B}$ is square and nonsingular, and the matrices $\mathcal{S}$ and $\mathcal{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 $Ax-s=b$. The number of superbasic variables (${n}_{S}$, say), therefore, indicates the number of degrees of freedom remaining after the constraints have been satisfied. In broad terms, ${n}_{S}$ is a measure of how nonlinear the problem is. In particular, ${n}_{S}$ will always be zero for LP problems.
If it appears that no improvement can be made with the current definition of $\mathcal{B}$, $\mathcal{S}$ and $\mathcal{N}$, a nonbasic variable is selected to be added to $\mathcal{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 $Ax-s=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 ${\mathcal{Q}}^{\prime }-{A}^{\mathrm{T}}\pi$, where ${\mathcal{Q}}^{\prime }$ is the gradient of the QP objective (6), 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({\stackrel{^}{x}}_{k},{\stackrel{^}{s}}_{k},{\stackrel{^}{\pi }}_{k}\right)$ is found by relaxing these conditions.

### 11.4The Merit Function

After a QP subproblem has been solved, new estimates of the solution are computed using a linesearch on the augmented Lagrangian merit function
 $M (x,s,π) = f(x) - πT (g(x)-sN) + 12 (g(x)-sN) T D (f(x)-sN) ,$ (8)
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({\stackrel{^}{x}}_{k},{\stackrel{^}{s}}_{k},{\stackrel{^}{\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
 $( xk+1 sk+1 πk+1 ) = ( xk sk πk ) + αk ( x^k - xk s^k - sk π^k - πk )$ (9)
gives a sufficient decrease in the merit function $\mathcal{M}$. When necessary, the penalties in $D$ are increased by the minimum-norm 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)).

### 11.5Treatment of Constraint Infeasibilities

e04srf makes explicit allowance for infeasible constraints. First, infeasible linear constraints are detected by solving the linear program
 $minimize x,v,w eT (v+w) subject to ​ l ≤ ( x Bx - v + w ) ≤ u , v≥0 , w≥0 ,$ (10)
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.)
The linear constraints are infeasible if the optimal solution of (10) has $v\ne 0$ or $w\ne 0$. e04srf 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 nonlinear functions can be safely evaluated.) e04srf proceeds to solve the nonlinear problem as given, using search directions obtained from the sequence of QP subproblems, see (6).
If a QP subproblem proves to be infeasible or unbounded (or if the dual variables $\pi$ for the nonlinear constraints become large), e04srf enters ‘elastic’ mode and thereafter solves the problem
 $minimize x,v,w f(x) + γ eT (v+w) subject to ​ l ≤ ( x g(x) - v + w Bx ) ≤ u , v≥0 , w≥0 ,$ (11)
where $\gamma$ is a non-negative argument (the elastic weight), and $f\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).
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.

## 12Optional Parameters

Several optional parameters in e04srf define choices in the problem specification or the algorithm logic. In order to reduce the number of formal arguments of e04srf 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.
The remainder of this section can be skipped if you wish to use the default values for all optional parameters.
The optional parameters can be changed by calling e04zmf anytime between the initialization of the handle by e04raf and the call to the solver. Modification of the arguments during intermediate monitoring steps is not allowed. Once the solver finishes, the optional parameters can be altered again for the next solve.
The option values may be retrieved by e04znf.
The following is a list of the optional parameters available. A full description of each optional parameter is provided in Section 12.1.

### 12.1Description of the Optional Parameters

For each option, we give a summary line, a description of the optional parameter and details of constraints.
The summary line contains:
• the keywords;
• a parameter value, where the letters $a$, $i$ and $r$ denote options that take character, integer and real values respectively.
• the default value, where the symbol $\epsilon$ is a generic notation for machine precision (see x02ajf).
All options accept the value $\mathrm{DEFAULT}$ to return single options to their default states.
Keywords and character values are case and white space insensitive.
 Defaults
This special keyword may be used to reset all optional parameters to their default values. Any value given with this keyword will be ignored.
 Infinite Bound Size $r$ Default $\text{}={10}^{20}$
This defines the ‘infinite’ bound $\mathit{bigbnd}$ in the definition of the problem constraints. Any upper bound greater than or equal to $\mathit{bigbnd}$ will be regarded as $+\infty$ (and similarly any lower bound less than or equal to $-\mathit{bigbnd}$ will be regarded as $-\infty$). Note that a modification of this optional parameter does not influence constraints which have already been defined; only the constraints formulated after the change will be affected.
Constraint: ${\mathbf{Infinite Bound Size}}\ge 1000$.
 Monitoring File $i$ Default $=-1$
If $i\ge 0$, the unit number for the secondary (monitoring) output. If set to $-1$, no secondary output is provided. The following information is output to the unit:
• a listing of the optional parameters;
• problem statistics, the iteration log and the final status as set by Monitoring Level;
• the solution if set by Print Solution.
Constraint: ${\mathbf{Monitoring File}}\ge -1$.
 Monitoring Level $i$ Default $=4$
This parameter sets the amount of information detail that will be printed by the solver to the secondary output. The meaning of the levels is the same as with Print Level.
Constraint: $0\le {\mathbf{Monitoring Level}}\le 5$.
 Print File $i$ Default $=\text{advisory message unit number}$
If $i\ge 0$, the unit number for the primary output of the solver. If ${\mathbf{Print File}}=-1$, the primary output is completely turned off independently of other settings. The default value is the advisory message unit number as defined by x04abf at the time of the optional parameters initialization, e.g., at the initialization of the handle. The following information is output to the unit:
• a listing of optional parameters if set by Print Options;
• problem statistics, the iteration log and the final status from the solver as set by Print Level;
• the solution if set by Print Solution.
Constraint: ${\mathbf{Print File}}\ge -1$.
 Print Level $i$ Default $=2$
This parameter defines how detailed information should be printed by the solver to the primary output.
$\mathbit{i}$ Output
$0$ No output from the solver
$1$ Only the final status and the objective value
$2$ Problem statistics, one line per iteration showing the progress of the solution with respect to the convergence measures, final status and statistics
$3$ As level $2$ and additionally inner QP iteration log is reported every SSQP Print Frequency
$4,5$ As level $3$ but further details of each minor and major iteration are presented
Constraint: $0\le {\mathbf{Print Level}}\le 5$.
 Print Options $a$ Default $=\mathrm{YES}$
If ${\mathbf{Print Options}}=\mathrm{YES}$, a listing of optional parameters will be printed to the primary output.
Constraint: ${\mathbf{Print Options}}=\mathrm{YES}$ or $\mathrm{NO}$.
 Print Solution $a$ Default $=\mathrm{NO}$
If ${\mathbf{Print Solution}}=\mathrm{YES}$, the final values of the solution vector are printed on the primary and secondary outputs.
Constraint: ${\mathbf{Print Solution}}=\mathrm{YES}$ or $\mathrm{NO}$.
 SSQP Crash Option $a$ Default $\text{}=\mathrm{TRIPLE-PASS}$
Except on restarts ${\mathbf{SSQP Start Type}}=\mathrm{COLD}$, a crash procedure is used to select an initial basis $\mathcal{B}$ (see (7) in Section 11.3) from certain rows and columns of the constraint matrix $\left(\begin{array}{cc}A& -I\end{array}\right)$. The SSQP Crash Option determines which rows and columns of $A$ (see (5) in Section 11.2) are eligible initially, and how many times crash is called. Columns of $-I$ are used to pad the basis where necessary.
$\mathbit{a}$ Meaning
$\mathrm{ONLY-SLACKS}$ The initial basis contains only slack variables: $\mathcal{B}=I$ (see (7) in Section 11.3).
$\mathrm{SINGLE-PASS}$ crash is called once, looking for a triangular basis in all rows and columns of $A$.
$\mathrm{DOUBLE-PASS}$ crash is called twice (if there are nonlinear constraints). The first call looks for a triangular basis in linear rows, and the iteration proceeds with simplex iterations until the linear constraints are satisfied. The Jacobian is then evaluated for the first major iteration and crash is called again to find a triangular basis in the nonlinear rows (retaining the current basis for linear rows).
$\mathrm{TRIPLE-PASS}$ crash is called up to three times (if there are nonlinear constraints). The first two calls treat linear equalities and linear inequalities separately. As before, the last call treats nonlinear rows before the first major iteration.
If ${\mathbf{SSQP Crash Option}}\ne \mathrm{ONLY-SLACKS}$, certain slacks on inequality rows are selected for the basis first, furthermore if ${\mathbf{SSQP Crash Option}}=\mathrm{DOUBLE-PASS}$ or $\mathrm{TRIPLE-PASS}$, numerical values are used to exclude slacks that are close to a bound. crash then makes several passes through the columns of $A$, searching for a basis matrix that is essentially triangular. A column is assigned to pivot on a particular row if the column contains a suitably large element in a row that has not yet been assigned (the pivot elements ultimately form the diagonals of the triangular basis. For remaining unassigned rows, slack variables are inserted to complete the basis.
Note: each crash pass is counted as a minor iteration..
Constraint: ${\mathbf{SSQP Crash Option}}=\mathrm{ONLY-SLACKS}$, $\mathrm{SINGLE-PASS}$, $\mathrm{DOUBLE-PASS}$ or $\mathrm{TRIPLE-PASS}$.
 SSQP Crash Tolerance $r$ Default $=0.1$
The SSQP Crash Tolerance $r$ allows the starting procedure crash to ignore certain small nonzeros in each column of $A$ (see (5) in Section 11.2). If ${a}_{\mathrm{max}}$ is the largest element in column $j$, other nonzeros ${a}_{ij}$ in the column are ignored if $|{a}_{ij}|\le {a}_{\mathrm{max}}×r$. To be meaningful, $r$ needs to be in the range $0\le r\le 1$.)
When $r>0$, the basis obtained by crash may not be strictly triangular, but it is likely to be nonsingular and almost triangular. The intention is to obtain a starting basis containing more columns of $A$ and fewer (arbitrary) slacks (columns of $-I$). A feasible solution may be reached sooner on some problems.
When $r=-1$ (or $r<0$) the solver will automatically choose an appropriate tolerance.
Constraint: SSQP Crash Tolerance$\ge 0$ or SSQP Crash Tolerance$=-1.0$.
 SSQP Estimate Derivatives $a$ Default $\text{}=\mathrm{NO}$
This option specifies which nonlinear function gradients (objective and nonlinear constraints) are known analytically and will be supplied to the solver by you subroutines objgrd and congrd. This option has the following alternatives:
$\mathbit{a}$ Meaning
$\mathrm{NO}$ (default) Indicates to the solver that you are providing full derivative information, that is, all objective and constraint gradients are known.
$\mathrm{ONLYOBJ}$ Indicates that some or all components of the objective gradient will be missing but all constraint gradients are known.
$\mathrm{ONLYCON}$ Conversely to $\mathrm{ONLYOBJ}$, this indicated that some or all of the constraint gradients will be missing but the objective gradient is known.
$\mathrm{YES}$ Indicates that the user is not providing derivatives.
The value ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{NO}$ should be used whenever possible. It is the most reliable and will usually be the most efficient.
If ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{YES}$ or $\mathrm{ONLYOBJ}$, e04srf will estimate the missing components of the objective gradient, using finite differences. This may simplify the coding of subroutine objgrd. However, it could increase the total run-time substantially (since a special call to objfun is required for each missing element), and there is less assurance that an acceptable solution will be located. If the nonlinear variables are not well scaled, it may be necessary to specify a nonstandard difference interval changing SSQP Finite Diff Ctrl Interval or SSQP Finite Diff Interval.
Conversely, if ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{YES}$ or $\mathrm{ONLYCON}$, e04srf will estimate missing elements of the Jacobian. For each column of the Jacobian, one call to confun is needed to estimate all missing elements in that column, if any.
At times, central differences are used rather than forward differences. Twice as many calls to objfun and confun are then needed. (This is not under your control.)
Note: if ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{NO}$, $\mathrm{ONLYCON}$ or $\mathrm{ONLYOBJ}$ and the user does not provide one or more of the required gradient (or Jacobian) elements then these will be assumed to be constant, this might not be the desired behaviour. Under this circumstance, e04srf will inform the user with a message.
Constraint: ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{NO}$, $\mathrm{ONLYOBJ}$, $\mathrm{ONLYCON}$ or $\mathrm{YES}$.
 SSQP Finite Diff Ctrl Interval $r$ Default $\approx 1.7×{\epsilon }^{0.27}$
When ${\mathbf{SSQP Estimate Derivatives}}=\mathrm{YES}$, $\mathrm{ONLYOBJ}$ or $\mathrm{ONLYCON}$, the central-difference interval $r$ is used near an optimal solution to obtain more accurate (but more expensive) estimates of gradients. Twice as many function evaluations are required compared to forward differencing. The interval used for the $j$th variable is ${h}_{j}=r\left(1+|{x}_{j}|\right)$. The resulting derivative estimates should be accurate to $\mathit{O}\left({r}^{2}\right)$, unless the functions are badly scaled.
Constraint: $0.
 SSQP Finite Diff Interval $r$ Default $\approx 1.7×{\epsilon }^{0.4}$
This alters the interval $r$ used to estimate gradients by forward differences. It does so in the following circumstances:
• in the ‘cheap’ phase of verifying the problem derivatives, when ${\mathbf{Verify Derivatives}}=\mathrm{CHEAP}$;
• for verifying the problem derivatives (${\mathbf{Verify Derivatives}}=\mathrm{YES}$);
• for estimating missing derivatives.
In all cases, a derivative with respect to ${x}_{j}$ is estimated by perturbing that component of $x$ to the value ${x}_{j}+r\left(1+|{x}_{j}|\right)$, and then evaluating the objective or constraint functions at the perturbed point. The resulting gradient estimates should be accurate to $\mathit{O}\left(r\right)$ unless the functions are badly scaled. Judicious alteration of $r$ may sometimes lead to greater accuracy.
Constraint: $0.
 SSQP Function Precision $r$ Default $\approx 1.7×{\epsilon }^{0.8}$
The relative function precision ${\epsilon }_{r}$ is intended to be a measure of the relative accuracy with which the nonlinear functions can be computed. For example, if $f\left(x\right)$ is computed as $1000.56789$ for some relevant $x$ and if the first $6$ significant digits are known to be correct, the appropriate value for ${\epsilon }_{r}$ would be $\text{1.0E−6}$.
Ideally, the functions $f\left(x\right)$ and ${g}_{i}\left(x\right)$ should have a magnitude of order $1$. If all functions are substantially less than $1$ in magnitude, ${\epsilon }_{r}$ should be the absolute precision. For example, if $f\left(x\right)=\text{1.23456789E−4}$ at some point and if the first $6$ significant digits are known to be correct, the appropriate value for ${\epsilon }_{r}$ would be $\text{1.0E−10}$.)
The default value of ${\epsilon }_{r}$ is appropriate for simple analytic functions.
In some cases the function values will be the result of extensive computation, possibly involving a costly iterative procedure that can provide few digits of precision. Specifying an appropriate SSQP Function Precision may lead to savings, by allowing the linesearch procedure to terminate when the difference between function values along the search direction becomes as small as the absolute error in the values.
Constraint: $0<{\epsilon }_{r}<0.1$.
 SSQP Hessian $a$ Default $\text{}=\mathrm{AUTO}$
This option selects the method for storing and updating the approximate Hessian. (e04srf uses a quasi-Newton approximation to the Hessian of the Lagrangian. A BFGS update is applied after each major iteration.)
If ${\mathbf{SSQP Hessian}}=\mathrm{FULL-MEMORY}$ is specified, the approximate Hessian is treated as a dense matrix and the BFGS updates are applied explicitly. This option is most efficient when the number of variables $n$ is not too large (say, less than $75$). In this case, the storage requirement is fixed and one can expect $1$-step Q-superlinear convergence to the solution.
${\mathbf{SSQP Hessian}}=\mathrm{LIMITED-MEMORY}$ should be used on problems where $n$ is very large. In this case, a limited-memory procedure is used to update a diagonal Hessian approximation ${H}_{r}$ a limited number of times. Updates are accumulated as a list of vector pairs. They are discarded at regular intervals after ${H}_{r}$ has been reset to their diagonal, see optional parameter SSQP Hessian Updates.
If ${\mathbf{SSQP Hessian}}=\mathrm{AUTO}$ is specified, then the quasi-Newton method used is chosen based on a heuristic involving the number of nonlinear variables (see Section 3.1) in the problem.
Constraint: ${\mathbf{SSQP Hessian}}=\mathrm{AUTO}$, $\mathrm{FULL-MEMORY}$ or $\mathrm{LIMITED-MEMORY}$.
 SSQP Hessian Updates $i$ Default $\text{}=10$
If optional parameter ${\mathbf{SSQP Hessian}}=\mathrm{LIMITED-MEMORY}$ is in effect and $i$ BFGS updates have already been carried out, all but the diagonal elements of the accumulated updates are discarded and the updating process starts again.
Broadly speaking, the more updates stored, the better the quality of the approximate Hessian. However, the more vectors stored, the greater the cost of each QP iteration. The default value is likely to give a robust algorithm without significant expense, but faster convergence can sometimes be obtained with significantly fewer updates (e.g., $i=5$).
Constraint: $0<{\mathbf{SSQP Hessian Updates}}$.
 SSQP Iteration Limit $i$ Default $\text{}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(10000,10\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{nvar}},m\right)\right)$
The value of $i$ specifies the maximum number of minor iterations allowed (i.e., iterations of the simplex method or the QP algorithm), summed over all major iterations. (See also the description of the optional parameter SSQP Minor Iteration Limit.) For the default value, $m$ is the number of linear, quadratic and nonlinear constraints defined in the model.
Constraint: $0<{\mathbf{SSQP Iteration Limit}}$.
 SSQP Major Feasibility Tol $r$ Default $\approx {\epsilon }^{0.38}$
This tolerance, $r$, specifies how accurately the nonlinear constraints should be satisfied. The default value is appropriate when the linear and nonlinear constraints contain data to about that accuracy.
Let ${v}_{\mathrm{max}}$ be the maximum nonlinear constraint violation, normalized by the size of the solution, which is required to satisfy
 $vmax=maxi vi ​ ​/​ ​ ‖x‖ ≤ r ,$ (12)
where ${v}_{\mathit{i}}$ is the violation of the $\mathit{i}$th nonlinear constraint, for $\mathit{i}=1,2,\dots ,{\mathbf{ncnln}}$.
In the major iteration log (see Section 9.1), ${v}_{\mathrm{max}}$ appears as the quantity labelled ‘Feasible’. If some of the problem functions are known to be of low accuracy, a larger SSQP Major Feasibility Tol may be appropriate.
Constraint: $0.
 SSQP Major Iteration Limit $i$ Default $\text{}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1000,3\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{nvar}},m\right)\right)$
This is the maximum number of major iterations allowed. It is intended to guard against an excessive number of linearizations of the constraints. For the default value, $m$ is the number of linear, quadratic and nonlinear constraints defined in the model.
Constraint: $0<{\mathbf{SSQP Major Iteration Limit}}$.
 SSQP Major Optimality Tol $r$ Default $\approx {\epsilon }^{0.38}$
This tolerance, $r$, specifies the final accuracy of the dual variables. On successful termination, e04srf will have computed a solution $\left(x,s,\pi \right)$ such that
 $cmax=maxj cj ​ ​/​ ​ ‖π‖ ≤ r ,$ (13)
where ${c}_{\mathit{j}}$ is an estimate of the complementarity slackness for variable $\mathit{j}$, for $\mathit{j}=1,2,\dots ,n+m$. Here $m$ represents the number of constraints. The values ${c}_{i}$ are computed from the final QP solution using the reduced gradients ${d}_{j}={g}_{j}-{\pi }^{\mathrm{T}}{a}_{j}$ (where ${g}_{j}$ is the $j$th component of the objective gradient, ${a}_{j}$ is the associated column of the constraint matrix $\left(\begin{array}{cc}A& -I\end{array}\right)$, and $\pi$ is the set of QP dual variables):
 $cj = { - dj min( xj - lj ,1) if ​ dj≥0 ; -dj min( uj - xj ,1) if ​ dj<0 . )$ (14)
In the major iteration log, ${c}_{\mathrm{max}}$ appears as the quantity labelled ‘Optimal’, See section Section 9.1.5.
Constraint: $0.
 SSQP Minor Feasibility Tol $r$ Default $\approx {\epsilon }^{0.38}$
e04srf tries to ensure that all variables eventually satisfy their upper and lower bounds to within this tolerance, $r$. This includes slack variables. Hence, general linear constraints should also be satisfied to within $r$.
Feasibility with respect to nonlinear constraints is judged by the optional parameter SSQP Major Feasibility Tol (not by SSQP Minor Feasibility Tol).
If the bounds and linear constraints cannot be satisfied to within $r$, the problem is declared infeasible. If rinfo$\left(12\right)$ is quite small, it may be appropriate to raise $r$ by a factor of $10$ or $100$. Otherwise, some error in the data should be suspected.
Nonlinear functions will be evaluated only at points that satisfy the bounds and linear constraints, except perhaps in presolve where the solver may require to evaluate confun, objgrd or congrd at a bound-feasible point. If there are regions where a function is undefined, every attempt should be made to eliminate these regions from the problem.
For example, if $f\left(x\right)=\sqrt{{x}_{1}}+\mathrm{log}\left({x}_{2}\right)$, it is essential to place lower bounds on both variables. If $r=\text{1.0E−6}$, the bounds ${x}_{1}\ge {10}^{-5}$ and ${x}_{2}\ge {10}^{-4}$ might be appropriate. (The log singularity is more serious. In general, keep $x$ as far away from singularities as possible.)
If ${\mathbf{SSQP Scale Option}}=\mathrm{ALL}$ or $\mathrm{VARS AND LC}$, feasibility is defined in terms of the scaled problem (since it is then more likely to be meaningful).
In reality, e04srf uses $r$ as a feasibility tolerance for satisfying the bounds on $x$ and $s$ in each QP subproblem. If the sum of infeasibilities cannot be reduced to zero, the QP subproblem is declared infeasible. e04srf then activates elastic mode for the rest of the solving process (with only the linearized nonlinear constraints defined to be elastic). For further details see Section 11.5.
Constraint: $0.
 SSQP Minor Iteration Limit $i$ Default $\text{}=500$
If the number of minor iterations for the optimality phase of the QP subproblem exceeds $i$, then all nonbasic QP variables that have not yet moved are frozen at their current values and the reduced QP is solved to optimality.
Note that more than $i$ minor iterations may be necessary to solve the reduced QP to optimality. These extra iterations are necessary to ensure that the terminated point gives a suitable direction for the linesearch.
In the major iteration log (see Section 9.1) the flag t at the end of a line indicates that the corresponding QP was artificially terminated using the limit $i$.
Compare with the optional parameter SSQP Iteration Limit, which defines an independent absolute limit on the total number of minor iterations (summed over all QP subproblems).
Constraint: $0<{\mathbf{SSQP Minor Iteration Limit}}$.
 SSQP Monitor Frequency $i$ Default $=0$
This parameter specifies the frequency on which to call the monitor routine monit. If zero, the monitor routine will not be called.
Constraint: ${\mathbf{SSQP Monitor Frequency}}\ge 0$.
 SSQP Penalty Parameter $r$ Default $\text{}=0.0$
Defines the initial value of the penalty parameter and sets $D=\text{diag}\left(r\right)$ in (8).
Constraint: $0.0\le {\mathbf{SSQP Penalty Parameter}}$.
 SSQP Print Frequency $i$ Default $\text{}=1$
Defines how frequently to print log information for the inner iterations. Only relevant if Print Level$\ge 3$ (with Print File$\ne -1$) or Monitoring Level$\ge 3$ (with Monitoring File$\ne -1$)
Constraint: $0\le {\mathbf{SSQP Print Frequency}}$.
 SSQP Scale Option $a$ Default $\text{}=\mathrm{NONE}$
Three scale options are available as follows:
$\mathbit{a}$ Meaning
$\mathrm{NONE}$ No scaling. This is recommended if it is known that $x$ and the constraint matrix never have very large elements (say, larger than $100$).
The constraints and variables are scaled by an iterative procedure that attempts to make the matrix coefficients as close as possible to $1.0$ (see Fourer (1982)). This will sometimes improve the performance of the solution procedures.
$\mathrm{ALL}$ The constraints and variables are scaled by the iterative procedure. Also, a certain additional scaling is performed that may be helpful if the right-hand side $b$ or the solution $x$ is large. This takes into account columns of $\left(\begin{array}{cc}A& -I\end{array}\right)$ that are fixed or have positive lower bounds or negative upper bounds.
Constraint: ${\mathbf{SSQP Scale Option}}=\mathrm{NONE}$, $\mathrm{VARS AND LC}$ or $\mathrm{ALL}$.
 SSQP Scale Print $a$ Default $\text{}=\mathrm{NO}$
Defines whether to print the row scales $r\left(i\right)$ and column scales $c\left(j\right)$. The scaled matrix coefficients are ${\overline{a}}_{ij}={a}_{ij}c\left(j\right)/r\left(i\right)$ and the scaled bounds on the variables and slacks are ${\overline{l}}_{j}={l}_{j}/c\left(j\right)$, ${\overline{u}}_{j}={u}_{j}/r\left(j-n\right)$ if $j>n$.
Constraint: ${\mathbf{SSQP Scale Print}}=\mathrm{YES}$ or $\mathrm{NO}$.
 SSQP Scale Tolerance $r$ Default $\text{}=0.9$
This parameter affects how many passes might be needed through the constraint matrix. On each pass, the scaling procedure computes the ratio of the largest and smallest nonzero coefficients in each column:
 $ρj=maxj |aij| / mini |aij| (aij≠0) .$
If $\underset{j}{\mathrm{max}}\phantom{\rule{0.25em}{0ex}}{\rho }_{j}$ is less than $r$ times its previous value, another scaling pass is performed to adjust the row and column scales. Raising $r$ from $0.9$ to $0.99$ (say) usually increases the number of scaling passes through $A$. At most $10$ passes are made. The value of $r$ should lie in the range $0.
Constraint: $\epsilon \le {\mathbf{SSQP Scale Tolerance}}$.
 SSQP Start Type $a$ Default $\text{}=\mathrm{COLD}$
Defines whether to perform a cold or warm start. If warm start data is not provided or is considered to have an unexpected size or content, then the solver will revert to perform a cold start on the problem. See Section 9.4 on how to correctly warm start a problem.
Constraint: ${\mathbf{SSQP Start Type}}=\mathrm{COLD}$ or $\mathrm{WARM}$.
 SSQP Superbasics Limit $i$ Default $\text{}=-1$
This option places a limit on the storage allocated for superbasic variables. Ideally, $i$ should be set slightly larger than the ‘number of degrees of freedom’ expected at an optimal solution. For nonlinear problems, the number of degrees of freedom is often called the ‘number of independent variables’. Normally, $i$ need not be greater than $n+1$.
If $i=-1$ (the default), then e04srf automatically sets the limit to the number of nonlinear variables. For many problems, $i$ may be considerably smaller than $n$. This will save storage if $n$ is very large.
If the limit set is too low then the solver will terminate with ${\mathbf{ifail}}={\mathbf{27}}$ error message.
Constraint: $-1\le {\mathbf{SSQP Superbasics Limit}}$.
 Stats Time $a$ Default $=\mathrm{NO}$
This parameter allows you to turn on the timings of various parts of the algorithm to give a better overview of where most of the time is spent. This might be helpful for a choice of different solving approaches. It is possible to choose between CPU and wall clock time. Choice $\mathrm{YES}$ is equivalent to $\mathrm{WALL CLOCK}$.
Constraint: ${\mathbf{Stats Time}}=\mathrm{YES}$, $\mathrm{NO}$, $\mathrm{CPU}$ or $\mathrm{WALL CLOCK}$.
 Task $a$ Default $=\mathrm{MINIMIZE}$
This parameter specifies the required direction of the optimization. If , the objective function (if set) is ignored and the algorithm stops as soon as a feasible point is found. If no objective function is set, Task reverts to automatically.
Constraint: ${\mathbf{Task}}=\mathrm{MINIMIZE}$, $\mathrm{MAXIMIZE}$ or $\mathrm{FEASIBLE POINT}$.
 Time Limit $r$ Default $\text{}={10}^{6}$
This parameter specifies a limit in seconds that the solver can use to solve one problem. If during the convergence check this limit is exceeded, the solver will terminate with ${\mathbf{ifail}}={\mathbf{23}}$ error message.
Constraint: ${\mathbf{Time Limit}}>0$.
 Verify Derivatives $a$ Default $=\mathrm{NO}$
This parameter specifies whether the routine should perform numerical checks on the consistency of the user-supplied gradient functions objgrd and congrd. If any discrepancies are found, ${\mathbf{ifail}}={\mathbf{26}}$ is returned. It is recommended that such checks are enabled when first developing the formulation of the problem, however, the verification process results in a significant increase in the number of the function evaluations and thus it shouldn't be used in production code.
a Meaning
$\mathrm{CHEAP}$ Only a ‘cheap’ test will be performed.
$\mathrm{YES}$ Individual columns of the problem Jacobian, as well as individual gradient elements, will be checked (with a more reliable test). A key of the form OK or Bad? indicates whether or not each component appears to be correct.
$\mathrm{NO}$ Derivative checking is disabled.
${\mathbf{Verify Derivatives}}=\mathrm{YES}$ should be specified whenever a new objgrd or congrd is being developed.
Constraint: ${\mathbf{Verify Derivatives}}=\mathrm{CHEAP}$, $\mathrm{YES}$ or $\mathrm{NO}$.