jrl::qp::test Namespace Reference

Classes

struct  FeasibilityConstraints
 
struct  LeastSquareProblem
 
struct  ProblemCharacteristics
 
struct  QPProblem
 
struct  RandomLeastSquare
 
struct  scalar_normal_random_op
 
struct  SeparatedFeasibilityConstraints
 

Functions

LeastSquareProblem JRLQP_DLLAPI generateBoxAndSingleConstraintProblem (int nbVar, bool act, double actLevel=0.5)
 
void JRLQP_DLLAPI checkDimensions (const MatrixConstRef &G, const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false)
 
void JRLQP_DLLAPI checkDimensions (const VectorConstRef &x, const VectorConstRef &u, const MatrixConstRef &G, const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false)
 
void JRLQP_DLLAPI checkDimensions (const VectorConstRef &x, const VectorConstRef &u, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false)
 
bool JRLQP_DLLAPI testKKT (const VectorConstRef &x, const VectorConstRef &u, const MatrixConstRef &G, const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false, double tau_p=1e-6, double tau_d=1e-6)
 
bool JRLQP_DLLAPI testKKT (const VectorConstRef &x, const VectorConstRef &u, const QPProblem<> &pb, double tau_p=1e-6, double tau_d=1e-6)
 
bool JRLQP_DLLAPI testKKTStationarity (const VectorConstRef &x, const VectorConstRef &u, const MatrixConstRef &G, const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false, double tau_d=1e-6)
 
bool JRLQP_DLLAPI testKKTStationarity (const VectorConstRef &x, const VectorConstRef &u, const QPProblem<> &pb, double tau_d=1e-6)
 
bool JRLQP_DLLAPI testKKTFeasibility (const VectorConstRef &x, const VectorConstRef &u, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC=false, double tau_p=1e-6, double tau_d=1e-6)
 
bool JRLQP_DLLAPI testKKTFeasibility (const VectorConstRef &x, const VectorConstRef &u, const FeasibilityConstraints &cstr, double tau_p=1e-6, double tau_d=1e-6)
 
auto randnVec (Eigen::Index size, double mean=0, double stddev=1)
 
auto randnMat (Eigen::Index rows, Eigen::Index cols, double mean=0, double stddev=1)
 
Eigen::VectorXd randUnitVec (Eigen::Index size)
 
Eigen::MatrixXd randOrtho (Eigen::Index size, bool special=false)
 
Eigen::MatrixXd randn (Eigen::Index rows, Eigen::Index cols, Eigen::Index rank=-1)
 
std::pair< Eigen::MatrixXd, Eigen::MatrixXd > randDependent (Eigen::Index cols, Eigen::Index rowsA, Eigen::Index rankA, Eigen::Index rowsB, Eigen::Index rankB, Eigen::Index rankAB)
 
RandomLeastSquare JRLQP_DLLAPI randomProblem (const ProblemCharacteristics &characs)
 
void checkDimensions ([[maybe_unused]] int n, [[maybe_unused]] const MatrixConstRef &C, [[maybe_unused]] const VectorConstRef &bl, [[maybe_unused]] const VectorConstRef &bu, [[maybe_unused]] const VectorConstRef &xl, [[maybe_unused]] const VectorConstRef &xu, [[maybe_unused]] bool transposedC)
 
void checkDimensions ([[maybe_unused]] const MatrixConstRef &G, [[maybe_unused]] const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC)
 
void checkDimensions ([[maybe_unused]] const VectorConstRef &x, [[maybe_unused]] const VectorConstRef &u, const MatrixConstRef &G, const VectorConstRef &a, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC)
 
void checkDimensions ([[maybe_unused]] const VectorConstRef &x, [[maybe_unused]] const VectorConstRef &u, const MatrixConstRef &C, const VectorConstRef &bl, const VectorConstRef &bu, const VectorConstRef &xl, const VectorConstRef &xu, bool transposedC)
 
template<typename Derived >
void disp (const std::string &name, const MatrixBase< Derived > &M)
 

Variables

struct JRLQP_DLLAPI SeparatedFeasibilityConstraints
 

Function Documentation

◆ checkDimensions() [1/7]

