constrained_nlp
Solves a general nonlinear programming problem using a sequential equality constrained quadratic programming method.
Synopsis
#include <imsl.h>
float *imsl_f_constrained_nlp (void fcn(), int m, int meq, int n, int ibtype, float xlb[], float xub[], …, 0)
The typedoublefunction is imsl_d_constrained_nlp.
Required Arguments
void fcn(int n, float x[], int iact, float *result, int *ierr)(Input)
User supplied function to evaluate the objective function and constraints at a given point.
int n (Input)
Number of variables.
float x[] (Input)
The point at which the objective function or a constraint is evaluated.
int iact (Input)
Integer indicating whether evaluation of the function is requested or evaluation of a constraint is requested. If iact is zero, then an objective function evaluation is requested. If iact is nonzero then the value of iact indicates the index of the constraint to evaluate. iact =1 to meq for equality constraints and iact =meq+1 to m for inequality constraints.
float result[] (Output)
If iact is zero, then result is the computed objective function at the point x. If iact is nonzero, then result is the requested constraint value at the point x.
int *ierr (Output)
Address of an integer. On input ierr is set to 0. If an error or other undesirable condition occurs during evaluation, then ierr should be set to 1. Setting ierr to 1 will result in the step size being reduced and the step being tried again. (If ierr is set to 1 for xguess, then an error is issued.)
int m (Input)
Total number of constraints.
int meq (Input)
Number of equality constraints.
int n (Input)
Number of variables.
int ibtype (Input)
Scalar indicating the types of bounds on variables.
ibtype | Action |
---|
0 | User will supply all the bounds. |
1 | All variables are nonnegative. |
2 | All variables are nonpositive. |
3 | User supplies only the bounds on first variable, all other variables will have the same bounds. |
float xlb[] (Input, Output, or Input/Output)
Array with n components containing the lower bounds on the variables. (Input, if ibtype = 0; output, if ibtype = 1 or 2; Input/Output, if ibtype = 3)
If there is no lower bound on a variable, then the corresponding xlb value should be set to imsl_f_machine(8).
float xub[] (Input, Output, or Input/Output)
Array with n components containing the upper bounds on the variables. (Input, if ibtype = 0; output, if ibtype 1 or 2; Input/Output, if ibtype = 3)
If there is no upper bound on a variable, then the corresponding xub value should be set to imsl_f_machine(7).
Return Value
A pointer to the solution x of the nonlinear programming problem. To release this space, use free. If no solution can be computed, then NULL is returned.
Synopsis with Optional Arugments
#include <imsl.h>
float *imsl_f_constrained_nlp (void fcn(), int m, int meq, int n, int ibtype, float xlb[], float xub[],
IMSL_GRADIENT, void grad(),
IMSL_PRINT, int iprint,
IMSL_XGUESS, float xguess[],
IMSL_ITMAX, int itmax,
IMSL_TAU0, float tau0,
IMSL_DEL0, float del0,
IMSL_SMALLW, float smallw,
IMSL_DELMIN, float delmin,
IMSL_SCFMAX, float scfmax,
IMSL_RETURN_USER, float x[],
IMSL_OBJ, float *obj,
IMSL_DIFFTYPE, int difftype,
IMSL_XSCALE, float xscale[],
IMSL_EPSDIF, float epsdif,
IMSL_EPSFCN, float epsfcn,
IMSL_TAUBND, float taubnd,
IMSL_FCN_W_DATA, void fcn(), void *data,
IMSL_GRADIENT_W_DATA, void grad(), void *data,
0)
Optional Arguments
IMSL_GRADIENT, void grad(int n, float x[], int iact, float result[]) (Input)
User-supplied function to evaluate the gradients at a given point where
Arguments
int n (Input)
Number of variables.
float x[] (Input)
The point at which the gradient of the objective function or gradient of a constraint is evaluated
int iact (Input)
Integer indicating whether evaluation of the function gradient is requested or evaluation of a constraint gradient is requested. If iact is zero, then an objective function gradient evaluation is requested. If iact is nonzero then the value of iact indicates the index of the constraint gradient to evaluate. iact =1 to meq for equality constraints and iact =meq+1 to m for inequality constraints.
float result[] (Output)
If iact is zero, then result is the computed gradient of the objective function at the point x. If iact is nonzero, then result is the computed gradient of the requested constraint value at the point x.
IMSL_PRINT, int iprint(Input)
Parameter indicating the desired output level. (Input)
iprint | Action |
---|
0 | No output printed. |
1 | One line of intermediate results is printed in each iteration. |
2 | Lines of intermediate results summarizing the most important data for each step are printed. |
3 | Lines of detailed intermediate results showing all primal and dual variables, the relevant values from the working set, progress in the backtracking and etc are printed |
4 | Lines of detailed intermediate results showing all primal and dual variables, the relevant values from the working set, progress in the backtracking, the gradients in the working set, the quasi-Newton updated and etc are printed. |
Default: iprint = 0.
IMSL_XGUESS, float xguess[](Input)
Array of length n containing an initial guess of the solution.
Default: xguess = X, (with the smallest value of ∥x∥2) that satisfies the bounds.
IMSL_ITMAX, int itmax(Input)
Maximum number of iterations allowed.
Default: itmax = 200.
IMSL_TAU0, float tau0(Input)
A universal bound describing how much the unscaled penalty-term may deviate from zero. imsl_f_constrained_nlp assumes that within the region described by
all functions may be evaluated safely. The initial guess, however, may violate these requirements. In that case an initial feasibility improvement phase is run by imsl_f_constrained_nlp until such a point is found. A small tau0diminishes the efficiency of imsl_f_constrained_nlp, because the iterates then will follow the boundary of the feasible set closely. Conversely, a large tau0 may degrade the reliability of the code.
Default tau0 = 1.0.
IMSL_DEL0, float del0(Input)
In the initial phase of minimization a constraint is considered binding if
Good values are between .01 and 1.0. If del0is chosen too small then identification of the correct set of binding constraints may be delayed. Contrary, if del0is too large, then the method will often escape to the full regularized SQP method, using individual slack variables for any active constraint, which is quite costly. For well-scaled problems del0 =1.0 is reasonable.
Default: del0 = .5* tau0
IMSL_SMALLW, float smallw(Input)
Scalar containing the error allowed in the multipliers. For example, a negative multiplier of an inequality constraint is accepted (as zero) if its absolute value is less than smallw.
Default: smallw = exp(2*log(eps/3)) where eps is the machine precision.
IMSL_DELMIN, float delmin (Input)
Scalar which defines allowable constraint violations of the final accepted result. Constraints are satisfied if ∣gi(x)∣≤ delmin for equality constraints, and gi(x) ≥(‑delmin) for equality constraints.
Default: delmin = min(.1*del0, max(epsdif, max(1.e-6*del0, smallw))
IMSL_SCFMAX, float scfmax(Input)
Scalar containing the bound for the internal automatic scaling of the objective function. (Input)
Default: scfmax = 1.0e4
IMSL_RETURN_USER, float x[] (Output)
A user allocated array of length n containing the solution x.
IMSL_OBJ, float *obj(Output)
Scalar containing the value of the objective function at the computed solution.
IMSL_LAGRANGE_MULTIPLIERS, float **lagrange(Output)
The address of a pointer, which on exit, points to an array containing the Lagrange multiplier estimates of the constraints.
IMSL_LAGRANGE_MULTIPLIERS_USER, float lagrange_user[](Output)
A user-supplied array of length ncon containing the Lagrange multiplier estimates of the constraints.
IMSL_CONSTRAINT_RESIDUALS, float **const_res(Output)
The address of a pointer, which on exit, points to an array containing the constraints residuals.
IMSL_CONSTRAINT_RESIDUALS_USER, float const_res_user[](Output)
A user-supplied array of length ncon containing the constraint residuals.
IMSL_FCN_W_DATA,
void fcn(
int n,
float x[],
int iact,
float *result,
int *ierr,
void *data),
void *data, (Input)
User supplied function to evaluate the objective function and constraints at a given point, which also accepts a pointer to data that is supplied by the user.
data is a pointer to the data to be passed to the user-supplied function. See
Passing Data to User-Supplied Functions in the introduction to this manual for more details.
IMSL_GRADIENT_W_DATA,
void grad(
int n,
float x[],
int iact,
float result[],
void *data),
void *data, (Input)
User-supplied function to evaluate the gradients at a given point, which also accepts a pointer to data that is supplied by the user.
data is a pointer to the data to be passed to the user-supplied function. See
Passing Data to User-Supplied Functions in the introduction to this manual for more details.
Note: The following optional arguments are valid only if IMSL_GRADIENT is not supplied. |
IMSL_DIFFTYPE, int difftype(Input)
Type of numerical differentiation to be used.
Default: difftype = 1
difftype | Action |
---|
1 | Use a forward difference quotient with discretization stepsize 0.1(epsfcn)1/2 componentwise relative. |
2 | Use the symmetric difference quotient with discretization stepsize 0.1(epsfcn)1/3 componentwise relative. |
3 | Use the sixth order approximation computing a Richardson extrapolation of three symmetric difference quotient values. This uses a discretization stepsize 0.01(epsfcn)1/7. |
IMSL_XSCALE, float xscale[](Input)
Vector of length n setting the internal scaling of the variables. The initial value given and the objective function and gradient evaluations however are always in the original unscaled variables. The first internal variable is obtained by dividing values x[i] by xscale[i]. In the absence of other information, set all entries to 1.0.
Default: xscale[] = 1.0.
IMSL_EPSDIF, float epsdif(Input)
Relative precision in gradients.
Default: epsdif = ɛ where ɛ is the machine precision.
IMSL_EPSFCN, float epsfcn(Input)
Relative precision of the function evaluation routine.
Default: epsfcn = ɛ where ɛ is the machine precision
IMSL_TAUBND, float taubnd(Input)
Amount by which bounds may be violated during numerical differentiation. Bounds are violated by taubnd(at most) only if a variable is on a bound and finite differences are taken taken for gradient evaluations.
Default: taubnd = 1.0
Description
The function constrained_nlp provides an interface to a licensed version of subroutine DONLP2, a code developed by Peter Spellucci (1998). It uses a sequential equality constrained quadratic programming method with an active set technique, and an alternative usage of a fully regularized mixed constrained subproblem in case of nonregular constraints (i.e. linear dependent gradients in the “working sets”). It uses a slightly modified version of the Pantoja-Mayne update for the Hessian of the Lagrangian, variable dual scaling and an improved Armjijo-type stepsize algorithm. Bounds on the variables are treated in a gradient-projection like fashion. Details may be found in the following two papers:
P. Spellucci: An SQP method for general nonlinear programs using only equality constrained subproblems. Math. Prog. 82, (1998), 413-448.
P. Spellucci: A new technique for inconsistent problems in the SQP method. Math. Meth. of Oper. Res. 47, (1998), 355-500. (published by Physica Verlag, Heidelberg, Germany).
The problem is stated as follows:
Although default values are provided for optional input arguments, it may be necessary to adjust these values for some problems. Through the use of optional arguments, imsl_f_constrained_nlp allows for several parameters of the algorithm to be adjusted to account for specific characteristics of problems. The DONLP2 Users Guide provides detailed descriptions of these parameters as well as strategies for maximizing the perfomance of the algorithm. The DONLP2 Users Guide is available in the “help” subdirectory of the main IMSL product installation directory. In addition, the following are a number of guidelines to consider when using imsl_f_constrained_nlp.
A good initial starting point is very problem specific and should be provided by the calling program whenever possible. See optional argument
IMSL_XGUESS.
Gradient approximation methods can have an effect on the success of
imsl_f_constrained_nlp. Selecting a higher order approximation method may be necessary for some problems. See optional argument
IMSL_DIFFTYPE.
If a two sided constraint
li ≤ gi(
x)
≤ ui is transformed into two constraints
g 2i(
x)
≥ 0 and
g 2i+1(
x)
≥ 0 , then choose del0
< ½(
ui-li) /max{1,
∥∇gi(
x)
∥}, or at least try to provide an estimate for that value. This will increase the efficiency of the algorithm. See optional argument
IMSL_DEL0.
The parameter
ierr provided in the interface to the user supplied function
fcn can be very useful in cases when evaluation is requested at a point that is not possible or reasonable. For example, if evaluation at the requested point would result in a floating point exception, then setting
ierr to
1 and returning without performing the evaluation will avoid the exception.
imsl_f_constrained_nlp will then reduce the stepsize and try the step again. Note, if
ierr is set to
1 for the initial guess, then an error is issued.
On some platforms,
constrained_nlp can evaluate the user-supplied functions
fcn and
grad in parallel. This is done only if the function
imsl_omp_options is called to flag user-defined functions as thread-safe. A function is thread-safe if there are no dependencies between calls. Such dependencies are usually the result of writing to global or static variables.
Example
The problem
is solved.
#include "imsl.h"
#define M 2
#define ME 1
#define N 2
void grad(int n, float x[], int iact, float result[]);
void fcn(int n, float x[], int iact, float *result, int *ierr);
int main()
{
int ibtype = 0;
float *x, ans[2];
static float xlb[N], xub[N];
imsl_omp_options(IMSL_SET_FUNCTIONS_THREAD_SAFE, 1, 0);
xlb[0] = xlb[1] = imsl_f_machine(8);
xub[0] = xub[1] = imsl_f_machine(7);
x = imsl_f_constrained_nlp(fcn, M, ME, N, ibtype, xlb, xub, 0);
imsl_f_write_matrix ("The solution is", 1, N, x, 0);
}
/* Himmelblau problem 1 */
void fcn(int n, float x[], int iact, float *result, int *ierr)
{
float tmp1, tmp2;
tmp1 = x[0] - 2.0e0;
tmp2 = x[1] - 1.0e0;
switch (iact) {
case 0:
*result = tmp1 * tmp1 + tmp2 * tmp2;
break;
case 1:
*result = x[0] - 2.0e0 * x[1] + 1.0e0;
break;
case 2:
*result = -(x[0]*x[0]) / 4.0e0 - x[1]*x[1] + 1.0e0;
break;
default: ;
break;
}
*ierr = 0;
return;
}
Output
The solution is
1 2
0.8229 0.9114
Fatal Errors
IMSL_STOP_USER_FCN | Request from user supplied function to stop algorithm. User flag = "#". |
IMSL_BAD_CONSTR_EVAL | Constraint evaluation returns an error with current point. |
IMSL_BAD_OBJ_EVAL | Objective evaluation returns an error with current point. |
IMSL_WORKING_SET_SINGULAR | Working set is singular in dual extended QP. |
IMSL_QP_INFEASIBLE | QP problem is seemingly infeasible. The solutionprocess is severely corrupted by roundoff, mostprobably a problem of bad scaling which was notovercome by the internal scaling techniques. |
IMSL_STATIONARY_PT_ERR_1 | A stationary point of the penalty function has been located which is not feasible. |
IMSL_STATIONARY_PT_ERR_2 | A stationary point of the penalty function has been located which is not feasible. Or, limiting accuracy reached for a singular problem, with termination criteria being too strong. |
IMSL_STATIONARY_PT_ERR_3 | A stationary point of the penalty function which is not feasible for the original problem has been located. Try some other initial guess. |
IMSL_ITMAX_EXCEEDED | Maximum number of iterations limit "itmax" = # exceeded. The best answer found is returned. |
IMSL_STEPSIZE_SELECTION | No acceptable stepsize in [sigsm, sigla]. This is often due to a programming error in the user supplied (analytic) gradients. It may also be due to termination criteria that are too stringent for the problem at hand, because of evaluation impreciseness of functions and/or gradients or because of ill conditioning of the (projected) Hessian matrix. |
IMSL_SLOW_PRIMAL_PROGRESS | Very slow primal progress. The problem is singular or ill-conditioned. |
IMSL_SLOW_PROGRESS_IN_X | Very slow progress in X, the problem is singular. |
IMSL_LIN_DEP_GRAD | The gradients in the working set are linearly dependent, such that a full regularized QP is solved. It may occur in a problem that the second order sufficiency condition is not satisfied and the matrix of gradients of binding constraints is singular or very ill-conditioned. |
IMSL_SMALL_CHANGE | For max(n, 10) consecutive steps, there were only small changes in the penalty function without the other termination criteria satisfied. |