e04st: FL CL CPP AD PY MB

NAG FL Interface
e04stf (handle_solve_ipopt)

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.

Keyword Search:

NAG Library Manual, Mark 30.2

Interfaces: FL CL CPP AD PY MB

NAG FL Interface Introduction

E04 (Opt) Chapter Contents

E04 (Opt) Chapter Introduction

e04st: FL CL CPP AD PY MB

▸▿ Contents

1 Purpose

2 Specification

▸▿ 3 Description

3.1 Structure of the Lagrange Multipliers

4 References

5 Arguments

6 Error Indicators and Warnings

7 Accuracy

8 Parallelism and Performance

▸▿ 9 Further Comments

9.1 Description of the Printed Output

9.2 Internal Changes

9.3 Additional Licensor

▸▿ 10 Example

10.1 Program Text

10.2 Program Results

▸▿ 11 Algorithmic Details

11.1 Stopping Criteria

11.2 Scaling the NLP

▸▿ 12 Optional Parameters

12.1 Description of the Optional Parameters

1 Purpose

e04stf is a solver from the NAG optimization modelling suite for constrained large-scale Nonlinear Programming (NLP) problems. It is an interior point method optimization solver based on the IPOPT software package.

2 Specification

Fortran Interface

Subroutine e04stf (

handle, objfun, objgrd, confun, congrd, hess, monit, nvar, x, nnzu, u, rinfo, stats, iuser, ruser, cpuser, ifail)

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

C Header Interface

#include <nag.h>

void