void jrl::qp::test::checkDimensions ( [[maybe_unused] ] const MatrixConstRef G,
[[maybe_unused] ] const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC 
)

◆ checkDimensions() [2/7]

void jrl::qp::test::checkDimensions ( [[maybe_unused] ] const VectorConstRef x,
[[maybe_unused] ] const VectorConstRef u,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC 
)

◆ checkDimensions() [3/7]

void jrl::qp::test::checkDimensions ( [[maybe_unused] ] const VectorConstRef x,
[[maybe_unused] ] const VectorConstRef u,
const MatrixConstRef G,
const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC 
)

◆ checkDimensions() [4/7]

void jrl::qp::test::checkDimensions ( [[maybe_unused] ] int  n,
[[maybe_unused] ] const MatrixConstRef C,
[[maybe_unused] ] const VectorConstRef bl,
[[maybe_unused] ] const VectorConstRef bu,
[[maybe_unused] ] const VectorConstRef xl,
[[maybe_unused] ] const VectorConstRef xu,
[[maybe_unused] ] bool  transposedC 
)

◆ checkDimensions() [5/7]

void JRLQP_DLLAPI jrl::qp::test::checkDimensions ( const MatrixConstRef G,
const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false 
)

Check that all the matrix and vector dimensions of the following QP problem are consistent: min. .5 x^T G x + a^T x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. When transposedC is true, the constraints read bl <= C^T x <= bu

◆ checkDimensions() [6/7]

void JRLQP_DLLAPI jrl::qp::test::checkDimensions ( const VectorConstRef x,
const VectorConstRef u,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false 
)

Check that all the matrix and vector dimensions of the following FP problem are consistent: find x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. Further check that the solution vector x and the vector of Lagrange multipliers have the correct dimensions. When transposedC is true, the constraints read bl <= C^T x <= bu

◆ checkDimensions() [7/7]

void JRLQP_DLLAPI jrl::qp::test::checkDimensions ( const VectorConstRef x,
const VectorConstRef u,
const MatrixConstRef G,
const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false 
)

Check that all the matrix and vector dimensions of the following QP problem are consistent: min. .5 x^T G x + a^T x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. Further check that the solution vector x and the vector of Lagrange multipliers have the correct dimensions. When transposedC is true, the constraints read bl <= C^T x <= bu

◆ disp()

template<typename Derived >
void jrl::qp::test::disp ( const std::string &  name,
const MatrixBase< Derived > &  M 
)

◆ generateBoxAndSingleConstraintProblem()

LeastSquareProblem jrl::qp::test::generateBoxAndSingleConstraintProblem ( int  nbVar,
bool  act,
double  actLevel = 0.5 
)

Generate a problem min. || x- x0 || s.t. c'x >= bl (1) xl<=x<=xu with nbVar variables and the constraint c'x >= bl active if act is true Noting xb the solution without the constraint (1) and su the vertex of the box [xl, xu] the furthest away in the direction given by c, actLevel indicates where the constraint plan intersects the segement [xb, su].

◆ randDependent()

std::pair<Eigen::MatrixXd, Eigen::MatrixXd> jrl::qp::test::randDependent ( Eigen::Index  cols,
Eigen::Index  rowsA,
Eigen::Index  rankA,
Eigen::Index  rowsB,
Eigen::Index  rankB,
Eigen::Index  rankAB 
)
inline

Returns two matrices A and B with specified sizes and ranks and such that [A;B] has a given rank

Parameters
colsColumn size of A and B.
rowsARow size of A.
rankARank of A.
rowsBRow size of B.
rankARank of B.
rankABRank of [AB].

◆ randn()

Eigen::MatrixXd jrl::qp::test::randn ( Eigen::Index  rows,
Eigen::Index  cols,
Eigen::Index  rank = -1 
)
inline

Generate a random matrix with a specified rank.

Coefficients of the matrix have empirically a normal distribution with a mean of 0 and a standard deviation of 1.

Parameters
rowsNumber of rows
colsNumber of cols
rankRank of the matrix. Must be lower or equal to both rows and cols. If let to is default value or set negative, the minimum of rows and cols is taken instead.

◆ randnMat()

auto jrl::qp::test::randnMat ( Eigen::Index  rows,
Eigen::Index  cols,
double  mean = 0,
double  stddev = 1 
)
inline

