qpmad
Eigen-based C++ QP solver.
solver.h
Go to the documentation of this file.
1 /**
2  @file
3  @author Alexander Sherikov
4 
5  @copyright 2017 Alexander Sherikov. Licensed under the Apache License,
6  Version 2.0. (see LICENSE or http://www.apache.org/licenses/LICENSE-2.0)
7 
8  @brief
9 */
10 
11 
12 #pragma once
13 
14 #include "common.h"
15 #include "givens.h"
16 #include "input_parser.h"
17 #include "inverse.h"
18 #include "solver_parameters.h"
19 #include "constraint_status.h"
20 #include "chosen_constraint.h"
21 #include "active_set.h"
22 #include "factorization_data.h"
23 
24 #ifdef QPMAD_ENABLE_TRACING
25 # include "testing.h"
26 #endif
27 
28 namespace qpmad
29 {
30  template <typename t_Scalar, int t_primal_size, int t_has_bounds, int t_general_ctr_number>
31  class SolverTemplate : public InputParser
32  {
33  public:
35  {
36  OK = 0,
38  };
39 
40  template <int t_rows>
41  using Vector = Eigen::Matrix<t_Scalar, t_rows, 1>;
42  template <int t_rows, int t_cols>
43  using Matrix = Eigen::Matrix<t_Scalar, t_rows, t_cols>;
44 
45 
46  protected:
48  Eigen::Dynamic == t_general_ctr_number ?
49  Eigen::Dynamic :
50  (0 == t_has_bounds ? t_general_ctr_number :
51  (Eigen::Dynamic == t_primal_size ? Eigen::Dynamic :
52  t_general_ctr_number + t_primal_size));
55 
58 
60 
63 
65 
66  Eigen::Array<uint8_t, num_constraints_compile_time_, 1> constraints_status_;
67 
69 
70  std::ptrdiff_t iter_counter_;
71 
73 
74 
75  public:
77  {
78  iter_counter_ = 0;
79  machinery_initialized_ = false;
81  }
82 
83 
84  /**
85  * @brief Returns type of the Hessian produced by the latest execution
86  * of `solve()`.
87  */
89  {
90  return (hessian_type_);
91  }
92 
93 
94  /**
95  * @brief Returns number of inequality iterations during the latest
96  * execution of `solve()`.
97  */
98  std::ptrdiff_t getNumberOfInequalityIterations() const
99  {
100  return (iter_counter_);
101  }
102 
103 
104  /**
105  * @brief Returns dual variables (Lagrange multipliers) corresponding
106  * to inequality constraints. Must be called after successful
107  * `solve()`, the result is undefined if previous call to `solve()`
108  * failed.
109  *
110  * @tparam t_status_size
111  * @tparam t_dual_size
112  * @tparam t_index_size
113  *
114  * @param[out] dual dual variables
115  * @param[out] indices constraint indices corresponding to the dual
116  * variables, index 0 corresponds to the first simple bound if present
117  * or to the first general constraint otherwise
118  * @param[out] is_lower flags indicating if lower or upper bound is
119  * active
120  */
121  template <int t_status_size, int t_dual_size, int t_index_size>
123  Vector<t_dual_size> &dual,
124  Eigen::Matrix<MatrixIndex, t_index_size, 1> &indices,
125  Eigen::Matrix<bool, t_status_size, 1> &is_lower) const
126  {
128  {
129  const MatrixIndex size = active_set_.size_ - active_set_.num_equalities_;
130 
131  dual.resize(size);
132  indices.resize(size);
133  is_lower.resize(size);
134 
135  for (MatrixIndex i = active_set_.num_equalities_; i < active_set_.size_; ++i)
136  {
137  const std::size_t output_index = i - active_set_.num_equalities_;
138 
139  dual(output_index) = dual_(i);
140  indices(output_index) = active_set_.getIndex(i);
141  is_lower(output_index) =
143  }
144  }
145  else
146  {
147  dual.resize(0);
148  indices.resize(0);
149  is_lower.resize(0);
150  }
151  }
152 
153 
154  template <
155  int t_rows_primal,
156  int t_rows_H,
157  int t_cols_H,
158  int t_rows_h,
159  int t_rows_A,
160  int t_cols_A,
161  int t_rows_Alb,
162  int t_rows_Aub>
164  Vector<t_rows_primal> &primal,
166  const Vector<t_rows_h> &h,
168  const Vector<t_rows_Alb> &Alb,
169  const Vector<t_rows_Aub> &Aub)
170  {
171  return (solve(primal, H, h, Eigen::VectorXd(), Eigen::VectorXd(), A, Alb, Aub, SolverParameters()));
172  }
173 
174 
175  template <
176  int t_rows_primal,
177  int t_rows_H,
178  int t_cols_H,
179  int t_rows_h,
180  int t_rows_lb,
181  int t_rows_ub,
182  int t_rows_A,
183  int t_cols_A,
184  int t_rows_Alb,
185  int t_rows_Aub>
187  Vector<t_rows_primal> &primal,
189  const Vector<t_rows_h> &h,
190  const Vector<t_rows_lb> &lb,
191  const Vector<t_rows_ub> &ub,
193  const Vector<t_rows_Alb> &Alb,
194  const Vector<t_rows_Aub> &Aub)
195  {
196  return (solve(primal, H, h, lb, ub, A, Alb, Aub, SolverParameters()));
197  }
198 
199 
200  template <int t_rows_primal, int t_rows_H, int t_cols_H, int t_rows_h, int t_rows_lb, int t_rows_ub>
202  Vector<t_rows_primal> &primal,
204  const Vector<t_rows_h> &h,
205  const Vector<t_rows_lb> &lb,
206  const Vector<t_rows_ub> &ub,
207  const SolverParameters &param)
208  {
209  return (solve(primal, H, h, lb, ub, Eigen::MatrixXd(), Eigen::VectorXd(), Eigen::VectorXd(), param));
210  }
211 
212 
213  template <int t_rows_primal, int t_rows_H, int t_cols_H, int t_rows_h, int t_rows_lb, int t_rows_ub>
215  Vector<t_rows_primal> &primal,
217  const Vector<t_rows_h> &h,
218  const Vector<t_rows_lb> &lb,
219  const Vector<t_rows_ub> &ub)
220  {
221  return (solve(
222  primal, H, h, lb, ub, Eigen::MatrixXd(), Eigen::VectorXd(), Eigen::VectorXd(), SolverParameters()));
223  }
224 
225 
226 
227  template <
228  int t_rows_primal,
229  int t_rows_H,
230  int t_cols_H,
231  int t_rows_h,
232  int t_rows_lb,
233  int t_rows_ub,
234  int t_rows_A,
235  int t_cols_A,
236  int t_rows_Alb,
237  int t_rows_Aub>
239  Vector<t_rows_primal> &primal,
241  const Vector<t_rows_h> &h,
242  const Vector<t_rows_lb> &lb,
243  const Vector<t_rows_ub> &ub,
245  const Vector<t_rows_Alb> &Alb,
246  const Vector<t_rows_Aub> &Aub,
247  const SolverParameters &param)
248  {
249  QPMAD_TRACE(std::setprecision(std::numeric_limits<double>::digits10));
250 
251  machinery_initialized_ = false;
252  iter_counter_ = 0;
253 
254  parseObjective(H, h);
255  parseSimpleBounds(lb, ub);
256  parseGeneralConstraints(A, Alb, Aub);
257 
258 
260 
261  if (0 == h_size_)
262  {
263  // trivial unconstrained optimum
264  primal.setZero(primal_size_);
265 
266  if (0 == num_constraints_)
267  {
268  // trivial solution
269  return (OK);
270  }
271  }
272 
273 
274  switch (param.hessian_type_)
275  {
277  {
278  const Eigen::LLT<Eigen::Ref<Eigen::MatrixXd>, Eigen::Lower> llt(H);
280  Eigen::Success == llt.info(), "Could not perform Cholesky decomposition of the Hessian.");
281  }
282  // no break here!
283  /* Falls through. */
286  // unconstrained optimum
287  if (h_size_ > 0)
288  {
289  primal = H.template triangularView<Eigen::Lower>().solve(-h);
290  H.transpose().template triangularView<Eigen::Upper>().solveInPlace(primal);
291  }
292  break;
293 
296  // unconstrained optimum
297  if (h_size_ > 0)
298  {
299  primal_step_direction_.noalias() = H.template triangularView<Eigen::Upper>().transpose() * -h;
300  primal.noalias() = H.template triangularView<Eigen::Upper>() * primal_step_direction_;
301  }
302  break;
303 
304  default:
305  QPMAD_UTILS_THROW("Malformed solver parameters!");
306  break;
307  }
308 
309 
310  if (0 == num_constraints_)
311  {
312  // return unconstrained optimum
313  return (OK);
314  }
315 
316 
317  // check consistency of general constraints and activate
318  // equality constraints
320  MatrixIndex num_equalities = 0;
321  for (MatrixIndex i = 0; i < num_constraints_; ++i)
322  {
324 
325  double lb_i;
326  double ub_i;
327 
328  if (true == chosen_ctr_.is_simple_)
329  {
330  lb_i = lb(i);
331  ub_i = ub(i);
332  }
333  else
334  {
338  }
339 
340 
341  if (lb_i - param.tolerance_ > ub_i)
342  {
344  QPMAD_UTILS_THROW("Inconsistent constraints (lb > ub).");
345  }
346 
347  if (std::abs(lb_i - ub_i) > param.tolerance_)
348  {
350  }
351  else
352  {
354  ++num_equalities;
355 
356 
357  if (true == chosen_ctr_.is_simple_)
358  {
359  chosen_ctr_.violation_ = lb_i - primal(i);
360  }
361  else
362  {
364  }
365 
367 
368 
369  // if 'primal_size_' constraints are already activated
370  // all other constraints are linearly dependent
371  if (active_set_.hasEmptySpace())
372  {
373  double ctr_i_dot_primal_step_direction;
374 
375  if (true == chosen_ctr_.is_simple_)
376  {
377  factorization_data_.computeEqualityPrimalStep(primal_step_direction_, i, active_set_.size_);
378 
379  ctr_i_dot_primal_step_direction = primal_step_direction_(i);
380  }
381  else
382  {
383  factorization_data_.computeEqualityPrimalStep(
386  active_set_.size_);
387 
388  ctr_i_dot_primal_step_direction =
390  }
391 
392  // if step direction is a zero vector, constraint is
393  // linearly dependent with previously added constraints
394  if (ctr_i_dot_primal_step_direction < -param.tolerance_)
395  {
396  double primal_step_length_ = chosen_ctr_.violation_ / ctr_i_dot_primal_step_direction;
397 
398  primal.noalias() += primal_step_length_ * primal_step_direction_;
399 
400  if (false
401  == factorization_data_.update(
403  {
404  QPMAD_UTILS_THROW("Failed to add an equality constraint -- is this possible?");
405  }
406  active_set_.addEquality(i);
407 
408  continue;
409  }
410  } // otherwise -- linear dependence
411 
412  // this point is reached if constraint is linearly dependent
413 
414  // check if this constraint is actually satisfied
415  if (std::abs(chosen_ctr_.violation_) > param.tolerance_)
416  {
417  // nope it is not
419  QPMAD_UTILS_THROW("Infeasible equality constraints");
420  }
421  // otherwise keep going
422  }
423  }
424 
425 
426  if (num_equalities == num_constraints_)
427  {
428  // exit early -- avoid unnecessary memory allocations
429  return (OK);
430  }
431 
432 
433  ReturnStatus return_status;
434  chooseConstraint(primal, lb, ub, A, Alb, Aub, param.tolerance_);
435 
436 
437  if (std::abs(chosen_ctr_.violation_) < param.tolerance_)
438  {
439  // all constraints are satisfied
440  return_status = OK;
441  }
442  else
443  {
444  return_status = MAXIMAL_NUMBER_OF_ITERATIONS;
445 
446 
448  dual_.resize(primal_size_);
450 
451 
452  double chosen_ctr_dot_primal_step_direction = 0.0;
453 
454  //
455  factorization_data_.computeInequalityDualStep(dual_step_direction_, chosen_ctr_, A, active_set_);
456  if (active_set_.hasEmptySpace())
457  {
458  // compute step direction in primal space
459  factorization_data_.computeInequalityPrimalStep(primal_step_direction_, active_set_);
460  chosen_ctr_dot_primal_step_direction =
462  }
463 
464 
465  // last iteration is not counted, so iter_counter_ starts with 1.
466  for (iter_counter_ = 1; (param.max_iter_ < 0) or (iter_counter_ <= param.max_iter_); ++iter_counter_)
467  {
468  QPMAD_TRACE(">>>>>>>>>" << iter_counter_ << "<<<<<<<<<");
469 #ifdef QPMAD_ENABLE_TRACING
470  testing::computeObjective(H, h, primal);
471 #endif
472  QPMAD_TRACE("||| Chosen ctr index = " << chosen_ctr_.index_);
473  QPMAD_TRACE("||| Chosen ctr dual = " << chosen_ctr_.dual_);
474  QPMAD_TRACE("||| Chosen ctr violation = " << chosen_ctr_.violation_);
475 
476 
477  // check dual feasibility
478  MatrixIndex dual_blocking_index = primal_size_;
479  double dual_step_length = std::numeric_limits<double>::infinity();
480  for (MatrixIndex i = active_set_.num_equalities_; i < active_set_.size_; ++i)
481  {
482  if (dual_step_direction_(i) < -param.tolerance_)
483  {
484  const double dual_step_length_i = -dual_(i) / dual_step_direction_(i);
485  if (dual_step_length_i < dual_step_length)
486  {
487  dual_step_length = dual_step_length_i;
488  dual_blocking_index = i;
489  }
490  }
491  }
492 
493 
494 #ifdef QPMAD_ENABLE_TRACING
496  H,
497  h,
498  primal,
499  A,
500  active_set_,
503  dual_,
505 #endif
506 
507 
508  if (active_set_.hasEmptySpace()
509  // if step direction is a zero vector, constraint is
510  // linearly dependent with previously added constraints
511  && (std::abs(chosen_ctr_dot_primal_step_direction) > param.tolerance_))
512  {
513  double step_length = -chosen_ctr_.violation_ / chosen_ctr_dot_primal_step_direction;
514 
515  QPMAD_TRACE("======================");
516  QPMAD_TRACE("||| Primal step length = " << step_length);
517  QPMAD_TRACE("||| Dual step length = " << dual_step_length);
518  QPMAD_TRACE("======================");
519 
520 
521  bool partial_step = false;
523  (step_length >= 0.0) && (dual_step_length >= 0.0),
524  "Non-negative step lengths expected.");
525  if (dual_step_length <= step_length)
526  {
527  step_length = dual_step_length;
528  partial_step = true;
529  }
530 
531 
532  primal.noalias() += step_length * primal_step_direction_;
533 
534  dual_.segment(active_set_.num_equalities_, active_set_.num_inequalities_).noalias() +=
535  step_length
536  * dual_step_direction_.segment(
537  active_set_.num_equalities_, active_set_.num_inequalities_);
538  chosen_ctr_.dual_ += step_length;
539  chosen_ctr_.violation_ += step_length * chosen_ctr_dot_primal_step_direction;
540 
541  QPMAD_TRACE("||| Chosen ctr dual = " << chosen_ctr_.dual_);
542  QPMAD_TRACE("||| Chosen ctr violation = " << chosen_ctr_.violation_);
543 
544 
545  if ((partial_step)
546  // if violation is almost zero -- assume that a full step is made
547  && (std::abs(chosen_ctr_.violation_) > param.tolerance_))
548  {
549  QPMAD_TRACE("||| PARTIAL STEP");
550  // deactivate blocking constraint
551  constraints_status_[active_set_.getIndex(dual_blocking_index)] = ConstraintStatus::INACTIVE;
552 
553  dropElementWithoutResize(dual_, dual_blocking_index, active_set_.size_);
554 
555  factorization_data_.downdate(dual_blocking_index, active_set_.size_);
556 
557  active_set_.removeInequality(dual_blocking_index);
558 
559  // compute step direction in primal & dual space
560  factorization_data_.updateStepsAfterPartialStep(
562  chosen_ctr_dot_primal_step_direction =
564  }
565  else
566  {
567  QPMAD_TRACE("||| FULL STEP");
568  // activate constraint
569  if (false
570  == factorization_data_.update(
572  {
573  QPMAD_UTILS_THROW("Failed to add an inequality constraint -- is this possible?");
574  }
575 
577  {
579  }
580  else
581  {
583  }
585  active_set_.addInequality(chosen_ctr_.index_);
586 
587  chooseConstraint(primal, lb, ub, A, Alb, Aub, param.tolerance_);
588 
589  if (std::abs(chosen_ctr_.violation_) < param.tolerance_)
590  {
591  // all constraints are satisfied
592  return_status = OK;
593  break;
594  }
595 
596  chosen_ctr_dot_primal_step_direction = 0.0;
597  factorization_data_.computeInequalityDualStep(
599  if (active_set_.hasEmptySpace())
600  {
601  // compute step direction in primal & dual space
602  factorization_data_.computeInequalityPrimalStep(primal_step_direction_, active_set_);
603  chosen_ctr_dot_primal_step_direction =
605  }
606  }
607  }
608  else
609  {
610  if (dual_blocking_index == primal_size_)
611  {
612  QPMAD_UTILS_THROW("Infeasible inequality constraints.");
613  }
614  else
615  {
616  QPMAD_TRACE("======================");
617  QPMAD_TRACE("||| Dual step length = " << dual_step_length);
618  QPMAD_TRACE("======================");
619 
620  // otherwise -- deactivate
621  dual_.segment(active_set_.num_equalities_, active_set_.num_inequalities_).noalias() +=
622  dual_step_length
623  * dual_step_direction_.segment(
624  active_set_.num_equalities_, active_set_.num_inequalities_);
625  chosen_ctr_.dual_ += dual_step_length;
626 
627  constraints_status_[active_set_.getIndex(dual_blocking_index)] = ConstraintStatus::INACTIVE;
628 
629  dropElementWithoutResize(dual_, dual_blocking_index, active_set_.size_);
630 
631  factorization_data_.downdate(dual_blocking_index, active_set_.size_);
632 
633  active_set_.removeInequality(dual_blocking_index);
634 
635  // compute step direction in primal & dual space
636  factorization_data_.updateStepsAfterPureDualStep(
638  chosen_ctr_dot_primal_step_direction =
640  }
641  }
642  }
643  }
644 
645 #ifdef QPMAD_ENABLE_TRACING
647  {
649 
652  }
653  else
654  {
655  QPMAD_TRACE("||| NO ACTIVE CONSTRAINTS");
656  }
657 #endif
658  return (return_status);
659  }
660 
661 
662  private:
663  template <class t_MatrixType>
664  void initializeMachineryLazy(t_MatrixType &H, const bool return_inverted_cholesky_factor)
665  {
666  if (not machinery_initialized_)
667  {
668  active_set_.initialize(primal_size_);
671 
672  factorization_data_.initialize(H, hessian_type_, primal_size_, return_inverted_cholesky_factor);
673  if (return_inverted_cholesky_factor)
674  {
676  }
677 
678  machinery_initialized_ = true;
679  }
680  }
681 
682 
683  template <
684  class t_Primal,
685  class t_LowerBounds,
686  class t_UpperBounds,
687  class t_Constraints,
688  class t_ConstraintsLowerBounds,
689  class t_ConstraintsUpperBounds>
691  const t_Primal &primal,
692  const t_LowerBounds &lb,
693  const t_UpperBounds &ub,
694  const t_Constraints &A,
695  const t_ConstraintsLowerBounds &Alb,
696  const t_ConstraintsUpperBounds &Aub,
697  const double tolerance)
698  {
699  chosen_ctr_.reset();
700 
701 
702  for (MatrixIndex i = 0; i < num_simple_bounds_; ++i)
703  {
705  {
706  checkConstraintViolation(i, lb(i), ub(i), primal(i));
707  }
708  }
709 
710  if ((std::abs(chosen_ctr_.violation_) < tolerance) && (num_general_constraints_ > 0))
711  {
712  general_ctr_dot_primal_.noalias() = A * primal;
714  {
716  {
718  i,
719  Alb(i - num_simple_bounds_),
720  Aub(i - num_simple_bounds_),
722  }
723  }
725  {
727  }
728  }
729 
732  }
733 
734 
736  const MatrixIndex i,
737  const double lb_i,
738  const double ub_i,
739  const double ctr_i_dot_primal)
740  {
741  double ctr_violation_i = ctr_i_dot_primal - lb_i;
742  if (ctr_violation_i < -std::abs(chosen_ctr_.violation_))
743  {
744  chosen_ctr_.violation_ = ctr_violation_i;
745  chosen_ctr_.index_ = i;
746  }
747  else
748  {
749  ctr_violation_i = ctr_i_dot_primal - ub_i;
750  if (ctr_violation_i > std::abs(chosen_ctr_.violation_))
751  {
752  chosen_ctr_.violation_ = ctr_violation_i;
753  chosen_ctr_.index_ = i;
754  }
755  }
756  }
757 
758 
759  template <class t_VectorType, class t_MatrixType>
760  double getConstraintDotPrimalStepDirection(const t_VectorType &primal_step_direction, const t_MatrixType &A)
761  const
762  {
764  {
765  return (primal_step_direction(chosen_ctr_.index_));
766  }
767  else
768  {
769  return (A.row(chosen_ctr_.general_constraint_index_) * primal_step_direction);
770  }
771  }
772  };
773 
774 
776 } // namespace qpmad
MatrixIndex general_constraint_index_
MatrixIndex num_simple_bounds_
Definition: input_parser.h:20
void parseObjective(const t_MatrixType &H, const t_VectorType &h)
Definition: input_parser.h:35
MatrixIndex h_size_
Definition: input_parser.h:19
MatrixIndex num_general_constraints_
Definition: input_parser.h:21
void parseSimpleBounds(const t_VectorTypelb &lb, const t_VectorTypeub &ub)
Definition: input_parser.h:48
void parseGeneralConstraints(const t_MatrixTypeA &A, const t_VectorTypelb &lb, const t_VectorTypeub &ub)
Definition: input_parser.h:68
MatrixIndex primal_size_
Definition: input_parser.h:18
static const MatrixIndex num_constraints_compile_time_
Definition: solver.h:47
void initializeMachineryLazy(t_MatrixType &H, const bool return_inverted_cholesky_factor)
Definition: solver.h:664
FactorizationData< t_Scalar, t_primal_size > factorization_data_
Definition: solver.h:57
Vector< t_primal_size > dual_
Definition: solver.h:59
void checkConstraintViolation(const MatrixIndex i, const double lb_i, const double ub_i, const double ctr_i_dot_primal)
Definition: solver.h:735
SolverParameters::HessianType getHessianType() const
Returns type of the Hessian produced by the latest execution of solve().
Definition: solver.h:88
SolverParameters::HessianType hessian_type_
Definition: solver.h:72
ReturnStatus solve(Vector< t_rows_primal > &primal, Matrix< t_rows_H, t_cols_H > &H, const Vector< t_rows_h > &h, const Vector< t_rows_lb > &lb, const Vector< t_rows_ub > &ub, const SolverParameters &param)
Definition: solver.h:201
ReturnStatus solve(Vector< t_rows_primal > &primal, Matrix< t_rows_H, t_cols_H > &H, const Vector< t_rows_h > &h, const Vector< t_rows_lb > &lb, const Vector< t_rows_ub > &ub, const Matrix< t_rows_A, t_cols_A > &A, const Vector< t_rows_Alb > &Alb, const Vector< t_rows_Aub > &Aub, const SolverParameters &param)
Definition: solver.h:238
Vector< t_primal_size > dual_step_direction_
Definition: solver.h:62
ReturnStatus solve(Vector< t_rows_primal > &primal, Matrix< t_rows_H, t_cols_H > &H, const Vector< t_rows_h > &h, const Vector< t_rows_lb > &lb, const Vector< t_rows_ub > &ub)
Definition: solver.h:214
ReturnStatus solve(Vector< t_rows_primal > &primal, Matrix< t_rows_H, t_cols_H > &H, const Vector< t_rows_h > &h, const Vector< t_rows_lb > &lb, const Vector< t_rows_ub > &ub, const Matrix< t_rows_A, t_cols_A > &A, const Vector< t_rows_Alb > &Alb, const Vector< t_rows_Aub > &Aub)
Definition: solver.h:186
@ MAXIMAL_NUMBER_OF_ITERATIONS
Definition: solver.h:37
void chooseConstraint(const t_Primal &primal, const t_LowerBounds &lb, const t_UpperBounds &ub, const t_Constraints &A, const t_ConstraintsLowerBounds &Alb, const t_ConstraintsUpperBounds &Aub, const double tolerance)
Definition: solver.h:690
Eigen::Matrix< t_Scalar, t_rows, 1 > Vector
Definition: solver.h:41
std::ptrdiff_t iter_counter_
Definition: solver.h:70
MatrixIndex num_constraints_
Definition: solver.h:53
ChosenConstraint chosen_ctr_
Definition: solver.h:68
ReturnStatus solve(Vector< t_rows_primal > &primal, Matrix< t_rows_H, t_cols_H > &H, const Vector< t_rows_h > &h, const Matrix< t_rows_A, t_cols_A > &A, const Vector< t_rows_Alb > &Alb, const Vector< t_rows_Aub > &Aub)
Definition: solver.h:163
std::ptrdiff_t getNumberOfInequalityIterations() const
Returns number of inequality iterations during the latest execution of solve().
Definition: solver.h:98
Vector< t_general_ctr_number > general_ctr_dot_primal_
Definition: solver.h:64
ActiveSet< t_primal_size > active_set_
Definition: solver.h:56
Eigen::Array< uint8_t, num_constraints_compile_time_, 1 > constraints_status_
Definition: solver.h:66
Vector< t_primal_size > primal_step_direction_
Definition: solver.h:61
bool machinery_initialized_
Definition: solver.h:54
void getInequalityDual(Vector< t_dual_size > &dual, Eigen::Matrix< MatrixIndex, t_index_size, 1 > &indices, Eigen::Matrix< bool, t_status_size, 1 > &is_lower) const
Returns dual variables (Lagrange multipliers) corresponding to inequality constraints....
Definition: solver.h:122
double getConstraintDotPrimalStepDirection(const t_VectorType &primal_step_direction, const t_MatrixType &A) const
Definition: solver.h:760
Eigen::Matrix< t_Scalar, t_rows, t_cols > Matrix
Definition: solver.h:43
#define QPMAD_TRACE(info)
Definition: common.h:26
#define QPMAD_UTILS_PERSISTENT_ASSERT(condition, message)
#define QPMAD_UTILS_THROW(s)
#define QPMAD_UTILS_ASSERT(condition, message)
double computeObjective(const Eigen::MatrixXd &H, const Eigen::VectorXd &h, const Eigen::VectorXd &primal)
Definition: testing.h:21
void printActiveSet(const t_ActiveSet &active_set, const t_ConstraintStatuses &constraints_status, const Eigen::VectorXd &dual)
Definition: testing.h:165
void checkLagrangeMultipliers(const Eigen::MatrixXd &H, const Eigen::VectorXd &h, const Eigen::VectorXd &primal, const Eigen::MatrixXd &A, const t_ActiveSet &active_set, const MatrixIndex &num_simple_bounds, const t_ConstraintStatuses &constraints_status, const Eigen::VectorXd &dual, const Eigen::VectorXd &dual_direction=Eigen::VectorXd())
Definition: testing.h:39
void dropElementWithoutResize(t_VectorType &vector, const MatrixIndex index, const MatrixIndex size)
Definition: common.h:39
EIGEN_DEFAULT_DENSE_INDEX_TYPE MatrixIndex
Definition: common.h:32