e04stf_ (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 e04stf or nagf_opt_handle_solve_ipopt.

3 Description

e04stf is typically used to solve the following nonlinear programming problem

\begin{array}{l} \underset{x \in ℝ^{n}}{minimize} & f (x) \\ subject to & l_{g} \leq g (x) \leq u_{g}, \\ \frac{1}{2} x^{T} Q_{i} x + r_{i}^{T} x + s_{i} \leq 0, 1 \leq i \leq m_{Q}, \\ l_{B} \leq B x \leq u_{B}, \\ l_{x} \leq x \leq u_{x}, \end{array}

where

$n$ is the number of the decision variables,
$m_{g}$ is the number of the nonlinear constraints and $g (x)$ , $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 a $m_{B} \times 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.

The objective

f (x)

can be specified in a number of ways: e04ref for a dense linear function, e04rff (and e04rsf or e04rtf) for a sparse linear or quadratic function, e04rtf is used to provide the quadratic function in factorized form, and e04rgf for a general nonlinear function. In the last case, objfun and objgrd will be used to compute values and gradients of the objective function. Variable box bounds

l_{x}, u_{x}

can be specified with e04rhf. The special case of linear constraints

l_{B}, B, u_{B}

is handled by e04rjf, quadratic constraints are handled by e04rsf and e04rtf while general nonlinear constraints

l_{g}, g (x), u_{g}

are specified by e04rkf (all can be specified). Again, in the last case, confun and congrd will be used to compute values and gradients of the nonlinear constraint functions.

Finally, if it is viable to calculate second derivatives, the sparsity structure of the second partial derivatives of a general nonlinear objective and/or of any general nonlinear constraints is specified by e04rlf and the values of these derivatives themselves will be computed by user-supplied hess. While there is an option (see Hessian Mode) that forces internal approximation of second derivatives, no such option exists for first derivatives which must be computed accurately. If e04rlf has been called and hess is used to calculate values for second derivatives, both the nonlinear objective and all the nonlinear constraints must be included; it is not possible to provide a subset of these. If the problem has only linear or quadratic objective and constraints, then hess is never called since the required Hessian information is already provided by the calls to e04ref, e04rff, e04rjf, e04rsf and e04rtf. If e04rlf is not called, then internal approximation of second derivatives will take place.

See Section 3.1 in the E04 Chapter Introduction for more details about the NAG optimization modelling suite.

3.1 Structure of the Lagrange Multipliers

For a problem consisting of

n

variable bounds,

m_{B}

linear constraints,

m_{Q}

quadratic constraints, and

m_{g}

nonlinear constraints, the number of Lagrange multipliers, and consequently the correct value for nnzu, will be

q = 2 * n + 2 * m_{B} + 2 * m_{Q} + 2 * m_{g}

. The order these will be found in the u array is

\begin{array}{l} z_{1_{L}}, z_{1_{U}}, z_{2_{L}}, z_{2_{U}} \dots z_{n_{L}}, z_{n_{U}}, \\ λ_{1_{L}}, λ_{1_{U}}, λ_{2_{L}}, λ_{2_{U}} \dots λ_{{m_{B}}_{L}}, λ_{{m_{B}}_{U}}, \\ λ_{{(m_{B} + 1)}_{L}}, λ_{{(m_{B} + 1)}_{U}}, λ_{{(m_{B} + 2)}_{L}}, λ_{{(m_{B} + 2)}_{U}} \dots λ_{{(m_{B} + m_{Q})}_{L}}, λ_{{(m_{B} + m_{Q})}_{U}}, \\ λ_{{(m_{B} + m_{Q} + 1)}_{L}}, λ_{{(m_{B} + m_{Q} + 1)}_{U}}, λ_{{(m_{B} + m_{Q} + 2)}_{L}}, λ_{{(m_{B} + m_{Q} + 2)}_{U}} \dots λ_{{(m_{B} + m_{Q} + m_{g})}_{L}}, λ_{{(m_{B} + m_{Q} + m_{g})}_{U}} . \end{array}

where the

L

and

U

subscripts refer to lower and upper bounds, respectively, and the variable bound constraint multipliers come first (if present, i.e., if e04rhf was called), followed by the linear constraint multipliers (if present, i.e., if e04rjf was called), followed by the quadratic constraint multipliers (if present, i.e., if e04rsf or e04rtf were called), and the nonlinear constraint multipliers (if present, i.e., if e04rkf was called).

Significantly nonzero values for any of these, after the solver has terminated, indicates that the corresponding constraint is active. Significance is judged in the first instance by the relative scale of any value compared to the smallest among them.

4 References

Byrd R H, Gilbert J Ch and Nocedal J (2000) A trust region method based on interior point techniques for nonlinear programming Mathematical Programming 89 149–185

Byrd R H, Liu G and Nocedal J (1997) On the local behavior of an interior point method for nonlinear programming Numerical Analysis (eds D F Griffiths and D J Higham) Addison–Wesley

Conn A R, Gould N I M, Orban D and Toint Ph L (2000) A primal-dual trust-region algorithm for non-convex nonlinear programming Mathematical Programming 87 (2) 215–249

Conn A R, Gould N I M and Toint Ph L (2000) Trust Region Methods SIAM, Philadephia

Fiacco A V and McCormick G P (1990) Nonlinear Programming: Sequential Unconstrained Minimization Techniques SIAM, Philadelphia

Gould N I M, Orban D, Sartenaer A and Toint Ph L (2001) Superlinear convergence of primal-dual interior point algorithms for nonlinear programming SIAM Journal on Optimization 11 (4) 974–1002

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

Hogg J D and Scott J A (2010) An indefinite sparse direct solver for large problems on multicore machines RAL Technical Report. RAL-TR-2010-011

Hogg J D and Scott J A (2011) HSL MA97: a bit-compatible multifrontal code for sparse symmetric systems RAL Technical Report. RAL-TR-2011-024

Wächter A and Biegler L T (2006) On the implementation of a primal-dual interior point filter line search algorithm for large-scale nonlinear programming Mathematical Programming 106(1) 25–57

Williams P and Lang B (2013) A framework for the

M R^{3}

Algorithm: theory and implementation SIAM J. Sci. Comput. 35 740–766

Yamashita H (1998) A globally convergent primal-dual interior-point method for constrained optimization Optimization Methods and Software 10 443–469

5 Arguments

1: $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 e04stf. It must not be changed between calls to the NAG optimization modelling suite.

2: $objfun$ – Subroutine, supplied by the NAG Library or the user. External Procedure

objfun must calculate the value of the nonlinear objective function

f (x)

at a specified value of the

n

-element vector of

x

variables. If there is no nonlinear objective (e.g., e04ref, e04rff, e04rsf or e04rtf was called to define a linear or quadratic objective function), objfun will never be called by e04stf and objfun may be the dummy routine e04stv. (e04stv is included in the NAG Library.)

The specification of objfun is:

Fortran Interface

Subroutine objfun (

nvar, x, fx, inform, iuser, ruser, cpuser)

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

C Header Interface

void	objfun (const Integer nvar, const double x[], double fx, Integer inform, Integer iuser[], double ruser[], void *cpuser)

1: $nvar$ – Integer Input: On entry: $n$ , the current number of decision variables $x$ in the model.
2: $x (nvar)$ – Real (Kind=nag_wp) array Input: On entry: the vector $x$ of variable values at which the objective function is to be evaluated.
3: $fx$ – Real (Kind=nag_wp) Output: On exit: the value of the objective function at $x$ .
4: $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 objfun. Specifically, if the value is negative, then the value of fx will be discarded and the solver will either attempt to find a different trial point or terminate immediately with $ifail = 25$ ; otherwise, the solver will proceed normally.
5: $iuser (*)$ – Integer array User Workspace
6: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace
7: $cpuser$ – Type (c_ptr) User Workspace: objfun is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: objfun should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

3: $objgrd$ – Subroutine, supplied by the NAG Library or the user. External Procedure

objgrd must calculate the values of the nonlinear objective function gradients

\frac{\partial f}{\partial x}

at a specified value of the

n

-element vector of

x

variables. If there is no nonlinear objective (e.g., e04ref, e04rff, e04rsf or e04rtf was called to define a linear or quadratic objective function), objgrd will never be called by e04stf and objgrd may be the dummy subroutine e04stw included in the NAG Library.

The specification of objgrd is:

Fortran Interface

Subroutine objgrd (

nvar, x, nnzfd, fdx, inform, iuser, ruser, cpuser)

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

C Header Interface

void	objgrd (const Integer nvar, const double x[], const Integer nnzfd, double fdx[], Integer inform, Integer iuser[], double ruser[], void *cpuser)

1: $nvar$ – Integer Input: On entry: $n$ , the current number of decision variables $x$ in the model.
2: $x (nvar)$ – 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: $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: $fdx (nnzfd)$ – 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. $fdx (i)$ will be the gradient $\frac{\partial f}{\partial x_{idxfd (i)}}$ .
5: $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 objgrd. Specifically, if the value is negative then the value of fdx will be discarded and the solver will either attempt to find a different trial point or will terminate immediately with $ifail = 25$ ; otherwise, computations will continue.
6: $iuser (*)$ – Integer array User Workspace
7: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace
8: $cpuser$ – Type (c_ptr) User Workspace: objgrd is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: objgrd should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

4: $confun$ – Subroutine, supplied by the NAG Library or the user. External Procedure

confun must calculate the values of the

m_{g}

-element vector

g_{i} (x)

of nonlinear constraint functions at a specified value of the

n

-element vector of

x

variables. If there are no nonlinear constraints then confun will never be called by e04stf and it may be the dummy subroutine e04stx included in the NAG Library.

The specification of confun is:

Fortran Interface

Subroutine confun (

nvar, x, ncnln, gx, inform, iuser, ruser, cpuser)

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

C Header Interface

void	confun (const Integer nvar, const double x[], const Integer ncnln, double gx[], Integer inform, Integer iuser[], double ruser[], void *cpuser)

1: $nvar$ – Integer Input: On entry: $n$ , the current number of decision variables $x$ in the model.
2: $x (nvar)$ – Real (Kind=nag_wp) array Input: On entry: the vector $x$ of variable values at which the constraint functions are to be evaluated.
3: $ncnln$ – Integer Input: On entry: $m_{g}$ , the number of nonlinear constraints, as specified in an earlier call to e04rkf.
4: $gx (ncnln)$ – Real (Kind=nag_wp) array Output: On exit: the $m_{g}$ values of the nonlinear constraint functions at $x$ .
5: $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 confun. Specifically, if the value is negative, then the value of gx will be discarded and the solver will either attempt to find a different trial point or terminate immediately with $ifail = 25$ ; otherwise, the solver will proceed normally.
6: $iuser (*)$ – Integer array User Workspace
7: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace
8: $cpuser$ – Type (c_ptr) User Workspace: confun is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: confun should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

5: $congrd$ – Subroutine, supplied by the NAG Library or the user. External Procedure

congrd must calculate the nonzero values of the sparse Jacobian of the nonlinear constraint functions

\frac{\partial g_{i}}{\partial x}

at a specified value of the

n

-element vector of

x

variables. If there are no nonlinear constraints, congrd will never be called by e04stf and congrd may be the dummy subroutine e04sty included in the NAG Library.)

The specification of congrd is:

Fortran Interface

Subroutine congrd (

nvar, x, nnzgd, gdx, inform, iuser, ruser, cpuser)

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

C Header Interface

void	congrd (const Integer nvar, const double x[], const Integer nnzgd, double gdx[], Integer inform, Integer iuser[], double ruser[], void *cpuser)

1: $nvar$ – Integer Input: On entry: $n$ , the current number of decision variables $x$ in the model.
2: $x (nvar)$ – 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: $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: $gdx (nnzgd)$ – 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. $gdx (i)$ will be the gradient $\frac{\partial g_{j}}{\partial x_{k}}$ , where $j = irowgd (i)$ and $k = icolgd (i)$ .
5: $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 congrd. Specifically, if the value is negative the solution of the current problem will terminate immediately with $ifail = 25$ ; otherwise, computations will continue.
6: $iuser (*)$ – Integer array User Workspace
7: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace
8: $cpuser$ – Type (c_ptr) User Workspace: congrd is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: congrd should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

