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  {
127  const MatrixIndex size = active_set_.size_ - active_set_.num_equalities_;
128 
129  dual.resize(size);
130  indices.resize(size);
131  is_lower.resize(size);
132 
133  for (MatrixIndex i = active_set_.num_equalities_; i < active_set_.size_; ++i)
134  {
135  const std::size_t output_index = i - active_set_.num_equalities_;
136 
137  dual(output_index) = dual_(i);
138  indices(output_index) = active_set_.getIndex(i);
139  is_lower(output_index) =
141  }
142  }
143 
144 
145  template <
146  int t_rows_primal,
147  int t_rows_H,
148  int t_cols_H,
149  int t_rows_h,
150  int t_rows_A,
151  int t_cols_A,
152  int t_rows_Alb,
153  int t_rows_Aub>
155  Vector<t_rows_primal> &primal,
157  const Vector<t_rows_h> &h,
159  const Vector<t_rows_Alb> &Alb,
160  const Vector<t_rows_Aub> &Aub)
161  {
162  return (solve(primal, H, h, Eigen::VectorXd(), Eigen::VectorXd(), A, Alb, Aub, SolverParameters()));
163  }
164 
165 
166  template <
167  int t_rows_primal,
168  int t_rows_H,
169  int t_cols_H,
170  int t_rows_h,
171  int t_rows_lb,
172  int t_rows_ub,
173  int t_rows_A,
174  int t_cols_A,
175  int t_rows_Alb,
176  int t_rows_Aub>
178  Vector<t_rows_primal> &primal,
180  const Vector<t_rows_h> &h,
181  const Vector<t_rows_lb> &lb,
182  const Vector<t_rows_ub> &ub,
184  const Vector<t_rows_Alb> &Alb,
185  const Vector<t_rows_Aub> &Aub)
186  {
187  return (solve(primal, H, h, lb, ub, A, Alb, Aub, SolverParameters()));
188  }
189 
190 
191  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>
193  Vector<t_rows_primal> &primal,
195  const Vector<t_rows_h> &h,
196  const Vector<t_rows_lb> &lb,
197  const Vector<t_rows_ub> &ub,
198  const SolverParameters &param)
199  {
200  return (solve(primal, H, h, lb, ub, Eigen::MatrixXd(), Eigen::VectorXd(), Eigen::VectorXd(), param));
201  }
202 
203 
204  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>
206  Vector<t_rows_primal> &primal,
208  const Vector<t_rows_h> &h,
209  const Vector<t_rows_lb> &lb,
210  const Vector<t_rows_ub> &ub)
211  {
212  return (solve(
213  primal, H, h, lb, ub, Eigen::MatrixXd(), Eigen::VectorXd(), Eigen::VectorXd(), SolverParameters()));
214  }
215 
216 
217 
218  template <
219  int t_rows_primal,
220  int t_rows_H,
221  int t_cols_H,
222  int t_rows_h,
223  int t_rows_lb,
224  int t_rows_ub,
225  int t_rows_A,
226  int t_cols_A,
227  int t_rows_Alb,
228  int t_rows_Aub>
230  Vector<t_rows_primal> &primal,
232  const Vector<t_rows_h> &h,
233  const Vector<t_rows_lb> &lb,
234  const Vector<t_rows_ub> &ub,
236  const Vector<t_rows_Alb> &Alb,
237  const Vector<t_rows_Aub> &Aub,
238  const SolverParameters &param)
239  {
240  QPMAD_TRACE(std::setprecision(std::numeric_limits<double>::digits10));
241 
242  machinery_initialized_ = false;
243  iter_counter_ = 0;
244 
245  parseObjective(H, h);
246  parseSimpleBounds(lb, ub);
247  parseGeneralConstraints(A, Alb, Aub);
248 
249 
251 
252  if (0 == h_size_)
253  {
254  // trivial unconstrained optimum
255  primal.setZero(primal_size_);
256 
257  if (0 == num_constraints_)
258  {
259  // trivial solution
260  return (OK);
261  }
262  }
263 
264 
265  switch (param.hessian_type_)
266  {
268  {
269  const Eigen::LLT<Eigen::Ref<Eigen::MatrixXd>, Eigen::Lower> llt(H);
271  Eigen::Success == llt.info(), "Could not perform Cholesky decomposition of the Hessian.");
272  }
273  // no break here!
274  /* Falls through. */
277  // unconstrained optimum
278  if (h_size_ > 0)
279  {
280  primal = H.template triangularView<Eigen::Lower>().solve(-h);
281  H.transpose().template triangularView<Eigen::Upper>().solveInPlace(primal);
282  }
283  break;
284 
287  // unconstrained optimum
288  if (h_size_ > 0)
289  {
290  primal_step_direction_.noalias() = H.template triangularView<Eigen::Upper>().transpose() * -h;
291  primal.noalias() = H.template triangularView<Eigen::Upper>() * primal_step_direction_;
292  }
293  break;
294 
295  default:
296  QPMAD_UTILS_THROW("Malformed solver parameters!");
297  break;
298  }
299 
300 
301  if (0 == num_constraints_)
302  {
303  // return unconstrained optimum
304  return (OK);
305  }
306 
307 
308  // check consistency of general constraints and activate
309  // equality constraints
311  MatrixIndex num_equalities = 0;
312  for (MatrixIndex i = 0; i < num_constraints_; ++i)
313  {
315 
316  double lb_i;
317  double ub_i;
318 
319  if (true == chosen_ctr_.is_simple_)
320  {
321  lb_i = lb(i);
322  ub_i = ub(i);
323  }
324  else
325  {
329  }
330 
331 
332  if (lb_i - param.tolerance_ > ub_i)
333  {
335  QPMAD_UTILS_THROW("Inconsistent constraints (lb > ub).");
336  }
337 
338  if (std::abs(lb_i - ub_i) > param.tolerance_)
339  {
341  }
342  else
343  {
345  ++num_equalities;
346 
347 
348  if (true == chosen_ctr_.is_simple_)
349  {
350  chosen_ctr_.violation_ = lb_i - primal(i);
351  }
352  else
353  {
355  }
356 
358 
359 
360  // if 'primal_size_' constraints are already activated
361  // all other constraints are linearly dependent
362  if (active_set_.hasEmptySpace())
363  {
364  double ctr_i_dot_primal_step_direction;
365 
366  if (true == chosen_ctr_.is_simple_)
367  {
368  factorization_data_.computeEqualityPrimalStep(primal_step_direction_, i, active_set_.size_);
369 
370  ctr_i_dot_primal_step_direction = primal_step_direction_(i);
371  }
372  else
373  {
374  factorization_data_.computeEqualityPrimalStep(
377  active_set_.size_);
378 
379  ctr_i_dot_primal_step_direction =
381  }
382 
383  // if step direction is a zero vector, constraint is
384  // linearly dependent with previously added constraints
385  if (ctr_i_dot_primal_step_direction < -param.tolerance_)
386  {
387  double primal_step_length_ = chosen_ctr_.violation_ / ctr_i_dot_primal_step_direction;
388 
389  primal.noalias() += primal_step_length_ * primal_step_direction_;
390 
391  if (false
392  == factorization_data_.update(
394  {
395  QPMAD_UTILS_THROW("Failed to add an equality constraint -- is this possible?");
396  }
397  active_set_.addEquality(i);
398 
399  continue;
400  }
401  } // otherwise -- linear dependence
402 
403  // this point is reached if constraint is linearly dependent
404 
405  // check if this constraint is actually satisfied
406  if (std::abs(chosen_ctr_.violation_) > param.tolerance_)
407  {
408  // nope it is not
410  QPMAD_UTILS_THROW("Infeasible equality constraints");
411  }
412  // otherwise keep going
413  }
414  }
415 
416 
417  if (num_equalities == num_constraints_)
418  {
419  // exit early -- avoid unnecessary memory allocations
420  return (OK);
421  }
422 
423 
424  ReturnStatus return_status;
425  chooseConstraint(primal, lb, ub, A, Alb, Aub, param.tolerance_);
426 
427 
428  if (std::abs(chosen_ctr_.violation_) < param.tolerance_)
429  {
430  // all constraints are satisfied
431  return_status = OK;
432  }
433  else
434  {
435  return_status = MAXIMAL_NUMBER_OF_ITERATIONS;
436 
437 
439  dual_.resize(primal_size_);
441 
442 
443  double chosen_ctr_dot_primal_step_direction = 0.0;
444 
445  //
446  factorization_data_.computeInequalityDualStep(dual_step_direction_, chosen_ctr_, A, active_set_);
447  if (active_set_.hasEmptySpace())
448  {
449  // compute step direction in primal space
450  factorization_data_.computeInequalityPrimalStep(primal_step_direction_, active_set_);
451  chosen_ctr_dot_primal_step_direction =
453  }
454 
455 
456  // last iteration is not counted, so iter_counter_ starts with 1.
457  for (iter_counter_ = 1; (param.max_iter_ < 0) or (iter_counter_ <= param.max_iter_); ++iter_counter_)
458  {
459  QPMAD_TRACE(">>>>>>>>>" << iter_counter_ << "<<<<<<<<<");
460 #ifdef QPMAD_ENABLE_TRACING
461  testing::computeObjective(H, h, primal);
462 #endif
463  QPMAD_TRACE("||| Chosen ctr index = " << chosen_ctr_.index_);
464  QPMAD_TRACE("||| Chosen ctr dual = " << chosen_ctr_.dual_);
465  QPMAD_TRACE("||| Chosen ctr violation = " << chosen_ctr_.violation_);
466 
467 
468  // check dual feasibility
469  MatrixIndex dual_blocking_index = primal_size_;
470  double dual_step_length = std::numeric_limits<double>::infinity();
471  for (MatrixIndex i = active_set_.num_equalities_; i < active_set_.size_; ++i)
472  {
473  if (dual_step_direction_(i) < -param.tolerance_)
474  {
475  const double dual_step_length_i = -dual_(i) / dual_step_direction_(i);
476  if (dual_step_length_i < dual_step_length)
477  {
478  dual_step_length = dual_step_length_i;
479  dual_blocking_index = i;
480  }
481  }
482  }
483 
484 
485 #ifdef QPMAD_ENABLE_TRACING
487  H,
488  h,
489  primal,
490  A,
491  active_set_,
494  dual_,
496 #endif
497 
498 
499  if (active_set_.hasEmptySpace()
500  // if step direction is a zero vector, constraint is
501  // linearly dependent with previously added constraints
502  && (std::abs(chosen_ctr_dot_primal_step_direction) > param.tolerance_))
503  {
504  double step_length = -chosen_ctr_.violation_ / chosen_ctr_dot_primal_step_direction;
505 
506  QPMAD_TRACE("======================");
507  QPMAD_TRACE("||| Primal step length = " << step_length);
508  QPMAD_TRACE("||| Dual step length = " << dual_step_length);
509  QPMAD_TRACE("======================");
510 
511 
512  bool partial_step = false;
514  (step_length >= 0.0) && (dual_step_length >= 0.0),
515  "Non-negative step lengths expected.");
516  if (dual_step_length <= step_length)
517  {
518  step_length = dual_step_length;
519  partial_step = true;
520  }
521 
522 
523  primal.noalias() += step_length * primal_step_direction_;
524 
525  dual_.segment(active_set_.num_equalities_, active_set_.num_inequalities_).noalias() +=
526  step_length
527  * dual_step_direction_.segment(
528  active_set_.num_equalities_, active_set_.num_inequalities_);
529  chosen_ctr_.dual_ += step_length;
530  chosen_ctr_.violation_ += step_length * chosen_ctr_dot_primal_step_direction;
531 
532  QPMAD_TRACE("||| Chosen ctr dual = " << chosen_ctr_.dual_);
533  QPMAD_TRACE("||| Chosen ctr violation = " << chosen_ctr_.violation_);
534 
535 
536  if ((partial_step)
537  // if violation is almost zero -- assume that a full step is made
538  && (std::abs(chosen_ctr_.violation_) > param.tolerance_))
539  {
540  QPMAD_TRACE("||| PARTIAL STEP");
541  // deactivate blocking constraint
542  constraints_status_[active_set_.getIndex(dual_blocking_index)] = ConstraintStatus::INACTIVE;
543 
544  dropElementWithoutResize(dual_, dual_blocking_index, active_set_.size_);
545 
546  factorization_data_.downdate(dual_blocking_index, active_set_.size_);
547 
548  active_set_.removeInequality(dual_blocking_index);
549 
550  // compute step direction in primal & dual space
551  factorization_data_.updateStepsAfterPartialStep(
553  chosen_ctr_dot_primal_step_direction =
555  }
556  else
557  {
558  QPMAD_TRACE("||| FULL STEP");
559  // activate constraint
560  if (false
561  == factorization_data_.update(
563  {
564  QPMAD_UTILS_THROW("Failed to add an inequality constraint -- is this possible?");
565  }
566 
568  {
570  }
571  else
572  {
574  }
576  active_set_.addInequality(chosen_ctr_.index_);
577 
578  chooseConstraint(primal, lb, ub, A, Alb, Aub, param.tolerance_);
579 
580  if (std::abs(chosen_ctr_.violation_) < param.tolerance_)
581  {
582  // all constraints are satisfied
583  return_status = OK;
584  break;
585  }
586 
587  chosen_ctr_dot_primal_step_direction = 0.0;
588  factorization_data_.computeInequalityDualStep(
590  if (active_set_.hasEmptySpace())
591  {
592  // compute step direction in primal & dual space
593  factorization_data_.computeInequalityPrimalStep(primal_step_direction_, active_set_);
594  chosen_ctr_dot_primal_step_direction =
596  }
597  }
598  }
599  else
600  {
601  if (dual_blocking_index == primal_size_)
602  {
603  QPMAD_UTILS_THROW("Infeasible inequality constraints.");
604  }
605  else
606  {
607  QPMAD_TRACE("======================");
608  QPMAD_TRACE("||| Dual step length = " << dual_step_length);
609  QPMAD_TRACE("======================");
610 
611  // otherwise -- deactivate
612  dual_.segment(active_set_.num_equalities_, active_set_.num_inequalities_).noalias() +=
613  dual_step_length
614  * dual_step_direction_.segment(
615  active_set_.num_equalities_, active_set_.num_inequalities_);
616  chosen_ctr_.dual_ += dual_step_length;
617 
618  constraints_status_[active_set_.getIndex(dual_blocking_index)] = ConstraintStatus::INACTIVE;
619 
620  dropElementWithoutResize(dual_, dual_blocking_index, active_set_.size_);
621 
622  factorization_data_.downdate(dual_blocking_index, active_set_.size_);
623 
624  active_set_.removeInequality(dual_blocking_index);
625 
626  // compute step direction in primal & dual space
627  factorization_data_.updateStepsAfterPureDualStep(
629  chosen_ctr_dot_primal_step_direction =
631  }
632  }
633  }
634  }
635 
636 #ifdef QPMAD_ENABLE_TRACING
638  {
640 
643  }
644  else
645  {
646  QPMAD_TRACE("||| NO ACTIVE CONSTRAINTS");
647  }
648 #endif
649  return (return_status);
650  }
651 
652 
653  private:
654  template <class t_MatrixType>
655  void initializeMachineryLazy(t_MatrixType &H, const bool return_inverted_cholesky_factor)
656  {
657  if (not machinery_initialized_)
658  {
659  active_set_.initialize(primal_size_);
662 
663  factorization_data_.initialize(H, hessian_type_, primal_size_, return_inverted_cholesky_factor);
664  if (return_inverted_cholesky_factor)
665  {
667  }
668 
669  machinery_initialized_ = true;
670  }
671  }
672 
673 
674  template <
675  class t_Primal,
676  class t_LowerBounds,
677  class t_UpperBounds,
678  class t_Constraints,
679  class t_ConstraintsLowerBounds,
680  class t_ConstraintsUpperBounds>
682  const t_Primal &primal,
683  const t_LowerBounds &lb,
684  const t_UpperBounds &ub,
685  const t_Constraints &A,
686  const t_ConstraintsLowerBounds &Alb,
687  const t_ConstraintsUpperBounds &Aub,
688  const double tolerance)
689  {
690  chosen_ctr_.reset();
691 
692 
693  for (MatrixIndex i = 0; i < num_simple_bounds_; ++i)
694  {
696  {
697  checkConstraintViolation(i, lb(i), ub(i), primal(i));
698  }
699  }
700 
701  if ((std::abs(chosen_ctr_.violation_) < tolerance) && (num_general_constraints_ > 0))
702  {
703  general_ctr_dot_primal_.noalias() = A * primal;
705  {
707  {
709  i,
710  Alb(i - num_simple_bounds_),
711  Aub(i - num_simple_bounds_),
713  }
714  }
716  {
718  }
719  }
720 
723  }
724 
725 
727  const MatrixIndex i,
728  const double lb_i,
729  const double ub_i,
730  const double ctr_i_dot_primal)
731  {
732  double ctr_violation_i = ctr_i_dot_primal - lb_i;
733  if (ctr_violation_i < -std::abs(chosen_ctr_.violation_))
734  {
735  chosen_ctr_.violation_ = ctr_violation_i;
736  chosen_ctr_.index_ = i;
737  }
738  else
739  {
740  ctr_violation_i = ctr_i_dot_primal - ub_i;
741  if (ctr_violation_i > std::abs(chosen_ctr_.violation_))
742  {
743  chosen_ctr_.violation_ = ctr_violation_i;
744  chosen_ctr_.index_ = i;
745  }
746  }
747  }
748 
749 
750  template <class t_VectorType, class t_MatrixType>
751  double getConstraintDotPrimalStepDirection(const t_VectorType &primal_step_direction, const t_MatrixType &A)
752  const
753  {
755  {
756  return (primal_step_direction(chosen_ctr_.index_));
757  }
758  else
759  {
760  return (A.row(chosen_ctr_.general_constraint_index_) * primal_step_direction);
761  }
762  }
763  };
764 
765 
767 } // namespace qpmad
#define QPMAD_UTILS_THROW(s)
Eigen::Array< uint8_t, num_constraints_compile_time_, 1 > constraints_status_
Definition: solver.h:66
double computeObjective(const Eigen::MatrixXd &H, const Eigen::VectorXd &h, const Eigen::VectorXd &primal)
Definition: testing.h:21
bool machinery_initialized_
Definition: solver.h:54
Vector< t_primal_size > primal_step_direction_
Definition: solver.h:61
#define QPMAD_TRACE(info)
Definition: common.h:26
MatrixIndex general_constraint_index_
double getConstraintDotPrimalStepDirection(const t_VectorType &primal_step_direction, const t_MatrixType &A) const
Definition: solver.h:751
void parseObjective(const t_MatrixType &H, const t_VectorType &h)
Definition: input_parser.h:35
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
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:205
MatrixIndex num_simple_bounds_
Definition: input_parser.h:20
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:229
#define QPMAD_UTILS_PERSISTENT_ASSERT(condition, message)
Eigen::Matrix< t_Scalar, t_rows, t_cols > Matrix
Definition: solver.h:43
EIGEN_DEFAULT_DENSE_INDEX_TYPE MatrixIndex
Definition: common.h:32
void initializeMachineryLazy(t_MatrixType &H, const bool return_inverted_cholesky_factor)
Definition: solver.h:655
MatrixIndex h_size_
Definition: input_parser.h:19
SolverParameters::HessianType hessian_type_
Definition: solver.h:72
Vector< t_primal_size > dual_step_direction_
Definition: solver.h:62
Vector< t_primal_size > dual_
Definition: solver.h:59
MatrixIndex num_general_constraints_
Definition: input_parser.h:21
static const MatrixIndex num_constraints_compile_time_
Definition: solver.h:47
Eigen::Matrix< t_Scalar, t_rows, 1 > Vector
Definition: solver.h:41
std::ptrdiff_t iter_counter_
Definition: solver.h:70
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:192
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:154
ActiveSet< t_primal_size > active_set_
Definition: solver.h:56
void checkConstraintViolation(const MatrixIndex i, const double lb_i, const double ub_i, const double ctr_i_dot_primal)
Definition: solver.h:726
Vector< t_general_ctr_number > general_ctr_dot_primal_
Definition: solver.h:64
std::ptrdiff_t getNumberOfInequalityIterations() const
Returns number of inequality iterations during the latest execution of solve().
Definition: solver.h:98
void printActiveSet(const t_ActiveSet &active_set, const t_ConstraintStatuses &constraints_status, const Eigen::VectorXd &dual)
Definition: testing.h:165
#define QPMAD_UTILS_ASSERT(condition, message)
MatrixIndex primal_size_
Definition: input_parser.h:18
MatrixIndex num_constraints_
Definition: solver.h:53
SolverParameters::HessianType getHessianType() const
Returns type of the Hessian produced by the latest execution of solve().
Definition: solver.h:88
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:681
void parseGeneralConstraints(const t_MatrixTypeA &A, const t_VectorTypelb &lb, const t_VectorTypeub &ub)
Definition: input_parser.h:68
FactorizationData< t_Scalar, t_primal_size > factorization_data_
Definition: solver.h:57
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:177
void dropElementWithoutResize(t_VectorType &vector, const MatrixIndex index, const MatrixIndex size)
Definition: common.h:39
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 parseSimpleBounds(const t_VectorTypelb &lb, const t_VectorTypeub &ub)
Definition: input_parser.h:48
ChosenConstraint chosen_ctr_
Definition: solver.h:68