PDF version (NAG web site
, 64-bit version, 64-bit version)
NAG Toolbox: nag_mip_sqp (h02da)
Purpose
nag_mip_sqp (h02da) solves general nonlinear programming problems with integer constraints on some of the variables.
Syntax
[
ax,
x,
c,
cjac,
objgrd,
objmip,
user,
ifail] = h02da(
ncnln,
a,
d,
bl,
bu,
varcon,
x,
confun,
objfun,
iopts,
opts, 'n',
n, 'nclin',
nclin, 'maxit',
maxit, 'acc',
acc, 'user',
user)
[
ax,
x,
c,
cjac,
objgrd,
objmip,
user,
ifail] = nag_mip_sqp(
ncnln,
a,
d,
bl,
bu,
varcon,
x,
confun,
objfun,
iopts,
opts, 'n',
n, 'nclin',
nclin, 'maxit',
maxit, 'acc',
acc, 'user',
user)
Before calling
nag_mip_sqp (h02da),
nag_mip_optset (h02zk) must be called with
optstr set to ‘’. Optional parameters may also be specified by calling
nag_mip_optset (h02zk) before the call to
nag_mip_sqp (h02da).
Description
nag_mip_sqp (h02da) solves mixed integer nonlinear programming problems using a modified sequential quadratic programming method. The problem is assumed to be stated in the following general form:
with
continuous variables and
binary and integer variables in a total of
variables;
equality constraints in a total of
constraint functions.
Partial derivatives of and are not required for the integer variables. Gradients with respect to integer variables are approximated by difference formulae.
No assumptions are made regarding except that it is twice continuously differentiable with respect to continuous elements of . It is not assumed that integer variables are relaxable. In other words, problem functions are evaluated only at integer points.
The method seeks to minimize the exact penalty function:
where
is adapted by the algorithm and
is given by:
Successive quadratic approximations are applied under the assumption that integer variables have a smooth influence on the model functions, that is function values do not change drastically when incrementing or decrementing an integer value. In practice this requires integer variables to be ordinal not categorical. The algorithm is stabilised by a trust region method including Yuan's second order corrections, see
Yuan and Sun (2006). The Hessian of the Lagrangian function is approximated by BFGS (see
The Quasi-Newton Update in
nag_opt_nlp1_solve (e04uc)) updates subject to the continuous and integer variables.
The mixed-integer quadratic programming subproblems of the SQP-trust region method are solved by a branch and cut method with continuous subproblem solutions obtained by the primal-dual method of Goldfarb and Idnani, see
Powell (1983). Different strategies are available for selecting a branching variable:
- Maximal fractional branching. Select an integer variable from the relaxed subproblem solution with largest distance from next integer value
- Minimal fractional branching. Select an integer variable from the relaxed subproblem solution with smallest distance from next integer value
and a node from where branching, that is the generation of two new subproblems, begins:
- Best of two. The optimal objective function values of the two child nodes are compared and the node with a lower value is chosen
- Best of all. Select an integer variable from the relaxed subproblem solution with the smallest distance from the next integer value
- Depth first. Select a child node whenever possible.
This implementation is based on the algorithm MISQP as described in
Exler et al. (2013).
Linear constraints may optionally be supplied by a matrix
and vector
rather than the constraint functions
such that
Partial derivatives with respect to of these constraint functions are not requested by nag_mip_sqp (h02da).
References
Exler O, Lehmann T and Schittkowski K (2013) A comparative study of SQP-type algorithms for nonlinear and nonconvex mixed-integer optimization Mathematical Programming Computation 4 383–412
Mann A (1986) GAMS/MINOS: Three examples Department of Operations Research Technical Report Stanford University
Powell M J D (1983) On the quadratic programming algorithm of Goldfarb and Idnani Report DAMTP 1983/Na 19 University of Cambridge, Cambridge
Yuan Y-x and Sun W (2006) Optimization Theory and Methods Springer–Verlag
Parameters
Compulsory Input Parameters
- 1:
– int64int32nag_int scalar
-
, the number of constraints supplied by .
Constraint:
.
- 2:
– double array
-
The first dimension of the array
a must be at least
.
The second dimension of the array
a must be at least
if
.
The
th row of
a must contain the coefficients of the
th general linear constraint, for
. Any equality constraints must be specified first.
If
, the array
a is not referenced.
- 3:
– double array
-
, the constant for the
th linear constraint.
If
, the array
d is not referenced.
- 4:
– double array
- 5:
– double array
-
must contain the lower bounds, , and the upper bounds, , for the variables; bounds on integer variables are rounded, bounds on binary variables need not be supplied.
Constraint:
, for .
- 6:
– int64int32nag_int array
-
varcon indicates the nature of each variable and constraint in the problem. The first
elements of the array must describe the nature of the variables, the next
elements the nature of the general linear constraints (if any) and the next
elements the general constraints (if any).
- A continuous variable.
- A binary variable.
- An integer variable.
- An equality constraint.
- An inequality constraint.
Constraints:
- , or , for ;
- or , for ;
- At least one variable must be either binary or integer;
- Any equality constraints must precede any inequality constraints.
- 7:
– double array
-
An initial estimate of the solution, which need not be feasible. Values corresponding to integer variables are rounded; if an initial value less than half is supplied for a binary variable the value zero is used, otherwise the value one is used.
- 8:
– function handle or string containing name of m-file
-
confun must calculate the constraint functions supplied by
and their Jacobian at
. If all constraints are supplied by
and
(i.e.,
),
confun will never be called by
nag_mip_sqp (h02da) and
confun may be the string
nag_mip_sqp_dummy_confun (h02ddm). (
nag_mip_sqp_dummy_confun (h02ddm) is included in the NAG Toolbox.) If
, the first call to
confun will occur after the first call to
objfun.
[mode, c, cjac, user] = confun(mode, ncnln, n, varcon, x, cjac, nstate, user)
Input Parameters
- 1:
– int64int32nag_int scalar
-
Indicates which values must be assigned during each call of
objfun. Only the following values need be assigned:
- Elements of c containing continuous variables.
- Elements of cjac containing continuous variables.
- 2:
– int64int32nag_int scalar
-
The dimension of the array
c and the first dimension of the array
cjac. the number of constraints supplied by
,
.
- 3:
– int64int32nag_int scalar
-
The second dimension of the array
cjac.
, the total number of variables,
.
- 4:
– int64int32nag_int array
-
The array
varcon as supplied to
nag_mip_sqp (h02da).
- 5:
– double array
-
The vector of variables at which the objective function and/or all continuous elements of its gradient are to be evaluated.
- 6:
– double array
-
Note: the derivative of the th constraint with respect to the th variable, , is stored in .
Continuous elements of
cjac are set to the value of NaN.
- 7:
– int64int32nag_int scalar
-
If
,
nag_mip_sqp (h02da) is calling
confun for the first time. This argument setting allows you to save computation time if certain data must be read or calculated only once.
- 8:
– Any MATLAB object
confun is called from
nag_mip_sqp (h02da) with the object supplied to
nag_mip_sqp (h02da).
Output Parameters
- 1:
– int64int32nag_int scalar
-
May be set to a negative value if you wish to terminate the solution to the current problem, and in this case
nag_mip_sqp (h02da) will terminate with
ifail set to
mode.
- 2:
– double array
-
Must contain
ncnln constraint values, with the value of the
th constraint
in
.
- 3:
– double array
-
Note: the derivative of the th constraint with respect to the th variable, , is stored in .
The
th row of
cjac must contain elements of
for each continuous variable
.
- 4:
– Any MATLAB object
- 9:
– function handle or string containing name of m-file
-
objfun must calculate the objective function
and its gradient for a specified
-element vector
.
[mode, objmip, objgrd, user] = objfun(mode, n, varcon, x, objgrd, nstate, user)
Input Parameters
- 1:
– int64int32nag_int scalar
-
Indicates which values must be assigned during each call of
objfun. Only the following values need be assigned:
- The objective function value, objmip.
- The continuous elements of objgrd.
- 2:
– int64int32nag_int scalar
-
, the total number of variables, .
- 3:
– int64int32nag_int array
-
The array
varcon as supplied to
nag_mip_sqp (h02da).
- 4:
– double array
-
The vector of variables at which the objective function and/or all continuous elements of its gradient are to be evaluated.
- 5:
– double array
-
Continuous elements of
objgrd are set to the value of NaN.
- 6:
– int64int32nag_int scalar
-
If
,
nag_mip_sqp (h02da) is calling
objfun for the first time. This argument setting allows you to save computation time if certain data must be read or calculated only once.
- 7:
– Any MATLAB object
objfun is called from
nag_mip_sqp (h02da) with the object supplied to
nag_mip_sqp (h02da).
Output Parameters
- 1:
– int64int32nag_int scalar
-
May be set to a negative value if you wish to terminate the solution to the current problem, and in this case
nag_mip_sqp (h02da) will terminate with
ifail set to
mode.
- 2:
– double scalar
-
Must be set to the objective function value,
, if
; otherwise
objmip is not referenced.
- 3:
– double array
-
Must contain the gradient vector of the objective function if
, with
containing the partial derivative of
with respect to continuous variable
; otherwise
objgrd is not referenced.
- 4:
– Any MATLAB object
- 10:
– int64int32nag_int array
-
Note: the dimension of this array is dictated by the requirements of associated functions that must have been previously called. This array
must be the same array passed as argument
iopts in the previous call to
nag_mip_optset (h02zk).
- 11:
– double array
-
Note: the dimension of this array is dictated by the requirements of associated functions that must have been previously called. This array
must be the same array passed as argument
opts in the previous call to
nag_mip_optset (h02zk).
The real option array as returned by
nag_mip_optset (h02zk).
Optional Input Parameters
- 1:
– int64int32nag_int scalar
-
Default:
the dimension of the arrays
bl,
bu,
x and the second dimension of the array
a. (An error is raised if these dimensions are not equal.)
, the total number of variables, .
Constraint:
.
- 2:
– int64int32nag_int scalar
-
Default:
the dimension of the array
d and the first dimension of the array
a. (An error is raised if these dimensions are not equal.)
, the number of general linear constraints defined by and .
Constraint:
.
- 3:
– int64int32nag_int scalar
Default:
The maximum number of iterations within which to find a solution. If
maxit is less than or equal to zero, the suggested value below is used.
- 4:
– double scalar
Default:
The requested accuracy for determining feasible points during iterations and for halting the method when the predicted improvement in objective function is less than
acc. If
acc is less than or equal to
(
being the
machine precision as given by
nag_machine_precision (x02aj)), the below suggested value is used.
- 5:
– Any MATLAB object
user is not used by
nag_mip_sqp (h02da), but is passed to
confun and
objfun. Note that for large objects it may be more efficient to use a global variable which is accessible from the m-files than to use
user.
Output Parameters
- 1:
– double array
-
The final values of the linear constraints
.
If
,
ax is not referenced.
- 2:
– double array
-
The final estimate of the solution.
- 3:
– double array
-
If
,
contains the value of the
th constraint function
at the final iterate, for
.
If
, the array
c is not referenced.
- 4:
– double array
-
Note: the derivative of the th constraint with respect to the th variable, , is stored in .
If
,
cjac contains the Jacobian matrix of the constraint functions at the final iterate, i.e.,
contains the partial derivative of the
th constraint function with respect to the
th variable, for
and
. (See the discussion of argument
cjac under
confun.)
If
, the array
cjac is not referenced.
- 5:
– double array
-
The objective function gradient at the solution.
- 6:
– double scalar
-
With
,
objmip contains the value of the objective function for the MINLP solution.
- 7:
– Any MATLAB object
- 8:
– int64int32nag_int scalar
unless the function detects an error (see
Error Indicators and Warnings).
Error Indicators and Warnings
Errors or warnings detected by the function:
-
-
Constraint: .
-
-
Constraint: .
-
-
Constraint: .
-
-
Constraint: .
-
-
Constraint: , for .
-
-
Constraint: , or , for .
-
-
Constraint: or , for .
-
-
The supplied
objfun returned a NaN value.
-
-
The supplied
confun returned a NaN value.
-
-
On entry, the optional parameter arrays
iopts and
opts have either not been initialized or been corrupted.
-
-
On entry, there are no binary or integer variables.
-
-
On entry, linear equality constraints do not precede linear inequality constraints.
-
-
On entry, nonlinear equality constraints do not precede nonlinear inequality constraints.
-
-
One or more objective gradients appear to be incorrect.
-
-
One or more constraint gradients appear to be incorrect.
-
-
On entry. Exceeded the maximum number of iterations.
-
-
More than the maximum number of feasible steps without improvement in the objective function.
-
-
Penalty parameter tends to infinity in an underlying mixed-integer quadratic program; the problem may be infeasible.
-
-
Termination at an infeasible iterate; if the problem is feasible, try a different starting value.
-
-
Termination with zero integer trust region for integer variables; try a different starting value.
-
-
The optimization failed due to numerical difficulties. Set optional parameter for more information.
-
-
The optimization halted because you set
mode negative in
objfun or
mode negative in
confun, to
.
-
An unexpected error has been triggered by this routine. Please
contact
NAG.
-
Your licence key may have expired or may not have been installed correctly.
-
Dynamic memory allocation failed.
Accuracy
Not applicable.
Further Comments
None.
Example
Select a portfolio of at most
assets from
available with expected return
, is fully invested and that minimizes
where
- is a vector of proportions of selected assets
- is an indicator variable that describes if an asset is in or out
- is a vector of mean returns
- is the covariance matrix of returns.
This example is taken from
Mann (1986) with
Linear constraints are supplied through both
and
, and
confun.
Open in the MATLAB editor:
h02da_example
function h02da_example
diary h02da_example.r
fprintf('h02da example results\n\n');
n = int64(8);
nclin = int64(5);
ncnln = int64(2);
varcon = zeros(n+nclin+ncnln, 1, 'int64');
varcon(1:4) = int64([0;0;0;0]);
varcon(5:8) = int64([1;1;1;1]);
bl = zeros(n,1);
bu = zeros(n,1);
bl(1:4) = 0;
bu(1:4) = 1000;
varcon(n+1) = 3;
varcon(n+2:n+nclin) = 4;
a = zeros(nclin,n);
a(1,1:4) = 1;
a(2,[1,5]) = [-1,1];
a(3,[2,6]) = [-1,1];
a(4,[3,7]) = [-1,1];
a(5,[4,8]) = [-1,1];
d = zeros(nclin,1);
d(1) = 1;
varcon(n+nclin+1) = 3;
varcon(n+nclin+2) = 4;
iopts = zeros(200, 1, 'int64');
opts = zeros(100, 1);
[iopts, opts, ifail] = h02zk( ...
'Initialize = h02da', iopts, opts);
maxit = int64(500);
acc = 1e-6;
x = zeros(n,1);
x(1:4) = 1;
user.p = 3;
user.rho = 10;
[ax, x, c, cjac, objgrd, objmip, user, ifail] = ...
h02da( ...
ncnln, a, d, bl, bu, varcon, x, @confun, @objfun, iopts, opts, ...
'user', user);
[ivalue, accqp, cvalue, optype, ifail] = ...
h02zl('QP Accuracy', iopts, opts);
[ifail] = x04ca('G', ' ', x, 'Final estimate');
fprintf('\nMixed integer QP Accuracy: %g', accqp);
fprintf('\nOptimised value: %g\n\n', objmip);
diary off
function [mode, c, cjac, user] = ...
confun(mode, ncnln, n, varcon, x, cjac, nstate, user)
c = zeros(ncnln, 1);
if (mode == 0)
c(1) = 8*x(1) + 9*x(2) + 12*x(3) + 7*x(4) - user.rho;
c(2) = user.p - x(5) - x(6) - x(7) - x(8);
else
cjac(1,1:4) = [8,9,12,7];
cjac(2,1:4) = 0;
end
function [mode, objmip, objgrd, user] = ...
objfun(mode, n, varcon, x, objgrd, nstate, user)
objmip = 0;
if (mode == 0)
objmip = x(1)*(4*x(1)+3*x(2)-x(3)) + ...
x(2)*(3*x(1)+6*x(2)+x(3)) + ...
x(3)*(x(2)-x(1)+10*x(3));
else
objgrd(1) = 8*x(1) + 6*x(2) - 2*x(3);
objgrd(2) = 6*x(1) + 12*x(2) + 2*x(3);
objgrd(3) = 2*(x(2)-x(1)) + 20*x(3);
objgrd(4) = 0;
end
h02da example results
Final estimate
1
1 0.3750
2 0.0000
3 0.5250
4 0.1000
5 1.0000
6 0.0000
7 1.0000
8 1.0000
Mixed integer QP Accuracy: 1e-10
Optimised value: 2.925
Optional Parameters
This section can be skipped if you wish to use the default values for all optional parameters, otherwise, the following is a list of the optional parameters available and a full description of each optional parameter is provided in
Description of the s.
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;
- a parameter value,
where the letters , and denote options that take character, integer and real values respectively.
All options accept the value in order to return single options to their default states.
Keywords and character values are case insensitive, however they must be separated by at least one space.
nag_mip_optset (h02zk) can be called to supply options, one call being necessary for each optional parameter. For example,
Call H02ZKF('Check Gradients = Yes', iopts, liopts, opts, lopts, ifail)
nag_mip_optset (h02zk) should be consulted for a full description of the method of supplying optional parameters.
For
nag_mip_sqp (h02da) the maximum length of the argument
cvalue used by
nag_mip_optget (h02zl) is
.
Branch Bound Steps Default
Maximum number of branch-and-bound steps for solving the mixed integer quadratic problems.
Constraint:
.
Branching Rule Default
Branching rule for branch and bound search.
- Maximum fractional branching.
- Minimum fractional branching.
Check Gradients Default
Perform an internal check of supplied objective and constraint gradients. It is advisable to set during code development to avoid difficulties associated with incorrect user-supplied gradients.
Continuous Trust Radius Default
Initial continuous trust region radius, ; the initial trial step for the SQP approximation must satisfy .
Constraint:
.
Descent Default
Initial descent parameter, , larger values of allow penalty parameter to increase faster which can lead to faster convergence.
Constraint:
.
Descent Factor Default
Factor for decreasing the internal descent parameter, , between iterations.
Constraint:
.
Feasible Steps Default
Maximum number of feasible steps without improvements, where feasibility is measured by .
Constraint:
.
Improved Bounds Default
Calculate improved bounds in case of ‘Best of all’ node selection strategy.
Integer Trust Radius Default
Initial integer trust region radius, ; the initial trial step for the SQP approximation must satisfy .
Constraint:
.
Maximum Restarts Default
Maximum number of restarts that allow the mixed integer SQP algorithm to return to a better solution. Setting a value higher than the default might lead to better results at the expense of function evaluations.
Constraint:
.
Minor Print Level Default
Print level of the subproblem solver. Active only if .
Constraint:
.
Modify Hessian Default
Modify the Hessian approximation in an attempt to get more accurate search directions. Calculation time is increased when the number of integer variables is large.
Node Selection Default
Node selection strategy for branch and bound.
- Large tree search; this method is the slowest as it solves all subproblem QPs independently.
- Uses warm starts and less memory.
- Uses more warm starts. If warm starts are applied, they can speed up the solution of mixed integer subproblems significantly when solving almost identical QPs.
Non Monotone Default
Maximum number of successive iterations considered for the non-monotone trust region algorithm, allowing the penalty function to increase.
Constraint:
.
Objective Scale Bound Default
When
internally scale absolute function values greater than
or
Objective Scale Bound.
Constraint:
.
Penalty Default
Initial penalty parameter, .
Constraint:
.
Penalty Factor Default
Factor for increasing penalty parameter when the trust regions cannot be enlarged at a trial step.
Constraint:
.
Print Level Default
Specifies the desired output level of printing.
- No output.
- Final convergence analysis.
- One line of intermediate results per iteration.
- Detailed information printed per iteration.
QP Accuracy Default
Termination tolerance of the relaxed quadratic program subproblems.
Constraint:
.
Scale Continuous Variables Default
Internally scale continuous variables values.
Scale Objective Function Default
Internally scale objective function values.
- No scaling.
- Scale absolute values greater than Objective Scale Bound.
Warm Starts Default
Maximum number of warm starts within the mixed integer QP solver, see
Node Selection.
Constraint:
.
PDF version (NAG web site
, 64-bit version, 64-bit version)
© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2015