6: $hess$ – Subroutine, supplied by the NAG Library or the user. External Procedure

hess must calculate the nonzero values of one of a set of second derivative quantities:

the Hessian of the Lagrangian function $σ \nabla^{2} f (x) + \sum_{i = 1}^{m_{g}} λ_{i} \nabla^{2} g_{i} (x)$ ,
the Hessian of the objective function $\nabla^{2} f (x)$ ,
the Hessian of the $i$ th constraint function $\nabla^{2} g_{i} (x)$ .

The value of argument idf determines which one of these is to be computed and this, in turn, is determined by earlier calls to e04rlf, when the nonzero sparsity structure of these Hessians was registered. Please note that it is not possible to only supply a subset of the Hessians (see

ifail = 6

). If there were no calls to e04rlf, hess will never be called by e04stf and hess may be the dummy subroutine e04stz (e04stz is included in the NAG Library). In this case, the Hessian of the Lagrangian will be approximated by a limited-memory quasi-Newton method (L-BFGS).

The specification of hess is:

Fortran Interface

Subroutine hess (

nvar, x, ncnln, idf, sigma, lambda, nnzh, hx, inform, iuser, ruser, cpuser)

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

C Header Interface

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: $nvar$ – Integer Input

On entry:

n

, the current number of decision variables

x

in the model.

2: $x (nvar)$ – Real (Kind=nag_wp) array Input

On entry: the vector

x

of variable values at which the Hessian functions are to be evaluated.

3: $ncnln$ – Integer Input

On entry:

m_{g}

, the number of nonlinear constraints, as specified in an earlier call to e04rkf.

4: $idf$ – Integer Input

On entry: specifies the quantities to be computed in hx.

$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.
$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.
$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: $sigma$ – Real (Kind=nag_wp) Input

On entry: if

idf = −1

, the value of the

σ

quantity in the definition of the Hessian of the Lagrangian. Otherwise, sigma should not be referenced.

6: $lambda (ncnln)$ – Real (Kind=nag_wp) array Input

On entry: if

idf = −1

, the values of the

λ_{i}

quantities in the definition of the Hessian of the Lagrangian. Otherwise, lambda should not be referenced.

7: $nnzh$ – Integer Input

On entry: the number of nonzero elements in the Hessian to be computed.

8: $hx (nnzh)$ – 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: $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

ifail = 25

; otherwise, computations will continue.

10: $iuser (*)$ – Integer array User Workspace

11: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace

12: $cpuser$ – Type (c_ptr) User Workspace

hess is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: hess should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

7: $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.

monit may be the dummy subroutine e04stu included in the NAG Library.

The specification of monit is:

Fortran Interface

Subroutine monit (

nvar, x, nnzu, u, inform, rinfo, stats, iuser, ruser, cpuser)

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

C Header Interface

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: $nvar$ – Integer Input: On entry: $n$ , the current number of decision variables $x$ in the model.
2: $x (nvar)$ – Real (Kind=nag_wp) array Input: On entry: $x^{i}$ , the values of the decision variables $x$ at the current iteration.
3: $nnzu$ – Integer Input: On entry: the dimension of array u.
4: $u (nnzu)$ – Real (Kind=nag_wp) array Input: On entry: if $nnzu > 0$ , u holds the values of Lagrange multipliers (dual variables) for the constraints at the current iteration. See Section 3.1 for layout information.
5: $inform$ – Integer Input/Output: On entry: a non-negative value.

On exit: may be used to request the solver to stop immediately. Specifically, if $inform < 0$ the solver will terminate immediately with $ifail = 20$ ; otherwise, the solver will proceed normally.
6: $rinfo (100)$ – Real (Kind=nag_wp) array Input: On entry: error measures and various indicators at the end of the current iteration as described in Section 9.1.
7: $stats (100)$ – Real (Kind=nag_wp) array Input: On entry: solver statistics at the end of the current iteration. It reports only the iteration count and the number of backtracking trial steps taken. See Section 9.1.
8: $iuser (*)$ – Integer array User Workspace
9: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace
10: $cpuser$ – Type (c_ptr) User Workspace: monit is called with the arguments iuser, ruser and cpuser as supplied to e04stf. 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 e04stf is called. Arguments denoted as Input must not be changed by this procedure.

Note: monit should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e04stf. If your code inadvertently does return any NaNs or infinities, e04stf is likely to produce unexpected results.

8: $nvar$ – Integer Input

On entry:

n

, the current number of decision variables

x

in the model.

9: $x (nvar)$ – 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: $nnzu$ – Integer Input

On entry: the number of Lagrange multipliers that are to be returned in array u.

nnzu = 0

, u will not be referenced; otherwise, it needs to match the dimension

q

as explained in Section 3.1.

Constraints:

$nnzu \geq 0$ ;
if $nnzu > 0$ , $nnzu = q$ .

11: $u (nnzu)$ – Real (Kind=nag_wp) array Input/Output

On entry: the input of u is reserved for future releases of the NAG Library and it is ignored at the moment.

Note: if

nnzu > 0

, u holds Lagrange multipliers (dual variables) for the constraints. See Section 3.1 for layout information. If

nnzu = 0

, u will not be referenced.

On exit: the final value of Lagrange multipliers

(z, λ)

12: $rinfo (100)$ – Real (Kind=nag_wp) array Output

On exit: error measures and various indicators at the end of the final iteration as given in the list below:

$1$	Objective function value $f (x)$ .
$2$	Constraint violation (primal infeasibility), see (7).
$3$	Dual infeasibility, see (6).
$4$	Complementarity.
$5$	Karush–Kuhn–Tucker error.
$6$	Regularization term for the Hessian of the Lagrangian. This value is only available in monit, see Iteration log in Section 9.1.
$7$	Step size for the dual variables. This value is only available in monit, see Iteration log in Section 9.1.
$8$	Step size for the primal variables. This value is only available in monit, see Iteration log in Section 9.1.
$9$ – $100$	Reserved for future use.

13: $stats (100)$ – Real (Kind=nag_wp) array Output

On exit: solver statistics at the end of the final iteration as given in the list below:

$1$	Number of the iterations.
$2$	Reserved for future use.
$3$	Number of backtracking trial steps.
$4$	Number of Hessian evaluations.
$5$	Number of objective gradient evaluations.
$6$ , $7$	Reserved for future use.
$8$	Total wall clock time elapsed.
$9$ – $18$	Reserved for future use.
$19$	Number of objective function evaluations.
$20$	Number of constraint function evaluations.
$21$	Number of constraint Jacobian evaluations.
$22$ – $100$	Reserved for future use.

