1 /* 2 * This file is part of qpOASES. 3 * 4 * qpOASES -- An Implementation of the Online Active Set Strategy. 5 * Copyright (C) 2007-2015 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/QProblemB.hpp 27 * \author Hans Joachim Ferreau, Andreas Potschka, Christian Kirches 28 * \version 3.2 29 * \date 2007-2015 30 * 31 * Declaration of the QProblemB class which is able to use the newly 32 * developed online active set strategy for parametric quadratic programming 33 * for problems with (simple) bounds only. 34 */ 35 36 37 38 #ifndef QPOASES_QPROBLEMB_HPP 39 #define QPOASES_QPROBLEMB_HPP 40 41 42 #include <qpOASES/Flipper.hpp> 43 #include <qpOASES/Options.hpp> 44 #include <qpOASES/Matrices.hpp> 45 46 47 BEGIN_NAMESPACE_QPOASES 48 49 50 class SolutionAnalysis; 51 52 /** 53 * \brief Implements the online active set strategy for box-constrained QPs. 54 * 55 * Class for setting up and solving quadratic programs with bounds (= box constraints) only. 56 * The main feature is the possibily to use the newly developed online active set strategy 57 * for parametric quadratic programming. 58 * 59 * \author Hans Joachim Ferreau, Andreas Potschka, Christian Kirches 60 * \version 3.2 61 * \date 2007-2015 62 */ 63 class QProblemB 64 { 65 /* allow SolutionAnalysis class to access private members */ 66 friend class SolutionAnalysis; 67 68 /* 69 * PUBLIC MEMBER FUNCTIONS 70 */ 71 public: 72 /** Default constructor. */ 73 QProblemB( ); 74 75 /** Constructor which takes the QP dimension and Hessian type 76 * information. If the Hessian is the zero (i.e. HST_ZERO) or the 77 * identity matrix (i.e. HST_IDENTITY), respectively, no memory 78 * is allocated for it and a NULL pointer can be passed for it 79 * to the init() functions. */ 80 QProblemB( int_t _nV, /**< Number of variables. */ 81 HessianType _hessianType = HST_UNKNOWN /**< Type of Hessian matrix. */ 82 ); 83 84 /** Copy constructor (deep copy). */ 85 QProblemB( const QProblemB& rhs /**< Rhs object. */ 86 ); 87 88 /** Destructor. */ 89 virtual ~QProblemB( ); 90 91 /** Assignment operator (deep copy). */ 92 virtual QProblemB& operator=( const QProblemB& 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 simply bounded 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 empty (or all implicit equality bounds), \n 105 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB by "clipping", \n 106 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB from yOpt != 0, \n 107 * 4. 0, 0, gB: starts with xOpt = 0, yOpt = 0 and gB, \n 108 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB from yOpt != 0, \n 109 * 6. xOpt, 0, gB: starts with xOpt, yOpt = 0 and gB, \n 110 * 7. xOpt, yOpt, gB: starts with xOpt, yOpt and gB (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_HOTSTART \n 118 RET_INIT_FAILED_INFEASIBILITY \n 119 RET_INIT_FAILED_UNBOUNDEDNESS \n 120 RET_MAX_NWSR_REACHED \n 121 RET_INVALID_ARGUMENTS */ 122 returnValue init( SymmetricMatrix *_H, /**< Hessian matrix (a shallow copy is made). */ 123 const real_t* const _g, /**< Gradient vector. */ 124 const real_t* const _lb, /**< Lower bounds (on variables). \n 125 If no lower bounds exist, a NULL pointer can be passed. */ 126 const real_t* const _ub, /**< Upper bounds (on variables). \n 127 If no upper bounds exist, a NULL pointer can be passed. */ 128 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. \n 129 Output: Number of performed working set recalculations. */ 130 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 131 Output: CPU time spent for QP initialisation (if pointer passed). */ 132 const real_t* const xOpt = 0, /**< Optimal primal solution vector. A NULL pointer can be passed. \n 133 (If a null pointer is passed, the old primal solution is kept!) */ 134 const real_t* const yOpt = 0, /**< Optimal dual solution vector. A NULL pointer can be passed. \n 135 (If a null pointer is passed, the old dual solution is kept!) */ 136 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 137 (If a null pointer is passed, all bounds are assumed inactive!) */ 138 const real_t* const _R = 0 /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. 139 The Cholesky factor must be stored in a real_t array of size nV*nV 140 in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! \n 141 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 142 ); 143 144 /** Initialises a simply bounded QP problem with given QP data and tries to solve it 145 * using at most nWSR iterations. Depending on the parameter constellation it: \n 146 * 1. 0, 0, 0 : starts with xOpt = 0, yOpt = 0 and gB empty (or all implicit equality bounds), \n 147 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB by "clipping", \n 148 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB from yOpt != 0, \n 149 * 4. 0, 0, gB: starts with xOpt = 0, yOpt = 0 and gB, \n 150 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB from yOpt != 0, \n 151 * 6. xOpt, 0, gB: starts with xOpt, yOpt = 0 and gB, \n 152 * 7. xOpt, yOpt, gB: starts with xOpt, yOpt and gB (assume them to be consistent!) 153 * 154 * Note: This function internally calls solveInitialQP for initialisation! 155 * 156 * \return SUCCESSFUL_RETURN \n 157 RET_INIT_FAILED \n 158 RET_INIT_FAILED_CHOLESKY \n 159 RET_INIT_FAILED_HOTSTART \n 160 RET_INIT_FAILED_INFEASIBILITY \n 161 RET_INIT_FAILED_UNBOUNDEDNESS \n 162 RET_MAX_NWSR_REACHED \n 163 RET_INVALID_ARGUMENTS */ 164 returnValue init( const real_t* const _H, /**< Hessian matrix (a shallow copy is made). \n 165 If Hessian matrix is trivial, a NULL pointer can be passed. */ 166 const real_t* const _g, /**< Gradient vector. */ 167 const real_t* const _lb, /**< Lower bounds (on variables). \n 168 If no lower bounds exist, a NULL pointer can be passed. */ 169 const real_t* const _ub, /**< Upper bounds (on variables). \n 170 If no upper bounds exist, a NULL pointer can be passed. */ 171 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. \n 172 Output: Number of performed working set recalculations. */ 173 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 174 Output: CPU time spent for QP initialisation (if pointer passed). */ 175 const real_t* const xOpt = 0, /**< Optimal primal solution vector. A NULL pointer can be passed. \n 176 (If a null pointer is passed, the old primal solution is kept!) */ 177 const real_t* const yOpt = 0, /**< Optimal dual solution vector. A NULL pointer can be passed. \n 178 (If a null pointer is passed, the old dual solution is kept!) */ 179 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 180 (If a null pointer is passed, all bounds are assumed inactive!) */ 181 const real_t* const _R = 0 /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. 182 The Cholesky factor must be stored in a real_t array of size nV*nV 183 in row-major format. Note: Only used if xOpt/yOpt and gB are NULL! \n 184 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 185 ); 186 187 /** Initialises a simply bounded QP problem with given QP data to be read from files and solves it 188 * using at most nWSR iterations. Depending on the parameter constellation it: \n 189 * 1. 0, 0, 0 : starts with xOpt = 0, yOpt = 0 and gB empty (or all implicit equality bounds), \n 190 * 2. xOpt, 0, 0 : starts with xOpt, yOpt = 0 and obtain gB by "clipping", \n 191 * 3. 0, yOpt, 0 : starts with xOpt = 0, yOpt and obtain gB from yOpt != 0, \n 192 * 4. 0, 0, gB: starts with xOpt = 0, yOpt = 0 and gB, \n 193 * 5. xOpt, yOpt, 0 : starts with xOpt, yOpt and obtain gB from yOpt != 0, \n 194 * 6. xOpt, 0, gB: starts with xOpt, yOpt = 0 and gB, \n 195 * 7. xOpt, yOpt, gB: starts with xOpt, yOpt and gB (assume them to be consistent!) 196 * 197 * Note: This function internally calls solveInitialQP for initialisation! 198 * 199 * \return SUCCESSFUL_RETURN \n 200 RET_INIT_FAILED \n 201 RET_INIT_FAILED_CHOLESKY \n 202 RET_INIT_FAILED_HOTSTART \n 203 RET_INIT_FAILED_INFEASIBILITY \n 204 RET_INIT_FAILED_UNBOUNDEDNESS \n 205 RET_MAX_NWSR_REACHED \n 206 RET_UNABLE_TO_READ_FILE */ 207 returnValue init( const char* const H_file, /**< Name of file where Hessian matrix is stored. \n 208 If Hessian matrix is trivial, a NULL pointer can be passed. */ 209 const char* const g_file, /**< Name of file where gradient vector is stored. */ 210 const char* const lb_file, /**< Name of file where lower bound vector. \n 211 If no lower bounds exist, a NULL pointer can be passed. */ 212 const char* const ub_file, /**< Name of file where upper bound vector. \n 213 If no upper bounds exist, a NULL pointer can be passed. */ 214 int_t& nWSR, /**< Input: Maximum number of working set recalculations when using initial homotopy. \n 215 Output: Number of performed working set recalculations. */ 216 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP initialisation. \n 217 Output: CPU time spent for QP initialisation (if pointer passed). */ 218 const real_t* const xOpt = 0, /**< Optimal primal solution vector. A NULL pointer can be passed. \n 219 (If a null pointer is passed, the old primal solution is kept!) */ 220 const real_t* const yOpt = 0, /**< Optimal dual solution vector. A NULL pointer can be passed. \n 221 (If a null pointer is passed, the old dual solution is kept!) */ 222 const Bounds* const guessedBounds = 0, /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 223 (If a null pointer is passed, all bounds are assumed inactive!) */ 224 const char* const R_file = 0 /**< Name of the file where a pre-computed (upper triangular) Cholesky factor 225 of the Hessian matrix is stored. \n 226 (If a null pointer is passed, Cholesky decomposition is computed internally!) */ 227 ); 228 229 230 /** Solves an initialised QP sequence using the online active set strategy. 231 * By default, QP solution is started from previous solution. If a guess 232 * for the working set is provided, an initialised homotopy is performed. 233 * 234 * Note: This function internally calls solveQP/solveRegularisedQP 235 * for solving an initialised QP! 236 * 237 * \return SUCCESSFUL_RETURN \n 238 RET_MAX_NWSR_REACHED \n 239 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 240 RET_HOTSTART_FAILED \n 241 RET_SHIFT_DETERMINATION_FAILED \n 242 RET_STEPDIRECTION_DETERMINATION_FAILED \n 243 RET_STEPLENGTH_DETERMINATION_FAILED \n 244 RET_HOMOTOPY_STEP_FAILED \n 245 RET_HOTSTART_STOPPED_INFEASIBILITY \n 246 RET_HOTSTART_STOPPED_UNBOUNDEDNESS \n 247 RET_SETUP_AUXILIARYQP_FAILED */ 248 returnValue hotstart( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 249 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 250 If no lower bounds exist, a NULL pointer can be passed. */ 251 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 252 If no upper bounds exist, a NULL pointer can be passed. */ 253 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 254 Output: Number of performed working set recalculations. */ 255 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP solution. \n 256 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 257 const Bounds* const guessedBounds = 0 /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 258 (If a null pointer is passed, the previous working set is kept!) */ 259 ); 260 261 /** Solves an initialised QP sequence using the online active set strategy, 262 * where QP data is read from files. 263 * By default, QP solution is started from previous solution. If a guess 264 * for the working set is provided, an initialised homotopy is performed. 265 * 266 * Note: This function internally calls solveQP/solveRegularisedQP 267 * for solving an initialised QP! 268 * 269 * \return SUCCESSFUL_RETURN \n 270 RET_MAX_NWSR_REACHED \n 271 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 272 RET_HOTSTART_FAILED \n 273 RET_SHIFT_DETERMINATION_FAILED \n 274 RET_STEPDIRECTION_DETERMINATION_FAILED \n 275 RET_STEPLENGTH_DETERMINATION_FAILED \n 276 RET_HOMOTOPY_STEP_FAILED \n 277 RET_HOTSTART_STOPPED_INFEASIBILITY \n 278 RET_HOTSTART_STOPPED_UNBOUNDEDNESS \n 279 RET_UNABLE_TO_READ_FILE \n 280 RET_SETUP_AUXILIARYQP_FAILED \n 281 RET_INVALID_ARGUMENTS */ 282 returnValue hotstart( const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 283 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 284 If no lower bounds exist, a NULL pointer can be passed. */ 285 const char* const ub_file, /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 286 If no upper bounds exist, a NULL pointer can be passed. */ 287 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 288 Output: Number of performed working set recalculations. */ 289 real_t* const cputime = 0, /**< Input: Maximum CPU time allowed for QP solution. \n 290 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 291 const Bounds* const guessedBounds = 0 /**< Optimal working set of bounds for solution (xOpt,yOpt). \n 292 (If a null pointer is passed, the previous working set is kept!) */ 293 ); 294 295 296 /** Writes a vector with the state of the working set 297 * \return SUCCESSFUL_RETURN \n 298 * RET_INVALID_ARGUMENTS */ 299 virtual returnValue getWorkingSet( real_t* workingSet /** Output: array containing state of the working set. */ 300 ); 301 302 /** Writes a vector with the state of the working set of bounds 303 * \return SUCCESSFUL_RETURN \n 304 * RET_INVALID_ARGUMENTS */ 305 virtual returnValue getWorkingSetBounds( real_t* workingSetB /** Output: array containing state of the working set of bounds. */ 306 ); 307 308 /** Writes a vector with the state of the working set of constraints 309 * \return SUCCESSFUL_RETURN \n 310 * RET_INVALID_ARGUMENTS */ 311 virtual returnValue getWorkingSetConstraints( real_t* workingSetC /** Output: array containing state of the working set of constraints. */ 312 ); 313 314 315 /** Returns current bounds object of the QP (deep copy). 316 * \return SUCCESSFUL_RETURN \n 317 RET_QPOBJECT_NOT_SETUP */ 318 inline returnValue getBounds( Bounds& _bounds /** Output: Bounds object. */ 319 ) const; 320 321 322 /** Returns the number of variables. 323 * \return Number of variables. */ 324 inline int_t getNV( ) const; 325 326 /** Returns the number of free variables. 327 * \return Number of free variables. */ 328 inline int_t getNFR( ) const; 329 330 /** Returns the number of fixed variables. 331 * \return Number of fixed variables. */ 332 inline int_t getNFX( ) const; 333 334 /** Returns the number of implicitly fixed variables. 335 * \return Number of implicitly fixed variables. */ 336 inline int_t getNFV( ) const; 337 338 /** Returns the dimension of null space. 339 * \return Dimension of null space. */ 340 virtual int_t getNZ( ) const; 341 342 343 /** Returns the optimal objective function value. 344 * \return finite value: Optimal objective function value (QP was solved) \n 345 +infinity: QP was not yet solved */ 346 real_t getObjVal( ) const; 347 348 /** Returns the objective function value at an arbitrary point x. 349 * \return Objective function value at point x */ 350 real_t getObjVal( const real_t* const _x /**< Point at which the objective function shall be evaluated. */ 351 ) const; 352 353 /** Returns the primal solution vector. 354 * \return SUCCESSFUL_RETURN \n 355 RET_QP_NOT_SOLVED */ 356 returnValue getPrimalSolution( real_t* const xOpt /**< Output: Primal solution vector (if QP has been solved). */ 357 ) const; 358 359 /** Returns the dual solution vector. 360 * \return SUCCESSFUL_RETURN \n 361 RET_QP_NOT_SOLVED */ 362 virtual returnValue getDualSolution( real_t* const yOpt /**< Output: Dual solution vector (if QP has been solved). */ 363 ) const; 364 365 366 /** Returns status of the solution process. 367 * \return Status of solution process. */ 368 inline QProblemStatus getStatus( ) const; 369 370 371 /** Returns if the QProblem object is initialised. 372 * \return BT_TRUE: QProblemB initialised \n 373 BT_FALSE: QProblemB not initialised */ 374 inline BooleanType isInitialised( ) const; 375 376 /** Returns if the QP has been solved. 377 * \return BT_TRUE: QProblemB solved \n 378 BT_FALSE: QProblemB not solved */ 379 inline BooleanType isSolved( ) const; 380 381 /** Returns if the QP is infeasible. 382 * \return BT_TRUE: QP infeasible \n 383 BT_FALSE: QP feasible (or not known to be infeasible!) */ 384 inline BooleanType isInfeasible( ) const; 385 386 /** Returns if the QP is unbounded. 387 * \return BT_TRUE: QP unbounded \n 388 BT_FALSE: QP unbounded (or not known to be unbounded!) */ 389 inline BooleanType isUnbounded( ) const; 390 391 392 /** Returns Hessian type flag (type is not determined due to this call!). 393 * \return Hessian type. */ 394 inline HessianType getHessianType( ) const; 395 396 /** Changes the print level. 397 * \return SUCCESSFUL_RETURN */ 398 inline returnValue setHessianType( HessianType _hessianType /**< New Hessian type. */ 399 ); 400 401 /** Returns if the QP has been internally regularised. 402 * \return BT_TRUE: Hessian is internally regularised for QP solution \n 403 BT_FALSE: No internal Hessian regularisation is used for QP solution */ 404 inline BooleanType usingRegularisation( ) const; 405 406 /** Returns current options struct. 407 * \return Current options struct. */ 408 inline Options getOptions( ) const; 409 410 /** Overrides current options with given ones. 411 * \return SUCCESSFUL_RETURN */ 412 inline returnValue setOptions( const Options& _options /**< New options. */ 413 ); 414 415 /** Returns the print level. 416 * \return Print level. */ 417 inline PrintLevel getPrintLevel( ) const; 418 419 /** Changes the print level. 420 * \return SUCCESSFUL_RETURN */ 421 returnValue setPrintLevel( PrintLevel _printlevel /**< New print level. */ 422 ); 423 424 425 /** Returns the current number of QP problems solved. 426 * \return Number of QP problems solved. */ 427 inline uint_t getCount( ) const; 428 429 /** Resets QP problem counter (to zero). 430 * \return SUCCESSFUL_RETURN. */ 431 inline returnValue resetCounter( ); 432 433 434 /** Prints concise list of properties of the current QP. 435 * \return SUCCESSFUL_RETURN \n */ 436 virtual returnValue printProperties( ); 437 438 /** Prints a list of all options and their current values. 439 * \return SUCCESSFUL_RETURN \n */ 440 returnValue printOptions( ) const; 441 442 443 /* 444 * PROTECTED MEMBER FUNCTIONS 445 */ 446 protected: 447 /** Frees all allocated memory. 448 * \return SUCCESSFUL_RETURN */ 449 returnValue clear( ); 450 451 /** Copies all members from given rhs object. 452 * \return SUCCESSFUL_RETURN */ 453 returnValue copy( const QProblemB& rhs /**< Rhs object. */ 454 ); 455 456 /** If Hessian type has been set by the user, nothing is done. 457 * Otherwise the Hessian type is set to HST_IDENTITY, HST_ZERO, or 458 * HST_POSDEF (default), respectively. 459 * \return SUCCESSFUL_RETURN \n 460 RET_HESSIAN_INDEFINITE */ 461 returnValue determineHessianType( ); 462 463 /** Determines type of existing constraints and bounds (i.e. implicitly fixed, unbounded etc.). 464 * \return SUCCESSFUL_RETURN \n 465 RET_SETUPSUBJECTTOTYPE_FAILED */ 466 virtual returnValue setupSubjectToType( ); 467 468 /** Determines type of new constraints and bounds (i.e. implicitly fixed, unbounded etc.). 469 * \return SUCCESSFUL_RETURN \n 470 RET_SETUPSUBJECTTOTYPE_FAILED */ 471 virtual returnValue setupSubjectToType( const real_t* const lb_new, /**< New lower bounds. */ 472 const real_t* const ub_new /**< New upper bounds. */ 473 ); 474 475 /** Computes the Cholesky decomposition of the (simply projected) Hessian 476 * (i.e. R^T*R = Z^T*H*Z). It only works in the case where Z is a simple 477 * projection matrix! 478 * Note: If Hessian turns out not to be positive definite, the Hessian type 479 * is set to HST_SEMIDEF accordingly. 480 * \return SUCCESSFUL_RETURN \n 481 * RET_HESSIAN_NOT_SPD \n 482 * RET_INDEXLIST_CORRUPTED */ 483 virtual returnValue computeCholesky( ); 484 485 486 /** Computes initial Cholesky decomposition of the (simply projected) Hessian 487 * making use of the function computeCholesky(). 488 * \return SUCCESSFUL_RETURN \n 489 * RET_HESSIAN_NOT_SPD \n 490 * RET_INDEXLIST_CORRUPTED */ 491 virtual returnValue setupInitialCholesky( ); 492 493 /** Obtains the desired working set for the auxiliary initial QP in 494 * accordance with the user specifications 495 * \return SUCCESSFUL_RETURN \n 496 RET_OBTAINING_WORKINGSET_FAILED \n 497 RET_INVALID_ARGUMENTS */ 498 returnValue obtainAuxiliaryWorkingSet( const real_t* const xOpt, /**< Optimal primal solution vector. 499 * If a NULL pointer is passed, all entries are assumed to be zero. */ 500 const real_t* const yOpt, /**< Optimal dual solution vector. 501 * If a NULL pointer is passed, all entries are assumed to be zero. */ 502 const Bounds* const guessedBounds, /**< Guessed working set for solution (xOpt,yOpt). */ 503 Bounds* auxiliaryBounds /**< Input: Allocated bound object. \n 504 * Output: Working set for auxiliary QP. */ 505 ) const; 506 507 /** Decides if lower bounds are smaller than upper bounds 508 * 509 * \return SUCCESSFUL_RETURN \n 510 * RET_QP_INFEASIBLE */ 511 512 returnValue areBoundsConsistent(const real_t* const lb, /**< Vector of lower bounds*/ 513 const real_t* const ub /**< Vector of upper bounds*/ 514 ) const; 515 516 /** Solves the system Ra = b or R^Ta = b where R is an upper triangular matrix. 517 * \return SUCCESSFUL_RETURN \n 518 RET_DIV_BY_ZERO */ 519 virtual returnValue backsolveR( const real_t* const b, /**< Right hand side vector. */ 520 BooleanType transposed, /**< Indicates if the transposed system shall be solved. */ 521 real_t* const a /**< Output: Solution vector */ 522 ) const; 523 524 /** Solves the system Ra = b or R^Ta = b where R is an upper triangular matrix. \n 525 * Special variant for the case that this function is called from within "removeBound()". 526 * \return SUCCESSFUL_RETURN \n 527 RET_DIV_BY_ZERO */ 528 virtual returnValue backsolveR( const real_t* const b, /**< Right hand side vector. */ 529 BooleanType transposed, /**< Indicates if the transposed system shall be solved. */ 530 BooleanType removingBound, /**< Indicates if function is called from "removeBound()". */ 531 real_t* const a /**< Output: Solution vector */ 532 ) const; 533 534 535 /** Determines step direction of the shift of the QP data. 536 * \return SUCCESSFUL_RETURN */ 537 returnValue determineDataShift( const real_t* const g_new, /**< New gradient vector. */ 538 const real_t* const lb_new, /**< New lower bounds. */ 539 const real_t* const ub_new, /**< New upper bounds. */ 540 real_t* const delta_g, /**< Output: Step direction of gradient vector. */ 541 real_t* const delta_lb, /**< Output: Step direction of lower bounds. */ 542 real_t* const delta_ub, /**< Output: Step direction of upper bounds. */ 543 BooleanType& Delta_bB_isZero/**< Output: Indicates if active bounds are to be shifted. */ 544 ); 545 546 547 /** Sets up internal QP data. 548 * \return SUCCESSFUL_RETURN \n 549 RET_INVALID_ARGUMENTS */ 550 returnValue setupQPdata( SymmetricMatrix *_H, /**< Hessian matrix.*/ 551 const real_t* const _g, /**< Gradient vector. */ 552 const real_t* const _lb, /**< Lower bounds (on variables). \n 553 If no lower bounds exist, a NULL pointer can be passed. */ 554 const real_t* const _ub /**< Upper bounds (on variables). \n 555 If no upper bounds exist, a NULL pointer can be passed. */ 556 ); 557 558 /** Sets up internal QP data. If the current Hessian is trivial 559 * (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, 560 * memory for Hessian is allocated and it is set to the given one. 561 * \return SUCCESSFUL_RETURN \n 562 RET_INVALID_ARGUMENTS \n 563 RET_NO_HESSIAN_SPECIFIED */ 564 returnValue setupQPdata( const real_t* const _H, /**< Hessian matrix. \n 565 If Hessian matrix is trivial,a NULL pointer can be passed. */ 566 const real_t* const _g, /**< Gradient vector. */ 567 const real_t* const _lb, /**< Lower bounds (on variables). \n 568 If no lower bounds exist, a NULL pointer can be passed. */ 569 const real_t* const _ub /**< Upper bounds (on variables). \n 570 If no upper bounds exist, a NULL pointer can be passed. */ 571 ); 572 573 /** Sets up internal QP data by loading it from files. If the current Hessian 574 * is trivial (i.e. HST_ZERO or HST_IDENTITY) but a non-trivial one is given, 575 * memory for Hessian is allocated and it is set to the given one. 576 * \return SUCCESSFUL_RETURN \n 577 RET_UNABLE_TO_OPEN_FILE \n 578 RET_UNABLE_TO_READ_FILE \n 579 RET_INVALID_ARGUMENTS \n 580 RET_NO_HESSIAN_SPECIFIED */ 581 returnValue setupQPdataFromFile( const char* const H_file, /**< Name of file where Hessian matrix, of neighbouring QP to be solved, is stored. \n 582 If Hessian matrix is trivial,a NULL pointer can be passed. */ 583 const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 584 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 585 If no lower bounds exist, a NULL pointer can be passed. */ 586 const char* const ub_file /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 587 If no upper bounds exist, a NULL pointer can be passed. */ 588 ); 589 590 /** Loads new QP vectors from files (internal members are not affected!). 591 * \return SUCCESSFUL_RETURN \n 592 RET_UNABLE_TO_OPEN_FILE \n 593 RET_UNABLE_TO_READ_FILE \n 594 RET_INVALID_ARGUMENTS */ 595 returnValue loadQPvectorsFromFile( const char* const g_file, /**< Name of file where gradient, of neighbouring QP to be solved, is stored. */ 596 const char* const lb_file, /**< Name of file where lower bounds, of neighbouring QP to be solved, is stored. \n 597 If no lower bounds exist, a NULL pointer can be passed. */ 598 const char* const ub_file, /**< Name of file where upper bounds, of neighbouring QP to be solved, is stored. \n 599 If no upper bounds exist, a NULL pointer can be passed. */ 600 real_t* const g_new, /**< Output: Gradient of neighbouring QP to be solved. */ 601 real_t* const lb_new, /**< Output: Lower bounds of neighbouring QP to be solved */ 602 real_t* const ub_new /**< Output: Upper bounds of neighbouring QP to be solved */ 603 ) const; 604 605 606 /** Sets internal infeasibility flag and throws given error in case the far bound 607 * strategy is not enabled (as QP might actually not be infeasible in this case). 608 * \return RET_HOTSTART_STOPPED_INFEASIBILITY \n 609 RET_ENSURELI_FAILED_CYCLING \n 610 RET_ENSURELI_FAILED_NOINDEX */ 611 returnValue setInfeasibilityFlag( returnValue returnvalue, /**< Returnvalue to be tunneled. */ 612 BooleanType doThrowError = BT_FALSE /**< Flag forcing to throw an error. */ 613 ); 614 615 616 /** Determines if next QP iteration can be performed within given CPU time limit. 617 * \return BT_TRUE: CPU time limit is exceeded, stop QP solution. \n 618 BT_FALSE: Sufficient CPU time for next QP iteration. */ 619 BooleanType isCPUtimeLimitExceeded( const real_t* const cputime, /**< Maximum CPU time allowed for QP solution. */ 620 real_t starttime, /**< Start time of current QP solution. */ 621 int_t nWSR /**< Number of working set recalculations performed so far. */ 622 ) const; 623 624 625 /** Regularise Hessian matrix by adding a scaled identity matrix to it. 626 * \return SUCCESSFUL_RETURN \n 627 RET_HESSIAN_ALREADY_REGULARISED */ 628 returnValue regulariseHessian( ); 629 630 631 /** Sets Hessian matrix of the QP. 632 * \return SUCCESSFUL_RETURN */ 633 inline returnValue setH( SymmetricMatrix* H_new /**< New Hessian matrix (a shallow copy is made). */ 634 ); 635 636 /** Sets dense Hessian matrix of the QP. 637 * If a null pointer is passed and 638 * a) hessianType is HST_IDENTITY, nothing is done, 639 * b) hessianType is not HST_IDENTITY, Hessian matrix is set to zero. 640 * \return SUCCESSFUL_RETURN */ 641 inline returnValue setH( const real_t* const H_new /**< New dense Hessian matrix (with correct dimension!), a shallow copy is made. */ 642 ); 643 644 /** Changes gradient vector of the QP. 645 * \return SUCCESSFUL_RETURN \n 646 * RET_INVALID_ARGUMENTS */ 647 inline returnValue setG( const real_t* const g_new /**< New gradient vector (with correct dimension!). */ 648 ); 649 650 /** Changes lower bound vector of the QP. 651 * \return SUCCESSFUL_RETURN \n 652 * RET_QPOBJECT_NOT_SETUP */ 653 inline returnValue setLB( const real_t* const lb_new /**< New lower bound vector (with correct dimension!). */ 654 ); 655 656 /** Changes single entry of lower bound vector of the QP. 657 * \return SUCCESSFUL_RETURN \n 658 * RET_QPOBJECT_NOT_SETUP \n 659 * RET_INDEX_OUT_OF_BOUNDS */ 660 inline returnValue setLB( int_t number, /**< Number of entry to be changed. */ 661 real_t value /**< New value for entry of lower bound vector. */ 662 ); 663 664 /** Changes upper bound vector of the QP. 665 * \return SUCCESSFUL_RETURN \n 666 * RET_QPOBJECT_NOT_SETUP */ 667 inline returnValue setUB( const real_t* const ub_new /**< New upper bound vector (with correct dimension!). */ 668 ); 669 670 /** Changes single entry of upper bound vector of the QP. 671 * \return SUCCESSFUL_RETURN \n 672 * RET_QPOBJECT_NOT_SETUP \n 673 * RET_INDEX_OUT_OF_BOUNDS */ 674 inline returnValue setUB( int_t number, /**< Number of entry to be changed. */ 675 real_t value /**< New value for entry of upper bound vector. */ 676 ); 677 678 679 /** Computes parameters for the Givens matrix G for which [x,y]*G = [z,0] 680 * \return SUCCESSFUL_RETURN */ 681 inline void computeGivens( real_t xold, /**< Matrix entry to be normalised. */ 682 real_t yold, /**< Matrix entry to be annihilated. */ 683 real_t& xnew, /**< Output: Normalised matrix entry. */ 684 real_t& ynew, /**< Output: Annihilated matrix entry. */ 685 real_t& c, /**< Output: Cosine entry of Givens matrix. */ 686 real_t& s /**< Output: Sine entry of Givens matrix. */ 687 ) const; 688 689 /** Applies Givens matrix determined by c and s (cf. computeGivens). 690 * \return SUCCESSFUL_RETURN */ 691 inline void applyGivens( real_t c, /**< Cosine entry of Givens matrix. */ 692 real_t s, /**< Sine entry of Givens matrix. */ 693 real_t nu, /**< Further factor: s/(1+c). */ 694 real_t xold, /**< Matrix entry to be transformed corresponding to 695 * the normalised entry of the original matrix. */ 696 real_t yold, /**< Matrix entry to be transformed corresponding to 697 * the annihilated entry of the original matrix. */ 698 real_t& xnew, /**< Output: Transformed matrix entry corresponding to 699 * the normalised entry of the original matrix. */ 700 real_t& ynew /**< Output: Transformed matrix entry corresponding to 701 * the annihilated entry of the original matrix. */ 702 ) const; 703 704 705 706 /** Compute relative length of homotopy in data space for termination 707 * criterion. 708 * \return Relative length in data space. */ 709 real_t getRelativeHomotopyLength( const real_t* const g_new, /**< Final gradient. */ 710 const real_t* const lb_new, /**< Final lower variable bounds. */ 711 const real_t* const ub_new /**< Final upper variable bounds. */ 712 ); 713 714 /** Ramping Strategy to avoid ties. Modifies homotopy start without 715 * changing current active set. 716 * \return SUCCESSFUL_RETURN */ 717 virtual returnValue performRamping( ); 718 719 720 /** ... */ 721 returnValue updateFarBounds( real_t curFarBound, /**< ... */ 722 int_t nRamp, /**< ... */ 723 const real_t* const lb_new, /**< ... */ 724 real_t* const lb_new_far, /**< ... */ 725 const real_t* const ub_new, /**< ... */ 726 real_t* const ub_new_far /**< ... */ 727 ) const; 728 729 730 /** Performs robustified ratio test yield the maximum possible step length 731 * along the homotopy path. 732 * \return SUCCESSFUL_RETURN */ 733 returnValue performRatioTest( int_t nIdx, /**< Number of ratios to be checked. */ 734 const int_t* const idxList, /**< Array containing the indices of all ratios to be checked. */ 735 const SubjectTo* const subjectTo, /**< Bound/Constraint object corresponding to ratios to be checked. */ 736 const real_t* const num, /**< Array containing all numerators for performing the ratio test. */ 737 const real_t* const den, /**< Array containing all denominators for performing the ratio test. */ 738 real_t epsNum, /**< Numerator tolerance. */ 739 real_t epsDen, /**< Denominator tolerance. */ 740 real_t& t, /**< Output: Maximum possible step length along the homotopy path. */ 741 int_t& BC_idx /**< Output: Index of blocking constraint. */ 742 ) const; 743 744 /** Checks whether given ratio is blocking, i.e. limits the maximum step length 745 * along the homotopy path to a value lower than given one. 746 * \return SUCCESSFUL_RETURN */ 747 inline BooleanType isBlocking( real_t num, /**< Numerator for performing the ratio test. */ 748 real_t den, /**< Denominator for performing the ratio test. */ 749 real_t epsNum, /**< Numerator tolerance. */ 750 real_t epsDen, /**< Denominator tolerance. */ 751 real_t& t /**< Input: Current maximum step length along the homotopy path, 752 * Output: Updated maximum possible step length along the homotopy path. */ 753 ) const; 754 755 756 /** Creates a sparse diagonal (square-)matrix which is a given 757 * multiple of the identity matrix. 758 * \return Diagonal matrix \n 759 */ 760 SymSparseMat* createDiagSparseMat( int_t n, /**< Row/column dimension of matrix to be created. */ 761 real_t diagVal = 1.0 /**< Value of all diagonal entries. */ 762 ); 763 764 765 /* 766 * PRIVATE MEMBER FUNCTIONS 767 */ 768 private: 769 /** Solves a QProblemB whose QP data is assumed to be stored in the member variables. 770 * A guess for its primal/dual optimal solution vectors and the corresponding 771 * optimal working set can be provided. 772 * Note: This function is internally called by all init functions! 773 * \return SUCCESSFUL_RETURN \n 774 RET_INIT_FAILED \n 775 RET_INIT_FAILED_CHOLESKY \n 776 RET_INIT_FAILED_HOTSTART \n 777 RET_INIT_FAILED_INFEASIBILITY \n 778 RET_INIT_FAILED_UNBOUNDEDNESS \n 779 RET_MAX_NWSR_REACHED */ 780 returnValue solveInitialQP( const real_t* const xOpt, /**< Optimal primal solution vector.*/ 781 const real_t* const yOpt, /**< Optimal dual solution vector. */ 782 const Bounds* const guessedBounds, /**< Optimal working set of bounds for solution (xOpt,yOpt). */ 783 const real_t* const _R, /**< Pre-computed (upper triangular) Cholesky factor of Hessian matrix. */ 784 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 785 * Output: Number of performed working set recalculations. */ 786 real_t* const cputime /**< Input: Maximum CPU time allowed for QP solution. \n 787 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 788 ); 789 790 /** Solves an initialised QProblemB using online active set strategy. 791 * Note: This function is internally called by all hotstart functions! 792 * \return SUCCESSFUL_RETURN \n 793 RET_MAX_NWSR_REACHED \n 794 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 795 RET_HOTSTART_FAILED \n 796 RET_SHIFT_DETERMINATION_FAILED \n 797 RET_STEPDIRECTION_DETERMINATION_FAILED \n 798 RET_STEPLENGTH_DETERMINATION_FAILED \n 799 RET_HOMOTOPY_STEP_FAILED \n 800 RET_HOTSTART_STOPPED_INFEASIBILITY \n 801 RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 802 returnValue solveQP( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 803 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 804 If no lower bounds exist, a NULL pointer can be passed. */ 805 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 806 If no upper bounds exist, a NULL pointer can be passed. */ 807 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 808 Output: Number of performed working set recalculations. */ 809 real_t* const cputime, /**< Input: Maximum CPU time allowed for QP solution. \n 810 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 811 int_t nWSRperformed = 0, /**< Number of working set recalculations already performed to solve 812 this QP within previous solveQP() calls. This number is 813 always zero, except for successive calls from solveRegularisedQP() 814 or when using the far bound strategy. */ 815 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 816 ); 817 818 819 /** Solves an initialised QProblemB using online active set strategy. 820 * Note: This function is internally called by all hotstart functions! 821 * \return SUCCESSFUL_RETURN \n 822 RET_MAX_NWSR_REACHED \n 823 RET_HOTSTART_FAILED_AS_QP_NOT_INITIALISED \n 824 RET_HOTSTART_FAILED \n 825 RET_SHIFT_DETERMINATION_FAILED \n 826 RET_STEPDIRECTION_DETERMINATION_FAILED \n 827 RET_STEPLENGTH_DETERMINATION_FAILED \n 828 RET_HOMOTOPY_STEP_FAILED \n 829 RET_HOTSTART_STOPPED_INFEASIBILITY \n 830 RET_HOTSTART_STOPPED_UNBOUNDEDNESS */ 831 returnValue solveRegularisedQP( const real_t* const g_new, /**< Gradient of neighbouring QP to be solved. */ 832 const real_t* const lb_new, /**< Lower bounds of neighbouring QP to be solved. \n 833 If no lower bounds exist, a NULL pointer can be passed. */ 834 const real_t* const ub_new, /**< Upper bounds of neighbouring QP to be solved. \n 835 If no upper bounds exist, a NULL pointer can be passed. */ 836 int_t& nWSR, /**< Input: Maximum number of working set recalculations; \n 837 Output: Number of performed working set recalculations. */ 838 real_t* const cputime, /**< Input: Maximum CPU time allowed for QP solution. \n 839 Output: CPU time spent for QP solution (or to perform nWSR iterations). */ 840 int_t nWSRperformed = 0, /**< Number of working set recalculations already performed to solve 841 this QP within previous solveRegularisedQP() calls. This number is 842 always zero, except for successive calls when using the far bound strategy. */ 843 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 844 ); 845 846 847 /** Sets up bound data structure according to auxiliaryBounds. 848 * (If the working set shall be setup afresh, make sure that 849 * bounds data structure has been resetted!) 850 * \return SUCCESSFUL_RETURN \n 851 RET_SETUP_WORKINGSET_FAILED \n 852 RET_INVALID_ARGUMENTS \n 853 RET_UNKNOWN_BUG */ 854 returnValue setupAuxiliaryWorkingSet( const Bounds* const auxiliaryBounds, /**< Working set for auxiliary QP. */ 855 BooleanType setupAfresh /**< Flag indicating if given working set shall be 856 * setup afresh or by updating the current one. */ 857 ); 858 859 /** Sets up the optimal primal/dual solution of the auxiliary initial QP. 860 * \return SUCCESSFUL_RETURN */ 861 returnValue setupAuxiliaryQPsolution( const real_t* const xOpt, /**< Optimal primal solution vector. 862 * If a NULL pointer is passed, all entries are set to zero. */ 863 const real_t* const yOpt /**< Optimal dual solution vector. 864 * If a NULL pointer is passed, all entries are set to zero. */ 865 ); 866 867 /** Sets up gradient of the auxiliary initial QP for given 868 * optimal primal/dual solution and given initial working set 869 * (assumes that members X, Y and BOUNDS have already been initialised!). 870 * \return SUCCESSFUL_RETURN */ 871 returnValue setupAuxiliaryQPgradient( ); 872 873 /** Sets up bounds of the auxiliary initial QP for given 874 * optimal primal/dual solution and given initial working set 875 * (assumes that members X, Y and BOUNDS have already been initialised!). 876 * \return SUCCESSFUL_RETURN \n 877 RET_UNKNOWN_BUG */ 878 returnValue setupAuxiliaryQPbounds( BooleanType useRelaxation /**< Flag indicating if inactive bounds shall be relaxed. */ 879 ); 880 881 882 protected: 883 /** Updates QP vectors, working sets and internal data structures in order to 884 start from an optimal solution corresponding to initial guesses of the working 885 set for bounds 886 * \return SUCCESSFUL_RETURN \n 887 * RET_SETUP_AUXILIARYQP_FAILED */ 888 virtual returnValue setupAuxiliaryQP( const Bounds* const guessedBounds /**< Initial guess for working set of bounds. */ 889 ); 890 891 private: 892 /** Determines step direction of the homotopy path. 893 * \return SUCCESSFUL_RETURN \n 894 RET_STEPDIRECTION_FAILED_CHOLESKY */ 895 returnValue determineStepDirection( const real_t* const delta_g, /**< Step direction of gradient vector. */ 896 const real_t* const delta_lb, /**< Step direction of lower bounds. */ 897 const real_t* const delta_ub, /**< Step direction of upper bounds. */ 898 BooleanType Delta_bB_isZero, /**< Indicates if active bounds are to be shifted. */ 899 real_t* const delta_xFX, /**< Output: Primal homotopy step direction of fixed variables. */ 900 real_t* const delta_xFR, /**< Output: Primal homotopy step direction of free variables. */ 901 real_t* const delta_yFX /**< Output: Dual homotopy step direction of fixed variables' multiplier. */ 902 ); 903 904 /** Determines the maximum possible step length along the homotopy path 905 * and performs this step (without changing working set). 906 * \return SUCCESSFUL_RETURN \n 907 * RET_QP_INFEASIBLE \n 908 */ 909 returnValue performStep( const real_t* const delta_g, /**< Step direction of gradient. */ 910 const real_t* const delta_lb, /**< Step direction of lower bounds. */ 911 const real_t* const delta_ub, /**< Step direction of upper bounds. */ 912 const real_t* const delta_xFX, /**< Primal homotopy step direction of fixed variables. */ 913 const real_t* const delta_xFR, /**< Primal homotopy step direction of free variables. */ 914 const real_t* const delta_yFX, /**< Dual homotopy step direction of fixed variables' multiplier. */ 915 int_t& BC_idx, /**< Output: Index of blocking constraint. */ 916 SubjectToStatus& BC_status /**< Output: Status of blocking constraint. */ 917 ); 918 919 /** Updates active set. 920 * \return SUCCESSFUL_RETURN \n 921 RET_REMOVE_FROM_ACTIVESET_FAILED \n 922 RET_ADD_TO_ACTIVESET_FAILED */ 923 returnValue changeActiveSet( int_t BC_idx, /**< Index of blocking constraint. */ 924 SubjectToStatus BC_status /**< Status of blocking constraint. */ 925 ); 926 927 /** Drift correction at end of each active set iteration 928 * \return SUCCESSFUL_RETURN */ 929 virtual returnValue performDriftCorrection( ); 930 931 /** Determines if it is more efficient to refactorise the matrices when 932 * hotstarting or not (i.e. better to update the existing factorisations). 933 * \return BT_TRUE iff matrices shall be refactorised afresh 934 */ 935 BooleanType shallRefactorise( const Bounds* const guessedBounds /**< Guessed new working set. */ 936 ) const; 937 938 939 /** Adds a bound to active set (specialised version for the case where no constraints exist). 940 * \return SUCCESSFUL_RETURN \n 941 RET_ADDBOUND_FAILED */ 942 returnValue addBound( int_t number, /**< Number of bound to be added to active set. */ 943 SubjectToStatus B_status, /**< Status of new active bound. */ 944 BooleanType updateCholesky /**< Flag indicating if Cholesky decomposition shall be updated. */ 945 ); 946 947 /** Removes a bounds from active set (specialised version for the case where no constraints exist). 948 * \return SUCCESSFUL_RETURN \n 949 RET_HESSIAN_NOT_SPD \n 950 RET_REMOVEBOUND_FAILED */ 951 returnValue removeBound( int_t number, /**< Number of bound to be removed from active set. */ 952 BooleanType updateCholesky /**< Flag indicating if Cholesky decomposition shall be updated. */ 953 ); 954 955 956 /** Prints concise information on the current iteration. 957 * \return SUCCESSFUL_RETURN \n */ 958 returnValue printIteration( int_t iter, /**< Number of current iteration. */ 959 int_t BC_idx, /**< Index of blocking bound. */ 960 SubjectToStatus BC_status, /**< Status of blocking bound. */ 961 real_t homotopyLength, /**< Current homotopy distance. */ 962 BooleanType isFirstCall = BT_TRUE /**< Indicating whether this is the first call for current QP. */ 963 ); 964 965 966 /* 967 * PROTECTED MEMBER VARIABLES 968 */ 969 protected: 970 BooleanType freeHessian; /**< Flag indicating whether the Hessian matrix needs to be de-allocated. */ 971 SymmetricMatrix* H; /**< Hessian matrix. */ 972 973 real_t* g; /**< Gradient. */ 974 real_t* lb; /**< Lower bound vector (on variables). */ 975 real_t* ub; /**< Upper bound vector (on variables). */ 976 977 Bounds bounds; /**< Data structure for problem's bounds. */ 978 979 real_t* R; /**< Cholesky factor of H (i.e. H = R^T*R). */ 980 BooleanType haveCholesky; /**< Flag indicating whether Cholesky decomposition has already been setup. */ 981 982 real_t* x; /**< Primal solution vector. */ 983 real_t* y; /**< Dual solution vector. */ 984 985 real_t tau; /**< Last homotopy step length. */ 986 987 QProblemStatus status; /**< Current status of the solution process. */ 988 989 BooleanType infeasible; /**< QP infeasible? */ 990 BooleanType unbounded; /**< QP unbounded? */ 991 992 HessianType hessianType; /**< Type of Hessian matrix. */ 993 real_t regVal; /**< Holds the offset used to regularise Hessian matrix (zero by default). */ 994 995 uint_t count; /**< Counts the number of hotstart function calls. */ 996 997 real_t *delta_xFR_TMP; /**< Temporary for determineStepDirection */ 998 999 real_t ramp0; /**< Start value for Ramping Strategy. */ 1000 real_t ramp1; /**< Final value for Ramping Strategy. */ 1001 int_t rampOffset; /**< Offset index for Ramping. */ 1002 1003 Options options; /**< Struct containing all user-defined options for solving QPs. */ 1004 1005 Flipper flipper; /**< Struct for making a temporary copy of the matrix factorisations. */ 1006 1007 TabularOutput tabularOutput; /**< Struct storing information for tabular output (printLevel == PL_TABULAR). */ 1008 }; 1009 1010 1011 END_NAMESPACE_QPOASES 1012 1013 #include <qpOASES/QProblemB.ipp> 1014 1015 #endif /* QPOASES_QPROBLEMB_HPP */ 1016 1017 1018 /* 1019 * end of file 1020 */ 1021