Generate a random matrix with elements following a normal distribution.

Parameters
rowsNumber of rows of the matrix.
colsNumber of columns of the matrix.
meanMean parameter of the normal distribution.
stddevStandard deviation parameter of the normal distribution.

◆ randnVec()

auto jrl::qp::test::randnVec ( Eigen::Index  size,
double  mean = 0,
double  stddev = 1 
)
inline

Generate a random vector with elements following a normal distribution.

Parameters
sizeSize of the vector.
meanMean parameter of the normal distribution.
stddevStandard deviation parameter of the normal distribution.

◆ randomProblem()

RandomLeastSquare jrl::qp::test::randomProblem ( const ProblemCharacteristics characs)

◆ randOrtho()

Eigen::MatrixXd jrl::qp::test::randOrtho ( Eigen::Index  size,
bool  special = false 
)
inline

Generate a random orthogonal matrix with Haar distribution.

Parameters
sizeNumber of rows and columns of the matrix.
specialIf true, the determinant of the generated matrix is 1 (i.e the matrix belongs to SO(size)). Otherwise it is 1 or -1 with 50% probability each (it belongs to O(size))

◆ randUnitVec()

Eigen::VectorXd jrl::qp::test::randUnitVec ( Eigen::Index  size)
inline

Generate a random unit vector with uniform distribution on the unit sphere.

Parameters
sizeSize of the vector.

◆ testKKT() [1/2]

bool jrl::qp::test::testKKT ( const VectorConstRef x,
const VectorConstRef u,
const MatrixConstRef G,
const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false,
double  tau_p = 1e-6,
double  tau_d = 1e-6 
)

Check the KKT conditions at the dual pair (x,u) for the problem min. .5 x^T G x + a^T x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. u is the vector of Lagrange multipliers with multipliers corresponding to the constraints first and multipliers corresponding to the bound last.

◆ testKKT() [2/2]

bool JRLQP_DLLAPI jrl::qp::test::testKKT ( const VectorConstRef x,
const VectorConstRef u,
const QPProblem<> &  pb,
double  tau_p = 1e-6,
double  tau_d = 1e-6 
)

Overload for a problem given by a QPProblem object.

◆ testKKTFeasibility() [1/2]

bool JRLQP_DLLAPI jrl::qp::test::testKKTFeasibility ( const VectorConstRef x,
const VectorConstRef u,
const FeasibilityConstraints cstr,
double  tau_p = 1e-6,
double  tau_d = 1e-6 
)

Overload for a problem given by a QPProblem object.

◆ testKKTFeasibility() [2/2]

bool jrl::qp::test::testKKTFeasibility ( const VectorConstRef x,
const VectorConstRef u,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false,
double  tau_p = 1e-6,
double  tau_d = 1e-6 
)

Check all KKT conditions other than the stationarity at the dual pair (x,u) for the problem min. .5 x^T G x + a^T x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. u is the vector of Lagrange multipliers with multipliers corresponding to the constraints first and multipliers corresponding to the bound last.

◆ testKKTStationarity() [1/2]

bool jrl::qp::test::testKKTStationarity ( const VectorConstRef x,
const VectorConstRef u,
const MatrixConstRef G,
const VectorConstRef a,
const MatrixConstRef C,
const VectorConstRef bl,
const VectorConstRef bu,
const VectorConstRef xl,
const VectorConstRef xu,
bool  transposedC = false,
double  tau_d = 1e-6 
)

Check the stationnarity part of the KKT conditions (i.e. dL/dx = 0, where L is the Lagragian), at the dual pair (x,u) for the problem min. .5 x^T G x + a^T x s.t. bl <= C x <= bu xl <= x <= xu where the variable bounds are optional. u is the vector of Lagrange multipliers with multipliers corresponding to the constraints first and multipliers corresponding to the bound last.

◆ testKKTStationarity() [2/2]

bool JRLQP_DLLAPI jrl::qp::test::testKKTStationarity ( const VectorConstRef x,
const VectorConstRef u,
const QPProblem<> &  pb,
double  tau_d = 1e-6 
)

Overload for a problem given by a QPProblem object.

Variable Documentation

◆ SeparatedFeasibilityConstraints