14: $iuser (*)$ – Integer array User Workspace

15: $ruser (*)$ – Real (Kind=nag_wp) array User Workspace

16: $cpuser$ – Type (c_ptr) User Workspace

iuser, ruser and cpuser are not used by e04stf, 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: $ifail$ – Integer Input/Output

On entry: ifail must be set to

0

−1

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

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

ifail \neq 0

on exit. When the value $- 1$ or $1$ is used it is essential to test the value of ifail on exit.

On exit:

ifail = 0

unless the routine detects an error or a warning has been flagged (see Section 6).

6 Error Indicators and Warnings

If on entry

ifail = 0

−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 e04stf may return useful information.

$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 initialized by e04raf or it has been corrupted.

$ifail = 2$: The problem is already being solved.

This solver does not support the model defined in the handle.

$ifail = 4$: On entry, $nvar = ⟨ value ⟩$ , expected $value = ⟨ value ⟩$ .
Constraint: nvar must match the current number of variables of the model in the handle.

$ifail = 5$: On entry, $nnzu = ⟨ value ⟩$ .
Constraint: $nnzu = ⟨ value ⟩$ or $0$ .

On entry, $nnzu = ⟨ value ⟩$ .
Constraint: no constraints present, so $nnzu$ must be $0$ .

$ifail = 6$: Either all of the constraint and objective Hessian structures must be defined or none (in which case, the Hessians will be approximated by a limited-memory quasi-Newton L-BFGS method).
On entry, a nonlinear objective function has been defined but no objective Hessian sparsity structure has been defined through e04rlf.

On entry, a nonlinear constraint function has been defined but no constraint Hessian sparsity structure has been defined through e04rlf, for constraint number $⟨ value ⟩$ .

$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. Either change the optional parameter Hessian Mode or 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.

$ifail = 20$: User requested termination during a monitoring step. inform was set to a negative value in monit.

$ifail = 22$: Maximum number of iterations exceeded.

$ifail = 23$: The solver terminated after an error in the step computation. This message is printed if the solver is unable to compute a search direction, despite several attempts to modify the iteration matrix. Usually, the value of the regularization parameter then becomes too large. One situation where this can happen is when values in the Hessian are invalid (NaN or Infinity). You can check whether this is true by using the Verify Derivatives option.

The solver terminated after failure in the restoration phase. This indicates that the restoration phase failed to find a feasible point that was acceptable to the filter line search for the original problem. This could happen if the problem is highly degenerate, does not satisfy the constraint qualification, or if your NLP code provides incorrect derivative information.

The solver terminated after the maximum time allowed was exceeded. Maximum number of seconds exceeded. Use optional parameter Time Limit to reset the limit.

The solver terminated due to an invalid option. Please contact NAG with details of the call to e04stf.

The solver terminated due to an invalid problem definition. Please contact NAG with details of the call to e04stf.

The solver terminated with not enough degrees of freedom. This indicates that your problem, as specified, has too few degrees of freedom. This can happen if you have too many equality constraints, or if you fix too many variables.

$ifail = 24$: The solver terminated after the search direction became too small. This indicates that the solver is calculating very small step sizes and is making very little progress. This could happen if the problem has been solved to the best numerical accuracy possible given the current NLP scaling.

$ifail = 25$: Invalid number detected in user function. 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.

$ifail = 50$: The solver reports NLP solved to acceptable level. This indicates that the algorithm did not converge to the desired tolerances, but that it was able to obtain a point satisfying the acceptable tolerance level. This may happen if the desired tolerances are too small for the current problem.

$ifail = 51$: The solver detected an infeasible problem. The restoration phase converged to a point that is a minimizer for the constraint violation (in the $ℓ_{1}$ -norm), but is not feasible for the original problem. This indicates that the problem may be infeasible (or at least that the algorithm is stuck at a locally infeasible point). The returned point (the minimizer of the constraint violation) might help you to find which constraint is causing the problem. If you believe that the NLP is feasible, it might help to start the optimization from a different point.

$ifail = 54$: The solver terminated due to diverging iterates. The max-norm of the iterates has become larger than a preset value. This can happen if the problem is unbounded below and the iterates are diverging.

$ifail = - 199$: e04stf is not available in this implementation.

$ifail = - 99$: An unexpected error has been triggered by this routine. Please contact NAG.
See Section 7 in the Introduction to the NAG Library FL Interface for further information.

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

$ifail = - 999$: Dynamic memory allocation failed.
See Section 9 in the Introduction to the NAG Library FL Interface for further information.

7 Accuracy

The accuracy of the solution is driven by optional parameter Stop Tolerance 1.

ifail = 0

on the final exit, the returned point satisfies Karush–Kuhn–Tucker (KKT) conditions to the requested accuracy (under the default settings close to

\sqrt{ε}

where

ε

is the machine precision) and thus it is a good estimate of a local solution. If

ifail = 50

, some of the convergence conditions were not fully satisfied but the point still seems to be a reasonable estimate and should be usable. Please refer to Section 11.1 and the description of the particular options.

8 Parallelism and Performance

Background information to multithreading can be found in the Multithreading documentation.

e04stf is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.

e04stf makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.

Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the Users' Note for your implementation for any additional implementation-specific information.

9 Further Comments

9.1 Description 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 and Monitoring Level determine the exposed level of detail. This allows, for example, the generation of a detailed log in a file while the condensed information is displayed on the screen. This section also describes what kind of information is made available to the monitoring routine monit via rinfo and stats.

There are four sections printed to the primary output with the default settings (level

2

): a derivative check, a header, an iteration log and a summary. At higher levels more information will be printed, including any internal IPOPT options that have been changed from their default values.

Header

Print Level \geq 1

, the header will contain option settings and statistics about the size of the problem how the solver sees it, i.e., it reflects any changes imposed by preprocessing and problem transformations. The header may look similar to:

Banner and optional parameters list

 ------------------------------------------------------------------------------
  E04ST, Interior point method for large-scale nonlinear optimization problems
 ------------------------------------------------------------------------------
 
 Begin of Options
     Print File                    =                   6     * d
     Print Level                   =                   2     * U
     Monitoring File               =                  67     * U
     Monitoring Level              =                   2     * U
 
     Infinite Bound Size           =         1.00000E+20     * d
     Task                          =            Minimize     * d
     Stats Time                    =                  No     * d
     Time Limit                    =         1.00000E+01     * U
     Verify Derivatives            =                  No     * d
 
     Hessian Mode                  =                Auto     * d
     Nlp Factorization Method      =                Ma97     * d
     Matrix Ordering               =                Auto     * d
     Outer Iteration Limit         =                  26     * U
     Stop Tolerance 1              =         2.50000E-08     * U
 End of Options

Summary of the problem

