qpOASES  3.2.1
An Implementation of the Online Active Set Strategy
SQProblemSchur.hpp
Go to the documentation of this file.
1 /*
2  * This file is part of qpOASES.
3  *
4  * qpOASES -- An Implementation of the Online Active Set Strategy.
5  * Copyright (C) 2007-2014 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 
38 #ifndef QPOASES_SQPROBLEMSCHUR_HPP
39 #define QPOASES_SQPROBLEMSCHUR_HPP
40 
41 
42 #include <qpOASES/SQProblem.hpp>
43 #include <qpOASES/SparseSolver.hpp>
45 
46 
48 
49 
61 class SQProblemSchur : public SQProblem
62 {
63  /* allow SolutionAnalysis class to access private members */
64  friend class SolutionAnalysis;
65 
66  /*
67  * PUBLIC MEMBER FUNCTIONS
68  */
69  public:
71  SQProblemSchur( );
72 
78  SQProblemSchur( int_t _nV,
79  int_t _nC,
80  HessianType _hessianType = HST_UNKNOWN,
81  int_t maxSchurUpdates = 75
82  );
83 
85  SQProblemSchur( const SQProblemSchur& rhs
86  );
87 
89  virtual ~SQProblemSchur( );
90 
92  virtual SQProblemSchur& operator=( const SQProblemSchur& rhs
93  );
94 
98  virtual returnValue reset( );
99 
102  returnValue resetSchurComplement( BooleanType allowInertiaCorrection );
103 
105  inline int_t getNumFactorizations( ) const;
106 
107  /*
108  * PROTECTED MEMBER FUNCTIONS
109  */
110  protected:
113  returnValue clear( );
114 
117  returnValue copy( const SQProblemSchur& rhs
118  );
119 
123 
128 
132 
144  Matrix *A_new,
146  const real_t *lb_new,
147  const real_t *ub_new,
148  const real_t *lbA_new,
149  const real_t *ubA_new
150  );
151 
158  virtual returnValue setupAuxiliaryWorkingSet( const Bounds* const auxiliaryBounds,
159  const Constraints* const auxiliaryConstraints,
160  BooleanType setupAfresh
161  );
162 
163 
169  virtual returnValue addConstraint( int_t number,
170  SubjectToStatus C_status,
171  BooleanType updateCholesky,
172  BooleanType ensureLI = BT_TRUE
173  );
174 
180  virtual returnValue addConstraint_checkLI( int_t number
181  );
182 
191  virtual returnValue addConstraint_ensureLI( int_t number,
192  SubjectToStatus C_status
193  );
194 
200  virtual returnValue addBound( int_t number,
201  SubjectToStatus B_status,
202  BooleanType updateCholesky,
203  BooleanType ensureLI = BT_TRUE
204  );
205 
210  virtual returnValue addBound_checkLI( int_t number
211  );
212 
221  virtual returnValue addBound_ensureLI( int_t number,
222  SubjectToStatus B_status
223  );
224 
230  virtual returnValue removeConstraint( int_t number,
231  BooleanType updateCholesky,
232  BooleanType allowFlipping = BT_FALSE,
233  BooleanType ensureNZC = BT_FALSE
234  );
235 
241  virtual returnValue removeBound( int_t number,
242  BooleanType updateCholesky,
243  BooleanType allowFlipping = BT_FALSE,
244  BooleanType ensureNZC = BT_FALSE
245  );
246 
249  virtual returnValue backsolveT( const real_t* const b,
250  BooleanType transposed,
251  real_t* const a
252  ) const;
253 
256  virtual returnValue backsolveR( const real_t* const b,
257  BooleanType transposed,
258  real_t* const a
259  ) const;
260 
263  virtual returnValue backsolveR( const real_t* const b,
264  BooleanType transposed,
265  BooleanType removingBound,
266  real_t* const a
267  ) const;
268 
269 
274  virtual returnValue determineStepDirection( const real_t* const delta_g,
275  const real_t* const delta_lbA,
276  const real_t* const delta_ubA,
277  const real_t* const delta_lb,
278  const real_t* const delta_ub,
279  BooleanType Delta_bC_isZero,
280  BooleanType Delta_bB_isZero,
281  real_t* const delta_xFX,
282  real_t* const delta_xFR,
283  real_t* const delta_yAC,
284  real_t* const delta_yFX
285  );
286 
287  virtual returnValue determineStepDirection2( const real_t* const delta_g,
288  const real_t* const delta_lbA,
289  const real_t* const delta_ubA,
290  const real_t* const delta_lb,
291  const real_t* const delta_ub,
292  BooleanType Delta_bC_isZero,
293  BooleanType Delta_bB_isZero,
294  real_t* const delta_xFX,
295  real_t* const delta_xFR,
296  real_t* const delta_yAC,
297  real_t* const delta_yFX
298  );
299 
300  /*
301  * PRIVATE MEMBER FUNCTION
302  */
303  private:
310  real_t* const xiC,
311  real_t* const xiX
312  );
313 
320  real_t* const xiC,
321  real_t* const xiX
322  );
323 
326  returnValue computeMTimes( real_t alpha, const real_t* const x, real_t beta, real_t* const y );
327 
330  returnValue computeMTransTimes( real_t alpha, const real_t* const x, real_t beta, real_t* const y );
331 
333  returnValue addToSchurComplement( int_t number, SchurUpdateType update, int_t numNonzerosM, const sparse_int_t* M_pos, const real_t* const M_vals, int_t numNonzerosN, const sparse_int_t* Npos, const real_t* const Nvals, real_t N_diag );
334 
337 
340 
342  real_t calcDetSchur( int_t idxDel );
343 
345  returnValue updateSchurQR( int_t idxDel );
346 
348  returnValue backsolveSchurQR( int_t dimS, const real_t* const rhs, int_t dimRhs, real_t* const sol );
349 
352 
355 
356  returnValue stepCalcRhs( int_t nFR, int_t nFX, int_t nAC, int_t* FR_idx, int_t* FX_idx, int_t* AC_idx, real_t& rhs_max, const real_t* const delta_g,
357  const real_t* const delta_lbA, const real_t* const delta_ubA,
358  const real_t* const delta_lb, const real_t* const delta_ub,
359  BooleanType Delta_bC_isZero, BooleanType Delta_bB_isZero,
360  real_t* const delta_xFX, real_t* const delta_xFR,
361  real_t* const delta_yAC, real_t* const delta_yFX
362  );
363 
364  returnValue stepCalcReorder(int_t nFR, int_t nAC, int_t* FR_idx, int_t* AC_idx, int_t nFRStart, int_t nACStart,
365  int_t* FR_idxStart, int_t* AC_idxStart, int_t* FR_iSort, int_t* FR_iSortStart,
366  int_t* AC_iSort, int_t* AC_iSortStart, real_t* rhs
367  );
368 
369  returnValue stepCalcBacksolveSchur( int_t nFR, int_t nFX, int_t nAC, int_t* FR_idx, int_t* FX_idx, int_t* AC_idx,
370  int_t dim, real_t* rhs, real_t* sol
371  );
372 
373  returnValue stepCalcReorder2( int_t nFR, int_t nAC, int_t* FR_idx, int_t* AC_idx, int_t nFRStart, int_t nACStart,
374  int_t* FR_idxStart, int_t* AC_idxStart, int_t* FR_iSort, int_t* FR_iSortStart,
375  int_t* AC_iSort, int_t* AC_iSortStart, real_t* sol, real_t* const delta_xFR, real_t* const delta_yAC
376  );
377 
378  returnValue stepCalcResid( int_t nFR, int_t nFX, int_t nAC, int_t* FR_idx, int_t* FX_idx, int_t* AC_idx,
379  BooleanType Delta_bC_isZero, real_t* const delta_xFX, real_t* const delta_xFR,
380  real_t* const delta_yAC, const real_t* const delta_g,
381  const real_t* const delta_lbA, const real_t* const delta_ubA, real_t& rnrm
382  );
383 
384  returnValue stepCalcDeltayFx( int_t nFR, int_t nFX, int_t nAC, int_t* FX_idx, const real_t* const delta_g,
385  real_t* const delta_xFX, real_t* const delta_xFR, real_t* const delta_yAC, real_t* const delta_yFX
386  );
387 
388  /*
389  * PROTECTED MEMBER VARIABLES
390  */
391  protected:
414 };
415 
416 
418 
420 
421 #endif /* QPOASES_QPROBLEMSCHUR_HPP */
422 
423 
424 /*
425  * end of file
426  */
#define BEGIN_NAMESPACE_QPOASES
Definition: Types.hpp:107
virtual returnValue addBound(int_t number, SubjectToStatus B_status, BooleanType updateCholesky, BooleanType ensureLI=BT_TRUE)
Definition: SQProblemSchur.cpp:1115
#define END_NAMESPACE_QPOASES
Definition: Types.hpp:110
int_t nSmax
Definition: SQProblemSchur.hpp:396
Manages working sets of constraints.
Definition: Constraints.hpp:56
virtual returnValue setupTQfactorisation()
Definition: SQProblemSchur.cpp:639
int_t M_physicallength
Definition: SQProblemSchur.hpp:407
returnValue deleteFromSchurComplement(int_t idx, BooleanType allowUndo=BT_FALSE)
Definition: SQProblemSchur.cpp:3244
int int_t
Definition: Types.hpp:180
SchurUpdateType
Definition: Types.hpp:283
returnValue undoDeleteFromSchurComplement(int_t idx)
Definition: SQProblemSchur.cpp:3343
Implements the online active set strategy for QPs with varying matrices.
Definition: SQProblem.hpp:59
HessianType
Definition: Types.hpp:248
sparse_int_t * M_jc
Definition: SQProblemSchur.hpp:410
virtual returnValue backsolveR(const real_t *const b, BooleanType transposed, real_t *const a) const
Definition: SQProblemSchur.cpp:2173
virtual returnValue addConstraint_checkLI(int_t number)
Definition: SQProblemSchur.cpp:834
virtual returnValue computeProjectedCholesky()
Definition: SQProblemSchur.cpp:621
returnValue resetSchurComplement(BooleanType allowInertiaCorrection)
Definition: SQProblemSchur.cpp:2984
returnValue computeMTimes(real_t alpha, const real_t *const x, real_t beta, real_t *const y)
Definition: SQProblemSchur.cpp:3116
SparseSolver * sparseSolver
Definition: SQProblemSchur.hpp:392
Definition: Types.hpp:206
SubjectToStatus
Definition: Types.hpp:272
SQProblemSchur()
Definition: SQProblemSchur.cpp:68
int_t * schurUpdateIndex
Definition: SQProblemSchur.hpp:404
real_t * Q_
Definition: SQProblemSchur.hpp:398
virtual returnValue determineStepDirection2(const real_t *const delta_g, const real_t *const delta_lbA, const real_t *const delta_ubA, const real_t *const delta_lb, const real_t *const delta_ub, BooleanType Delta_bC_isZero, BooleanType Delta_bB_isZero, real_t *const delta_xFX, real_t *const delta_xFR, real_t *const delta_yAC, real_t *const delta_yFX)
Definition: SQProblemSchur.cpp:2841
real_t calcDetSchur(int_t idxDel)
Definition: SQProblemSchur.cpp:2191
sparse_int_t * M_ir
Definition: SQProblemSchur.hpp:409
virtual returnValue reset()
Definition: SQProblemSchur.cpp:191
returnValue addToSchurComplement(int_t number, SchurUpdateType update, int_t numNonzerosM, const sparse_int_t *M_pos, const real_t *const M_vals, int_t numNonzerosN, const sparse_int_t *Npos, const real_t *const Nvals, real_t N_diag)
Definition: SQProblemSchur.cpp:3162
returnValue correctInertia()
Definition: SQProblemSchur.cpp:3428
virtual returnValue addBound_ensureLI(int_t number, SubjectToStatus B_status)
Definition: SQProblemSchur.cpp:1357
virtual SQProblemSchur & operator=(const SQProblemSchur &rhs)
Definition: SQProblemSchur.cpp:176
virtual returnValue addBound_checkLI(int_t number)
Definition: SQProblemSchur.cpp:1250
returnValue updateSchurQR(int_t idxDel)
Definition: SQProblemSchur.cpp:2284
Indexlist constraintsActiveStart
Definition: SQProblemSchur.hpp:413
int_t numFactorizations
Definition: SQProblemSchur.hpp:402
virtual ~SQProblemSchur()
Definition: SQProblemSchur.cpp:165
Definition: Types.hpp:205
real_t * M_vals
Definition: SQProblemSchur.hpp:408
real_t * y
Definition: QProblemB.hpp:984
returnValue repairSingularWorkingSet()
Definition: SQProblemSchur.cpp:3525
real_t * S
Definition: SQProblemSchur.hpp:394
virtual returnValue addConstraint_ensureLI(int_t number, SubjectToStatus C_status)
Definition: SQProblemSchur.cpp:947
Provides additional tools for analysing QP solutions.
Definition: SolutionAnalysis.hpp:57
real_t * x
Definition: QProblemB.hpp:983
virtual returnValue addConstraint(int_t number, SubjectToStatus C_status, BooleanType updateCholesky, BooleanType ensureLI=BT_TRUE)
Definition: SQProblemSchur.cpp:648
returnValue
Defines all symbols for global return values.
Definition: MessageHandling.hpp:64
SchurUpdateType * schurUpdate
Definition: SQProblemSchur.hpp:405
returnValue addConstraint_checkLISchur(int_t number, real_t *const xiC, real_t *const xiX)
Definition: SQProblemSchur.cpp:855
Abstract base class for interfacing tailored matrix-vector operations.
Definition: Matrices.hpp:59
Definition: Types.hpp:256
virtual returnValue setupAuxiliaryWorkingSet(const Bounds *const auxiliaryBounds, const Constraints *const auxiliaryConstraints, BooleanType setupAfresh)
Definition: SQProblemSchur.cpp:512
int_t nS
Definition: SQProblemSchur.hpp:395
virtual returnValue determineStepDirection(const real_t *const delta_g, const real_t *const delta_lbA, const real_t *const delta_ubA, const real_t *const delta_lb, const real_t *const delta_ub, BooleanType Delta_bC_isZero, BooleanType Delta_bB_isZero, real_t *const delta_xFX, real_t *const delta_xFR, real_t *const delta_yAC, real_t *const delta_yFX)
Definition: SQProblemSchur.cpp:2810
virtual returnValue removeConstraint(int_t number, BooleanType updateCholesky, BooleanType allowFlipping=BT_FALSE, BooleanType ensureNZC=BT_FALSE)
Definition: SQProblemSchur.cpp:1523
returnValue computeMTransTimes(real_t alpha, const real_t *const x, real_t beta, real_t *const y)
Definition: SQProblemSchur.cpp:3133
real_t * R_
Definition: SQProblemSchur.hpp:399
real_t rcondS
Definition: SQProblemSchur.hpp:401
returnValue copy(const SQProblemSchur &rhs)
Definition: SQProblemSchur.cpp:235
Base class for linear solvers that are used in a Schur-complement implementation in qpOASES...
Definition: SparseSolver.hpp:52
Stores and manages index lists.
Definition: Indexlist.hpp:55
virtual returnValue setupAuxiliaryQP(SymmetricMatrix *H_new, Matrix *A_new, const real_t *lb_new, const real_t *ub_new, const real_t *lbA_new, const real_t *ubA_new)
Definition: SQProblemSchur.cpp:312
Indexlist boundsFreeStart
Definition: SQProblemSchur.hpp:412
int_t getNumFactorizations() const
Definition: SQProblemSchur.ipp:47
real_t detS
Definition: SQProblemSchur.hpp:400
virtual returnValue removeBound(int_t number, BooleanType updateCholesky, BooleanType allowFlipping=BT_FALSE, BooleanType ensureNZC=BT_FALSE)
Definition: SQProblemSchur.cpp:1798
virtual returnValue computeInitialCholesky()
Definition: SQProblemSchur.cpp:630
returnValue addBound_checkLISchur(int_t number, real_t *const xiC, real_t *const xiX)
Definition: SQProblemSchur.cpp:1270
Manages working sets of bounds (i.e. box constraints).
Definition: Bounds.hpp:56
Implements the online active set strategy for QPs with varying, sparse matrices.
Definition: SQProblemSchur.hpp:61
returnValue clear()
Definition: SQProblemSchur.cpp:211
virtual returnValue backsolveT(const real_t *const b, BooleanType transposed, real_t *const a) const
Definition: SQProblemSchur.cpp:2164
returnValue backsolveSchurQR(int_t dimS, const real_t *const rhs, int_t dimRhs, real_t *const sol)
Definition: SQProblemSchur.cpp:2423
BooleanType
Definition: Types.hpp:203
BEGIN_NAMESPACE_QPOASES typedef double real_t
Definition: Types.hpp:171
int_t sparse_int_t
Definition: Types.hpp:199
Abstract base class for interfacing matrix-vector operations tailored to symmetric matrices...
Definition: Matrices.hpp:292