#include <qpOASES_e/Bounds.h>
#include <qpOASES_e/Options.h>
#include <qpOASES_e/Matrices.h>
#include <qpOASES_e/Constraints.h>
#include <qpOASES_e/Flipper.h>
#include <qpOASES_e/ConstraintProduct.h>
Go to the source code of this file.
Classes | |
struct | QProblem |
Implements the online active set strategy for QPs with general constraints. More... | |
Functions | |
returnValue | QProblem_addBound (QProblem *_THIS, int number, SubjectToStatus B_status, BooleanType updateCholesky, BooleanType ensureLI) |
returnValue | QProblem_addBound_checkLI (QProblem *_THIS, int number) |
returnValue | QProblem_addBound_ensureLI (QProblem *_THIS, int number, SubjectToStatus B_status) |
returnValue | QProblem_addConstraint (QProblem *_THIS, int number, SubjectToStatus C_status, BooleanType updateCholesky, BooleanType ensureLI) |
returnValue | QProblem_addConstraint_checkLI (QProblem *_THIS, int number) |
returnValue | QProblem_addConstraint_ensureLI (QProblem *_THIS, int number, SubjectToStatus C_status) |
returnValue | QProblem_areBoundsConsistent (QProblem *_THIS, const real_t *const lb, const real_t *const ub, const real_t *const lbA, const real_t *const ubA) |
returnValue | QProblem_backsolveR (QProblem *_THIS, const real_t *const b, BooleanType transposed, real_t *const a) |
returnValue | QProblem_backsolveRrem (QProblem *_THIS, const real_t *const b, BooleanType transposed, BooleanType removingBound, real_t *const a) |
returnValue | QProblem_backsolveT (QProblem *_THIS, const real_t *const b, BooleanType transposed, real_t *const a) |
returnValue | QProblem_changeActiveSet (QProblem *_THIS, int BC_idx, SubjectToStatus BC_status, BooleanType BC_isBound) |
returnValue | QProblem_computeProjectedCholesky (QProblem *_THIS) |
returnValue | QProblem_determineDataShift (QProblem *_THIS, const real_t *const g_new, const real_t *const lbA_new, const real_t *const ubA_new, const real_t *const lb_new, const real_t *const ub_new, real_t *const delta_g, real_t *const delta_lbA, real_t *const delta_ubA, real_t *const delta_lb, real_t *const delta_ub, BooleanType *Delta_bC_isZero, BooleanType *Delta_bB_isZero) |
returnValue | QProblem_determineHessianType (QProblem *_THIS) |
returnValue | QProblem_determineStepDirection (QProblem *_THIS, const real_t *const delta_g, const real_t *const delta_lbA, const real_t *const delta_ubA, const real_t *const delta_lb, const real_t *const delta_ub, BooleanType Delta_bC_isZero, BooleanType Delta_bB_isZero, real_t *const delta_xFX, real_t *const delta_xFR, real_t *const delta_yAC, real_t *const delta_yFX) |
returnValue | QProblem_dropInfeasibles (QProblem *_THIS, int BC_number, SubjectToStatus BC_status, BooleanType BC_isBound, real_t *xiB, real_t *xiC) |
returnValue | QProblem_ensureNonzeroCurvature (QProblem *_THIS, BooleanType removeBoundNotConstraint, int remIdx, BooleanType *exchangeHappened, BooleanType *addBoundNotConstraint, int *addIdx, SubjectToStatus *addStatus) |
static returnValue | QProblem_getBounds (QProblem *_THIS, Bounds *_bounds) |
static returnValue | QProblem_getConstraints (QProblem *_THIS, Constraints *_constraints) |
static unsigned int | QProblem_getCount (QProblem *_THIS) |
returnValue | QProblem_getDualSolution (QProblem *_THIS, real_t *const yOpt) |
static HessianType | QProblem_getHessianType (QProblem *_THIS) |
static int | QProblem_getNAC (QProblem *_THIS) |
static int | QProblem_getNC (QProblem *_THIS) |
static int | QProblem_getNEC (QProblem *_THIS) |
static int | QProblem_getNFR (QProblem *_THIS) |
static int | QProblem_getNFV (QProblem *_THIS) |
static int | QProblem_getNFX (QProblem *_THIS) |
static int | QProblem_getNIAC (QProblem *_THIS) |
static int | QProblem_getNV (QProblem *_THIS) |
int | QProblem_getNZ (QProblem *_THIS) |
real_t | QProblem_getObjVal (QProblem *_THIS) |
real_t | QProblem_getObjValX (QProblem *_THIS, const real_t *const _x) |
static Options | QProblem_getOptions (QProblem *_THIS) |
returnValue | QProblem_getPrimalSolution (QProblem *_THIS, real_t *const xOpt) |
static PrintLevel | QProblem_getPrintLevel (QProblem *_THIS) |
real_t | QProblem_getRelativeHomotopyLength (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new) |
static QProblemStatus | QProblem_getStatus (QProblem *_THIS) |
returnValue | QProblem_getWorkingSet (QProblem *_THIS, real_t *workingSet) |
returnValue | QProblem_getWorkingSetBounds (QProblem *_THIS, real_t *workingSetB) |
returnValue | QProblem_getWorkingSetConstraints (QProblem *_THIS, real_t *workingSetC) |
returnValue | QProblem_hotstart (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new, int *nWSR, real_t *const cputime) |
returnValue | QProblem_hotstartF (QProblem *_THIS, const char *const g_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file, int *nWSR, real_t *const cputime) |
returnValue | QProblem_hotstartFW (QProblem *_THIS, const char *const g_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file, int *nWSR, real_t *const cputime, Bounds *const guessedBounds, Constraints *const guessedConstraints) |
returnValue | QProblem_hotstartW (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new, int *nWSR, real_t *const cputime, Bounds *const guessedBounds, Constraints *const guessedConstraints) |
returnValue | QProblem_init (QProblem *_THIS, real_t *const _H, const real_t *const _g, real_t *const _A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA, int *nWSR, real_t *const cputime) |
returnValue | QProblem_initF (QProblem *_THIS, const char *const H_file, const char *const g_file, const char *const A_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file, int *nWSR, real_t *const cputime) |
returnValue | QProblem_initFW (QProblem *_THIS, const char *const H_file, const char *const g_file, const char *const A_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file, int *nWSR, real_t *const cputime, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Constraints *const guessedConstraints, const char *const R_file) |
returnValue | QProblem_initM (QProblem *_THIS, DenseMatrix *_H, const real_t *const _g, DenseMatrix *_A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA, int *nWSR, real_t *const cputime) |
returnValue | QProblem_initMW (QProblem *_THIS, DenseMatrix *_H, const real_t *const _g, DenseMatrix *_A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA, int *nWSR, real_t *const cputime, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Constraints *const guessedConstraints, const real_t *const _R) |
returnValue | QProblem_initW (QProblem *_THIS, real_t *const _H, const real_t *const _g, real_t *const _A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA, int *nWSR, real_t *const cputime, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Constraints *const guessedConstraints, const real_t *const _R) |
static BooleanType | QProblem_isBlocking (QProblem *_THIS, real_t num, real_t den, real_t epsNum, real_t epsDen, real_t *t) |
BooleanType | QProblem_isCPUtimeLimitExceeded (QProblem *_THIS, const real_t *const cputime, real_t starttime, int nWSR) |
static BooleanType | QProblem_isInfeasible (QProblem *_THIS) |
static BooleanType | QProblem_isInitialised (QProblem *_THIS) |
static BooleanType | QProblem_isSolved (QProblem *_THIS) |
static BooleanType | QProblem_isUnbounded (QProblem *_THIS) |
returnValue | QProblem_loadQPvectorsFromFile (QProblem *_THIS, const char *const g_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file, real_t *const g_new, real_t *const lb_new, real_t *const ub_new, real_t *const lbA_new, real_t *const ubA_new) |
returnValue | QProblem_obtainAuxiliaryWorkingSet (QProblem *_THIS, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Constraints *const guessedConstraints, Bounds *auxiliaryBounds, Constraints *auxiliaryConstraints) |
returnValue | QProblem_performDriftCorrection (QProblem *_THIS) |
returnValue | QProblem_performPlainRatioTest (QProblem *_THIS, int nIdx, const int *const idxList, const real_t *const num, const real_t *const den, real_t epsNum, real_t epsDen, real_t *t, int *BC_idx) |
returnValue | QProblem_performRamping (QProblem *_THIS) |
returnValue | QProblem_performRatioTestB (QProblem *_THIS, int nIdx, const int *const idxList, Bounds *const subjectTo, const real_t *const num, const real_t *const den, real_t epsNum, real_t epsDen, real_t *t, int *BC_idx) |
returnValue | QProblem_performRatioTestC (QProblem *_THIS, int nIdx, const int *const idxList, Constraints *const subjectTo, const real_t *const num, const real_t *const den, real_t epsNum, real_t epsDen, real_t *t, int *BC_idx) |
returnValue | QProblem_performStep (QProblem *_THIS, const real_t *const delta_g, const real_t *const delta_lbA, const real_t *const delta_ubA, const real_t *const delta_lb, const real_t *const delta_ub, const real_t *const delta_xFX, const real_t *const delta_xFR, const real_t *const delta_yAC, const real_t *const delta_yFX, int *BC_idx, SubjectToStatus *BC_status, BooleanType *BC_isBound) |
returnValue | QProblem_printIteration (QProblem *_THIS, int iter, int BC_idx, SubjectToStatus BC_status, BooleanType BC_isBound, real_t homotopyLength, BooleanType isFirstCall) |
returnValue | QProblem_printOptions (QProblem *_THIS) |
returnValue | QProblem_printProperties (QProblem *_THIS) |
returnValue | QProblem_regulariseHessian (QProblem *_THIS) |
returnValue | QProblem_removeBound (QProblem *_THIS, int number, BooleanType updateCholesky, BooleanType allowFlipping, BooleanType ensureNZC) |
returnValue | QProblem_removeConstraint (QProblem *_THIS, int number, BooleanType updateCholesky, BooleanType allowFlipping, BooleanType ensureNZC) |
returnValue | QProblem_reset (QProblem *_THIS) |
static returnValue | QProblem_resetCounter (QProblem *_THIS) |
static returnValue | QProblem_setA (QProblem *_THIS, real_t *const A_new) |
static returnValue | QProblem_setAM (QProblem *_THIS, DenseMatrix *A_new) |
returnValue | QProblem_setConstraintProduct (QProblem *_THIS, ConstraintProduct _constraintProduct) |
static returnValue | QProblem_setG (QProblem *_THIS, const real_t *const g_new) |
static returnValue | QProblem_setH (QProblem *_THIS, real_t *const H_new) |
static returnValue | QProblem_setHessianType (QProblem *_THIS, HessianType _hessianType) |
static returnValue | QProblem_setHM (QProblem *_THIS, DenseMatrix *H_new) |
returnValue | QProblem_setInfeasibilityFlag (QProblem *_THIS, returnValue returnvalue, BooleanType doThrowError) |
static returnValue | QProblem_setLB (QProblem *_THIS, const real_t *const lb_new) |
static returnValue | QProblem_setLBA (QProblem *_THIS, const real_t *const lbA_new) |
static returnValue | QProblem_setLBAn (QProblem *_THIS, int number, real_t value) |
static returnValue | QProblem_setLBn (QProblem *_THIS, int number, real_t value) |
static returnValue | QProblem_setOptions (QProblem *_THIS, Options _options) |
returnValue | QProblem_setPrintLevel (QProblem *_THIS, PrintLevel _printlevel) |
static returnValue | QProblem_setUB (QProblem *_THIS, const real_t *const ub_new) |
static returnValue | QProblem_setUBA (QProblem *_THIS, const real_t *const ubA_new) |
static returnValue | QProblem_setUBAn (QProblem *_THIS, int number, real_t value) |
static returnValue | QProblem_setUBn (QProblem *_THIS, int number, real_t value) |
returnValue | QProblem_setupAuxiliaryQP (QProblem *_THIS, Bounds *const guessedBounds, Constraints *const guessedConstraints) |
returnValue | QProblem_setupAuxiliaryQPbounds (QProblem *_THIS, Bounds *const auxiliaryBounds, Constraints *const auxiliaryConstraints, BooleanType useRelaxation) |
returnValue | QProblem_setupAuxiliaryQPgradient (QProblem *_THIS) |
returnValue | QProblem_setupAuxiliaryQPsolution (QProblem *_THIS, const real_t *const xOpt, const real_t *const yOpt) |
returnValue | QProblem_setupAuxiliaryWorkingSet (QProblem *_THIS, Bounds *const auxiliaryBounds, Constraints *const auxiliaryConstraints, BooleanType setupAfresh) |
returnValue | QProblem_setupInitialCholesky (QProblem *_THIS) |
returnValue | QProblem_setupQPdata (QProblem *_THIS, real_t *const _H, const real_t *const _g, real_t *const _A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA) |
returnValue | QProblem_setupQPdataFromFile (QProblem *_THIS, const char *const H_file, const char *const g_file, const char *const A_file, const char *const lb_file, const char *const ub_file, const char *const lbA_file, const char *const ubA_file) |
returnValue | QProblem_setupQPdataM (QProblem *_THIS, DenseMatrix *_H, const real_t *const _g, DenseMatrix *_A, const real_t *const _lb, const real_t *const _ub, const real_t *const _lbA, const real_t *const _ubA) |
returnValue | QProblem_setupSubjectToType (QProblem *_THIS) |
returnValue | QProblem_setupSubjectToTypeNew (QProblem *_THIS, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new) |
returnValue | QProblem_setupTQfactorisation (QProblem *_THIS) |
BooleanType | QProblem_shallRefactorise (QProblem *_THIS, Bounds *const guessedBounds, Constraints *const guessedConstraints) |
returnValue | QProblem_solveCurrentEQP (QProblem *_THIS, const int n_rhs, const real_t *g_in, const real_t *lb_in, const real_t *ub_in, const real_t *lbA_in, const real_t *ubA_in, real_t *x_out, real_t *y_out) |
returnValue | QProblem_solveInitialQP (QProblem *_THIS, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Constraints *const guessedConstraints, const real_t *const _R, int *nWSR, real_t *const cputime) |
returnValue | QProblem_solveQP (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new, int *nWSR, real_t *const cputime, int nWSRperformed, BooleanType isFirstCall) |
returnValue | QProblem_solveRegularisedQP (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, const real_t *const lbA_new, const real_t *const ubA_new, int *nWSR, real_t *const cputime, int nWSRperformed, BooleanType isFirstCall) |
returnValue | QProblem_updateFarBounds (QProblem *_THIS, real_t curFarBound, int nRamp, const real_t *const lb_new, real_t *const lb_new_far, const real_t *const ub_new, real_t *const ub_new_far, const real_t *const lbA_new, real_t *const lbA_new_far, const real_t *const ubA_new, real_t *const ubA_new_far) |
static BooleanType | QProblem_usingRegularisation (QProblem *_THIS) |
returnValue | QProblem_writeQpDataIntoMatFile (QProblem *_THIS, const char *const filename) |
returnValue | QProblem_writeQpWorkspaceIntoMatFile (QProblem *_THIS, const char *const filename) |
returnValue | QProblemBCPY_computeCholesky (QProblem *_THIS) |
returnValue | QProblemBCPY_determineDataShift (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new, real_t *const delta_g, real_t *const delta_lb, real_t *const delta_ub, BooleanType *Delta_bB_isZero) |
real_t | QProblemBCPY_getRelativeHomotopyLength (QProblem *_THIS, const real_t *const g_new, const real_t *const lb_new, const real_t *const ub_new) |
returnValue | QProblemBCPY_loadQPvectorsFromFile (QProblem *_THIS, const char *const g_file, const char *const lb_file, const char *const ub_file, real_t *const g_new, real_t *const lb_new, real_t *const ub_new) |
returnValue | QProblemBCPY_obtainAuxiliaryWorkingSet (QProblem *_THIS, const real_t *const xOpt, const real_t *const yOpt, Bounds *const guessedBounds, Bounds *auxiliaryBounds) |
returnValue | QProblemBCPY_setupQPdata (QProblem *_THIS, real_t *const _H, const real_t *const _g, const real_t *const _lb, const real_t *const _ub) |
returnValue | QProblemBCPY_setupQPdataFromFile (QProblem *_THIS, const char *const H_file, const char *const g_file, const char *const lb_file, const char *const ub_file) |
returnValue | QProblemBCPY_setupQPdataM (QProblem *_THIS, DenseMatrix *_H, const real_t *const _g, const real_t *const _lb, const real_t *const _ub) |
returnValue | QProblemBCPY_updateFarBounds (QProblem *_THIS, real_t curFarBound, int nRamp, const real_t *const lb_new, real_t *const lb_new_far, const real_t *const ub_new, real_t *const ub_new_far) |
void | QProblemCON (QProblem *_THIS, int _nV, int _nC, HessianType _hessianType) |
void | QProblemCPY (QProblem *FROM, QProblem *TO) |
Declaration of the QProblem class which is able to use the newly developed online active set strategy for parametric quadratic programming.
Definition in file QProblem.h.
returnValue QProblem_addBound | ( | QProblem * | _THIS, |
int | number, | ||
SubjectToStatus | B_status, | ||
BooleanType | updateCholesky, | ||
BooleanType | ensureLI | ||
) |
Adds a bound to active set.
number | Number of bound to be added to active set. |
B_status | Status of new active bound. |
updateCholesky | Flag indicating if Cholesky decomposition shall be updated. |
ensureLI | Ensure linear independence by exchange rules by default. |
Definition at line 4107 of file QProblem.c.
returnValue QProblem_addBound_checkLI | ( | QProblem * | _THIS, |
int | number | ||
) |
Checks if new active bound to be added is linearly dependent from from row of the active constraints matrix.
number | Number of bound to be added to active set. |
Definition at line 4283 of file QProblem.c.
returnValue QProblem_addBound_ensureLI | ( | QProblem * | _THIS, |
int | number, | ||
SubjectToStatus | B_status | ||
) |
Ensures linear independence of constraint matrix when a new bound is added. To _THIS end a bound or constraint is removed simultaneously if necessary.
number | Number of bound to be added to active set. |
B_status | Status of new active bound. |
Definition at line 4388 of file QProblem.c.
returnValue QProblem_addConstraint | ( | QProblem * | _THIS, |
int | number, | ||
SubjectToStatus | C_status, | ||
BooleanType | updateCholesky, | ||
BooleanType | ensureLI | ||
) |
Adds a constraint to active set.
number | Number of constraint to be added to active set. |
C_status | Status of new active constraint. |
updateCholesky | Flag indicating if Cholesky decomposition shall be updated. |
ensureLI | Ensure linear independence by exchange rules by default. |
Definition at line 3639 of file QProblem.c.
returnValue QProblem_addConstraint_checkLI | ( | QProblem * | _THIS, |
int | number | ||
) |
Checks if new active constraint to be added is linearly dependent from from row of the active constraints matrix.
number | Number of constraint to be added to active set. |
Definition at line 3803 of file QProblem.c.
returnValue QProblem_addConstraint_ensureLI | ( | QProblem * | _THIS, |
int | number, | ||
SubjectToStatus | C_status | ||
) |
Ensures linear independence of constraint matrix when a new constraint is added. To _THIS end a bound or constraint is removed simultaneously if necessary.
number | Number of constraint to be added to active set. |
C_status | Status of new active bound. |
Definition at line 3925 of file QProblem.c.
returnValue QProblem_areBoundsConsistent | ( | QProblem * | _THIS, |
const real_t *const | lb, | ||
const real_t *const | ub, | ||
const real_t *const | lbA, | ||
const real_t *const | ubA | ||
) |
Decides if lower bounds are smaller than upper bounds
lb | Vector of lower bounds |
ub | Vector of upper bounds |
lbA | Vector of lower constraints |
ubA | Vector of upper constraints |
Definition at line 7090 of file QProblem.c.
returnValue QProblem_backsolveR | ( | QProblem * | _THIS, |
const real_t *const | b, | ||
BooleanType | transposed, | ||
real_t *const | a | ||
) |
Solves the system Ra = b or R^Ta = b where R is an upper triangular matrix.
b | Right hand side vector. |
transposed | Indicates if the transposed system shall be solved. |
a | Output: Solution vector |
Definition at line 1748 of file QProblem.c.
returnValue QProblem_backsolveRrem | ( | QProblem * | _THIS, |
const real_t *const | b, | ||
BooleanType | transposed, | ||
BooleanType | removingBound, | ||
real_t *const | a | ||
) |
Solves the system Ra = b or R^Ta = b where R is an upper triangular matrix.
Special variant for the case that _THIS function is called from within "removeBound()".
b | Right hand side vector. |
transposed | Indicates if the transposed system shall be solved. |
removingBound | Indicates if function is called from "removeBound()". |
a | Output: Solution vector |
Definition at line 1760 of file QProblem.c.
returnValue QProblem_backsolveT | ( | QProblem * | _THIS, |
const real_t *const | b, | ||
BooleanType | transposed, | ||
real_t *const | a | ||
) |
Solves the system Ta = b or T^Ta = b where T is a reverse upper triangular matrix.
b | Right hand side vector. |
transposed | Indicates if the transposed system shall be solved. |
a | Output: Solution vector |
Definition at line 5280 of file QProblem.c.
returnValue QProblem_changeActiveSet | ( | QProblem * | _THIS, |
int | BC_idx, | ||
SubjectToStatus | BC_status, | ||
BooleanType | BC_isBound | ||
) |
Updates the active set.
BC_idx | Index of blocking constraint. |
BC_status | Status of blocking constraint. |
BC_isBound | Indicates if blocking constraint is a bound. |
Definition at line 6055 of file QProblem.c.
returnValue QProblem_computeProjectedCholesky | ( | QProblem * | _THIS | ) |
Computes the Cholesky decomposition of the projected Hessian (i.e. R^T*R = Z^T*H*Z). Note: If Hessian turns out not to be positive definite, the Hessian type is set to HST_SEMIDEF accordingly.
Definition at line 2917 of file QProblem.c.
returnValue QProblem_determineDataShift | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
real_t *const | delta_g, | ||
real_t *const | delta_lbA, | ||
real_t *const | delta_ubA, | ||
real_t *const | delta_lb, | ||
real_t *const | delta_ub, | ||
BooleanType * | Delta_bC_isZero, | ||
BooleanType * | Delta_bB_isZero | ||
) |
Determines step direction of the shift of the QP data.
g_new | New gradient vector. |
lbA_new | New lower constraints' bounds. |
ubA_new | New upper constraints' bounds. |
lb_new | New lower bounds. |
ub_new | New upper bounds. |
delta_g | Output: Step direction of gradient vector. |
delta_lbA | Output: Step direction of lower constraints' bounds. |
delta_ubA | Output: Step direction of upper constraints' bounds. |
delta_lb | Output: Step direction of lower bounds. |
delta_ub | Output: Step direction of upper bounds. |
Delta_bC_isZero | Output: Indicates if active constraints' bounds are to be shifted. |
Delta_bB_isZero | Output: Indicates if active bounds are to be shifted. |
Definition at line 5333 of file QProblem.c.
returnValue QProblem_determineHessianType | ( | QProblem * | _THIS | ) |
If Hessian type has been set by the user, nothing is done. Otherwise the Hessian type is set to HST_IDENTITY, HST_ZERO, or HST_POSDEF (default), respectively.
Definition at line 1421 of file QProblem.c.
returnValue QProblem_determineStepDirection | ( | QProblem * | _THIS, |
const real_t *const | delta_g, | ||
const real_t *const | delta_lbA, | ||
const real_t *const | delta_ubA, | ||
const real_t *const | delta_lb, | ||
const real_t *const | delta_ub, | ||
BooleanType | Delta_bC_isZero, | ||
BooleanType | Delta_bB_isZero, | ||
real_t *const | delta_xFX, | ||
real_t *const | delta_xFR, | ||
real_t *const | delta_yAC, | ||
real_t *const | delta_yFX | ||
) |
Determines step direction of the homotopy path.
delta_g | Step direction of gradient vector. |
delta_lbA | Step direction of lower constraints' bounds. |
delta_ubA | Step direction of upper constraints' bounds. |
delta_lb | Step direction of lower bounds. |
delta_ub | Step direction of upper bounds. |
Delta_bC_isZero | Indicates if active constraints' bounds are to be shifted. |
Delta_bB_isZero | Indicates if active bounds are to be shifted. |
delta_xFX | Output: Primal homotopy step direction of fixed variables. |
delta_xFR | Output: Primal homotopy step direction of free variables. |
delta_yAC | Output: Dual homotopy step direction of active constraints' multiplier. |
delta_yFX | Output: Dual homotopy step direction of fixed variables' multiplier. |
Definition at line 5397 of file QProblem.c.
returnValue QProblem_dropInfeasibles | ( | QProblem * | _THIS, |
int | BC_number, | ||
SubjectToStatus | BC_status, | ||
BooleanType | BC_isBound, | ||
real_t * | xiB, | ||
real_t * | xiC | ||
) |
Drops the blocking bound/constraint that led to infeasibility, or finds another bound/constraint to drop according to drop priorities.
BC_number | Number of the bound or constraint to be added |
BC_status | New status of the bound or constraint to be added |
BC_isBound | Whether a bound or a constraint is to be added |
Definition at line 7121 of file QProblem.c.
returnValue QProblem_ensureNonzeroCurvature | ( | QProblem * | _THIS, |
BooleanType | removeBoundNotConstraint, | ||
int | remIdx, | ||
BooleanType * | exchangeHappened, | ||
BooleanType * | addBoundNotConstraint, | ||
int * | addIdx, | ||
SubjectToStatus * | addStatus | ||
) |
Ensure non-zero curvature by primal jump.
removeBoundNotConstraint | SubjectTo to be removed is a bound. |
remIdx | Index of bound/constraint to be removed. |
exchangeHappened | Output: Exchange was necessary to ensure. |
addBoundNotConstraint | SubjectTo to be added is a bound. |
addIdx | Index of bound/constraint to be added. |
addStatus | Status of bound/constraint to be added. |
Definition at line 5050 of file QProblem.c.
|
inlinestatic |
Returns current bounds object of the QP (deep copy).
Definition at line 1675 of file QProblem.h.
|
inlinestatic |
Returns current constraints object of the QP (deep copy).
_constraints | Output: Constraints object. |
Definition at line 2044 of file QProblem.h.
|
inlinestatic |
Returns the current number of QP problems solved.
Definition at line 1843 of file QProblem.h.
returnValue QProblem_getDualSolution | ( | QProblem * | _THIS, |
real_t *const | yOpt | ||
) |
Returns the dual solution vector (deep copy).
yOpt | Output: Dual solution vector (if QP has been solved). |
Definition at line 1058 of file QProblem.c.
|
inlinestatic |
Returns Hessian type flag (type is not determined due to _THIS call!).
Definition at line 1778 of file QProblem.h.
|
inlinestatic |
Returns the number of active constraints.
Definition at line 2079 of file QProblem.h.
|
inlinestatic |
Returns the number of constraints.
Definition at line 2061 of file QProblem.h.
|
inlinestatic |
Returns the number of (implicitly defined) equality constraints.
Definition at line 2070 of file QProblem.h.
|
inlinestatic |
Returns the number of free variables.
Definition at line 1700 of file QProblem.h.
|
inlinestatic |
Returns the number of implicitly fixed variables.
Definition at line 1718 of file QProblem.h.
|
inlinestatic |
Returns the number of fixed variables.
Definition at line 1709 of file QProblem.h.
|
inlinestatic |
Returns the number of inactive constraints.
Definition at line 2088 of file QProblem.h.
|
inlinestatic |
Returns the number of variables.
Definition at line 1691 of file QProblem.h.
int QProblem_getNZ | ( | QProblem * | _THIS | ) |
Returns the dimension of null space.
Definition at line 1048 of file QProblem.c.
Returns the optimal objective function value.
Definition at line 1095 of file QProblem.c.
Returns the objective function value at an arbitrary point x.
_x | Point at which the objective function shall be evaluated. |
Definition at line 1119 of file QProblem.c.
Returns current options struct.
Definition at line 1809 of file QProblem.h.
returnValue QProblem_getPrimalSolution | ( | QProblem * | _THIS, |
real_t *const | xOpt | ||
) |
Returns the primal solution vector.
xOpt | Output: Primal solution vector (if QP has been solved). |
Definition at line 1168 of file QProblem.c.
|
inlinestatic |
real_t QProblem_getRelativeHomotopyLength | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new | ||
) |
Compute relative length of homotopy in data space for termination criterion.
g_new | Final gradient. |
lb_new | Final lower variable bounds. |
ub_new | Final upper variable bounds. |
lbA_new | Final lower constraint bounds. |
ubA_new | Final upper constraint bounds. |
Definition at line 6142 of file QProblem.c.
|
inlinestatic |
Returns status of the solution process.
Definition at line 1727 of file QProblem.h.
returnValue QProblem_getWorkingSet | ( | QProblem * | _THIS, |
real_t * | workingSet | ||
) |
Writes a vector with the state of the working set
workingSet | Output: array containing state of the working set. |
Definition at line 981 of file QProblem.c.
returnValue QProblem_getWorkingSetBounds | ( | QProblem * | _THIS, |
real_t * | workingSetB | ||
) |
Writes a vector with the state of the working set of bounds
workingSetB | Output: array containing state of the working set of bounds. |
Definition at line 1001 of file QProblem.c.
returnValue QProblem_getWorkingSetConstraints | ( | QProblem * | _THIS, |
real_t * | workingSetC | ||
) |
Writes a vector with the state of the working set of constraints
workingSetC | Output: array containing state of the working set of constraints. |
Definition at line 1022 of file QProblem.c.
returnValue QProblem_hotstart | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Solves an initialised QP sequence using the online active set strategy. QP solution is started from previous solution.
Note: This function internally calls solveQP/solveRegularisedQP for solving an initialised QP!
g_new | Gradient of neighbouring QP to be solved. |
lb_new | Lower bounds of neighbouring QP to be solved. If no lower bounds exist, a NULL pointer can be passed. |
ub_new | Upper bounds of neighbouring QP to be solved. If no upper bounds exist, a NULL pointer can be passed. |
lbA_new | Lower constraints' bounds of neighbouring QP to be solved. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_new | Upper constraints' bounds of neighbouring QP to be solved. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
Definition at line 578 of file QProblem.c.
returnValue QProblem_hotstartF | ( | QProblem * | _THIS, |
const char *const | g_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Solves an initialised QP sequence using the online active set strategy, where QP data is read from files. QP solution is started from previous solution.
Note: This function internally calls solveQP/solveRegularisedQP for solving an initialised QP!
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
Definition at line 759 of file QProblem.c.
returnValue QProblem_hotstartFW | ( | QProblem * | _THIS, |
const char *const | g_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints | ||
) |
Solves an initialised QP sequence using the online active set strategy, where QP data is read from files. By default, QP solution is started from previous solution. If a guess for the working set is provided, an initialised homotopy is performed.
Note: This function internally calls solveQP/solveRegularisedQP for solving an initialised QP!
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). (If a null pointer is passed, the previous working set of bounds is kept!) |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). (If a null pointer is passed, the previous working set of constraints is kept!) |
Definition at line 864 of file QProblem.c.
returnValue QProblem_hotstartW | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints | ||
) |
Solves an initialised QP sequence using the online active set strategy. By default, QP solution is started from previous solution. If a guess for the working set is provided, an initialised homotopy is performed.
Note: This function internally calls solveQP/solveRegularisedQP for solving an initialised QP!
g_new | Gradient of neighbouring QP to be solved. |
lb_new | Lower bounds of neighbouring QP to be solved. If no lower bounds exist, a NULL pointer can be passed. |
ub_new | Upper bounds of neighbouring QP to be solved. If no upper bounds exist, a NULL pointer can be passed. |
lbA_new | Lower constraints' bounds of neighbouring QP to be solved. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_new | Upper constraints' bounds of neighbouring QP to be solved. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). (If a null pointer is passed, the previous working set of bounds is kept!) |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). (If a null pointer is passed, the previous working set of constraints is kept!) |
Definition at line 801 of file QProblem.c.
returnValue QProblem_init | ( | QProblem * | _THIS, |
real_t *const | _H, | ||
const real_t *const | _g, | ||
real_t *const | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Initialises a QP problem with given QP data and tries to solve it using at most nWSR iterations.
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_INVALID_ARGUMENTS
_H | Hessian matrix. If Hessian matrix is trivial, a NULL pointer can be passed. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation (if pointer passed). |
Definition at line 299 of file QProblem.c.
returnValue QProblem_initF | ( | QProblem * | _THIS, |
const char *const | H_file, | ||
const char *const | g_file, | ||
const char *const | A_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Initialises a QP problem with given QP data to be read from files and tries to solve it using at most nWSR iterations.
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_UNABLE_TO_READ_FILE
H_file | Name of file where Hessian matrix is stored. If Hessian matrix is trivial, a NULL pointer can be passed. |
g_file | Name of file where gradient vector is stored. |
A_file | Name of file where constraint matrix is stored. |
lb_file | Name of file where lower bound vector. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bound vector. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bound vector. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation (if pointer passed). |
Definition at line 327 of file QProblem.c.
returnValue QProblem_initFW | ( | QProblem * | _THIS, |
const char *const | H_file, | ||
const char *const | g_file, | ||
const char *const | A_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints, | ||
const char *const | R_file | ||
) |
Initialises a QP problem with given QP data to be ream from files and tries to solve it using at most nWSR iterations. Depending on the parameter constellation it:
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_UNABLE_TO_READ_FILE
H_file | Name of file where Hessian matrix is stored. If Hessian matrix is trivial, a NULL pointer can be passed. |
g_file | Name of file where gradient vector is stored. |
A_file | Name of file where constraint matrix is stored. |
lb_file | Name of file where lower bound vector. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bound vector. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bound vector. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation. |
xOpt | Optimal primal solution vector. (If a null pointer is passed, the old primal solution is kept!) |
yOpt | Optimal dual solution vector. (If a null pointer is passed, the old dual solution is kept!) |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). |
R_file | Pre-computed (upper triangular) Cholesky factor of Hessian matrix. The Cholesky factor must be stored in a real_t array of size nV*nV in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! (If a null pointer is passed, Cholesky decomposition is computed internally!) |
Definition at line 472 of file QProblem.c.
returnValue QProblem_initM | ( | QProblem * | _THIS, |
DenseMatrix * | _H, | ||
const real_t *const | _g, | ||
DenseMatrix * | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Initialises a QP problem with given QP data and tries to solve it using at most nWSR iterations.
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_INVALID_ARGUMENTS
_H | Hessian matrix. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation (if pointer passed). |
Definition at line 256 of file QProblem.c.
returnValue QProblem_initMW | ( | QProblem * | _THIS, |
DenseMatrix * | _H, | ||
const real_t *const | _g, | ||
DenseMatrix * | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints, | ||
const real_t *const | _R | ||
) |
Initialises a QP problem with given QP data and tries to solve it using at most nWSR iterations. Depending on the parameter constellation it:
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_INVALID_ARGUMENTS
_H | Hessian matrix. If Hessian matrix is trivial, a NULL pointer can be passed. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation. |
xOpt | Optimal primal solution vector. (If a null pointer is passed, the old primal solution is kept!) |
yOpt | Optimal dual solution vector. (If a null pointer is passed, the old dual solution is kept!) |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). |
_R | Pre-computed (upper triangular) Cholesky factor of Hessian matrix. The Cholesky factor must be stored in a real_t array of size nV*nV in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! (If a null pointer is passed, Cholesky decomposition is computed internally!) |
Definition at line 355 of file QProblem.c.
returnValue QProblem_initW | ( | QProblem * | _THIS, |
real_t *const | _H, | ||
const real_t *const | _g, | ||
real_t *const | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints, | ||
const real_t *const | _R | ||
) |
Initialises a QP problem with given QP data and tries to solve it using at most nWSR iterations. Depending on the parameter constellation it:
Note: This function internally calls solveInitialQP for initialisation!
\return SUCCESSFUL_RETURN \n RET_INIT_FAILED \n RET_INIT_FAILED_CHOLESKY \n RET_INIT_FAILED_TQ \n RET_INIT_FAILED_HOTSTART \n RET_INIT_FAILED_INFEASIBILITY \n RET_INIT_FAILED_UNBOUNDEDNESS \n RET_MAX_NWSR_REACHED \n RET_INVALID_ARGUMENTS
_H | Hessian matrix. If Hessian matrix is trivial, a NULL pointer can be passed. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations when using initial homotopy. Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP initialisation. Output: CPU time spent for QP initialisation. |
xOpt | Optimal primal solution vector. (If a null pointer is passed, the old primal solution is kept!) |
yOpt | Optimal dual solution vector. (If a null pointer is passed, the old dual solution is kept!) |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). |
_R | Pre-computed (upper triangular) Cholesky factor of Hessian matrix. The Cholesky factor must be stored in a real_t array of size nV*nV in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! (If a null pointer is passed, Cholesky decomposition is computed internally!) |
Definition at line 414 of file QProblem.c.
|
inlinestatic |
Checks whether given ratio is blocking, i.e. limits the maximum step length along the homotopy path to a value lower than given one.
num | Numerator for performing the ratio test. |
den | Denominator for performing the ratio test. |
epsNum | Numerator tolerance. |
epsDen | Denominator tolerance. |
t | Input: Current maximum step length along the homotopy path, Output: Updated maximum possible step length along the homotopy path. |
Definition at line 2022 of file QProblem.h.
BooleanType QProblem_isCPUtimeLimitExceeded | ( | QProblem * | _THIS, |
const real_t *const | cputime, | ||
real_t | starttime, | ||
int | nWSR | ||
) |
Determines if next QP iteration can be performed within given CPU time limit.
cputime | Maximum CPU time allowed for QP solution. |
starttime | Start time of current QP solution. |
nWSR | Number of working set recalculations performed so far. |
Definition at line 2064 of file QProblem.c.
|
inlinestatic |
Returns if the QP is infeasible.
Definition at line 1760 of file QProblem.h.
|
inlinestatic |
Returns if the QProblem object is initialised.
Definition at line 1736 of file QProblem.h.
|
inlinestatic |
Returns if the QP has been solved.
Definition at line 1748 of file QProblem.h.
|
inlinestatic |
Returns if the QP is unbounded.
Definition at line 1769 of file QProblem.h.
returnValue QProblem_loadQPvectorsFromFile | ( | QProblem * | _THIS, |
const char *const | g_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file, | ||
real_t *const | g_new, | ||
real_t *const | lb_new, | ||
real_t *const | ub_new, | ||
real_t *const | lbA_new, | ||
real_t *const | ubA_new | ||
) |
Loads new QP vectors from files (internal members are not affected!).
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. If no upper constraints' bounds exist, a NULL pointer can be passed. |
g_new | Output: Gradient of neighbouring QP to be solved. |
lb_new | Output: Lower bounds of neighbouring QP to be solved |
ub_new | Output: Upper bounds of neighbouring QP to be solved |
lbA_new | Output: Lower constraints' bounds of neighbouring QP to be solved |
ubA_new | Output: Upper constraints' bounds of neighbouring QP to be solved |
Definition at line 6785 of file QProblem.c.
returnValue QProblem_obtainAuxiliaryWorkingSet | ( | QProblem * | _THIS, |
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints, | ||
Bounds * | auxiliaryBounds, | ||
Constraints * | auxiliaryConstraints | ||
) |
Obtains the desired working set for the auxiliary initial QP in accordance with the user specifications (assumes that member AX has already been initialised!)
xOpt | Optimal primal solution vector. If a NULL pointer is passed, all entries are assumed to be zero. |
yOpt | Optimal dual solution vector. If a NULL pointer is passed, all entries are assumed to be zero. |
guessedBounds | Guessed working set of bounds for solution (xOpt,yOpt). |
guessedConstraints | Guessed working set for solution (xOpt,yOpt). |
auxiliaryBounds | Input: Allocated bound object. Ouput: Working set of constraints for auxiliary QP. |
auxiliaryConstraints | Input: Allocated bound object. Ouput: Working set for auxiliary QP. |
Definition at line 3040 of file QProblem.c.
returnValue QProblem_performDriftCorrection | ( | QProblem * | _THIS | ) |
Drift correction at end of each active set iteration
Definition at line 6429 of file QProblem.c.
returnValue QProblem_performPlainRatioTest | ( | QProblem * | _THIS, |
int | nIdx, | ||
const int *const | idxList, | ||
const real_t *const | num, | ||
const real_t *const | den, | ||
real_t | epsNum, | ||
real_t | epsDen, | ||
real_t * | t, | ||
int * | BC_idx | ||
) |
Performs robustified ratio test yield the maximum possible step length along the homotopy path.
nIdx | Number of ratios to be checked. |
idxList | Array containing the indices of all ratios to be checked. |
num | Array containing all numerators for performing the ratio test. |
den | Array containing all denominators for performing the ratio test. |
epsNum | Numerator tolerance. |
epsDen | Denominator tolerance. |
t | Output: Maximum possible step length along the homotopy path. |
BC_idx | Output: Index of blocking constraint. |
Definition at line 5027 of file QProblem.c.
returnValue QProblem_performRamping | ( | QProblem * | _THIS | ) |
Ramping Strategy to avoid ties. Modifies homotopy start without changing current active set.
Definition at line 6296 of file QProblem.c.
returnValue QProblem_performRatioTestB | ( | QProblem * | _THIS, |
int | nIdx, | ||
const int *const | idxList, | ||
Bounds *const | subjectTo, | ||
const real_t *const | num, | ||
const real_t *const | den, | ||
real_t | epsNum, | ||
real_t | epsDen, | ||
real_t * | t, | ||
int * | BC_idx | ||
) |
Performs robustified ratio test yield the maximum possible step length along the homotopy path.
nIdx | Number of ratios to be checked. |
idxList | Array containing the indices of all ratios to be checked. |
subjectTo | Bound object corresponding to ratios to be checked. |
num | Array containing all numerators for performing the ratio test. |
den | Array containing all denominators for performing the ratio test. |
epsNum | Numerator tolerance. |
epsDen | Denominator tolerance. |
t | Output: Maximum possible step length along the homotopy path. |
BC_idx | Output: Index of blocking constraint. |
Definition at line 2133 of file QProblem.c.
returnValue QProblem_performRatioTestC | ( | QProblem * | _THIS, |
int | nIdx, | ||
const int *const | idxList, | ||
Constraints *const | subjectTo, | ||
const real_t *const | num, | ||
const real_t *const | den, | ||
real_t | epsNum, | ||
real_t | epsDen, | ||
real_t * | t, | ||
int * | BC_idx | ||
) |
Performs robustified ratio test yield the maximum possible step length along the homotopy path.
nIdx | Number of ratios to be checked. |
idxList | Array containing the indices of all ratios to be checked. |
subjectTo | Constraint object corresponding to ratios to be checked. |
num | Array containing all numerators for performing the ratio test. |
den | Array containing all denominators for performing the ratio test. |
epsNum | Numerator tolerance. |
epsDen | Denominator tolerance. |
t | Output: Maximum possible step length along the homotopy path. |
BC_idx | Output: Index of blocking constraint. |
Definition at line 6379 of file QProblem.c.
returnValue QProblem_performStep | ( | QProblem * | _THIS, |
const real_t *const | delta_g, | ||
const real_t *const | delta_lbA, | ||
const real_t *const | delta_ubA, | ||
const real_t *const | delta_lb, | ||
const real_t *const | delta_ub, | ||
const real_t *const | delta_xFX, | ||
const real_t *const | delta_xFR, | ||
const real_t *const | delta_yAC, | ||
const real_t *const | delta_yFX, | ||
int * | BC_idx, | ||
SubjectToStatus * | BC_status, | ||
BooleanType * | BC_isBound | ||
) |
Determines the maximum possible step length along the homotopy path and performs _THIS step (without changing working set).
delta_g | Step direction of gradient. |
delta_lbA | Step direction of lower constraints' bounds. |
delta_ubA | Step direction of upper constraints' bounds. |
delta_lb | Step direction of lower bounds. |
delta_ub | Step direction of upper bounds. |
delta_xFX | Primal homotopy step direction of fixed variables. |
delta_xFR | Primal homotopy step direction of free variables. |
delta_yAC | Dual homotopy step direction of active constraints' multiplier. |
delta_yFX | Dual homotopy step direction of fixed variables' multiplier. |
BC_idx | Output: Index of blocking constraint. |
BC_status | Output: Status of blocking constraint. |
BC_isBound | Output: Indicates if blocking constraint is a bound. |
Definition at line 5741 of file QProblem.c.
returnValue QProblem_printIteration | ( | QProblem * | _THIS, |
int | iter, | ||
int | BC_idx, | ||
SubjectToStatus | BC_status, | ||
BooleanType | BC_isBound, | ||
real_t | homotopyLength, | ||
BooleanType | isFirstCall | ||
) |
Prints concise information on the current iteration.
iter | Number of current iteration. |
BC_idx | Index of blocking constraint. |
BC_status | Status of blocking constraint. |
BC_isBound | Indicates if blocking constraint is a bound. |
homotopyLength | Current homotopy distance. |
isFirstCall | Indicating whether this is the first call for current QP. |
Definition at line 6844 of file QProblem.c.
returnValue QProblem_printOptions | ( | QProblem * | _THIS | ) |
Prints a list of all options and their current values.
Definition at line 1240 of file QProblem.c.
returnValue QProblem_printProperties | ( | QProblem * | _THIS | ) |
Prints concise list of properties of the current QP.
Definition at line 1249 of file QProblem.c.
returnValue QProblem_regulariseHessian | ( | QProblem * | _THIS | ) |
Regularise Hessian matrix by adding a scaled identity matrix to it.
Definition at line 2094 of file QProblem.c.
returnValue QProblem_removeBound | ( | QProblem * | _THIS, |
int | number, | ||
BooleanType | updateCholesky, | ||
BooleanType | allowFlipping, | ||
BooleanType | ensureNZC | ||
) |
Removes a bounds from active set.
number | Number of bound to be removed from active set. |
updateCholesky | Flag indicating if Cholesky decomposition shall be updated. |
allowFlipping | Flag indicating if flipping bounds are allowed. |
ensureNZC | Flag indicating if non-zero curvature is ensured by exchange rules. |
Definition at line 4796 of file QProblem.c.
returnValue QProblem_removeConstraint | ( | QProblem * | _THIS, |
int | number, | ||
BooleanType | updateCholesky, | ||
BooleanType | allowFlipping, | ||
BooleanType | ensureNZC | ||
) |
Removes a constraint from active set.
number | Number of constraint to be removed from active set. |
updateCholesky | Flag indicating if Cholesky decomposition shall be updated. |
allowFlipping | Flag indicating if flipping bounds are allowed. |
ensureNZC | Flag indicating if non-zero curvature is ensured by exchange rules. |
Definition at line 4565 of file QProblem.c.
returnValue QProblem_reset | ( | QProblem * | _THIS | ) |
Clears all data structures of QProblem except for QP data.
Definition at line 202 of file QProblem.c.
|
inlinestatic |
Resets QP problem counter (to zero).
Definition at line 1852 of file QProblem.h.
|
inlinestatic |
Sets dense constraint matrix of the QP.
Note: Also internal vector Ax is recomputed!
A_new | New dense constraint matrix (with correct dimension!). |
Definition at line 2115 of file QProblem.h.
|
inlinestatic |
Sets constraint matrix of the QP.
Note: Also internal vector Ax is recomputed!
A_new | New constraint matrix. |
Definition at line 2103 of file QProblem.h.
returnValue QProblem_setConstraintProduct | ( | QProblem * | _THIS, |
ConstraintProduct | _constraintProduct | ||
) |
Defines user-defined routine for calculating the constraint product A*x
Definition at line 1084 of file QProblem.c.
|
inlinestatic |
Changes gradient vector of the QP.
g_new | New gradient vector (with correct dimension!). |
Definition at line 1905 of file QProblem.h.
|
inlinestatic |
Sets dense Hessian matrix of the QP. If a null pointer is passed and a) hessianType is HST_IDENTITY, nothing is done, b) hessianType is not HST_IDENTITY, Hessian matrix is set to zero.
H_new | New dense Hessian matrix (with correct dimension!). |
Definition at line 1880 of file QProblem.h.
|
inlinestatic |
Changes the print level.
_hessianType | New Hessian type. |
Definition at line 1787 of file QProblem.h.
|
inlinestatic |
Sets Hessian matrix of the QP.
H_new | New Hessian matrix. |
Definition at line 1868 of file QProblem.h.
returnValue QProblem_setInfeasibilityFlag | ( | QProblem * | _THIS, |
returnValue | returnvalue, | ||
BooleanType | doThrowError | ||
) |
Sets internal infeasibility flag and throws given error in case the far bound strategy is not enabled (as QP might actually not be infeasible in _THIS case).
returnvalue | Returnvalue to be tunneled. |
doThrowError | Flag forcing to throw an error. |
Definition at line 2048 of file QProblem.c.
|
inlinestatic |
Changes lower bound vector of the QP.
lb_new | New lower bound vector (with correct dimension!). |
Definition at line 1924 of file QProblem.h.
|
inlinestatic |
Sets constraints' lower bound vector of the QP.
lbA_new | New constraints' lower bound vector (with correct dimension!). |
Definition at line 2148 of file QProblem.h.
|
inlinestatic |
Changes single entry of lower constraints' bound vector of the QP.
number | Number of entry to be changed. |
value | New value for entry of lower constraints' bound vector (with correct dimension!). |
Definition at line 2175 of file QProblem.h.
|
inlinestatic |
Changes single entry of lower bound vector of the QP.
number | Number of entry to be changed. |
value | New value for entry of lower bound vector. |
Definition at line 1950 of file QProblem.h.
|
inlinestatic |
Overrides current options with given ones.
_options | New options. |
Definition at line 1818 of file QProblem.h.
returnValue QProblem_setPrintLevel | ( | QProblem * | _THIS, |
PrintLevel | _printlevel | ||
) |
Changes the print level.
_printlevel | New print level. |
Definition at line 1193 of file QProblem.c.
|
inlinestatic |
Changes upper bound vector of the QP.
ub_new | New upper bound vector (with correct dimension!). |
Definition at line 1972 of file QProblem.h.
|
inlinestatic |
Sets constraints' upper bound vector of the QP.
ubA_new | New constraints' upper bound vector (with correct dimension!). |
Definition at line 2196 of file QProblem.h.
|
inlinestatic |
Changes single entry of upper constraints' bound vector of the QP.
number | Number of entry to be changed. |
value | New value for entry of upper constraints' bound vector (with correct dimension!). |
Definition at line 2223 of file QProblem.h.
|
inlinestatic |
Changes single entry of upper bound vector of the QP.
number | Number of entry to be changed. |
value | New value for entry of upper bound vector. |
Definition at line 1998 of file QProblem.h.
returnValue QProblem_setupAuxiliaryQP | ( | QProblem * | _THIS, |
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints | ||
) |
Updates QP vectors, working sets and internal data structures in order to start from an optimal solution corresponding to initial guesses of the working set for bounds and constraints.
guessedBounds | Initial guess for working set of bounds. |
guessedConstraints | Initial guess for working set of constraints. |
Definition at line 6528 of file QProblem.c.
returnValue QProblem_setupAuxiliaryQPbounds | ( | QProblem * | _THIS, |
Bounds *const | auxiliaryBounds, | ||
Constraints *const | auxiliaryConstraints, | ||
BooleanType | useRelaxation | ||
) |
Sets up (constraints') bounds of the auxiliary initial QP for given optimal primal/dual solution and given initial working set (assumes that members X, Y and BOUNDS, CONSTRAINTS have already been initialised!).
auxiliaryBounds | Working set of bounds for auxiliary QP. |
auxiliaryConstraints | Working set of constraints for auxiliary QP. |
useRelaxation | Flag indicating if inactive (constraints') bounds shall be relaxed. |
Definition at line 3489 of file QProblem.c.
returnValue QProblem_setupAuxiliaryQPgradient | ( | QProblem * | _THIS | ) |
Sets up gradient of the auxiliary initial QP for given optimal primal/dual solution and given initial working set (assumes that members X, Y and BOUNDS, CONSTRAINTS have already been initialised!).
Definition at line 3444 of file QProblem.c.
returnValue QProblem_setupAuxiliaryQPsolution | ( | QProblem * | _THIS, |
const real_t *const | xOpt, | ||
const real_t *const | yOpt | ||
) |
Sets up the optimal primal/dual solution of the auxiliary initial QP.
xOpt | Optimal primal solution vector. If a NULL pointer is passed, all entries are set to zero. |
yOpt | Optimal dual solution vector. If a NULL pointer is passed, all entries are set to zero. |
Definition at line 3386 of file QProblem.c.
returnValue QProblem_setupAuxiliaryWorkingSet | ( | QProblem * | _THIS, |
Bounds *const | auxiliaryBounds, | ||
Constraints *const | auxiliaryConstraints, | ||
BooleanType | setupAfresh | ||
) |
Sets up bound and constraints data structures according to auxiliaryBounds/Constraints. (If the working set shall be setup afresh, make sure that bounds and constraints data structure have been resetted and the TQ factorisation has been initialised!)
auxiliaryBounds | Working set of bounds for auxiliary QP. |
auxiliaryConstraints | Working set of constraints for auxiliary QP. |
setupAfresh | Flag indicating if given working set shall be setup afresh or by updating the current one. |
Definition at line 3192 of file QProblem.c.
returnValue QProblem_setupInitialCholesky | ( | QProblem * | _THIS | ) |
Computes initial Cholesky decomposition of the projected Hessian making use of the function setupCholeskyDecomposition() or setupCholeskyDecompositionProjected().
Definition at line 544 of file QProblem.c.
returnValue QProblem_setupQPdata | ( | QProblem * | _THIS, |
real_t *const | _H, | ||
const real_t *const | _g, | ||
real_t *const | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA | ||
) |
Sets up dense internal QP data. If the current Hessian is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, memory for Hessian is allocated and it is set to the given one.
_H | Hessian matrix. If Hessian matrix is trivial,a NULL pointer can be passed. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
Definition at line 6682 of file QProblem.c.
returnValue QProblem_setupQPdataFromFile | ( | QProblem * | _THIS, |
const char *const | H_file, | ||
const char *const | g_file, | ||
const char *const | A_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
const char *const | lbA_file, | ||
const char *const | ubA_file | ||
) |
Sets up internal QP data by loading it from files. If the current Hessian is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, memory for Hessian is allocated and it is set to the given one.
H_file | Name of file where Hessian matrix, of neighbouring QP to be solved, is stored. If Hessian matrix is trivial,a NULL pointer can be passed. |
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
A_file | Name of file where constraint matrix, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
lbA_file | Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_file | Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. If no upper constraints' bounds exist, a NULL pointer can be passed. |
Definition at line 6718 of file QProblem.c.
returnValue QProblem_setupQPdataM | ( | QProblem * | _THIS, |
DenseMatrix * | _H, | ||
const real_t *const | _g, | ||
DenseMatrix * | _A, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub, | ||
const real_t *const | _lbA, | ||
const real_t *const | _ubA | ||
) |
Setups internal QP data.
_H | Hessian matrix. If Hessian matrix is trivial,a NULL pointer can be passed. |
_g | Gradient vector. |
_A | Constraint matrix. |
_lb | Lower bound vector (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bound vector (on variables). If no upper bounds exist, a NULL pointer can be passed. |
_lbA | Lower constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
_ubA | Upper constraints' bound vector. If no lower constraints' bounds exist, a NULL pointer can be passed. |
Definition at line 6656 of file QProblem.c.
returnValue QProblem_setupSubjectToType | ( | QProblem * | _THIS | ) |
Determines type of existing constraints and bounds (i.e. implicitly fixed, unbounded etc.).
Definition at line 2760 of file QProblem.c.
returnValue QProblem_setupSubjectToTypeNew | ( | QProblem * | _THIS, |
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new | ||
) |
Determines type of new constraints and bounds (i.e. implicitly fixed, unbounded etc.).
lb_new | New lower bounds. |
ub_new | New upper bounds. |
lbA_new | New lower constraints' bounds. |
ubA_new | New upper constraints' bounds. |
Definition at line 2769 of file QProblem.c.
returnValue QProblem_setupTQfactorisation | ( | QProblem * | _THIS | ) |
Initialises TQ factorisation of A (i.e. A*Q = [0 T]) if NO constraint is active.
Definition at line 3011 of file QProblem.c.
BooleanType QProblem_shallRefactorise | ( | QProblem * | _THIS, |
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints | ||
) |
Determines if it is more efficient to refactorise the matrices when hotstarting or not (i.e. better to update the existing factorisations).
guessedBounds | Guessed new working set of bounds. |
guessedConstraints | Guessed new working set of constraints. |
Definition at line 6617 of file QProblem.c.
returnValue QProblem_solveCurrentEQP | ( | QProblem * | _THIS, |
const int | n_rhs, | ||
const real_t * | g_in, | ||
const real_t * | lb_in, | ||
const real_t * | ub_in, | ||
const real_t * | lbA_in, | ||
const real_t * | ubA_in, | ||
real_t * | x_out, | ||
real_t * | y_out | ||
) |
Solves using the current working set
n_rhs | Number of consecutive right hand sides |
g_in | Gradient of neighbouring QP to be solved. |
lb_in | Lower bounds of neighbouring QP to be solved. If no lower bounds exist, a NULL pointer can be passed. |
ub_in | Upper bounds of neighbouring QP to be solved. If no upper bounds exist, a NULL pointer can be passed. |
lbA_in | Lower constraints' bounds of neighbouring QP to be solved. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_in | Upper constraints' bounds of neighbouring QP to be solved. |
x_out | Output: Primal solution |
y_out | Output: Dual solution |
Definition at line 910 of file QProblem.c.
returnValue QProblem_solveInitialQP | ( | QProblem * | _THIS, |
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Constraints *const | guessedConstraints, | ||
const real_t *const | _R, | ||
int * | nWSR, | ||
real_t *const | cputime | ||
) |
Solves a QProblem whose QP data is assumed to be stored in the member variables. A guess for its primal/dual optimal solution vectors and the corresponding working sets of bounds and constraints can be provided. Note: This function is internally called by all init functions!
xOpt | Optimal primal solution vector. |
yOpt | Optimal dual solution vector. |
guessedBounds | Optimal working set of bounds for solution (xOpt,yOpt). |
guessedConstraints | Optimal working set of constraints for solution (xOpt,yOpt). |
_R | Pre-computed (upper triangular) Cholesky factor of Hessian matrix. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
Definition at line 2232 of file QProblem.c.
returnValue QProblem_solveQP | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
int | nWSRperformed, | ||
BooleanType | isFirstCall | ||
) |
Solves QProblem using online active set strategy. Note: This function is internally called by all hotstart functions!
g_new | Gradient of neighbouring QP to be solved. |
lb_new | Lower bounds of neighbouring QP to be solved. If no lower bounds exist, a NULL pointer can be passed. |
ub_new | Upper bounds of neighbouring QP to be solved. If no upper bounds exist, a NULL pointer can be passed. |
lbA_new | Lower constraints' bounds of neighbouring QP to be solved. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_new | Upper constraints' bounds of neighbouring QP to be solved. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
nWSRperformed | Number of working set recalculations already performed to solve this QP within previous solveQP() calls. This number is always zero, except for successive calls from solveRegularisedQP() or when using the far bound strategy. |
isFirstCall | Indicating whether this is the first call for current QP. |
Definition at line 2414 of file QProblem.c.
returnValue QProblem_solveRegularisedQP | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
const real_t *const | lbA_new, | ||
const real_t *const | ubA_new, | ||
int * | nWSR, | ||
real_t *const | cputime, | ||
int | nWSRperformed, | ||
BooleanType | isFirstCall | ||
) |
Solves QProblem using online active set strategy. Note: This function is internally called by all hotstart functions!
g_new | Gradient of neighbouring QP to be solved. |
lb_new | Lower bounds of neighbouring QP to be solved. If no lower bounds exist, a NULL pointer can be passed. |
ub_new | Upper bounds of neighbouring QP to be solved. If no upper bounds exist, a NULL pointer can be passed. |
lbA_new | Lower constraints' bounds of neighbouring QP to be solved. If no lower constraints' bounds exist, a NULL pointer can be passed. |
ubA_new | Upper constraints' bounds of neighbouring QP to be solved. If no upper constraints' bounds exist, a NULL pointer can be passed. |
nWSR | Input: Maximum number of working set recalculations; Output: Number of performed working set recalculations. |
cputime | Input: Maximum CPU time allowed for QP solution. Output: CPU time spent for QP solution (or to perform nWSR iterations). |
nWSRperformed | Number of working set recalculations already performed to solve this QP within previous solveRegularisedQP() calls. This number is always zero, except for successive calls when using the far bound strategy. |
isFirstCall | Indicating whether this is the first call for current QP. |
Definition at line 2655 of file QProblem.c.
returnValue QProblem_updateFarBounds | ( | QProblem * | _THIS, |
real_t | curFarBound, | ||
int | nRamp, | ||
const real_t *const | lb_new, | ||
real_t *const | lb_new_far, | ||
const real_t *const | ub_new, | ||
real_t *const | ub_new_far, | ||
const real_t *const | lbA_new, | ||
real_t *const | lbA_new_far, | ||
const real_t *const | ubA_new, | ||
real_t *const | ubA_new_far | ||
) |
...
curFarBound | ... |
nRamp | ... |
lb_new | ... |
lb_new_far | ... |
ub_new | ... |
ub_new_far | ... |
lbA_new | ... |
lbA_new_far | ... |
ubA_new | ... |
ubA_new_far | ... |
Definition at line 6185 of file QProblem.c.
|
inlinestatic |
Returns if the QP has been internally regularised.
Definition at line 1797 of file QProblem.h.
returnValue QProblem_writeQpDataIntoMatFile | ( | QProblem * | _THIS, |
const char *const | filename | ||
) |
...
filename | Mat file name. |
returnValue QProblem_writeQpWorkspaceIntoMatFile | ( | QProblem * | _THIS, |
const char *const | filename | ||
) |
...
filename | Mat file name. |
returnValue QProblemBCPY_computeCholesky | ( | QProblem * | _THIS | ) |
Computes the Cholesky decomposition of the (simply projected) Hessian (i.e. R^T*R = Z^T*H*Z). It only works in the case where Z is a simple projection matrix! Note: If Hessian turns out not to be positive definite, the Hessian type is set to HST_SEMIDEF accordingly.
Definition at line 1525 of file QProblem.c.
returnValue QProblemBCPY_determineDataShift | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new, | ||
real_t *const | delta_g, | ||
real_t *const | delta_lb, | ||
real_t *const | delta_ub, | ||
BooleanType * | Delta_bB_isZero | ||
) |
Determines step direction of the shift of the QP data.
g_new | New gradient vector. |
lb_new | New lower bounds. |
ub_new | New upper bounds. |
delta_g | Output: Step direction of gradient vector. |
delta_lb | Output: Step direction of lower bounds. |
delta_ub | Output: Step direction of upper bounds. |
Delta_bB_isZero | Output: Indicates if active bounds are to be shifted. |
Definition at line 1818 of file QProblem.c.
real_t QProblemBCPY_getRelativeHomotopyLength | ( | QProblem * | _THIS, |
const real_t *const | g_new, | ||
const real_t *const | lb_new, | ||
const real_t *const | ub_new | ||
) |
Compute relative length of homotopy in data space for termination criterion.
g_new | Final gradient. |
lb_new | Final lower variable bounds. |
ub_new | Final upper variable bounds. |
Definition at line 2183 of file QProblem.c.
returnValue QProblemBCPY_loadQPvectorsFromFile | ( | QProblem * | _THIS, |
const char *const | g_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file, | ||
real_t *const | g_new, | ||
real_t *const | lb_new, | ||
real_t *const | ub_new | ||
) |
Loads new QP vectors from files (internal members are not affected!).
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
g_new | Output: Gradient of neighbouring QP to be solved. |
lb_new | Output: Lower bounds of neighbouring QP to be solved |
ub_new | Output: Upper bounds of neighbouring QP to be solved |
Definition at line 1987 of file QProblem.c.
returnValue QProblemBCPY_obtainAuxiliaryWorkingSet | ( | QProblem * | _THIS, |
const real_t *const | xOpt, | ||
const real_t *const | yOpt, | ||
Bounds *const | guessedBounds, | ||
Bounds * | auxiliaryBounds | ||
) |
Obtains the desired working set for the auxiliary initial QP in accordance with the user specifications
xOpt | Optimal primal solution vector. If a NULL pointer is passed, all entries are assumed to be zero. |
yOpt | Optimal dual solution vector. If a NULL pointer is passed, all entries are assumed to be zero. |
guessedBounds | Guessed working set for solution (xOpt,yOpt). |
auxiliaryBounds | Input: Allocated bound object. Output: Working set for auxiliary QP. |
Definition at line 1604 of file QProblem.c.
returnValue QProblemBCPY_setupQPdata | ( | QProblem * | _THIS, |
real_t *const | _H, | ||
const real_t *const | _g, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub | ||
) |
Sets up internal QP data. If the current Hessian is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, memory for Hessian is allocated and it is set to the given one.
_H | Hessian matrix. If Hessian matrix is trivial,a NULL pointer can be passed. |
_g | Gradient vector. |
_lb | Lower bounds (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bounds (on variables). If no upper bounds exist, a NULL pointer can be passed. |
Definition at line 1895 of file QProblem.c.
returnValue QProblemBCPY_setupQPdataFromFile | ( | QProblem * | _THIS, |
const char *const | H_file, | ||
const char *const | g_file, | ||
const char *const | lb_file, | ||
const char *const | ub_file | ||
) |
Sets up internal QP data by loading it from files. If the current Hessian is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, memory for Hessian is allocated and it is set to the given one.
H_file | Name of file where Hessian matrix, of neighbouring QP to be solved, is stored. If Hessian matrix is trivial,a NULL pointer can be passed. |
g_file | Name of file where gradient, of neighbouring QP to be solved, is stored. |
lb_file | Name of file where lower bounds, of neighbouring QP to be solved, is stored. If no lower bounds exist, a NULL pointer can be passed. |
ub_file | Name of file where upper bounds, of neighbouring QP to be solved, is stored. If no upper bounds exist, a NULL pointer can be passed. |
Definition at line 1919 of file QProblem.c.
returnValue QProblemBCPY_setupQPdataM | ( | QProblem * | _THIS, |
DenseMatrix * | _H, | ||
const real_t *const | _g, | ||
const real_t *const | _lb, | ||
const real_t *const | _ub | ||
) |
Sets up internal QP data.
_H | Hessian matrix. |
_g | Gradient vector. |
_lb | Lower bounds (on variables). If no lower bounds exist, a NULL pointer can be passed. |
_ub | Upper bounds (on variables). If no upper bounds exist, a NULL pointer can be passed. |
Definition at line 1881 of file QProblem.c.
returnValue QProblemBCPY_updateFarBounds | ( | QProblem * | _THIS, |
real_t | curFarBound, | ||
int | nRamp, | ||
const real_t *const | lb_new, | ||
real_t *const | lb_new_far, | ||
const real_t *const | ub_new, | ||
real_t *const | ub_new_far | ||
) |
...
curFarBound | ... |
nRamp | ... |
lb_new | ... |
lb_new_far | ... |
ub_new | ... |
ub_new_far | ... |
Definition at line 6245 of file QProblem.c.
void QProblemCON | ( | QProblem * | _THIS, |
int | _nV, | ||
int | _nC, | ||
HessianType | _hessianType | ||
) |
Constructor which takes the QP dimension and Hessian type information. If the Hessian is the zero (i.e. HST_ZERO) or the identity matrix (i.e. HST_IDENTITY), respectively, no memory is allocated for it and a NULL pointer can be passed for it to the init() functions.
_nV | Number of variables. |
_nC | Number of constraints. |
_hessianType | Type of Hessian matrix. |
Definition at line 51 of file QProblem.c.
Copies all members from given rhs object.
Definition at line 132 of file QProblem.c.