naginterfaces.library.opt.qp_dense_solve¶
- naginterfaces.library.opt.qp_dense_solve(bl, bu, h, x, comm, a=None, cvec=None, qphess=None, istate=None, data=None, io_manager=None)[source]¶
qp_dense_solve
solves general quadratic programming problems. It is not intended for large sparse problems.Note: this function uses optional algorithmic parameters, see also:
qp_dense_option_file()
,qp_dense_option_string()
,nlp1_init()
.For full information please refer to the NAG Library document for e04nf
https://www.nag.com/numeric/nl/nagdoc_29.2/flhtml/e04/e04nff.html
- Parameters
- blfloat, array-like, shape
must contain the lower bounds for all the constraints
- bufloat, array-like, shape
must contain the upper bounds for all the constraints
- hfloat, array-like, shape
May be used to store the quadratic term of the QP objective function if desired. In some cases, you need not use to store explicitly (see the specification of function ). The elements of are referenced only by function . The number of rows of is denoted by , whose default value is . (The option ‘Hessian Rows’ may be used to specify a value of .)
If the default version of is used and the problem is of type QP1 or QP2 (the default), the first rows and columns of must contain the leading rows and columns of the symmetric Hessian matrix .
Only the diagonal and upper triangular elements of the leading rows and columns of are referenced.
The remaining elements need not be assigned.
If the default version of is used and the problem is of type QP3 or QP4, the first rows of must contain an upper trapezoidal factor of the symmetric Hessian matrix .
The factor need not be of full rank, i.e., some of the diagonal elements may be zero.
However, as a general rule, the larger the dimension of the leading nonsingular sub-matrix of , the fewer iterations will be required.
Elements outside the upper trapezoidal part of the first rows of need not be assigned.
If a non-default version of is supplied, then in some cases it may be desirable to use a one-dimensional array to transmit data to . is then declared as an by array, where .
In other situations, it may be desirable to compute or without accessing – for example, if or is sparse or has special structure.
The arguments and may then refer to any convenient array.
If the problem is of type FP or LP, is not referenced.
- xfloat, array-like, shape
An initial estimate of the solution.
- commdict, communication object, modified in place
Communication structure.
This argument must have been initialized by a prior call to
nlp1_init()
.- aNone or float, array-like, shape , optional
Note: the required extent for this argument in dimension 2 is determined as follows: if : ; if : ; otherwise: .
The th row of must contain the coefficients of the th general linear constraint, for .
- cvecNone or float, array-like, shape , optional
The coefficients of the explicit linear term of the objective function when the problem is of type LP, QP2 (the default) and QP4.
If the problem is of type FP, QP1, or QP3, is not referenced.
- qphessNone or callable hx = qphess(jthcol, h, x, data=None), optional
Note: if this argument is None then a NAG-supplied facility will be used.
In general, you need not provide a version of .
However, the algorithm of
qp_dense_solve
requires only the product of or and a vector ; and in some cases you may obtain increased efficiency by providing a version of that avoids the need to define the elements of the matrices or explicitly.is not referenced if the problem is of type FP or LP, in which case may be None.
- Parameters
- jthcolint
Specifies whether or not the vector is a column of the identity matrix.
The vector is the th column of the identity matrix, and hence or is the th column of or , respectively. This may in some cases require very little computation and may be coded to take advantage of this. However special code is not necessary because is always stored explicitly in the array .
has no special form.
- hfloat, ndarray, shape
See:
- xfloat, ndarray, shape
The vector .
- dataarbitrary, optional, modifiable in place
User-communication data for callback functions.
- Returns
- hxfloat, array-like, shape
The product if the problem is of type QP1 or QP2 (the default), or the product if the problem is of type QP3 or QP4.
- istateNone or int, array-like, shape , optional
Need not be set if the (default) option ‘Cold Start’ is used.
If the option ‘Warm Start’ has been chosen, specifies the desired status of the constraints at the start of the feasibility phase.
More precisely, the first elements of refer to the upper and lower bounds on the variables, and the next elements refer to the general linear constraints (if any).
Possible values for are as follows:
Meaning
0
The corresponding constraint should not be in the initial working set.
1
The constraint should be in the initial working set at its lower bound.
2
The constraint should be in the initial working set at its upper bound.
3
The constraint should be in the initial working set as an equality. This value must not be specified unless .
The values , and are also acceptable but will be reset to zero by the function.
If
qp_dense_solve
has been called previously with the same values of and , already contains satisfactory information. (See also the description of the option ‘Warm Start’.) The function also adjusts (if necessary) the values supplied in to be consistent with .- dataarbitrary, optional
User-communication data for callback functions.
- io_managerFileObjManager, optional
Manager for I/O in this routine.
- Returns
- istateint, ndarray, shape
The status of the constraints in the working set at the point returned in . The significance of each possible value of is as follows:
Meaning
The constraint violates its lower bound by more than the feasibility tolerance.
The constraint violates its upper bound by more than the feasibility tolerance.
The constraint is satisfied to within the feasibility tolerance, but is not in the working set.
This inequality constraint is included in the working set at its lower bound.
This inequality constraint is included in the working set at its upper bound.
This constraint is included in the working set as an equality. This value of can occur only when .
This corresponds to optimality being declared with being temporarily fixed at its current value. This value of can occur only when = 1 on exit.
- xfloat, ndarray, shape
The point at which
qp_dense_solve
terminated. If the function exits successfully or = 1 or 4, contains an estimate of the solution.- iteraint
The total number of iterations performed.
- objfloat
The value of the objective function at if is feasible, or the sum of infeasibilities at otherwise. If the problem is of type FP and is feasible, is set to zero.
- axNone or float, ndarray, shape
The final values of the linear constraints .
- clamdafloat, ndarray, shape
The values of the Lagrange multipliers for each constraint with respect to the current working set. The first elements contain the multipliers for the bound constraints on the variables, and the next elements contain the multipliers for the general linear constraints (if any). If (i.e., constraint is not in the working set), is zero. If is optimal, should be non-negative if , non-positive if and zero if .
- Other Parameters
- ‘Check Frequency’float
Default
Every th iteration, a numerical test is made to see if the current solution satisfies the constraints in the working set. If the largest residual of the constraints in the working set is judged to be too large, the current working set is refactorized and the variables are recomputed to satisfy the constraints more accurately. If , the default value is used.
- ‘Cold Start’valueless
Default
This option specifies how the initial working set is chosen. With a ‘Cold Start’,
qp_dense_solve
chooses the initial working set based on the values of the variables and constraints at the initial point. Broadly speaking, the initial working set will include equality constraints and bounds or inequality constraints that violate or ‘nearly’ satisfy their bounds (to within ‘Crash Tolerance’).With a ‘Warm Start’, you must provide a valid definition of every element of the array .
qp_dense_solve
will override your specification of if necessary, so that a poor choice of the working set will not cause a fatal error. For instance, any elements of which are set to , or will be reset to zero, as will any elements which are set to when the corresponding elements of and are not equal. A warm start will be advantageous if a good estimate of the initial working set is available – for example, whenqp_dense_solve
is called repeatedly to solve related problems.- ‘Warm Start’valueless
This option specifies how the initial working set is chosen. With a ‘Cold Start’,
qp_dense_solve
chooses the initial working set based on the values of the variables and constraints at the initial point. Broadly speaking, the initial working set will include equality constraints and bounds or inequality constraints that violate or ‘nearly’ satisfy their bounds (to within ‘Crash Tolerance’).With a ‘Warm Start’, you must provide a valid definition of every element of the array .
qp_dense_solve
will override your specification of if necessary, so that a poor choice of the working set will not cause a fatal error. For instance, any elements of which are set to , or will be reset to zero, as will any elements which are set to when the corresponding elements of and are not equal. A warm start will be advantageous if a good estimate of the initial working set is available – for example, whenqp_dense_solve
is called repeatedly to solve related problems.- ‘Crash Tolerance’float
Default
This value is used in conjunction with the option ‘Cold Start’ (the default value) when
qp_dense_solve
selects an initial working set. If , the initial working set will include (if possible) bounds or general inequality constraints that lie within of their bounds. In particular, a constraint of the form will be included in the initial working set if . If or , the default value is used.- ‘Defaults’valueless
This special keyword may be used to reset all options to their default values.
- ‘Expand Frequency’int
Default
This option is part of an anti-cycling procedure designed to guarantee progress even on highly degenerate problems.
The strategy is to force a positive step at every iteration, at the expense of violating the constraints by a small amount. Suppose that the value of the option ‘Feasibility Tolerance’ is . Over a period of iterations, the feasibility tolerance actually used by
qp_dense_solve
(i.e., the working feasibility tolerance) increases from to (in steps of ).At certain stages the following ‘resetting procedure’ is used to remove constraint infeasibilities. First, all variables whose upper or lower bounds are in the working set are moved exactly onto their bounds. A count is kept of the number of nontrivial adjustments made. If the count is positive, iterative refinement is used to give variables that satisfy the working set to (essentially) machine precision. Finally, the working feasibility tolerance is reinitialized to .
If a problem requires more than iterations, the resetting procedure is invoked and a new cycle of iterations is started with incremented by . (The decision to resume the feasibility phase or optimality phase is based on comparing any constraint infeasibilities with .)
The resetting procedure is also invoked when
qp_dense_solve
reaches an apparently optimal, infeasible or unbounded solution, unless this situation has already occurred twice. If any nontrivial adjustments are made, iterations are continued.If , the default value is used. If , no anti-cycling procedure is invoked.
- ‘Feasibility Phase Iteration Limit’int
Default
For problems of type FP, the scalar specifies the maximum number of iterations allowed before temination. Setting and means that the workspace needed will be computed and printed, but no iterations will be performed.
For problems of type LP, the maximum number of iterations allowed before temination is taken as . Setting , and means that the workspace needed will be computed and printed, but no iterations will be performed.
For problems of type QP, the scalars and specify the maximum number of iterations allowed in the feasibility and optimality phases. ‘Optimality Phase Iteration Limit’ is equivalent to ‘Iteration Limit’. Setting and means that the workspace needed will be computed and printed, but no iterations will be performed.
If or , the default value is used.
- ‘Optimality Phase Iteration Limit’int
Default
For problems of type FP, the scalar specifies the maximum number of iterations allowed before temination. Setting and means that the workspace needed will be computed and printed, but no iterations will be performed.
For problems of type LP, the maximum number of iterations allowed before temination is taken as . Setting , and means that the workspace needed will be computed and printed, but no iterations will be performed.
For problems of type QP, the scalars and specify the maximum number of iterations allowed in the feasibility and optimality phases. ‘Optimality Phase Iteration Limit’ is equivalent to ‘Iteration Limit’. Setting and means that the workspace needed will be computed and printed, but no iterations will be performed.
If or , the default value is used.
- ‘Feasibility Tolerance’float
Default
If , defines the maximum acceptable absolute violation in each constraint at a ‘feasible’ point. For example, if the variables and the coefficients in the general constraints are of order unity, and the latter are correct to about decimal digits, it would be appropriate to specify as . If , the default value is used.
qp_dense_solve
attempts to find a feasible solution before optimizing the objective function. If the sum of infeasibilities cannot be reduced to zero, the option ‘Minimum Sum of Infeasibilities’ can be used to find the minimum value of the sum. Let Sinf be the corresponding sum of infeasibilities. If Sinf is quite small, it may be appropriate to raise by a factor of or . Otherwise, some error in the data should be suspected.Note that a ‘feasible solution’ is a solution that satisfies the current constraints to within the tolerance .
- ‘Hessian Rows’int
Default
Note that this option does not apply to problems of type FP or LP.
This specifies , the number of rows of the Hessian matrix . The default value of is , the number of variables of the problem.
If the problem is of type QP, then will usually be , the number of variables. However, a value of less than is appropriate for QP3 or QP4 if is an upper trapezoidal matrix with rows. Similarly, may be used to define the dimension of a leading block of nonzeros in the Hessian matrices of QP1 or QP2. In this case the last rows and columns of are assumed to be zero. In the QP case should not be greater than ; if it is, the last rows of are ignored.
If or , the default value is used.
- ‘Infinite Bound Size’float
Default
If , defines the ‘infinite’ bound in the definition of the problem constraints. Any upper bound greater than or equal to will be regarded as (and similarly any lower bound less than or equal to will be regarded as ). If , the default value is used.
- ‘Infinite Step Size’float
Default
If , specifies the magnitude of the change in variables that will be considered a step to an unbounded solution. (Note that an unbounded solution can occur only when the Hessian is not positive definite.) If the change in during an iteration would exceed the value of , then the objective function is considered to be unbounded below in the feasible region. If , the default value is used.
- ‘Iteration Limit’int
Default
See option ‘Feasibility Phase Iteration Limit’.
- ‘Iters’int
Default
See option ‘Feasibility Phase Iteration Limit’.
- ‘Itns’int
Default
See option ‘Feasibility Phase Iteration Limit’.
- ‘List’valueless
Option ‘List’ enables printing of each option specification as it is supplied. ‘Nolist’ suppresses this printing.
- ‘Nolist’valueless
Default
Option ‘List’ enables printing of each option specification as it is supplied. ‘Nolist’ suppresses this printing.
- ‘Maximum Degrees of Freedom’int
Default
Note that this option does not apply to problems of type FP or LP.
This places a limit on the storage allocated for the triangular factor of the reduced Hessian . Ideally, should be set slightly larger than the value of expected at the solution. It need not be larger than , where is the number of variables that appear nonlinearly in the quadratic objective function. For many problems it can be much smaller than .
For quadratic problems, a minimizer may lie on any number of constraints, so that may vary between and . The default value of is, therefore, the number of variables . If ‘Hessian Rows’ is specified, the default value of is the same number, .
- ‘Minimum Sum of Infeasibilities’str
Default
If no feasible point exists for the constraints, then this option is used to control whether or not
qp_dense_solve
will calculate a point that minimizes the constraint violations. If ,qp_dense_solve
will terminate as soon as it is evident that no feasible point exists for the constraints. The final point will generally not be the point at which the sum of infeasibilities is minimized. If ,qp_dense_solve
will continue until the sum of infeasibilities is minimized.- ‘Monitoring File’int
Default
If and , monitoring information produced by
qp_dense_solve
at every iteration is sent to a file with logical unit number .If and/or , no monitoring information is produced.
- ‘Optimality Tolerance’float
Default
If , defines the tolerance used to determine if the bounds and general constraints have the right ‘sign’ for the solution to be judged to be optimal.
If , the default value is used.
- ‘Print Level’int
Default
The value of controls the amount of printout produced by
qp_dense_solve
, as indicated below. A detailed description of the printed output is given in Further Comments (summary output at each iteration and the final solution) and Monitoring Information (monitoring information at each iteration). If , the default value is used.The following printout is sent to the file object associated with the advisory I/O unit (see
FileObjManager
):Output
No output.
The final solution only.
One line of summary output ( characters; see Further Comments) for each iteration (no printout of the final solution).
The final solution and one line of summary output for each iteration.
The following printout is sent to the unit number given by the option ‘Monitoring File’:
Output
No output.
One long line of output ( characters; see Monitoring Information) for each iteration (no printout of the final solution).
At each iteration: the Lagrange multipliers, the variables , the constraint values and the constraint status (see ).
At each iteration: the diagonal elements of the upper triangular matrix associated with the factorization (3) (see Definition of Search Direction) of the working set and the diagonal elements of the upper triangular matrix .
If and the unit number defined by the option ‘Monitoring File’ is the advisory unit number, the summary output is suppressed.
- ‘Problem Type’str
Default QP2
This option specifies the type of objective function to be minimized during the optimality phase. The following are the five optional keywords and the dimensions of the arrays that must be specified in order to define the objective function:
LP
not referenced, length- required;
QP1
rank- symmetric, not referenced;
QP2
rank- symmetric, length- required;
QP3
rank- upper trapezoidal, not referenced;
QP4
rank- upper trapezoidal, length- required.
FP
the objective function is omitted and neither nor are referenced.
For problems of type FP the objective function is omitted and neither nor are referenced.
The following keywords are also acceptable.
Option
Quadratic
QP2
Linear
LP
Feasible
FP
In addition, the keyword QP is equivalent to the default option QP2.
If (i.e., the objective function is purely linear), the efficiency of
qp_dense_solve
may be increased by specifying as LP.- ‘Rank Tolerance’float
Default
Note that this option does not apply to problems of type FP or LP.
This option enables you to control the condition number of the triangular factor (see Algorithmic Details). If denotes the function , the dimension of is defined to be smallest index such that . If , the default value is used.
- Raises
- NagValueError
- (errno )
Too many iterations.
- (errno )
Reduced Hessian exceeds assigned dimension.
- (errno )
On entry, .
Constraint: .
- (errno )
On entry, .
Constraint: .
- (errno )
On entry, the equal bounds on are infinite, because and , but : and .
- (errno )
On entry, the bounds on are inconsistent: and .
- (errno )
On entry with a Warm Start, .
- (errno )
On entry with a Cold Start, .
- (errno )
On entry, the equal bounds on variable are infinite, because and , but : and .
- (errno )
On entry, the equal bounds on linear constraint are infinite, because and , but : and .
- (errno )
On entry, the equal bounds on nonlinear constraint are infinite, because and , but : and .
- (errno )
On entry, the bounds on variable are inconsistent: and .
- (errno )
On entry, the bounds on linear constraint are inconsistent: and .
- (errno )
On entry, the bounds on nonlinear constraint are inconsistent: and .
- (errno )
‘Problem Type’ not recognized. Problem abandoned.
- Warns
- NagAlgorithmicWarning
- (errno )
Iterations terminated at a dead point.
- (errno )
Optimal solution is not unique.
- (errno )
solution is unbounded.
- (errno )
No feasible point for the linear constraints.
- (errno )
Cannot satisfy the working set constraints to the accuracy requested.
- Notes
In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.
qp_dense_solve
is designed to solve a class of quadratic programming problems that are assumed to be stated in the following general form:where is an matrix and may be specified in a variety of ways depending upon the particular problem to be solved. The available forms for are listed in Table [label omitted], in which the prefixes FP, LP and QP stand for ‘feasible point’, ‘linear programming’ and ‘quadratic programming’ respectively and is an -element vector.
Problem type
Matrix
FP
Not applicable
Not applicable
LP
Not applicable
QP1
symmetric
QP2
symmetric
QP3
upper trapezoidal
QP4
upper trapezoidal
There is no restriction on or apart from symmetry. If the quadratic function is convex, a global minimum is found; otherwise, a local minimum is found. The default problem type is QP2 and other objective functions are selected by using the option ‘Problem Type’. For problems of type FP, the objective function is omitted and the function attempts to find a feasible point for the set of constraints.
The constraints involving are called the general constraints. Note that upper and lower bounds are specified for all the variables and for all the general constraints. An equality constraint can be specified by setting . If certain bounds are not present, the associated elements of or can be set to special values that will be treated as or . (See the description of the option ‘Infinite Bound Size’.)
The defining feature of a quadratic function is that the second-derivative matrix (the Hessian matrix) is constant. For QP1 and QP2 (the default), ; for QP3 and QP4, ; and for the LP case, . If is positive semidefinite, it is usually more efficient to use
lsq_lincon_solve()
. If is defined as the zero matrix,qp_dense_solve
will still attempt to solve the resulting linear programming problem; however, this can be accomplished more efficiently by setting the option the problem is of type ‘LP’, or by usinglp_solve()
instead.You must supply an initial estimate of the solution.
In the QP case, you may supply either explicitly as an matrix, or implicitly in a function that computes the product or for any given vector .
In general, a successful run of
qp_dense_solve
will indicate one of three situations:a minimizer has been found;
the algorithm has terminated at a so-called dead-point; or
the problem has no bounded solution.
If a minimizer is found, and is positive definite or positive semidefinite,
qp_dense_solve
will obtain a global minimizer; otherwise, the solution will be a local minimizer (which may or may not be a global minimizer). A dead-point is a point at which the necessary conditions for optimality are satisfied but the sufficient conditions are not. At such a point, a feasible direction of decrease may or may not exist, so that the point is not necessarily a local solution of the problem. Verification of optimality in such instances requires further information, and is in general an NP-hard problem (see Pardalos and Schnitger (1988)). Termination at a dead-point can occur only if is not positive definite. If is positive semidefinite, the dead-point will be a weak minimizer (i.e., with a unique optimal objective value, but an infinite set of optimal ).The method used by
qp_dense_solve
(see Algorithmic Details) is most efficient when many constraints or bounds are active at the solution.
- References
Gill, P E, Hammarling, S, Murray, W, Saunders, M A and Wright, M H, 1986, Users’ guide for LSSOL (Version 1.0), Report SOL 86-1, Department of Operations Research, Stanford University
Gill, P E and Murray, W, 1978, Numerically stable methods for quadratic programming, Math. Programming (14), 349–372
Gill, P E, Murray, W, Saunders, M A and Wright, M H, 1984, Procedures for optimization problems with a mixture of bounds and general linear constraints, ACM Trans. Math. Software (10), 282–298
Gill, P E, Murray, W, Saunders, M A and Wright, M H, 1989, A practical anti-cycling procedure for linearly constrained optimization, Math. Programming (45), 437–474
Gill, P E, Murray, W, Saunders, M A and Wright, M H, 1991, Inertia-controlling methods for general quadratic programming, SIAM Rev. (33), 1–36
Gill, P E, Murray, W and Wright, M H, 1981, Practical Optimization, Academic Press
Pardalos, P M and Schnitger, G, 1988, Checking local optimality in constrained quadratic programming is NP-hard, Operations Research Letters (7), 33–35