Number of nonzeros in equality constraint Jacobian...:        4
Number of nonzeros in inequality constraint Jacobian.:        8
Number of nonzeros in Lagrangian Hessian.............:       10

Total number of variables............................:        4
                     variables with only lower bounds:        4
                variables with lower and upper bounds:        0
                     variables with only upper bounds:        0
Total number of equality constraints.................:        1
Total number of inequality constraints...............:        2
        inequality constraints with only lower bounds:        2
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:        0

Derivative Check

If Verify Derivatives is set, then information will appear about any errors detected in the user-supplied derivative routines objgrd, congrd or hess. It may look like this:

Starting derivative checker for first derivatives.

* grad_f[          1] = -2.000000E+00    ~  2.455000E+01  [ 1.081E+00]
* jac_g [    1,    4] =  4.700969E+01 v  ~  5.200968E+01  [ 9.614E-02]
Starting derivative checker for second derivatives.

*             obj_hess[    1,    1] =  1.881000E+03 v  ~  1.882000E+03  [ 5.314E-04]
*     1-th constr_hess[    1,    3] =  2.988964E+00 v  ~ -1.103543E-02  [ 3.000E+00]

Derivative checker detected 3 error(s).

The first line indicates that the value for the partial derivative of the objective with respect to the first variable as returned by objgrd (the first one printed) differs sufficiently from a finite difference estimation derived from objfun (the second one printed). The number in square brackets is the relative difference between these two numbers.

The second line reports on a discrepancy for the partial derivative of the first constraint with respect to the fourth variable. If the indicator v is absent, the discrepancy refers to a component that had not been included in the sparsity structure, in which case the nonzero structure of the derivatives should be corrected. Mistakes in the first derivatives should be corrected before attempting to correct mistakes in the second derivatives.

The third line reports on a discrepancy in a second derivative of the objective function, differentiated with respect to the first variable, twice.

The fourth line reports on a discrepancy in a second derivative of the first constraint, differentiated with respect to the first and third variables.

Iteration log

Print Level = 2

, the status of each iteration is condensed to one line. The line shows:

iter

The current iteration count. This includes regular iterations and iterations during the restoration phase. If the algorithm is in the restoration phase, the letter r will be appended to the iteration number. The iteration number

0

represents the starting point. This quantity is also available as

stats (1)

of monit.

objective

The unscaled objective value at the current point (given the current NLP scaling). During the restoration phase, this value remains the unscaled objective value for the original problem. This quantity is also available as

rinfo (1)

of monit.

inf_pr

The unscaled constraint violation at the current point (given the current NLP scaling). This quantity is the infinity-norm (max) of the (unscaled) constraints

g_{i}

. During the restoration phase, this value remains the constraint violation of the original problem at the current point. This quantity is also available as

rinfo (2)

of monit.

inf_du

The scaled dual infeasibility at the current point (given the current NLP scaling). This quantity measure the infinity-norm (max) of the internal dual infeasibility,

λ_{i}

of Eq. (4a) in the implementation paper Wächter and Biegler (2006), including inequality constraints reformulated using slack variables and NLP scaling. During the restoration phase, this is the value of the dual infeasibility for the restoration phase problem. This quantity is also available as

rinfo (3)

of monit.

lg(mu)

log10

of the value of the barrier parameter

μ

μ

itself is also available as

rinfo (4)

of monit.

||d||

The infinity norm (max) of the primal step (for the original variables x and the internal slack variables s). During the restoration phase, this value includes the values of additional variables,

\bar{p}

and

\bar{n}

(see Eq. (30) in Wächter and Biegler (2006)). This quantity is also available as

rinfo (5)

of monit.

lg(rg)

log10

of the value of the regularization term for the Hessian of the Lagrangian in the augmented system (

δ_{w}

of Eq. (26) and Section 3.1 in Wächter and Biegler (2006)). A dash (–) indicates that no regularization was done. The regularization term itself is also available as

rinfo (6)

of monit.

alpha_du

The step size for the dual variables (

α_{k}^{z}

of Eq. (14c) in Wächter and Biegler (2006)). This quantity is also available as

rinfo (7)

of monit.

alpha_pr

The step size for the primal variables (

α_{k}

of Eq. (14a) in Wächter and Biegler (2006)). This quantity is also available as

rinfo (8)

of monit. The number is usually followed by a character for additional diagnostic information regarding the step acceptance criterion.

f	f-type iteration in the filter method without second-order correction
F	f-type iteration in the filter method with second-order correction
h	h-type iteration in the filter method without second-order correction
H	h-type iteration in the filter method with second-order correction
k	penalty value unchanged in merit function method without second-order correction
K	penalty value unchanged in merit function method with second-order correction
n	penalty value updated in merit function method without second-order correction
N	penalty value updated in merit function method with second-order correction
R	Restoration phase just started
w	in watchdog procedure
s	step accepted in soft restoration phase
t/T	tiny step accepted without line search
r	some previous iterate restored

The number of backtracking line search steps (does not include second-order correction steps). This quantity is also available as

stats (3)

of monit.

Note that the step acceptance mechanisms in IPOPT consider the barrier objective function (4) which is usually different from the value reported in the objective column. Similarly, for the purposes of the step acceptance, the constraint violation is measured for the internal problem formulation, which includes slack variables for inequality constraints and potentially NLP scaling of the constraint functions. This value, too, is usually different from the value reported in inf_pr. As a consequence, a new iterate might have worse values both for the objective function and the constraint violation as reported in the iteration output, seemingly contradicting globalization procedure.

Note that all these values are also available in

rinfo (1), \dots, rinfo (8)

stats (1)

, and

stats (3)

of the monitoring routine monit.

