1 /* 2 * This file is part of qpOASES. 3 * 4 * qpOASES -- An Implementation of the Online Active Set Strategy. 5 * Copyright (C) 2007-2017 by Hans Joachim Ferreau, Andreas Potschka, 6 * Christian Kirches et al. All rights reserved. 7 * 8 * qpOASES is free software; you can redistribute it and/or 9 * modify it under the terms of the GNU Lesser General Public 10 * License as published by the Free Software Foundation; either 11 * version 2.1 of the License, or (at your option) any later version. 12 * 13 * qpOASES is distributed in the hope that it will be useful, 14 * but WITHOUT ANY WARRANTY; without even the implied warranty of 15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. 16 * See the GNU Lesser General Public License for more details. 17 * 18 * You should have received a copy of the GNU Lesser General Public 19 * License along with qpOASES; if not, write to the Free Software 20 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 21 * 22 */ 23 24 25 /** 26 * \file include/qpOASES/QProblem.hpp 27 * \author Hans Joachim Ferreau, Andreas Potschka, Christian Kirches 28 * \version 3.2 29 * \date 2007-2017 30 * 31 * Declaration of the QProblem class which is able to use the newly 32 * developed online active set strategy for parametric quadratic programming. 33 */ 34 35 36 37 #ifndef QPOASES_QPROBLEM_HPP 38 #define QPOASES_QPROBLEM_HPP 39 40 41 #include <qpOASES/QProblemB.hpp> 42 #include <qpOASES/Constraints.hpp> 43 #include <qpOASES/ConstraintProduct.hpp> 44 #include <qpOASES/Matrices.hpp> 45 46 47 BEGIN_NAMESPACE_QPOASES 48 49 50 /** 51 * \brief Implements the online active set strategy for QPs with general constraints. 52 * 53 * A class for setting up and solving quadratic programs. The main feature is 54 * the possibily to use the newly developed online active set strategy for 55 * parametric quadratic programming. 56 * 57 * \author Hans Joachim Ferreau, Andreas Potschka, Christian Kirches 58 * \version 3.2 59 * \date 2007-2017 60 */ 61 class QProblem : public QProblemB 62 { 63 /* allow SolutionAnalysis class to access private members */ 64 friend class SolutionAnalysis; 65 66 /* 67 * PUBLIC MEMBER FUNCTIONS 68 */ 69 public: 70 /** Default constructor. */ 71 QProblem( ); 72 73 /** Constructor which takes the QP dimension and Hessian type 74 * information. If the Hessian is the zero (i.e. HST_ZERO) or the 75 * identity matrix (i.e. HST_IDENTITY), respectively, no memory 76 * is allocated for it and a NULL pointer can be passed for it 77 * to the init() functions. */ 78 QProblem( int_t _nV, /**< Number of variables. */ 79 int_t _nC, /**< Number of constraints. */ 80 HessianType _hessianType = HST_UNKNOWN, /**< Type of Hessian matrix. */ 81 BooleanType allocDenseMats = BT_TRUE /**< Enable allocation of dense matrices. */ 82 ); 83 84 /** Copy constructor (deep copy). */ 85 QProblem( const QProblem& rhs /**< Rhs object. */ 86 ); 87 88 /** Destructor. */ 89 virtual ~QProblem( ); 90 91 /** Assignment operator (deep copy). */ 92 virtual QProblem& operator=( const QProblem& rhs /**< Rhs object. */ 93 ); 94 95 96 /** Clears all data structures of QProblemB except for QP data. 97 * \return SUCCESSFUL_RETURN \n 98 RET_RESET_FAILED */ 99 virtual returnValue reset( ); 100 101 102 /** Initialises a QP problem with given QP data and tries to solve it 103 * using at most nWSR iterations. Depending on the parameter constellation it: \n 104 * 1. 0, 0, 0 : starts with xOpt = 0, yOpt = 0 and gB/gC empty (or all implicit equality bounds), \n 105 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB/gC by "clipping", \n 106 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB/gC from yOpt != 0, \n 107 * 4. 0, 0, gB/gC: starts with xOpt = 0, yOpt = 0 and gB/gC, \n 108 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB/gC from yOpt != 0, \n 109 * 6. xOpt, 0, gB/gC: starts with xOpt, yOpt = 0 and gB/gC, \n 110 * 7. xOpt, yOpt, gB/gC: starts with xOpt, yOpt and gB/gC (assume them to be consistent!) 111 * 112 * Note: This function internally calls solveInitialQP for initialisation! 113 * 114 * \return SUCCESSFUL_RETURN \n 115 RET_INIT_FAILED \n 116 RET_INIT_FAILED_CHOLESKY \n 117 RET_INIT_FAILED_TQ \n 118 RET_INIT_FAILED_HOTSTART \n 119 RET_INIT_FAILED_INFEASIBILITY \n 120 RET_INIT_FAILED_UNBOUNDEDNESS \n 121 RET_MAX_NWSR_REACHED \n 122 RET_INVALID_ARGUMENTS */ 123 returnValue init( SymmetricMatrix *_H, /**< Hessian matrix (a shallow copy is made). */ 124 const real_t* const _g, /**< Gradient vector. */ 125 Matrix *_A, /**< Constraint matrix (a shallow copy is made). */ 126 const real_t* const _lb, /**< Lower bound vector (on variables). \n 127 If no lower bounds exist, a NULL pointer can be passed. */ 128 const real_t* const _ub, /**< Upper bound vector (on variables). \n 129 If no upper bounds exist, a NULL pointer can be passed. */ 130 const real_t* const _lbA, /**< Lower constraints' bound vector. \n 131 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 132 const real_t* const _ubA, /**< Upper constraints' bound vector. \n 133 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 134 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. 135 Output: Number of performed working set recalculations. */ 136 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 137 Output: CPU time spent for QP initialisation (if pointer passed). */ 138 const real_t* const xOpt = 0, /**< Optimal primal solution vector. \n 139 (If a null pointer is passed, the old primal solution is kept!) */ 140 const real_t* const yOpt = 0, /**< Optimal dual solution vector. \n 141 (If a null pointer is passed, the old dual solution is kept!) */ 142 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 143 (If a null pointer is passed, all bounds are assumed inactive!) */ 144 const Constraints* const guessedConstraints = 0,/**< Optimal working set of constraints for solution (xOpt,yOpt). \n 145 (If a null pointer is passed, all constraints are assumed inactive!) */ 146 const real_t* const _R = 0 /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. 147 The Cholesky factor must be stored in a real_t array of size nV*nV 148 in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! \n 149 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 150 ); 151 152 153 /** Initialises a QP problem with given QP data and tries to solve it 154 * using at most nWSR iterations. Depending on the parameter constellation it: \n 155 * 1. 0, 0, 0 : starts with xOpt = 0, yOpt = 0 and gB/gC empty (or all implicit equality bounds), \n 156 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB/gC by "clipping", \n 157 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB/gC from yOpt != 0, \n 158 * 4. 0, 0, gB/gC: starts with xOpt = 0, yOpt = 0 and gB/gC, \n 159 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB/gC from yOpt != 0, \n 160 * 6. xOpt, 0, gB/gC: starts with xOpt, yOpt = 0 and gB/gC, \n 161 * 7. xOpt, yOpt, gB/gC: starts with xOpt, yOpt and gB/gC (assume them to be consistent!) 162 * 163 * Note: This function internally calls solveInitialQP for initialisation! 164 * 165 * \return SUCCESSFUL_RETURN \n 166 RET_INIT_FAILED \n 167 RET_INIT_FAILED_CHOLESKY \n 168 RET_INIT_FAILED_TQ \n 169 RET_INIT_FAILED_HOTSTART \n 170 RET_INIT_FAILED_INFEASIBILITY \n 171 RET_INIT_FAILED_UNBOUNDEDNESS \n 172 RET_MAX_NWSR_REACHED \n 173 RET_INVALID_ARGUMENTS */ 174 returnValue init( const real_t* const _H, /**< Hessian matrix (a shallow copy is made). \n 175 If Hessian matrix is trivial, a NULL pointer can be passed. */ 176 const real_t* const _g, /**< Gradient vector. */ 177 const real_t* const _A, /**< Constraint matrix (a shallow copy is made). */ 178 const real_t* const _lb, /**< Lower bound vector (on variables). \n 179 If no lower bounds exist, a NULL pointer can be passed. */ 180 const real_t* const _ub, /**< Upper bound vector (on variables). \n 181 If no upper bounds exist, a NULL pointer can be passed. */ 182 const real_t* const _lbA, /**< Lower constraints' bound vector. \n 183 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 184 const real_t* const _ubA, /**< Upper constraints' bound vector. \n 185 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 186 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. 187 Output: Number of performed working set recalculations. */ 188 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 189 Output: CPU time spent for QP initialisation (if pointer passed). */ 190 const real_t* const xOpt = 0, /**< Optimal primal solution vector. \n 191 (If a null pointer is passed, the old primal solution is kept!) */ 192 const real_t* const yOpt = 0, /**< Optimal dual solution vector. \n 193 (If a null pointer is passed, the old dual solution is kept!) */ 194 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 195 (If a null pointer is passed, all bounds are assumed inactive!) */ 196 const Constraints* const guessedConstraints = 0,/**< Optimal working set of constraints for solution (xOpt,yOpt). \n 197 (If a null pointer is passed, all constraints are assumed inactive!) */ 198 const real_t* const _R = 0 /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. 199 The Cholesky factor must be stored in a real_t array of size nV*nV 200 in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! \n 201 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 202 ); 203 204 /** Initialises a QP problem with given data to be read from files and solves it 205 * using at most nWSR iterations. Depending on the parameter constellation it: \n 206 * 1. 0, 0, 0 : starts with xOpt = 0, yOpt = 0 and gB/gC empty (or all implicit equality bounds), \n 207 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB/gC by "clipping", \n 208 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB/gC from yOpt != 0, \n 209 * 4. 0, 0, gB/gC: starts with xOpt = 0, yOpt = 0 and gB/gC, \n 210 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB/gC from yOpt != 0, \n 211 * 6. xOpt, 0, gB/gC: starts with xOpt, yOpt = 0 and gB/gC, \n 212 * 7. xOpt, yOpt, gB/gC: starts with xOpt, yOpt and gB/gC (assume them to be consistent!) 213 * 214 * Note: This function internally calls solveInitialQP for initialisation! 215 * 216 * \return SUCCESSFUL_RETURN \n 217 RET_INIT_FAILED \n 218 RET_INIT_FAILED_CHOLESKY \n 219 RET_INIT_FAILED_TQ \n 220 RET_INIT_FAILED_HOTSTART \n 221 RET_INIT_FAILED_INFEASIBILITY \n 222 RET_INIT_FAILED_UNBOUNDEDNESS \n 223 RET_MAX_NWSR_REACHED \n 224 RET_UNABLE_TO_READ_FILE \n 225 RET_INVALID_ARGUMENTS */ 226 returnValue init( const char* const H_file, /**< Name of file where Hessian matrix is stored. \n 227 If Hessian matrix is trivial, a NULL pointer can be passed. */ 228 const char* const g_file, /**< Name of file where gradient vector is stored. */ 229 const char* const A_file, /**< Name of file where constraint matrix is stored. */ 230 const char* const lb_file, /**< Name of file where lower bound vector. \n 231 If no lower bounds exist, a NULL pointer can be passed. */ 232 const char* const ub_file, /**< Name of file where upper bound vector. \n 233 If no upper bounds exist, a NULL pointer can be passed. */ 234 const char* const lbA_file, /**< Name of file where lower constraints' bound vector. \n 235 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 236 const char* const ubA_file, /**< Name of file where upper constraints' bound vector. \n 237 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 238 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. 239 Output: Number of performed working set recalculations. */ 240 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 241 Output: CPU time spent for QP initialisation (if pointer passed). */ 242 const real_t* const xOpt = 0, /**< Optimal primal solution vector. \n 243 (If a null pointer is passed, the old primal solution is kept!) */ 244 const real_t* const yOpt = 0, /**< Optimal dual solution vector. \n 245 (If a null pointer is passed, the old dual solution is kept!) */ 246 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 247 (If a null pointer is passed, all bounds are assumed inactive!) */ 248 const Constraints* const guessedConstraints = 0,/**< Optimal working set of constraints for solution (xOpt,yOpt). \n 249 (If a null pointer is passed, all constraints are assumed inactive!) */ 250 const char* const R_file = 0 /**< Name of the file where a pre-computed (upper triangular) Cholesky factor 251 of the Hessian matrix is stored. \n 252 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 253 ); 254 255 256 /** Solves an initialised QP sequence using the online active set strategy. 257 * By default, QP solution is started from previous solution. If a guess 258 * for the working set is provided, an initialised homotopy is performed. 259 * 260 * Note: This function internally calls solveQP/solveRegularisedQP 261 * for solving an initialised QP! 262 * 263 * \return SUCCESSFUL_RETURN \n 264 RET_MAX_NWSR_REACHED \n 265 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 266 RET_HOTSTART_FAILED \n 267 RET_SHIFT_DETERMINATION_FAILED \n 268 RET_STEPDIRECTION_DETERMINATION_FAILED \n 269 RET_STEPLENGTH_DETERMINATION_FAILED \n 270 RET_HOMOTOPY_STEP_FAILED \n 271 RET_HOTSTART_STOPPED_INFEASIBILITY \n 272 RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 273 returnValue hotstart( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 274 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 275 If no lower bounds exist, a NULL pointer can be passed. */ 276 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 277 If no upper bounds exist, a NULL pointer can be passed. */ 278 const real_t* const lbA_new, /**< Lower constraints' bounds of neighbouring QP to be solved. \n 279 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 280 const real_t* const ubA_new, /**< Upper constraints' bounds of neighbouring QP to be solved. \n 281 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 282 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 283 Output: Number of performed working set recalculations. */ 284 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP solution. \n 285 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 286 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 287 (If a null pointer is passed, the previous working set of bounds is kept!) */ 288 const Constraints* const guessedConstraints = 0 /**< Optimal working set of constraints for solution (xOpt,yOpt). \n 289 (If a null pointer is passed, the previous working set of constraints is kept!) */ 290 ); 291 292 /** Solves an initialised QP sequence using the online active set strategy, 293 * where QP data is read from files. 294 * By default, QP solution is started from previous solution. If a guess 295 * for the working set is provided, an initialised homotopy is performed. 296 * 297 * Note: This function internally calls solveQP/solveRegularisedQP 298 * for solving an initialised QP! 299 * 300 * \return SUCCESSFUL_RETURN \n 301 RET_MAX_NWSR_REACHED \n 302 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 303 RET_HOTSTART_FAILED \n 304 RET_SHIFT_DETERMINATION_FAILED \n 305 RET_STEPDIRECTION_DETERMINATION_FAILED \n 306 RET_STEPLENGTH_DETERMINATION_FAILED \n 307 RET_HOMOTOPY_STEP_FAILED \n 308 RET_HOTSTART_STOPPED_INFEASIBILITY \n 309 RET_HOTSTART_STOPPED_UNBOUNDEDNESS \n 310 RET_UNABLE_TO_READ_FILE \n 311 RET_INVALID_ARGUMENTS */ 312 returnValue hotstart( const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 313 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 314 If no lower bounds exist, a NULL pointer can be passed. */ 315 const char* const ub_file, /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 316 If no upper bounds exist, a NULL pointer can be passed. */ 317 const char* const lbA_file, /**< Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. \n 318 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 319 const char* const ubA_file, /**< Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. \n 320 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 321 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 322 Output: Number of performed working set recalculations. */ 323 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP solution. \n 324 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 325 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 326 (If a null pointer is passed, the previous working set of bounds is kept!) */ 327 const Constraints* const guessedConstraints = 0 /**< Optimal working set of constraints for solution (xOpt,yOpt). \n 328 (If a null pointer is passed, the previous working set of constraints is kept!) */ 329 ); 330 331 332 /** Solves an equality-constrained QP problem resulting from the current working set. 333 * \return SUCCESSFUL_RETURN \n 334 * RET_STEPDIRECTION_FAILED_TQ \n 335 * RET_STEPDIRECTION_FAILED_CHOLESKY \n 336 * RET_INVALID_ARGUMENTS */ 337 returnValue solveCurrentEQP ( const int_t n_rhs, /**< Number of consecutive right hand sides */ 338 const real_t* g_in, /**< Gradient of neighbouring QP to be solved. */ 339 const real_t* lb_in, /**< Lower bounds of neighbouring QP to be solved. \n 340 If no lower bounds exist, a NULL pointer can be passed. */ 341 const real_t* ub_in, /**< Upper bounds of neighbouring QP to be solved. \n 342 If no upper bounds exist, a NULL pointer can be passed. */ 343 const real_t* lbA_in, /**< Lower constraints' bounds of neighbouring QP to be solved. \n 344 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 345 const real_t* ubA_in, /**< Upper constraints' bounds of neighbouring QP to be solved. \n */ 346 real_t* x_out, /**< Output: Primal solution */ 347 real_t* y_out /**< Output: Dual solution */ 348 ); 349 350 /** Writes a vector with the state of the working set 351 * \return SUCCESSFUL_RETURN \n 352 * RET_INVALID_ARGUMENTS */ 353 virtual returnValue getWorkingSet( real_t* workingSet /** Output: array containing state of the working set. */ 354 ); 355 356 /** Writes a vector with the state of the working set of bounds 357 * \return SUCCESSFUL_RETURN \n 358 * RET_INVALID_ARGUMENTS */ 359 virtual returnValue getWorkingSetBounds( real_t* workingSetB /** Output: array containing state of the working set of bounds. */ 360 ); 361 362 /** Writes a vector with the state of the working set of constraints 363 * \return SUCCESSFUL_RETURN \n 364 * RET_INVALID_ARGUMENTS */ 365 virtual returnValue getWorkingSetConstraints( real_t* workingSetC /** Output: array containing state of the working set of constraints. */ 366 ); 367 368 369 /** Returns current constraints object of the QP (deep copy). 370 * \return SUCCESSFUL_RETURN \n 371 RET_QPOBJECT_NOT_SETUP */ 372 inline returnValue getConstraints( Constraints& _constraints /** Output: Constraints object. */ 373 ) const; 374 375 376 /** Returns the number of constraints. 377 * \return Number of constraints. */ 378 inline int_t getNC( ) const; 379 380 /** Returns the number of (implicitly defined) equality constraints. 381 * \return Number of (implicitly defined) equality constraints. */ 382 inline int_t getNEC( ) const; 383 384 /** Returns the number of active constraints. 385 * \return Number of active constraints. */ 386 inline int_t getNAC( ) const; 387 388 /** Returns the number of inactive constraints. 389 * \return Number of inactive constraints. */ 390 inline int_t getNIAC( ) const; 391 392 /** Returns the dimension of null space. 393 * \return Dimension of null space. */ 394 virtual int_t getNZ( ) const; 395 396 397 /** Returns the dual solution vector (deep copy). 398 * \return SUCCESSFUL_RETURN \n 399 RET_QP_NOT_SOLVED */ 400 virtual returnValue getDualSolution( real_t* const yOpt /**< Output: Dual solution vector (if QP has been solved). */ 401 ) const; 402 403 404 /** Defines user-defined routine for calculating the constraint product A*x 405 * \return SUCCESSFUL_RETURN \n */ 406 returnValue setConstraintProduct( ConstraintProduct* const _constraintProduct 407 ); 408 409 410 /** Prints concise list of properties of the current QP. 411 * \return SUCCESSFUL_RETURN \n */ 412 virtual returnValue printProperties( ); 413 414 /** Set the incoming array to true for each variable entry that is 415 in the set of free variables */ 416 returnValue getFreeVariablesFlags( BooleanType* varIsFree ); 417 418 419 /* 420 * PROTECTED MEMBER FUNCTIONS 421 */ 422 protected: 423 /** Frees all allocated memory. 424 * \return SUCCESSFUL_RETURN */ 425 returnValue clear( ); 426 427 /** Copies all members from given rhs object. 428 * \return SUCCESSFUL_RETURN */ 429 returnValue copy( const QProblem& rhs /**< Rhs object. */ 430 ); 431 432 /** Solves a QProblem whose QP data is assumed to be stored in the member variables. 433 * A guess for its primal/dual optimal solution vectors and the corresponding 434 * working sets of bounds and constraints can be provided. 435 * Note: This function is internally called by all init functions! 436 * \return SUCCESSFUL_RETURN \n 437 RET_INIT_FAILED \n 438 RET_INIT_FAILED_CHOLESKY \n 439 RET_INIT_FAILED_TQ \n 440 RET_INIT_FAILED_HOTSTART \n 441 RET_INIT_FAILED_INFEASIBILITY \n 442 RET_INIT_FAILED_UNBOUNDEDNESS \n 443 RET_MAX_NWSR_REACHED */ 444 returnValue solveInitialQP( const real_t* const xOpt, /**< Optimal primal solution vector.*/ 445 const real_t* const yOpt, /**< Optimal dual solution vector. */ 446 const Bounds* const guessedBounds, /**< Optimal working set of bounds for solution (xOpt,yOpt). */ 447 const Constraints* const guessedConstraints, /**< Optimal working set of constraints for solution (xOpt,yOpt). */ 448 const real_t* const _R, /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. */ 449 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 450 Output: Number of performed working set recalculations. */ 451 real_t* const cputime /**< Input: Maximum CPU time allowed for QP solution. \n 452 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 453 ); 454 455 /** Solves QProblem using online active set strategy. 456 * Note: This function is internally called by all hotstart functions! 457 * \return SUCCESSFUL_RETURN \n 458 RET_MAX_NWSR_REACHED \n 459 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 460 RET_HOTSTART_FAILED \n 461 RET_SHIFT_DETERMINATION_FAILED \n 462 RET_STEPDIRECTION_DETERMINATION_FAILED \n 463 RET_STEPLENGTH_DETERMINATION_FAILED \n 464 RET_HOMOTOPY_STEP_FAILED \n 465 RET_HOTSTART_STOPPED_INFEASIBILITY \n 466 RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 467 returnValue solveQP( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 468 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 469 If no lower bounds exist, a NULL pointer can be passed. */ 470 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 471 If no upper bounds exist, a NULL pointer can be passed. */ 472 const real_t* const lbA_new, /**< Lower constraints' bounds of neighbouring QP to be solved. \n 473 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 474 const real_t* const ubA_new, /**< Upper constraints' bounds of neighbouring QP to be solved. \n 475 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 476 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 477 Output: Number of performed working set recalculations. */ 478 real_t* const cputime, /**< Input: Maximum CPU time allowed for QP solution. \n 479 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 480 int_t nWSRperformed = 0, /**< Number of working set recalculations already performed to solve 481 this QP within previous solveQP() calls. This number is 482 always zero, except for successive calls from solveRegularisedQP() 483 or when using the far bound strategy. */ 484 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 485 ); 486 487 488 /** Solves QProblem using online active set strategy. 489 * Note: This function is internally called by all hotstart functions! 490 * \return SUCCESSFUL_RETURN \n 491 RET_MAX_NWSR_REACHED \n 492 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 493 RET_HOTSTART_FAILED \n 494 RET_SHIFT_DETERMINATION_FAILED \n 495 RET_STEPDIRECTION_DETERMINATION_FAILED \n 496 RET_STEPLENGTH_DETERMINATION_FAILED \n 497 RET_HOMOTOPY_STEP_FAILED \n 498 RET_HOTSTART_STOPPED_INFEASIBILITY \n 499 RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 500 returnValue solveRegularisedQP( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 501 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 502 If no lower bounds exist, a NULL pointer can be passed. */ 503 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 504 If no upper bounds exist, a NULL pointer can be passed. */ 505 const real_t* const lbA_new, /**< Lower constraints' bounds of neighbouring QP to be solved. \n 506 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 507 const real_t* const ubA_new, /**< Upper constraints' bounds of neighbouring QP to be solved. \n 508 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 509 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 510 Output: Number of performed working set recalculations. */ 511 real_t* const cputime, /**< Input: Maximum CPU time allowed for QP solution. \n 512 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 513 int_t nWSRperformed = 0, /**< Number of working set recalculations already performed to solve 514 this QP within previous solveRegularisedQP() calls. This number is 515 always zero, except for successive calls when using the far bound strategy. */ 516 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 517 ); 518 519 520 /** Update activities in a hot start if some of the bounds have 521 become infinity or if variables have become fixed. */ 522 /* \return SUCCESSFUL_RETURN \n 523 RET_HOTSTART_FAILED */ 524 virtual returnValue updateActivitiesForHotstart( const real_t* const lb_new, /**< New lower bounds. */ 525 const real_t* const ub_new, /**< New upper bounds. */ 526 const real_t* const lbA_new, /**< New lower constraints' bounds. */ 527 const real_t* const ubA_new /**< New upper constraints' bounds. */ 528 ); 529 530 531 /** Determines type of existing constraints and bounds (i.e. implicitly fixed, unbounded etc.). 532 * \return SUCCESSFUL_RETURN \n 533 RET_SETUPSUBJECTTOTYPE_FAILED */ 534 virtual returnValue setupSubjectToType( ); 535 536 /** Determines type of new constraints and bounds (i.e. implicitly fixed, unbounded etc.). 537 * \return SUCCESSFUL_RETURN \n 538 RET_SETUPSUBJECTTOTYPE_FAILED */ 539 using QProblemB::setupSubjectToType; 540 virtual returnValue setupSubjectToType( const real_t* const lb_new, /**< New lower bounds. */ 541 const real_t* const ub_new, /**< New upper bounds. */ 542 const real_t* const lbA_new, /**< New lower constraints' bounds. */ 543 const real_t* const ubA_new /**< New upper constraints' bounds. */ 544 ); 545 546 /** Computes the Cholesky decomposition of the projected Hessian (i.e. R^T*R = Z^T*H*Z). 547 * Note: If Hessian turns out not to be positive definite, the Hessian type 548 * is set to HST_SEMIDEF accordingly. 549 * \return SUCCESSFUL_RETURN \n 550 * RET_HESSIAN_NOT_SPD \n 551 * RET_INDEXLIST_CORRUPTED */ 552 virtual returnValue computeProjectedCholesky( ); 553 554 /** Computes initial Cholesky decomposition of the projected Hessian making 555 * use of the function computeCholesky() or computeProjectedCholesky(). 556 * \return SUCCESSFUL_RETURN \n 557 * RET_HESSIAN_NOT_SPD \n 558 * RET_INDEXLIST_CORRUPTED */ 559 virtual returnValue setupInitialCholesky( ); 560 561 /** Initialises TQ factorisation of A (i.e. A*Q = [0 T]) if NO constraint is active. 562 * \return SUCCESSFUL_RETURN \n 563 RET_INDEXLIST_CORRUPTED */ 564 virtual returnValue setupTQfactorisation( ); 565 566 567 /** Obtains the desired working set for the auxiliary initial QP in 568 * accordance with the user specifications 569 * (assumes that member AX has already been initialised!) 570 * \return SUCCESSFUL_RETURN \n 571 RET_OBTAINING_WORKINGSET_FAILED \n 572 RET_INVALID_ARGUMENTS */ 573 returnValue obtainAuxiliaryWorkingSet( const real_t* const xOpt, /**< Optimal primal solution vector. 574 * If a NULL pointer is passed, all entries are assumed to be zero. */ 575 const real_t* const yOpt, /**< Optimal dual solution vector. 576 * If a NULL pointer is passed, all entries are assumed to be zero. */ 577 const Bounds* const guessedBounds, /**< Guessed working set of bounds for solution (xOpt,yOpt). */ 578 const Constraints* const guessedConstraints, /**< Guessed working set for solution (xOpt,yOpt). */ 579 Bounds* auxiliaryBounds, /**< Input: Allocated bound object. \n 580 * Ouput: Working set of constraints for auxiliary QP. */ 581 Constraints* auxiliaryConstraints /**< Input: Allocated bound object. \n 582 * Ouput: Working set for auxiliary QP. */ 583 ) const; 584 585 /** Sets up bound and constraints data structures according to auxiliaryBounds/Constraints. 586 * (If the working set shall be setup afresh, make sure that 587 * bounds and constraints data structure have been resetted 588 * and the TQ factorisation has been initialised!) 589 * \return SUCCESSFUL_RETURN \n 590 RET_SETUP_WORKINGSET_FAILED \n 591 RET_INVALID_ARGUMENTS \n 592 RET_UNKNOWN_BUG */ 593 virtual returnValue setupAuxiliaryWorkingSet( const Bounds* const auxiliaryBounds, /**< Working set of bounds for auxiliary QP. */ 594 const Constraints* const auxiliaryConstraints, /**< Working set of constraints for auxiliary QP. */ 595 BooleanType setupAfresh /**< Flag indicating if given working set shall be 596 * setup afresh or by updating the current one. */ 597 ); 598 599 /** Sets up the optimal primal/dual solution of the auxiliary initial QP. 600 * \return SUCCESSFUL_RETURN */ 601 returnValue setupAuxiliaryQPsolution( const real_t* const xOpt, /**< Optimal primal solution vector. 602 * If a NULL pointer is passed, all entries are set to zero. */ 603 const real_t* const yOpt /**< Optimal dual solution vector. 604 * If a NULL pointer is passed, all entries are set to zero. */ 605 ); 606 607 /** Sets up gradient of the auxiliary initial QP for given 608 * optimal primal/dual solution and given initial working set 609 * (assumes that members X, Y and BOUNDS, CONSTRAINTS have already been initialised!). 610 * \return SUCCESSFUL_RETURN */ 611 returnValue setupAuxiliaryQPgradient( ); 612 613 /** Sets up (constraints') bounds of the auxiliary initial QP for given 614 * optimal primal/dual solution and given initial working set 615 * (assumes that members X, Y and BOUNDS, CONSTRAINTS have already been initialised!). 616 * \return SUCCESSFUL_RETURN \n 617 RET_UNKNOWN_BUG */ 618 returnValue setupAuxiliaryQPbounds( const Bounds* const auxiliaryBounds, /**< Working set of bounds for auxiliary QP. */ 619 const Constraints* const auxiliaryConstraints, /**< Working set of constraints for auxiliary QP. */ 620 BooleanType useRelaxation /**< Flag indicating if inactive (constraints') bounds shall be relaxed. */ 621 ); 622 623 624 /** Adds a constraint to active set. 625 * \return SUCCESSFUL_RETURN \n 626 RET_ADDCONSTRAINT_FAILED \n 627 RET_ADDCONSTRAINT_FAILED_INFEASIBILITY \n 628 RET_ENSURELI_FAILED */ 629 virtual returnValue addConstraint( int_t number, /**< Number of constraint to be added to active set. */ 630 SubjectToStatus C_status, /**< Status of new active constraint. */ 631 BooleanType updateCholesky, /**< Flag indicating if Cholesky decomposition shall be updated. */ 632 BooleanType ensureLI = BT_TRUE /**< Ensure linear independence by exchange rules by default. */ 633 ); 634 635 /** Checks if new active constraint to be added is linearly dependent from 636 * from row of the active constraints matrix. 637 * \return RET_LINEARLY_DEPENDENT \n 638 RET_LINEARLY_INDEPENDENT \n 639 RET_INDEXLIST_CORRUPTED */ 640 virtual returnValue addConstraint_checkLI( int_t number /**< Number of constraint to be added to active set. */ 641 ); 642 643 /** Ensures linear independence of constraint matrix when a new constraint is added. 644 * To this end a bound or constraint is removed simultaneously if necessary. 645 * \return SUCCESSFUL_RETURN \n 646 RET_LI_RESOLVED \n 647 RET_ENSURELI_FAILED \n 648 RET_ENSURELI_FAILED_TQ \n 649 RET_ENSURELI_FAILED_NOINDEX \n 650 RET_REMOVE_FROM_ACTIVESET */ 651 virtual returnValue addConstraint_ensureLI( int_t number, /**< Number of constraint to be added to active set. */ 652 SubjectToStatus C_status /**< Status of new active bound. */ 653 ); 654 655 /** Adds a bound to active set. 656 * \return SUCCESSFUL_RETURN \n 657 RET_ADDBOUND_FAILED \n 658 RET_ADDBOUND_FAILED_INFEASIBILITY \n 659 RET_ENSURELI_FAILED */ 660 virtual returnValue addBound( int_t number, /**< Number of bound to be added to active set. */ 661 SubjectToStatus B_status, /**< Status of new active bound. */ 662 BooleanType updateCholesky, /**< Flag indicating if Cholesky decomposition shall be updated. */ 663 BooleanType ensureLI = BT_TRUE /**< Ensure linear independence by exchange rules by default. */ 664 ); 665 666 /** Checks if new active bound to be added is linearly dependent from 667 * from row of the active constraints matrix. 668 * \return RET_LINEARLY_DEPENDENT \n 669 RET_LINEARLY_INDEPENDENT */ 670 virtual returnValue addBound_checkLI( int_t number /**< Number of bound to be added to active set. */ 671 ); 672 673 /** Ensures linear independence of constraint matrix when a new bound is added. 674 * To this end a bound or constraint is removed simultaneously if necessary. 675 * \return SUCCESSFUL_RETURN \n 676 RET_LI_RESOLVED \n 677 RET_ENSURELI_FAILED \n 678 RET_ENSURELI_FAILED_TQ \n 679 RET_ENSURELI_FAILED_NOINDEX \n 680 RET_REMOVE_FROM_ACTIVESET */ 681 virtual returnValue addBound_ensureLI( int_t number, /**< Number of bound to be added to active set. */ 682 SubjectToStatus B_status /**< Status of new active bound. */ 683 ); 684 685 /** Removes a constraint from active set. 686 * \return SUCCESSFUL_RETURN \n 687 RET_CONSTRAINT_NOT_ACTIVE \n 688 RET_REMOVECONSTRAINT_FAILED \n 689 RET_HESSIAN_NOT_SPD */ 690 virtual returnValue removeConstraint( int_t number, /**< Number of constraint to be removed from active set. */ 691 BooleanType updateCholesky, /**< Flag indicating if Cholesky decomposition shall be updated. */ 692 BooleanType allowFlipping = BT_FALSE, /**< Flag indicating if flipping bounds are allowed. */ 693 BooleanType ensureNZC = BT_FALSE /**< Flag indicating if non-zero curvature is ensured by exchange rules. */ 694 ); 695 696 /** Removes a bounds from active set. 697 * \return SUCCESSFUL_RETURN \n 698 RET_BOUND_NOT_ACTIVE \n 699 RET_HESSIAN_NOT_SPD \n 700 RET_REMOVEBOUND_FAILED */ 701 virtual returnValue removeBound( int_t number, /**< Number of bound to be removed from active set. */ 702 BooleanType updateCholesky, /**< Flag indicating if Cholesky decomposition shall be updated. */ 703 BooleanType allowFlipping = BT_FALSE, /**< Flag indicating if flipping bounds are allowed. */ 704 BooleanType ensureNZC = BT_FALSE /**< Flag indicating if non-zero curvature is ensured by exchange rules. */ 705 ); 706 707 708 /** Performs robustified ratio test yield the maximum possible step length 709 * along the homotopy path. 710 * \return SUCCESSFUL_RETURN */ 711 returnValue performPlainRatioTest( int_t nIdx, /**< Number of ratios to be checked. */ 712 const int_t* const idxList, /**< Array containing the indices of all ratios to be checked. */ 713 const real_t* const num, /**< Array containing all numerators for performing the ratio test. */ 714 const real_t* const den, /**< Array containing all denominators for performing the ratio test. */ 715 real_t epsNum, /**< Numerator tolerance. */ 716 real_t epsDen, /**< Denominator tolerance. */ 717 real_t& t, /**< Output: Maximum possible step length along the homotopy path. */ 718 int_t& BC_idx /**< Output: Index of blocking constraint. */ 719 ) const; 720 721 722 /** Ensure non-zero curvature by primal jump. 723 * \return SUCCESSFUL_RETURN \n 724 * RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 725 returnValue ensureNonzeroCurvature( 726 BooleanType removeBoundNotConstraint, /**< SubjectTo to be removed is a bound. */ 727 int_t remIdx, /**< Index of bound/constraint to be removed. */ 728 BooleanType &exchangeHappened, /**< Output: Exchange was necessary to ensure. */ 729 BooleanType &addBoundNotConstraint, /**< SubjectTo to be added is a bound. */ 730 int_t &addIdx, /**< Index of bound/constraint to be added. */ 731 SubjectToStatus &addStatus /**< Status of bound/constraint to be added. */ 732 ); 733 734 735 /** Solves the system Ta = b or T^Ta = b where T is a reverse upper triangular matrix. 736 * \return SUCCESSFUL_RETURN \n 737 RET_DIV_BY_ZERO */ 738 virtual returnValue backsolveT( const real_t* const b, /**< Right hand side vector. */ 739 BooleanType transposed, /**< Indicates if the transposed system shall be solved. */ 740 real_t* const a /**< Output: Solution vector */ 741 ) const; 742 743 744 /** Determines step direction of the shift of the QP data. 745 * \return SUCCESSFUL_RETURN */ 746 returnValue determineDataShift( const real_t* const g_new, /**< New gradient vector. */ 747 const real_t* const lbA_new,/**< New lower constraints' bounds. */ 748 const real_t* const ubA_new,/**< New upper constraints' bounds. */ 749 const real_t* const lb_new, /**< New lower bounds. */ 750 const real_t* const ub_new, /**< New upper bounds. */ 751 real_t* const delta_g, /**< Output: Step direction of gradient vector. */ 752 real_t* const delta_lbA, /**< Output: Step direction of lower constraints' bounds. */ 753 real_t* const delta_ubA, /**< Output: Step direction of upper constraints' bounds. */ 754 real_t* const delta_lb, /**< Output: Step direction of lower bounds. */ 755 real_t* const delta_ub, /**< Output: Step direction of upper bounds. */ 756 BooleanType& Delta_bC_isZero,/**< Output: Indicates if active constraints' bounds are to be shifted. */ 757 BooleanType& Delta_bB_isZero/**< Output: Indicates if active bounds are to be shifted. */ 758 ); 759 760 /** Determines step direction of the homotopy path. 761 * \return SUCCESSFUL_RETURN \n 762 RET_STEPDIRECTION_FAILED_TQ \n 763 RET_STEPDIRECTION_FAILED_CHOLESKY */ 764 virtual returnValue determineStepDirection( const real_t* const delta_g, /**< Step direction of gradient vector. */ 765 const real_t* const delta_lbA, /**< Step direction of lower constraints' bounds. */ 766 const real_t* const delta_ubA, /**< Step direction of upper constraints' bounds. */ 767 const real_t* const delta_lb, /**< Step direction of lower bounds. */ 768 const real_t* const delta_ub, /**< Step direction of upper bounds. */ 769 BooleanType Delta_bC_isZero, /**< Indicates if active constraints' bounds are to be shifted. */ 770 BooleanType Delta_bB_isZero, /**< Indicates if active bounds are to be shifted. */ 771 real_t* const delta_xFX, /**< Output: Primal homotopy step direction of fixed variables. */ 772 real_t* const delta_xFR, /**< Output: Primal homotopy step direction of free variables. */ 773 real_t* const delta_yAC, /**< Output: Dual homotopy step direction of active constraints' multiplier. */ 774 real_t* const delta_yFX /**< Output: Dual homotopy step direction of fixed variables' multiplier. */ 775 ); 776 777 /** Determines the maximum possible step length along the homotopy path 778 * and performs this step (without changing working set). 779 * \return SUCCESSFUL_RETURN \n 780 * RET_ERROR_IN_CONSTRAINTPRODUCT \n 781 * RET_QP_INFEASIBLE */ 782 returnValue performStep( const real_t* const delta_g, /**< Step direction of gradient. */ 783 const real_t* const delta_lbA, /**< Step direction of lower constraints' bounds. */ 784 const real_t* const delta_ubA, /**< Step direction of upper constraints' bounds. */ 785 const real_t* const delta_lb, /**< Step direction of lower bounds. */ 786 const real_t* const delta_ub, /**< Step direction of upper bounds. */ 787 const real_t* const delta_xFX, /**< Primal homotopy step direction of fixed variables. */ 788 const real_t* const delta_xFR, /**< Primal homotopy step direction of free variables. */ 789 const real_t* const delta_yAC, /**< Dual homotopy step direction of active constraints' multiplier. */ 790 const real_t* const delta_yFX, /**< Dual homotopy step direction of fixed variables' multiplier. */ 791 int_t& BC_idx, /**< Output: Index of blocking constraint. */ 792 SubjectToStatus& BC_status, /**< Output: Status of blocking constraint. */ 793 BooleanType& BC_isBound /**< Output: Indicates if blocking constraint is a bound. */ 794 ); 795 796 /** Updates the active set. 797 * \return SUCCESSFUL_RETURN \n 798 RET_REMOVE_FROM_ACTIVESET_FAILED \n 799 RET_ADD_TO_ACTIVESET_FAILED */ 800 returnValue changeActiveSet( int_t BC_idx, /**< Index of blocking constraint. */ 801 SubjectToStatus BC_status, /**< Status of blocking constraint. */ 802 BooleanType BC_isBound /**< Indicates if blocking constraint is a bound. */ 803 ); 804 805 806 /** Compute relative length of homotopy in data space for termination 807 * criterion. 808 * \return Relative length in data space. */ 809 real_t getRelativeHomotopyLength( const real_t* const g_new, /**< Final gradient. */ 810 const real_t* const lb_new, /**< Final lower variable bounds. */ 811 const real_t* const ub_new, /**< Final upper variable bounds. */ 812 const real_t* const lbA_new, /**< Final lower constraint bounds. */ 813 const real_t* const ubA_new /**< Final upper constraint bounds. */ 814 ); 815 816 817 /** Ramping Strategy to avoid ties. Modifies homotopy start without 818 * changing current active set. 819 * \return SUCCESSFUL_RETURN */ 820 virtual returnValue performRamping( ); 821 822 823 /** ... */ 824 returnValue updateFarBounds( real_t curFarBound, /**< ... */ 825 int_t nRamp, /**< ... */ 826 const real_t* const lb_new, /**< ... */ 827 real_t* const lb_new_far, /**< ... */ 828 const real_t* const ub_new, /**< ... */ 829 real_t* const ub_new_far, /**< ... */ 830 const real_t* const lbA_new, /**< ... */ 831 real_t* const lbA_new_far, /**< ... */ 832 const real_t* const ubA_new, /**< ... */ 833 real_t* const ubA_new_far /**< ... */ 834 ) const; 835 836 837 /** Drift correction at end of each active set iteration 838 * \return SUCCESSFUL_RETURN */ 839 virtual returnValue performDriftCorrection( ); 840 841 842 /** Updates QP vectors, working sets and internal data structures in order to 843 start from an optimal solution corresponding to initial guesses of the working 844 set for bounds and constraints. 845 * \return SUCCESSFUL_RETURN \n 846 * RET_SETUP_AUXILIARYQP_FAILED \n 847 RET_INVALID_ARGUMENTS */ 848 using QProblemB::setupAuxiliaryQP; 849 virtual returnValue setupAuxiliaryQP( const Bounds* const guessedBounds, /**< Initial guess for working set of bounds. */ 850 const Constraints* const guessedConstraints /**< Initial guess for working set of constraints. */ 851 ); 852 853 /** Determines if it is more efficient to refactorise the matrices when 854 * hotstarting or not (i.e. better to update the existing factorisations). 855 * \return BT_TRUE iff matrices shall be refactorised afresh 856 */ 857 BooleanType shallRefactorise( const Bounds* const guessedBounds, /**< Guessed new working set of bounds. */ 858 const Constraints* const guessedConstraints /**< Guessed new working set of constraints. */ 859 ) const; 860 861 /** Sets up internal QP data. 862 * \return SUCCESSFUL_RETURN \n 863 RET_INVALID_ARGUMENTS \n 864 RET_UNKNONW_BUG */ 865 returnValue setupQPdata( SymmetricMatrix *_H, /**< Hessian matrix. \n 866 If Hessian matrix is trivial,a NULL pointer can be passed. */ 867 const real_t* const _g, /**< Gradient vector. */ 868 Matrix *_A, /**< Constraint matrix. */ 869 const real_t* const _lb, /**< Lower bound vector (on variables). \n 870 If no lower bounds exist, a NULL pointer can be passed. */ 871 const real_t* const _ub, /**< Upper bound vector (on variables). \n 872 If no upper bounds exist, a NULL pointer can be passed. */ 873 const real_t* const _lbA, /**< Lower constraints' bound vector. \n 874 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 875 const real_t* const _ubA /**< Upper constraints' bound vector. \n 876 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 877 ); 878 879 880 /** Sets up dense internal QP data. If the current Hessian is trivial 881 * (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, 882 * memory for Hessian is allocated and it is set to the given one. 883 * \return SUCCESSFUL_RETURN \n 884 RET_INVALID_ARGUMENTS \n 885 RET_UNKNONW_BUG */ 886 returnValue setupQPdata( const real_t* const _H, /**< Hessian matrix. \n 887 If Hessian matrix is trivial,a NULL pointer can be passed. */ 888 const real_t* const _g, /**< Gradient vector. */ 889 const real_t* const _A, /**< Constraint matrix. */ 890 const real_t* const _lb, /**< Lower bound vector (on variables). \n 891 If no lower bounds exist, a NULL pointer can be passed. */ 892 const real_t* const _ub, /**< Upper bound vector (on variables). \n 893 If no upper bounds exist, a NULL pointer can be passed. */ 894 const real_t* const _lbA, /**< Lower constraints' bound vector. \n 895 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 896 const real_t* const _ubA /**< Upper constraints' bound vector. \n 897 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 898 ); 899 900 /** Sets up internal QP data by loading it from files. If the current Hessian 901 * is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, 902 * memory for Hessian is allocated and it is set to the given one. 903 * \return SUCCESSFUL_RETURN \n 904 RET_UNABLE_TO_OPEN_FILE \n 905 RET_UNABLE_TO_READ_FILE \n 906 RET_INVALID_ARGUMENTS \n 907 RET_UNKNONW_BUG */ 908 returnValue setupQPdataFromFile( const char* const H_file, /**< Name of file where Hessian matrix, of neighbouring QP to be solved, is stored. \n 909 If Hessian matrix is trivial,a NULL pointer can be passed. */ 910 const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 911 const char* const A_file, /**< Name of file where constraint matrix, of neighbouring QP to be solved, is stored. */ 912 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 913 If no lower bounds exist, a NULL pointer can be passed. */ 914 const char* const ub_file, /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 915 If no upper bounds exist, a NULL pointer can be passed. */ 916 const char* const lbA_file, /**< Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. \n 917 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 918 const char* const ubA_file /**< Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. \n 919 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 920 ); 921 922 /** Loads new QP vectors from files (internal members are not affected!). 923 * \return SUCCESSFUL_RETURN \n 924 RET_UNABLE_TO_OPEN_FILE \n 925 RET_UNABLE_TO_READ_FILE \n 926 RET_INVALID_ARGUMENTS */ 927 returnValue loadQPvectorsFromFile( const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 928 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 929 If no lower bounds exist, a NULL pointer can be passed. */ 930 const char* const ub_file, /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 931 If no upper bounds exist, a NULL pointer can be passed. */ 932 const char* const lbA_file, /**< Name of file where lower constraints' bounds, of neighbouring QP to be solved, is stored. \n 933 If no lower constraints' bounds exist, a NULL pointer can be passed. */ 934 const char* const ubA_file, /**< Name of file where upper constraints' bounds, of neighbouring QP to be solved, is stored. \n 935 If no upper constraints' bounds exist, a NULL pointer can be passed. */ 936 real_t* const g_new, /**< Output: Gradient of neighbouring QP to be solved. */ 937 real_t* const lb_new, /**< Output: Lower bounds of neighbouring QP to be solved */ 938 real_t* const ub_new, /**< Output: Upper bounds of neighbouring QP to be solved */ 939 real_t* const lbA_new, /**< Output: Lower constraints' bounds of neighbouring QP to be solved */ 940 real_t* const ubA_new /**< Output: Upper constraints' bounds of neighbouring QP to be solved */ 941 ) const; 942 943 944 /** Prints concise information on the current iteration. 945 * \return SUCCESSFUL_RETURN \n */ 946 returnValue printIteration( int_t iter, /**< Number of current iteration. */ 947 int_t BC_idx, /**< Index of blocking constraint. */ 948 SubjectToStatus BC_status, /**< Status of blocking constraint. */ 949 BooleanType BC_isBound, /**< Indicates if blocking constraint is a bound. */ 950 real_t homotopyLength, /**< Current homotopy distance. */ 951 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 952 ); 953 954 955 /** Sets constraint matrix of the QP. \n 956 Note: Also internal vector Ax is recomputed! 957 * \return SUCCESSFUL_RETURN \n 958 * RET_INVALID_ARGUMENTS */ 959 inline returnValue setA( Matrix *A_new /**< New constraint matrix (a shallow copy is made). */ 960 ); 961 962 /** Sets dense constraint matrix of the QP. \n 963 Note: Also internal vector Ax is recomputed! 964 * \return SUCCESSFUL_RETURN \n 965 * RET_INVALID_ARGUMENTS */ 966 inline returnValue setA( const real_t* const A_new /**< New dense constraint matrix (with correct dimension!), a shallow copy is made. */ 967 ); 968 969 970 /** Sets constraints' lower bound vector of the QP. 971 * \return SUCCESSFUL_RETURN \n 972 * RET_QPOBJECT_NOT_SETUP */ 973 inline returnValue setLBA( const real_t* const lbA_new /**< New constraints' lower bound vector (with correct dimension!). */ 974 ); 975 976 /** Changes single entry of lower constraints' bound vector of the QP. 977 * \return SUCCESSFUL_RETURN \n 978 * RET_QPOBJECT_NOT_SETUP \n 979 * RET_INDEX_OUT_OF_BOUNDS */ 980 inline returnValue setLBA( int_t number, /**< Number of entry to be changed. */ 981 real_t value /**< New value for entry of lower constraints' bound vector (with correct dimension!). */ 982 ); 983 984 /** Sets constraints' upper bound vector of the QP. 985 * \return SUCCESSFUL_RETURN \n 986 * RET_QPOBJECT_NOT_SETUP */ 987 inline returnValue setUBA( const real_t* const ubA_new /**< New constraints' upper bound vector (with correct dimension!). */ 988 ); 989 990 /** Changes single entry of upper constraints' bound vector of the QP. 991 * \return SUCCESSFUL_RETURN \n 992 * RET_QPOBJECT_NOT_SETUP \n 993 * RET_INDEX_OUT_OF_BOUNDS */ 994 inline returnValue setUBA( int_t number, /**< Number of entry to be changed. */ 995 real_t value /**< New value for entry of upper constraints' bound vector (with correct dimension!). */ 996 ); 997 998 999 /** Drops the blocking bound/constraint that led to infeasibility, or finds another 1000 * bound/constraint to drop according to drop priorities. 1001 * \return SUCCESSFUL_RETURN \n 1002 */ 1003 returnValue dropInfeasibles ( int_t BC_number, /**< Number of the bound or constraint to be added. */ 1004 SubjectToStatus BC_status, /**< New status of the bound or constraint to be added. */ 1005 BooleanType BC_isBound, /**< Whether a bound or a constraint is to be added. */ 1006 real_t *xiB, /**< (not yet documented) */ 1007 real_t *xiC /**< (not yet documented) */ 1008 ); 1009 1010 /** Decides if lower bounds are smaller than upper bounds 1011 * 1012 * \return SUCCESSFUL_RETURN \n 1013 * RET_QP_INFEASIBLE */ 1014 1015 returnValue areBoundsConsistent(const real_t* const lb, /**< Vector of lower bounds*/ 1016 const real_t* const ub, /**< Vector of upper bounds*/ 1017 const real_t* const lbA, /**< Vector of lower constraints*/ 1018 const real_t* const ubA /**< Vector of upper constraints*/ 1019 ) const; 1020 1021 1022 public: 1023 /** ... 1024 * \return SUCCESSFUL_RETURN \n 1025 RET_UNABLE_TO_OPEN_FILE */ 1026 returnValue writeQpDataIntoMatFile( const char* const filename /**< Mat file name. */ 1027 ) const; 1028 1029 /** ... 1030 * \return SUCCESSFUL_RETURN \n 1031 RET_UNABLE_TO_OPEN_FILE */ 1032 returnValue writeQpWorkspaceIntoMatFile( const char* const filename /**< Mat file name. */ 1033 ); 1034 1035 1036 1037 /* 1038 * PROTECTED MEMBER VARIABLES 1039 */ 1040 protected: 1041 BooleanType freeConstraintMatrix; /**< Flag indicating whether the constraint matrix needs to be de-allocated. */ 1042 Matrix* A; /**< Constraint matrix. */ 1043 1044 real_t* lbA; /**< Lower constraints' bound vector. */ 1045 real_t* ubA; /**< Upper constraints' bound vector. */ 1046 1047 Constraints constraints; /**< Data structure for problem's constraints. */ 1048 1049 real_t* T; /**< Reverse triangular matrix, A = [0 T]*Q'. */ 1050 real_t* Q; /**< Orthonormal quadratic matrix, A = [0 T]*Q'. */ 1051 int_t sizeT; /**< Matrix T is stored in a (sizeT x sizeT) array. */ 1052 1053 real_t* Ax; /**< Stores the current A*x \n 1054 * (for increased efficiency only). */ 1055 real_t* Ax_l; /**< Stores the current distance to lower constraints' bounds A*x-lbA \n 1056 * (for increased efficiency only). */ 1057 real_t* Ax_u; /**< Stores the current distance to lower constraints' bounds ubA-A*x \n 1058 * (for increased efficiency only). */ 1059 1060 ConstraintProduct* constraintProduct; /**< Pointer to user-defined constraint product function. */ 1061 1062 real_t* tempA; /**< Temporary for determineStepDirection. */ 1063 real_t* tempB; /**< Temporary for determineStepDirection. */ 1064 real_t* ZFR_delta_xFRz; /**< Temporary for determineStepDirection. */ 1065 real_t* delta_xFRy; /**< Temporary for determineStepDirection. */ 1066 real_t* delta_xFRz; /**< Temporary for determineStepDirection. */ 1067 real_t* delta_yAC_TMP; /**< Temporary for determineStepDirection. */ 1068 1069 real_t* tempC; /**< Temporary for constraint types. */ 1070 }; 1071 1072 1073 END_NAMESPACE_QPOASES 1074 1075 #include <qpOASES/QProblem.ipp> 1076 1077 #endif /* QPOASES_QPROBLEM_HPP */ 1078 1079 1080 /* 1081 * end of file 1082 */ 1083