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