The output might look as follows:

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0  2.6603500E+05 1.55E+02 3.21E+01  -1.0 0.00E+00    -  0.00E+00 0.00E+00   0
   1  1.5053889E+05 7.95E+01 1.43E+01  -1.0 1.16E+00    -  4.55E-01 1.00E+00f  1
   2  8.9745785E+04 3.91E+01 6.45E+00  -1.0 3.07E+01    -  5.78E-03 1.00E+00f  1
   3  3.9878595E+04 1.63E+01 3.47E+00  -1.0 5.19E+00   0.0 2.43E-01 1.00E+00f  1
   4  2.7780042E+04 1.08E+01 1.64E+00  -1.0 3.66E+01    -  7.24E-01 8.39E-01f  1
   5  2.6194274E+04 1.01E+01 1.49E+00  -1.0 1.07E+01    -  1.00E+00 1.05E-01f  1
   6  1.5422960E+04 4.75E+00 6.82E-01  -1.0 1.74E+01    -  1.00E+00 1.00E+00f  1
   7  1.1975453E+04 3.14E+00 7.26E-01  -1.0 2.83E+01    -  1.00E+00 5.06E-01f  1
   8  8.3508421E+03 1.34E+00 2.04E-01  -1.0 3.96E+01    -  9.27E-01 1.00E+00f  1
   9  7.0657495E+03 4.85E-01 9.22E-02  -1.0 5.32E+01    -  1.00E+00 1.00E+00f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  10  6.8359393E+03 1.17E-01 1.28E-01  -1.7 4.69E+01    -  8.21E-01 1.00E+00h  1
  11  6.6508917E+03 1.52E-02 1.52E-02  -2.5 1.87E+01    -  1.00E+00 1.00E+00h  1
  12  6.4123213E+03 8.77E-03 1.49E-01  -3.8 1.85E+01    -  7.49E-01 1.00E+00f  1
  13  6.3157361E+03 4.33E-03 1.90E-03  -3.8 2.07E+01    -  1.00E+00 1.00E+00f  1
  14  6.2989280E+03 1.12E-03 4.06E-04  -3.8 1.54E+01    -  1.00E+00 1.00E+00h  1
  15  6.2996264E+03 9.90E-05 2.05E-04  -5.7 5.35E+00    -  9.63E-01 1.00E+00h  1
  16  6.2998436E+03 0.00E+00 1.86E-07  -5.7 4.55E-01    -  1.00E+00 1.00E+00h  1
  17  6.2998424E+03 0.00E+00 6.18E-12  -8.2 2.62E-03    -  1.00E+00 1.00E+00h  1

Print Level > 2

, each iteration produces significantly more detailed output comprising detailed error measures and output from internal operations. The output is reasonably self-explanatory so it is not featured here in detail.

Summary

Once the solver finishes, a detailed summary is produced if

Print Level \geq 1

. An example is shown below:

Number of Iterations....: 6

                                   (scaled)                 (unscaled)
Objective...............:   7.8692659500479623E-01    6.2324586324379867E+00
Dual infeasibility......:   7.9744615766675617E-10    6.3157735687207093E-09
Constraint violation....:   8.3555384833289281E-12    8.3555384833289281E-12
Complementarity.........:   0.0000000000000000E+00    0.0000000000000000E+00
Overall NLP error.......:   7.9744615766675617E-10    6.3157735687207093E-09


Number of objective function evaluations             = 7
Number of objective gradient evaluations             = 7
Number of equality constraint evaluations            = 7
Number of inequality constraint evaluations          = 0
Number of equality constraint Jacobian evaluations   = 7
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations             = 6
Total CPU secs in IPOPT (w/o function evaluations)   =      0.724
Total CPU secs in NLP function evaluations           =      0.343

EXIT: Optimal Solution Found.

It starts with the total number of iterations the algorithm went through. Then, five quantities are printed, all evaluated at the termination point: the value of the objective function, the dual infeasibility, the constraint violation, the complementarity and the NLP error.

This is followed by some statistics on the number of calls to user-supplied functions and CPU time taken in user-supplied functions and the main algorithm. Lastly, status at exit is indicated by a short message. Detailed timings of the algorithm are displayed only if Stats Time is set.

9.2 Internal Changes

Internal changes have been made to this routine as follows:

At Mark 26.1:
The default for the optional parameter Verify Derivatives was changed from $AUTO$ to $NO$ meaning that the derivatives will not be checked unless you explicitly request them to be and the description of the option $AUTO$ was removed.

$Print Level = 0$ and $Monitoring Level = 0$ no longer produce output, a banner was printed in the previous release.

A new option Task was introduced. It allows you to easily switch between minimization, maximization and feasible point. The previous release assumed minimization which is now the default choice.

A new option Matrix Ordering was introduced. It allows you to choose the fill-reducing ordering for the internal sparse linear algebra solver. Originally, at Mark 26, only AMD ordering was implemented. METIS ordering has now been introduced which is especially efficient for large-scale problems. A heuristic to automatically choose between the two orderings was also added and is now the default choice.
At Mark 27:
The name of the argument ‘mon’ was updated to monit to be consistent with the rest of the NAG optimization modelling suite routines.
At Mark 27.1:
The arguments fdx of objgrd and gdx of congrd are now Intent (Inout) instead of Intent (Out) to stay consistent with the other solvers of the NAG optimization modelling suite.

At Mark 28.3:

Solver e04stf interface was updated. The size of the information arrays rinfo and stats have been increased from

32

100

. The intent of the Lagrange multipliers, array u, changed from Intent(Out) to Intent(Inout); this array is only assigned to and not referenced. Also, the intent of the argument hx in user call-back hess was changed from Intent (Out) to Intent (Inout). In what follows, the changes to the interfaces are shown.

Old interface

Subroutine e04stf (handle, objfun, objgrd, confun, congrd, hess, monit, nvar, &
                   x, nnzu, u, rinfo, stats, iuser, ruser, cpuser, ifail)
Integer, Intent (In)               :: nvar, nnzu
Integer, Intent (Inout)            :: iuser(*), ifail
Real (Kind=nag_wp), Intent (Inout) :: x(nvar), ruser(*)
Real (Kind=nag_wp), Intent (Out)   :: u(nnzu), rinfo(32), stats(32)
Type (c_ptr), Intent (In)          :: handle, cpuser
External                           :: objfun, objgrd, confun, congrd, hess,  &
                                      monit

with Hessian user call-back

Subroutine hess (nvar, x, ncnln, idf, sigma, lambda, nnzh, hx, inform, iuser, &
                 ruser, cpuser)
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 (Out)   :: hx(nnzh), ruser(*)
Type (c_ptr), Intent (In)          :: cpuser

and monitor user call-back

Subroutine monit ( nvar, x, nnzu, u, inform, rinfo, stats, iuser, ruser, cpuser)
Integer, Intent (In)               :: nvar, nnzu
Integer, Intent (Inout)            :: inform, iuser(*)
Real (Kind=nag_wp), Intent (In)    :: x(nvar), u(nnzu), rinfo(32), stats(32)
Real (Kind=nag_wp), Intent (Inout) :: ruser(*)
Type (c_ptr), Intent (In)          :: cpuser

New Interface

Subroutine e04stf (handle, objfun, objgrd, confun, congrd, hess, monit, nvar, &
                   x, nnzu, u, rinfo, stats, iuser, ruser, cpuser, ifail)
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

with Hessian user call-back

Subroutine hess (nvar, x, ncnln, idf, sigma, lambda, nnzh, hx, inform, iuser, &
                 ruser, cpuser)
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

and monitor user call-back

Subroutine monit (nvar, x, nnzu, u, inform, rinfo, stats, iuser, ruser, cpuser)
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

At Mark 28.4:
Calls to vendor linear algebra routines have been rationalized in order to behave more predictably and consistently in processes which depend on external vendor libraries.
At Mark 30.2:
Another underlying sparse linear algebra solver, Harwell package MA86, was introduced. The option NLP Factorization Method can be used to choose between Harwell package MA86 or MA97.

