Processing math: 100%
KASKADE 7 development version
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Public Types | Public Member Functions | Static Public Attributes | List of all members
Kaskade::QPSolver< d, Real > Class Template Reference

An iterative solver for small instances of a particular class of quadratic programs. More...

#include <qp.hh>

Detailed Description

template<int d, class Real = double>
class Kaskade::QPSolver< d, Real >

An iterative solver for small instances of a particular class of quadratic programs.

The solver addresses problems of the form

\min_x \frac{1}{2} x^T Ax + c^T x \quad \text{s.t. } Bx \le b,

where A is positive definite. The implementation works with dense matrices and is intended to solve small problems of just a couple of variables, typically less than 30, and moreover well conditioned problems.

Template Parameters
Realthe real field type of matrix/vector entries. Explicit instantiations are provided for float and double.
dthe dimension of (vector valued) optimization variables

Solution algorithm

The solver forms the Lagrange dual problem

\max_\lambda - \frac{1}{2} \lambda^TBA^{-1}B^T\lambda + (BA^{-1}c+b)^T \lambda \quad\text{s.t. } \lambda \le 0,

which is bound constrained and can be solved by a projected gradient method. As the primal problem is convex, strong duality holds, and the unique primal minimizer x^* can be computed from a dual maximizer \lambda^* by

x^* = A^{-1}(B^T\lambda^*-c).

Definition at line 423 of file qp.hh.

Public Types

using VectorX = Dune::BlockVector< Dune::FieldVector< Real, d > >
 
using VectorB = Dune::BlockVector< Dune::FieldVector< Real, 1 > >
 
using MatrixA = DynamicMatrix< Dune::FieldMatrix< Real, d, d > >
 
using MatrixB = DynamicMatrix< Dune::FieldMatrix< Real, 1, d > >
 
using Self = QPSolver< d, Real >
 

Public Member Functions

 QPSolver (MatrixA const &A, MatrixB const &B, QPConvexificationStrategy convex=QPConvexificationStrategy::DONOTHING)
 Constructs the solver, providing the matrices A and B . More...
 
 QPSolver (Self const &qp)
 
Selfoperator= (Self const &qp)
 
std::tuple< VectorX, VectorB, int > solve (VectorX const &c, VectorB const &b, double tol) const
 Solves the QP for the given right hand side vectors. More...
 
MatrixB const & matrixB () const
 Provides access to the internal copy of B. More...
 

Static Public Attributes

static constexpr QPStructure sparse = QPStructure::DENSE
 Given info on whether dense or sparse arithmetic is used. More...
 

Member Typedef Documentation

◆ MatrixA

template<int d, class Real = double>
using Kaskade::QPSolver< d, Real >::MatrixA = DynamicMatrix<Dune::FieldMatrix<Real,d,d> >

Definition at line 430 of file qp.hh.

◆ MatrixB

template<int d, class Real = double>
using Kaskade::QPSolver< d, Real >::MatrixB = DynamicMatrix<Dune::FieldMatrix<Real,1,d> >

Definition at line 431 of file qp.hh.

◆ Self

template<int d, class Real = double>
using Kaskade::QPSolver< d, Real >::Self = QPSolver<d,Real>

Definition at line 432 of file qp.hh.

◆ VectorB

template<int d, class Real = double>
using Kaskade::QPSolver< d, Real >::VectorB = Dune::BlockVector<Dune::FieldVector<Real,1> >

Definition at line 429 of file qp.hh.

◆ VectorX

template<int d, class Real = double>
using Kaskade::QPSolver< d, Real >::VectorX = Dune::BlockVector<Dune::FieldVector<Real,d> >

Definition at line 428 of file qp.hh.

Constructor & Destructor Documentation

◆ QPSolver() [1/2]

template<int d, class Real = double>
Kaskade::QPSolver< d, Real >::QPSolver ( MatrixA const &  A,
MatrixB const &  B,
QPConvexificationStrategy  convex = QPConvexificationStrategy::DONOTHING 
)

Constructs the solver, providing the matrices A and B .

The argument matrices are copied internally, and are not referenced later.

Throws a NonpositiveMatrixException if A is not positive definite and convex is QPConvexificationStrategy::DONOTHING.

Todo:
write a move constructor

◆ QPSolver() [2/2]

template<int d, class Real = double>
Kaskade::QPSolver< d, Real >::QPSolver ( Self const &  qp)

Member Function Documentation

◆ matrixB()

template<int d, class Real = double>
MatrixB const & Kaskade::QPSolver< d, Real >::matrixB ( ) const
inline

Provides access to the internal copy of B.

Definition at line 474 of file qp.hh.

◆ operator=()

template<int d, class Real = double>
Self & Kaskade::QPSolver< d, Real >::operator= ( Self const &  qp)

◆ solve()

template<int d, class Real = double>
std::tuple< VectorX, VectorB, int > Kaskade::QPSolver< d, Real >::solve ( VectorX const &  c,
VectorB const &  b,
double  tol 
) const

Solves the QP for the given right hand side vectors.

Solves the QP by solving the dual QP with a projected gradient method starting at the dually feasible point \lambda = 0.

Todo:
add option to provide an initial guess for \lambda .
Parameters
cthe linear term of the objective
bthe constraints right hand side
tolthe tolerance for termination

The iteration is terminated as soon as the projected gradient s of the dual problem is small, i.e. \|s\|_2 \le \mathsf{tol}.

Returns
a tuple (x,\lambda,iter) of solution vector, multiplier, and iteration count

If the QP is infeasible, an InfeasibleProblemException is thrown.

Member Data Documentation

◆ sparse

template<int d, class Real = double>
constexpr QPStructure Kaskade::QPSolver< d, Real >::sparse = QPStructure::DENSE
staticconstexpr

Given info on whether dense or sparse arithmetic is used.

Definition at line 483 of file qp.hh.


The documentation for this class was generated from the following file: