ConstrainedOptPack: C++ Tools for Constrained (and Unconstrained) Optimization
Version of the Day
|
Represents the QP to be solved by QPSchur {abstract}. More...
#include <ConstrainedOptPack_QPSchur.hpp>
Public Types | |
typedef vector_one_based_checked < EBounds > | x_init_t |
typedef vector_one_based_checked < size_type > | l_x_X_map_t |
typedef vector_one_based_checked < size_type > | i_x_X_map_t |
typedef QPSchurPack::Constraints | Constraints |
Public Member Functions | |
virtual | ~QP () |
virtual size_type | n () const =0 |
virtual size_type | m () const =0 |
virtual const DVectorSlice | g () const =0 |
virtual const MatrixSymOp & | G () const =0 |
virtual const MatrixOp & | A () const =0 |
If m == 0 then don't call this, it may throw an exception or worse. More... | |
virtual size_type | n_R () const =0 |
virtual const x_init_t & | x_init () const =0 |
Return the status of a variable initially. More... | |
virtual const l_x_X_map_t & | l_x_X_map () const =0 |
Map from full x(i) to initially fixed x_X(l). More... | |
virtual const i_x_X_map_t & | i_x_X_map () const =0 |
Map from initially fixed x_X(l) to full x(i). More... | |
virtual const DVectorSlice | b_X () const =0 |
virtual const GenPermMatrixSlice & | Q_R () const =0 |
(Q_R().ordered_by() == BY_ROW) More... | |
virtual const GenPermMatrixSlice & | Q_X () const =0 |
(Q_X().ordered_by() == BY_ROW) More... | |
virtual const MatrixSymOpNonsing & | Ko () const =0 |
virtual const DVectorSlice | fo () const =0 |
virtual Constraints & | constraints ()=0 |
virtual const Constraints & | constraints () const =0 |
virtual void | dump_qp (std::ostream &out) |
Dump the definition of the QP to a stream. More... | |
Represents the QP to be solved by QPSchur {abstract}.
In order to solve a QP, clients must define subclasses for this interface and the Constraints
interface defined later. This is where the specialized properties of the QP are exploited. For certain types of QPs, standard implementation classes can be defined.
Here the QP is of the form:
(1.a) min g'*x + 1/2*x'*G*x (1.b) s.t. A'*x = c (1.c) cl_bar <= A_bar'*x <= cu_bar where: x <: R^n g <: R^n G <: R^(n x n) A <: R^(n x m) c <: R^m A_bar <: R^(n x m_bar) cl_bar, cu_bar <: R^m_bar
Above, cl_bar <= A_bar'*x <= cu_bar
may represent variable bounds, general inequalities and equality constraints and these constraints are represented by the class Constraints
.
When solving the above QP using the schur complement QP solver we start out with a KKT system with a subset of variables initially fixed at a bound:
[ G_RR G_RX A_F 0 ] [ x_R ] [ - g_R ] [ G_RX' G_XX A_X I ] [ x_X ] [ - g_X ] [ A_R' A_X' 0 0 ]*[ lambda ] = [ c ] [ I 0 0 ] [ mu_X ] [ b_X ]
We can simplify the above system by solving for the initially fixed variables and removing them from the initial KKT system to give:
x_X = b_X [ G_RR A_R ] [ x_R ] [ - g_R - G_RX*b_X ] [ A_R' 0 ] * [ lambda ] = [ c - A_X'*b_X ] \_______________/ \________/ \__________________/ Ko vo fo mu_X = - g_X - G_RX'*x_R - G_X*b_X - A_X*lambda where: n_X = n - n_R x_R = Q_R'*x <: R^n_R g_R = Q_R'*g <: R^n_R G_RR = Q_R'*G*Q_R <: R^(n_R x n_R) G_RX = Q_R'*G*Q_X <: R^(n_R x n_X) G_XX = Q_X'*G*Q_X <: R^(n_X x n_X) A_R = Q_R'*A <: R^(n_R x m) A_X = Q_X'*A <: R^(n_X x m) Q_R <: R^(n x n_R) Q_X <: R^(n x n_X)
This class is an interface for encapsulating the QP. Operations are available for accessing g
, G
, A
, Ko
, vo
, and fo
as well as the mapping matrices Q_R
and Q_X
(both ordered by row). Also, operations are available for accessing data structures that describe the set of initially fixed and free variables. See the Constraints
interface for how to access the constraints in (1.c) and the matrix A_bar
.
Definition at line 166 of file ConstrainedOptPack_QPSchur.hpp.
typedef vector_one_based_checked<EBounds> ConstrainedOptPack::QPSchurPack::QP::x_init_t |
Definition at line 173 of file ConstrainedOptPack_QPSchur.hpp.
Definition at line 175 of file ConstrainedOptPack_QPSchur.hpp.
Definition at line 177 of file ConstrainedOptPack_QPSchur.hpp.
Definition at line 180 of file ConstrainedOptPack_QPSchur.hpp.
|
inlinevirtual |
Definition at line 186 of file ConstrainedOptPack_QPSchur.hpp.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
If m == 0 then don't call this, it may throw an exception or worse.
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Return the status of a variable initially.
For 1 <= i <= n:
/ FREE : x(i) is initially free | LOWER : x(i) is initially fixed at xl(i) x_init(i) = | UPPER : x(i) is initially fixed at xu(i) \ EQUALITY : x(i) fixed at xl(i) = xu(i) always
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Map from full x(i) to initially fixed x_X(l).
For 1 <= i <= n:
/ l : x(i) = x_X(l) = b_X(l) initially (1 <= l <= n_X) l_x_X_map(i) = | \ 0 : otherwise
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Map from initially fixed x_X(l) to full x(i).
For 1 <= l <= n_X:
i_x_X_map(l) = i : x(i) = x_X(l) = b_X(l) initially (1 <= i <= n)
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
(Q_R().ordered_by() == BY_ROW)
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
(Q_X().ordered_by() == BY_ROW)
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
pure virtual |
Implemented in ConstrainedOptPack::QPSchurPack::QPInitFixedFreeStd.
|
virtual |
Dump the definition of the QP to a stream.
This function is only to be used for debugging small problems.
Definition at line 5022 of file ConstrainedOptPack_QPSchur.cpp.