For details of all known issues which have been reported for the NAG Library please refer to the Known Issues.

9.3 Additional Licensor

Parts of the code for e04stf are distributed according to terms imposed by the Eclipse Public License. Please refer to Library Licensors for further details.

10 Example

This example is based on Problem 73 in Hock and Schittkowski (1981) and involves the minimization of the linear function

f (x) = 24.55 x_{1} + 26.75 x_{2} + 39.00 x_{3} + 40.50 x_{4},

subject to the bounds

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

to the nonlinear constraint

12 x_{1} + 11.9 x_{2} + 41.8 x_{3} + 52.1 x_{4} - 21 - 1.645 \sqrt{0.28 x_{1}^{2} + 0.19 x_{2}^{2} + 20.5 x_{3}^{2} + 0.62 x_{4}^{2}} \geq 0,

and the linear constraints

\begin{matrix} 2.3 x_{1} + 5.6 x_{2} + 11.1 x_{3} + 1.3 x_{4} \geq 5, \\ x_{1} + x_{2} + x_{3} + x_{4} - 1 = 0 . \end{matrix}

The initial point, which is infeasible, is

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

and

f (x_{0}) = 130.8

. The optimal solution (to five significant figures) is

x^{*} = {(0.63552, 0.0, 0.31270, 0.051777)}^{T} .

10.1 Program Text

Program Text (e04stfe.f90)

10.2 Program Results

Program Results (e04stfe.r)

11 Algorithmic Details

e04stf is an implementation of IPOPT (see Wächter and Biegler (2006)) that is fully supported and maintained by NAG. It uses Harwell packages MA97 or MA86 for the underlying sparse linear algebra factorization. MA86 uses MC68 approximate minimum degree for the ordering and MA97 uses either MC68 approximate minimum degree or METIS algorithm. Any issues relating to e04stf should be directed to NAG who assume all responsibility for the e04stf routine and its implementation.

In the remainder of this section, we repeat part of Section 2.1 of Wächter and Biegler (2006).

To simplify notation, we describe the method for the problem formulation

\begin{array}{l} \underset{x \in ℝ^{n}}{minimize} & f (x) \end{array}

(1)

\begin{array}{l} subject to & g (x) = 0 \end{array}

(2)

\begin{array}{l} x \geq 0 . \end{array}

(3)

Range constraints of the form

l \leq c (x) \leq u

can be expressed in this formulation by introducing slack variables

x_{s} \geq 0

x_{t} \geq 0

(increasing

n

2

) and defining new equality constraints

g (x, x_{s}) \equiv c (x) - l - x_{s} = 0

and

g (x, x_{t}) \equiv u - c (x) - x_{t} = 0

e04stf, like the methods discussed in Williams and Lang (2013), Byrd et al. (2000), Conn et al. (2000) and Fiacco and McCormick (1990), computes (approximate) solutions for a sequence of barrier problems

\begin{array}{l} \underset{x \in ℝ^{n}}{minimize} & φ_{μ} (x) ≔ f (x) - μ \sum_{i = 1}^{n} \ln (x^{(i)}) \end{array}

(4)

\begin{array}{l} subject to & g (x) = 0 \end{array}

(5)

for a decreasing sequence of barrier parameters

μ

converging to zero.

The algorithm may be interpreted as a homotopy method to the primal-dual equations,

\begin{aligned} \nabla f (x) + \nabla g (x) λ - z = 0 \end{aligned}

(6)

\begin{aligned} g (x) = 0 \end{aligned}

(7)

\begin{aligned} X Z e - μ e = 0 \end{aligned}

(8)

with the homotopy parameter

μ

, which is driven to zero (see e.g., Byrd et al. (1997) and Gould et al. (2001)). Here,

X ≔ diag (x)

for a vector

x

, similarly

Z ≔ diag (z)

, and

e

stands for the vector of all ones for appropriate dimension, while

λ \in ℝ^{m}

and

z \in ℝ^{n}

correspond to the Lagrange multipliers for the equality constraints (2) and the bound constraints (3), respectively.

Note, that the equations (6), (7) and (8) for

μ = 0

together with ‘

x

z \geq 0

’ are the Karush–Kuhn–Tucker (KKT) conditions for the original problem (1), (2) and (3). Those are the first-order optimality conditions for (1), (2) and (3) if constraint qualifications are satisfied (Conn et al. (2000)).

Starting from an initial point supplied in x, e04stf computes an approximate solution to the barrier problem (4) and (5) for a fixed value of

μ

(by default,

0.1

), then decreases the barrier parameter, and continues the solution of the next barrier problem from the approximate solution of the previous one.

A sophisticated overall termination criterion for the algorithm is used to overcome potential difficulties when the Lagrange multipliers become large. This can happen, for example, when the gradients of the active constraints are nearly linear dependent. The termination criterion is described in detail by Wächter and Biegler (2006) (also see below Section 11.1).

11.1 Stopping Criteria

Using the individual parts of the primal-dual equations (6), (7) and (8), we define the optimality error for the barrier problem as

E_{μ} (x, λ, z) ≔ \max {\frac{{‖ \nabla f (x) + \nabla g (x) λ - z ‖}_{\infty}}{s_{d}}, {‖ g (x) ‖}_{\infty}, \frac{{‖ X Z e - μ e ‖}_{\infty}}{s_{c}}}

(9)

with scaling parameters

s_{d}

s_{c} \geq 1

defined below (not to be confused with NLP scaling factors described in Section 11.2). By

E_{0} (x, λ, z)

we denote (9) with

μ = 0

; this measures the optimality error for the original problem (1), (2) and (3). The overall algorithm terminates if an approximate solution

({\tilde{x}}_{*}, {\tilde{λ}}_{*}, {\tilde{z}}_{*})

(including multiplier estimates) satisfying

E_{0} ({\tilde{x}}_{*}, {\tilde{λ}}_{*}, {\tilde{z}}_{*}) \leq ε_{tol}

(10)

is found, where

ε_{tol} > 0

is the user-supplied error tolerance in optional parameter Stop Tolerance 1.

Even if the original problem is well scaled, the multipliers

λ

and

z

might become very large, for example, when the gradients of the active constraints are (nearly) linearly dependent at a solution of (1), (2) and (3). In this case, the algorithm might encounter numerical difficulties satisfying the unscaled primal-dual equations (6), (7) and (8) to a tight tolerance. In order to adapt the termination criteria to handle such circumstances, we choose the scaling factors

s_{d} ≔ \max {s_{\max}, \frac{{‖ λ ‖}_{1} + {‖ z ‖}_{1}}{(m + n)}} / s_{\max} s_{c} ≔ \max {s_{\max}, \frac{{‖ z ‖}_{1}}{n}} / s_{\max}

in (9). In this way, a component of the optimality error is scaled, whenever the average value of the multipliers becomes larger than a fixed number

