Quadratic Programming

From NEOS

The quadratic programming problem involves minimization of a quadratic function subject to linear constraints. Most codes use the formulation

minimize	$\frac{1}{2} x^T Q x + c^T x$
subject to	$a_i^T x = b_i$	for $i \in E\qquad$ (QP)
	$a_i^T x \geq b_i$	for $i \in I$ ,

where $Q \in R^{n\times n}$ is symmetric, and the index sets I and E specify the inequality and equality constraints, respectively.

The difficulty of solving the quadratic programming problem depends largely on the nature of the matrix Q. In convex quadratic programs, which are relatively easy to solve, the matrix Q is positive semidefinite (on the feasible set). If Q has negative eigenvalues-nonconvex quadratic programming-then the objective function may have more than one local minimizer. An extreme example is the problem

$\min \; - x^T x \; : \; -1 \leq x_i \leq 1, \; i=1,\ldots,n$ ,

which has a minimizer at any x with

| x i | = 1

for i = 1,..., n - a total of $2 n$ local minimizers.

The codes in the BQPD, LINDO, LSSOL, PORT 3, and QPOPT packages are based on active set methods. After finding a feasible point during an initial phase, these methods search for a solution along the edges and faces of the feasible set by solving a sequence of equality-constrained quadratic programming problems. Active set methods differ from the simplex method for linear programming in that neither the iterates nor the solution need be vertices of the feasible set. When the quadratic programming problem is nonconvex, these methods usually find a local minimizer. Finding a global minimizer is a more difficult task that is not addressed by the software currently available.

Equality-constrained quadratic programs arise, not only as subproblems in solving the general problem, but also in structural analysis and other areas of application. Null-space methods for solving

minimize	$\frac{1}{2} x^T Q x + c^T x$
subject to	$A x = b \qquad$ (EQP)

find a full-rank matrix $Z\in R^{n\times m}$ , such that Z spans the null space of A. This matrix can be computed with orthogonal factorizations or, in the case of sparse problems, by LU factorization of a submatrix of A, just as in the simplex method for linear programming. Given a feasible vector x₀, we can express any other feasible vector x in the form

x = x 0 + Z w

for some $w \in R^m$ . Direct computation shows that the equality-constrained subproblem (EQP) is equivalent to the unconstrained subproblem

$\min_w \; \frac{1}{2} w^T (Z^T Q Z) w + (Q x_0 + c)^T Z w$

If the reduced Hessian matrix $Z T Q Z$ is positive definite, then the unique solution $w *$ of this subproblem can be obtained by solving the linear system

(Z T Q Z) w = - Z T (Q x 0 + c)

The solution $x *$ of the equality-constrained subproblem (EQP) is then recovered by using $x = x 0 + Z w$ . Lagrange multipliers can be computed from $x *$ by noting that the first-order condition for optimality in (EQP) is that there exists a multiplier vector $y *$ such that

Q x * + c + A T y * = 0

If A has full rank, then

y * = - (A T A) - 1 A (Q x * + c)

is the unique set of multipliers. Most codes uses null-space methods. Range-space methods for (EQP) can be used when Q is positive definite and easy to invert, for example, diagonal or block-diagonal. In this approach, the solution and the multiplier vectors are calculated from the formulae

y * = - (A Q - 1 A) - 1 (b + A Q - 1 c)

and

x * = - Q - 1 (c + A T y *)

Although this approach works only for a subclass of problems, there are many applications in which it is useful.

Active set methods for the inequality-constrained problem (QP) solve a sequence of equality-constrained problems. Given a feasible $x k$ , these methods find a direction $d k$ by solving the subproblem

minimize	$q (x k + d)$
subject to	$a_i^T(x_k + d) = b_i\qquad i \in W_k$

where q is the objective function

$q(x) = \frac{1}{2} x^T Q x + c^T x$

and $W k$ is a working set of constraints. In all cases $W k$ is a subset of

$A(x_k) = \{ i \in I : a_i^T x_k = b_i \} \cup E$

the set of constraints that are active at $x k$ . Typically, $W k$ either is equal to $A (x k)$ or else has one fewer index than $A (x k)$ .

The working set $W k$ is updated at each iteration with the aim of determining the set $A *$ of active constraints at a solution $x *$ . When $W k$ is equal to $A *$ , a local minimizer of the original problem can be obtained as a solution of the equality-constrained subproblem (EQP). The updating of $W k$ depends on the solution of the direction-finding subproblem.

Subproblem (1.4) has a solution if the reduced Hessian matrix $Z_k^T Q Z_k$ is positive definite. This is always the case if $Q$ is positive definite. If subproblem [1.4) has a solution $d k$ , we compute the largest possible step

$\mu_k = \max\{ \frac{b_i - a_i^T x_k}{a_i^T d_k}: \; a_i^T d_k > 0, i \not \in W_k \}$

that does not violate any constraints, and we set $x k + 1 = x k + α k d k$ , where $α k = min{1,μ k}$ .

The step $α k = 1$ would take us to the minimizer of the objective function on the subspace defined by the current working set, but it may be necessary to truncate this step if a new constraint is encountered. The working set is updated by including in $W k + 1$ all constraints active at $x k + 1$ .

If the solution to subproblem (1.4) is $d k = 0$ , then $x k$ is the minimizer of the objective function on the subspace defined by $W k$ . First-order optimality conditions for (1.4) imply that there are multipliers $y_i^{(k)}$ such that $Q x_k + c + \sum_{i \in W_k} y_i^{(k)} a_i = 0$ .

If $y_i^{(k)} \geq 0$ for $i \in W_k$ , then $x k$ is a local minimizer of problem (1.1). Otherwise, we obtain $W k + 1$ by deleting one of the indices i for which $y_i^{(k)} < 0$ . As in the case of linear programming, various pricing schemes for making this choice can be implemented.

If the reduced Hessian matrix $Z_k^T Q Z_k$ is indefinite, then subproblem (1.4) is unbounded below. In this case we need to determine a direction $d k$ such that $q (x k + α d k$ is unbounded below, using techniques based on factorizations of the reduced Hessian matrix. Given $d k$ , we compute $μ k$ as in (1.5) and define $x k + 1 = x k + μ d k$ .

The new working set $W k + 1$ is obtained by adding to $W k$ all constraints active at $x k + 1$ .

A key to the efficient implementation of active set methods is the reuse of information from solving the equality-constrained subproblem at the next iteration. The only difference between consecutive subproblems is that the working set grows or shrinks by a single component. Efficient codes perform updates of the matrix factorizations obtained at the previous iteration, rather than calculating them from scratch each time.

The LSSOL package (duplicated in NAG) is specifically designed for convex quadratic programs and linearly constrained linear least squares problems. It is not aimed at large-scale problems; the constraint matrices and the Hessian Q are all specified in dense storage format. The quadratic programming routine in NLPQL has the same properties. IMSL contains codes for dense quadratic programs. If the matrix Q is not positive definite, it is replaced by

Q + γ I

where $\gamma \geq 0$ is chosen large enough to force convexity.

BQPD uses a null-space method to solve quadratic programs that are not necessarily convex. The linear algebra operations are performed in a modular way; the user is allowed to choose between sparse and dense matrix algebra. The reduced Hessian matrix is, however, processed as a dense matrix, even when sparse techniques are used to handle Q and the constraints. The code is efficient for large-scale problems when the size of the working set is close to n. LINDO also takes account of sparsity, while MATLAB, and QPOPT (also available in the NAG library) are designed for dense quadratic programs that are not necessarily convex.

Linear least squares problems are special instances of convex quadratic programs that arise frequently in data-fitting applications. The linear least squares problem

minimize	$\frac{1}{2} \\| C x - d \\|_2^2$
subject to	$a_i^T x = b_i$	for $i \in E\qquad$ (QP)
	$a_i^T x \geq b_i$	for $i \in I$ ,

where $C \in R^{m\times n}$ and $d \in R^m$ , is a special case of problem (1.1); we can see this by replacing Q by $C T C$ and c by $C T d$ in (1.1). In general, it is preferable to solve a least squares problem with a code that takes advantage of the special structure of the least squares problem (for example, LSSOL).

Algorithms for solving linear least squares problems tend to rely on null-space active set methods. For a least squares problem the null-space matrix $Z$ can be obtained from the $Q R$ factorization of $C$ explicit formation of $C T C$ is avoided, since $C T C$ is usually less well conditioned than $C$ .

Retrieved from "http://wiki.mcs.anl.gov/NEOS/index.php/Quadratic_Programming"

Quadratic Programming

From NEOS

Views

Personal tools

Navigation

Search

Toolbox