s_{\max} \geq 1

(

s_{\max} = 100

in our implementation). Also note, in the case that the multipliers diverge,

E_{0} (x, λ, z)

can only become small, if a Fritz John point for (1), (2) and (3) is approached, or if the primal variables diverge as well.

11.2 Scaling the NLP

Ideally, the formulated problem should be scaled so that, near the solution, all function gradients (objective and constraints), when nonzero, are of a similar order of a magnitude. e04stf will compute automatic NLP scaling factors for the objective and constraint functions (but not the decision variables) and apply them if large imbalances of scale are detected. This rescaling is only computed at the starting point. References to scaled or unscaled objective or constraints in Section 9.1 and Section 11 should be understood in this context.

12 Optional Parameters

Several optional parameters in e04stf define choices in the problem specification or the algorithm logic. In order to reduce the number of formal arguments of e04stf 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 and the call to the solver. Modification of the optional parameters during intermediate monitoring stops is not allowed. Once the solver finishes, the optional parameters can be altered again for the next solve.

If any options are set by the solver (typically those with the choice of

AUTO

), their value can be retrieved by e04znf. If the solver is called again, any such arguments are reset to their default values and the decision is made again.

The following is a list of the optional parameters available. A full description of each optional parameter is provided in Section 12.1.

Defaults
Hessian Mode
Infinite Bound Size
Matrix Ordering
Monitoring File
Monitoring Level
NLP Factorization Method
Outer Iteration Limit
Print File
Print Level
Print Options
Print Solution
Stats Time
Stop Tolerance 1
Task
Time Limit
Verify Derivatives

12.1 Description 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, where the minimum abbreviation of each keyword is underlined;
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 $ε$ is a generic notation for machine precision (see x02ajf).

All options accept the value

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.

Hessian Mode

a

Default

= AUTO

This parameter specifies whether the Hessian will be user-supplied (in hx) or approximated by e04stf using a limited-memory quasi-Newton L-BFGS method. In the

AUTO

setting, if no Hessian structure has been registered in the problem with a call to e04rlf and there are general nonlinear objective or constraints, then the Hessian will be approximated. Otherwise hess will be called if and only if any of e04rgf or e04rkf have been used to define the problem. Approximating the Hessian is likely to require more iterations to achieve convergence but will reduce the time spent in user-supplied functions.

Constraint:

Hessian Mode = AUTO

EXACT

APPROXIMATE

Infinite Bound Size

r

Default

= 10^{20}

This defines the ‘infinite’ bound

bigbnd

in the definition of the problem constraints. Any upper bound greater than or equal to

bigbnd

will be regarded as

+ \infty

(and similarly any lower bound less than or equal to

- 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.

It also serves as a limit for the objective function to be considered unbounded (

ifail = 54

Constraint:

Infinite Bound Size \geq 1000

Monitoring File

i

Default

= −1

i \geq 0

, the unit number for the secondary (monitoring) output. If set to

−1

, no secondary output is provided. The information output to this unit is controlled by Monitoring Level.

Constraint:

Monitoring File \geq −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 \leq Monitoring Level \leq 5

Matrix Ordering

a

Default

= AUTO

If NLP Factorization Method is

MA97

, this parameter specifies the ordering to be used by the internal sparse linear algebra solver. It affects the number of nonzeros in the factorized matrix and thus influences the cost per iteration.

$Matrix Ordering = AUTO$: A heuristic is used to choose automatically between METIS and AMD orderings.
$Matrix Ordering = BEST$: Both AMD and METIS orderings are computed at the beginning of the solve and the one with the fewest nonzeros in the factorized matrix is selected.
$Matrix Ordering = AMD$: An approximate minimum degree (AMD) ordering is used.
$Matrix Ordering = METIS$: METIS ordering is used.

Note: if NLP Factorization Method is

MA86

, this parameter is ignored and AMD ordering will be used.

Constraint:

Matrix Ordering = AUTO

BEST

AMD

METIS

NLP Factorization Method

a

Default

= MA97

This parameter controls whether Harwell package

MA97

MA86

is used for the sparse linear algebra factorization. To determine which best suits your application, it is recommended to try both

MA97

and

MA86

Constraint:

NLP Factorization Method = MA97

MA86

Outer Iteration Limit

i

Default

= 100

The maximum number of iterations to be performed by e04stf. Setting the option too low might lead to

ifail = 22

Constraint:

Outer Iteration Limit \geq 0

Print File

i

Default

= advisory message unit number

i \geq 0

, the unit number for the primary output of the solver. If

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 information output to this unit is controlled by Print Level.

Constraint:

Print File \geq −1

Print Level

i

Default

= 2

This parameter defines how detailed information should be printed by the solver to the primary output.

$i$	Output
$0$	No output from the solver
$1$	Additionally, derivative check information, the Header and Summary.
$2$	Additionally, the Iteration log.
$3$ , $4$	Additionally, details of each iteration with scalar quantities printed.
$5$	Additionally, individual components of arrays are printed resulting in large output.

Constraint:

0 \leq Print Level \leq 5

Print Options

a

Default

= YES

Print Options = YES

, a listing of optional parameters will be printed to the primary output.

Constraint:

Print Options = YES

NO

Print Solution

a

Default

= NO

Print Solution = X

, the final values of the primal variables are printed on the primary and secondary outputs.

Print Solution = YES

ALL

, in addition to the primal variables, the final values of the dual variables are printed on the primary and secondary outputs.

Constraint:

Print Solution = YES

NO

X

ALL

Stats Time

a

Default

= NO

This parameter allows you to turn on 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.

Constraint:

Stats Time = YES

NO

Stop Tolerance 1

r

Default

= \max (10^{−6}, \sqrt{ε})

This option sets the value

ε_{tol}

of (10) which is used for optimality and complementarity tests from KKT conditions. See Section 11.1.

Constraint:

Stop Tolerance 1 > ε

Task

a

Default

= MINIMIZE

This parameter specifies the required direction of the optimization. If

Task = FEASIBLE POINT

, the objective function (if set) is ignored and the algorithm stops as soon as a feasible point is found with respect to the given tolerance. If no objective function is set, Task reverts to

FEASIBLE POINT

automatically.

Constraint:

Task = MINIMIZE

MAXIMIZE

FEASIBLE POINT

Time Limit

r

Default

= 10^{6}

A limit to the number of seconds that the solver can use to solve one problem. If during the convergence check this limit is exceeded, the solver will terminate with a corresponding error message.

Constraint:

Time Limit > 0

Verify Derivatives

a

Default

= NO

This parameter specifies whether the routine should perform numerical checks on the consistency of the user-supplied functions. It is recommended that such checks are enabled when first developing the formulation of the problem, however, the derivative check results in a significant increase of the number of the function evaluations and thus it shouldn't be used in production code.

Constraint: