proxsuite 0.6.7
The Advanced Proximal Optimization Toolbox
Loading...
Searching...
No Matches
wrapper.hpp
Go to the documentation of this file.
1//
2// Copyright (c) 2022 INRIA
3//
8#ifndef PROXSUITE_PROXQP_DENSE_WRAPPER_HPP
9#define PROXSUITE_PROXQP_DENSE_WRAPPER_HPP
14#include <chrono>
15
16namespace proxsuite {
17namespace proxqp {
18namespace dense {
22
81template<typename T>
84 isize dim,
85 isize n_eq,
86 isize n_in,
87 bool box_constraints)
88{
89 if (_dense_backend == DenseBackend::Automatic) {
90 isize n_constraints(n_in);
91 if (box_constraints) {
92 n_constraints += dim;
93 }
94 T threshold(1.5);
95 T frequence(0.2);
96 T PrimalDualLDLTCost =
97 0.5 * std::pow(T(n_eq) / T(dim), 2) +
98 0.17 * (std::pow(T(n_eq) / T(dim), 3) +
99 std::pow(T(n_constraints) / T(dim), 3)) +
100 frequence * std::pow(T(n_eq + n_constraints) / T(dim), 2) / T(dim);
101 T PrimalLDLTCost =
102 threshold *
103 ((0.5 * T(n_eq) + T(n_constraints)) / T(dim) + frequence / T(dim));
104 bool choice = PrimalDualLDLTCost > PrimalLDLTCost;
105 if (choice) {
107 } else {
109 }
110 } else {
111 return _dense_backend;
112 }
113}
114template<typename T>
115struct QP
116{
117private:
118 // structure of the problem
119 // not supposed to change
120 DenseBackend dense_backend;
121 bool box_constraints;
122 HessianType hessian_type;
123
124public:
130
140 QP(isize _dim,
141 isize _n_eq,
142 isize _n_in,
143 bool _box_constraints,
144 proxsuite::proxqp::HessianType _hessian_type,
145 DenseBackend _dense_backend)
146 : dense_backend(dense_backend_choice<T>(_dense_backend,
147 _dim,
148 _n_eq,
149 _n_in,
150 _box_constraints))
151 , box_constraints(_box_constraints)
152 , hessian_type(_hessian_type)
153 , results(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
154 , settings(dense_backend)
155 , model(_dim, _n_eq, _n_in, _box_constraints)
156 , work(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
157 , ruiz(preconditioner::RuizEquilibration<T>{ _dim,
158 _n_eq,
159 _n_in,
160 _box_constraints })
161 {
162 work.timer.stop();
163 }
173 QP(isize _dim,
174 isize _n_eq,
175 isize _n_in,
176 bool _box_constraints,
177 DenseBackend _dense_backend,
178 proxsuite::proxqp::HessianType _hessian_type)
179 : dense_backend(dense_backend_choice<T>(_dense_backend,
180 _dim,
181 _n_eq,
182 _n_in,
183 _box_constraints))
184 , box_constraints(_box_constraints)
185 , hessian_type(_hessian_type)
186 , results(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
187 , settings(dense_backend)
188 , model(_dim, _n_eq, _n_in, _box_constraints)
189 , work(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
190 , ruiz(preconditioner::RuizEquilibration<T>{ _dim,
191 _n_eq,
192 _n_in,
193 _box_constraints })
194 {
195 work.timer.stop();
196 }
205 QP(isize _dim,
206 isize _n_eq,
207 isize _n_in,
208 bool _box_constraints,
209 proxsuite::proxqp::HessianType _hessian_type)
210 : dense_backend(dense_backend_choice<T>(DenseBackend::Automatic,
211 _dim,
212 _n_eq,
213 _n_in,
214 _box_constraints))
215 , box_constraints(_box_constraints)
216 , hessian_type(_hessian_type)
217 , results(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
218 , settings(dense_backend)
219 , model(_dim, _n_eq, _n_in, _box_constraints)
220 , work(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
221 , ruiz(preconditioner::RuizEquilibration<T>{ _dim,
222 _n_eq,
223 _n_in,
224 _box_constraints })
225 {
226 work.timer.stop();
227 }
237 QP(isize _dim,
238 isize _n_eq,
239 isize _n_in,
240 bool _box_constraints,
241 DenseBackend _dense_backend)
242 : dense_backend(dense_backend_choice<T>(_dense_backend,
243 _dim,
244 _n_eq,
245 _n_in,
246 _box_constraints))
247 , box_constraints(_box_constraints)
248 , hessian_type(HessianType::Dense)
249 , results(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
250 , settings(dense_backend)
251 , model(_dim, _n_eq, _n_in, _box_constraints)
252 , work(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
253 , ruiz(preconditioner::RuizEquilibration<T>{ _dim,
254 _n_eq,
255 _n_in,
256 _box_constraints })
257 {
258 work.timer.stop();
259 }
267 QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints)
268 : dense_backend(dense_backend_choice<T>(DenseBackend::Automatic,
269 _dim,
270 _n_eq,
271 _n_in,
272 _box_constraints))
273 , box_constraints(_box_constraints)
274 , hessian_type(proxsuite::proxqp::HessianType::Dense)
275 , results(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
276 , settings(dense_backend)
277 , model(_dim, _n_eq, _n_in, _box_constraints)
278 , work(_dim, _n_eq, _n_in, _box_constraints, dense_backend)
279 , ruiz(preconditioner::RuizEquilibration<T>{ _dim,
280 _n_eq,
281 _n_in,
282 _box_constraints })
283 {
284 work.timer.stop();
285 }
293 QP(isize _dim,
294 isize _n_eq,
295 isize _n_in,
296 proxsuite::proxqp::HessianType _hessian_type)
297 : dense_backend(dense_backend_choice<T>(DenseBackend::Automatic,
298 _dim,
299 _n_eq,
300 _n_in,
301 false))
302 , box_constraints(false)
303 , hessian_type(_hessian_type)
304 , results(_dim, _n_eq, _n_in, false, dense_backend)
305 , settings(dense_backend)
306 , model(_dim, _n_eq, _n_in, false)
307 , work(_dim, _n_eq, _n_in, false, dense_backend)
308 , ruiz(preconditioner::RuizEquilibration<T>{ _dim, _n_eq, _n_in, false })
309 {
310 work.timer.stop();
311 }
318 QP(isize _dim, isize _n_eq, isize _n_in)
319 : dense_backend(dense_backend_choice<T>(DenseBackend::Automatic,
320 _dim,
321 _n_eq,
322 _n_in,
323 false))
324 , box_constraints(false)
325 , hessian_type(proxsuite::proxqp::HessianType::Dense)
326 , results(_dim, _n_eq, _n_in, false, dense_backend)
327 , settings(dense_backend)
328 , model(_dim, _n_eq, _n_in, false)
329 , work(_dim, _n_eq, _n_in, false, dense_backend)
330 , ruiz(preconditioner::RuizEquilibration<T>{ _dim, _n_eq, _n_in, false })
331 {
332 work.timer.stop();
333 }
334 bool is_box_constrained() const { return box_constraints; };
335 DenseBackend which_dense_backend() const { return dense_backend; };
336 HessianType which_hessian_type() const { return hessian_type; };
361 bool compute_preconditioner = true,
362 optional<T> rho = nullopt,
363 optional<T> mu_eq = nullopt,
364 optional<T> mu_in = nullopt,
365 optional<T> manual_minimal_H_eigenvalue = nullopt)
366 {
368 box_constraints == true,
369 std::invalid_argument,
370 "wrong model setup: the QP object is designed with box "
371 "constraints, but is initialized without lower or upper box "
372 "inequalities.");
373 // dense case
374 if (settings.compute_timings) {
375 work.timer.stop();
376 work.timer.start();
377 }
378 settings.compute_preconditioner = compute_preconditioner;
379 // check the model is valid
380 if (g != nullopt && g.value().size() != 0) {
382 g.value().size(),
383 model.dim,
384 "the dimension wrt the primal variable x variable for initializing g "
385 "is not valid.");
386 } else {
387 g.reset();
388 }
389 if (b != nullopt && b.value().size() != 0) {
391 b.value().size(),
392 model.n_eq,
393 "the dimension wrt equality constrained variables for initializing b "
394 "is not valid.");
395 } else {
396 b.reset();
397 }
398 if (u != nullopt && u.value().size() != 0) {
400 u.value().size(),
401 model.n_in,
402 "the dimension wrt inequality constrained variables for initializing u "
403 "is not valid.");
404 } else {
405 u.reset();
406 }
407 if (l != nullopt && l.value().size() != 0) {
409 l.value().size(),
410 model.n_in,
411 "the dimension wrt inequality constrained variables for initializing l "
412 "is not valid.");
413 } else {
414 l.reset();
415 }
416 if (H != nullopt && H.value().size() != 0) {
418 H.value().rows(),
419 model.dim,
420 "the row dimension for initializing H is not valid.");
422 H.value().cols(),
423 model.dim,
424 "the column dimension for initializing H is not valid.");
425 } else {
426 H.reset();
427 }
428 if (A != nullopt && A.value().size() != 0) {
430 A.value().rows(),
431 model.n_eq,
432 "the row dimension for initializing A is not valid.");
434 A.value().cols(),
435 model.dim,
436 "the column dimension for initializing A is not valid.");
437 } else {
438 A.reset();
439 }
440 if (C != nullopt && C.value().size() != 0) {
442 C.value().rows(),
443 model.n_in,
444 "the row dimension for initializing C is not valid.");
446 C.value().cols(),
447 model.dim,
448 "the column dimension for initializing C is not valid.");
449 } else {
450 C.reset();
451 }
452 if (settings.initial_guess ==
454 work.refactorize =
455 true; // necessary for the first solve (then refactorize only if there
456 // is an update of the matrices)
457 } else {
458 work.refactorize = false;
459 }
460 work.proximal_parameter_update = false;
461 if (settings.compute_timings) {
462 work.timer.stop();
463 work.timer.start();
464 }
465 PreconditionerStatus preconditioner_status;
466 if (compute_preconditioner) {
468 } else {
470 }
472 settings, results, work, rho, mu_eq, mu_in);
475 manual_minimal_H_eigenvalue, results, settings);
476 typedef optional<VecRef<T>> optional_VecRef;
478 g,
479 A,
480 b,
481 C,
482 l,
483 u,
484 optional_VecRef(nullopt),
485 optional_VecRef(nullopt),
486 settings,
487 model,
488 work,
489 results,
490 box_constraints,
491 ruiz,
492 preconditioner_status,
493 hessian_type);
494 work.is_initialized = true;
495 if (settings.compute_timings) {
496 results.info.setup_time = work.timer.elapsed().user; // in microseconds
497 }
498 };
527 optional<VecRef<T>> l_box,
528 optional<VecRef<T>> u_box,
529 bool compute_preconditioner = true,
530 optional<T> rho = nullopt,
531 optional<T> mu_eq = nullopt,
532 optional<T> mu_in = nullopt,
533 optional<T> manual_minimal_H_eigenvalue = nullopt)
534 {
535
536 // dense case
537 if (settings.compute_timings) {
538 work.timer.stop();
539 work.timer.start();
540 }
541 settings.compute_preconditioner = compute_preconditioner;
543 box_constraints == false && (l_box != nullopt || u_box != nullopt),
544 std::invalid_argument,
545 "wrong model setup: the QP object is designed without box "
546 "constraints, but is initialized with lower or upper box inequalities.");
547 if (l_box != nullopt && l_box.value().size() != 0) {
548 PROXSUITE_CHECK_ARGUMENT_SIZE(l_box.value().size(),
549 model.dim,
550 "the dimension wrt the primal variable x "
551 "variable for initializing l_box "
552 "is not valid.");
553 } else {
554 l_box.reset();
555 }
556 if (u_box != nullopt && u_box.value().size() != 0) {
557 PROXSUITE_CHECK_ARGUMENT_SIZE(u_box.value().size(),
558 model.dim,
559 "the dimension wrt the primal variable x "
560 "variable for initializing u_box "
561 "is not valid.");
562 } else {
563 l_box.reset();
564 }
565 // check the model is valid
566 if (g != nullopt && g.value().size() != 0) {
568 g.value().size(),
569 model.dim,
570 "the dimension wrt the primal variable x variable for initializing g "
571 "is not valid.");
572 } else {
573 g.reset();
574 }
575 if (b != nullopt && b.value().size() != 0) {
577 b.value().size(),
578 model.n_eq,
579 "the dimension wrt equality constrained variables for initializing b "
580 "is not valid.");
581 } else {
582 b.reset();
583 }
584 if (u != nullopt && u.value().size() != 0) {
586 u.value().size(),
587 model.n_in,
588 "the dimension wrt inequality constrained variables for initializing u "
589 "is not valid.");
590 } else {
591 u.reset();
592 }
593 if (u_box != nullopt && u_box.value().size() != 0) {
595 u_box.value().size(),
596 model.dim,
597 "the dimension wrt box inequality constrained variables for "
598 "initializing u_box "
599 "is not valid.");
600 } else {
601 u_box.reset();
602 }
603 if (l != nullopt && l.value().size() != 0) {
605 l.value().size(),
606 model.n_in,
607 "the dimension wrt inequality constrained variables for initializing l "
608 "is not valid.");
609 } else {
610 l.reset();
611 }
612 if (l_box != nullopt && l_box.value().size() != 0) {
614 l_box.value().size(),
615 model.dim,
616 "the dimension wrt box inequality constrained variables for "
617 "initializing l_box "
618 "is not valid.");
619 } else {
620 l_box.reset();
621 }
622 if (H != nullopt && H.value().size() != 0) {
624 H.value().rows(),
625 model.dim,
626 "the row dimension for initializing H is not valid.");
628 H.value().cols(),
629 model.dim,
630 "the column dimension for initializing H is not valid.");
631 } else {
632 H.reset();
633 }
634 if (A != nullopt && A.value().size() != 0) {
636 A.value().rows(),
637 model.n_eq,
638 "the row dimension for initializing A is not valid.");
640 A.value().cols(),
641 model.dim,
642 "the column dimension for initializing A is not valid.");
643 } else {
644 A.reset();
645 }
646 if (C != nullopt && C.value().size() != 0) {
648 C.value().rows(),
649 model.n_in,
650 "the row dimension for initializing C is not valid.");
652 C.value().cols(),
653 model.dim,
654 "the column dimension for initializing C is not valid.");
655 } else {
656 C.reset();
657 }
658 if (settings.initial_guess ==
660 work.refactorize =
661 true; // necessary for the first solve (then refactorize only if there
662 // is an update of the matrices)
663 } else {
664 work.refactorize = false;
665 }
666 work.proximal_parameter_update = false;
667 if (settings.compute_timings) {
668 work.timer.stop();
669 work.timer.start();
670 }
671 PreconditionerStatus preconditioner_status;
672 if (compute_preconditioner) {
674 } else {
676 }
678 settings, results, work, rho, mu_eq, mu_in);
681 manual_minimal_H_eigenvalue, results, settings);
683 g,
684 A,
685 b,
686 C,
687 l,
688 u,
689 l_box,
690 u_box,
691 settings,
692 model,
693 work,
694 results,
695 box_constraints,
696 ruiz,
697 preconditioner_status,
698 hessian_type);
699 work.is_initialized = true;
700 if (settings.compute_timings) {
701 results.info.setup_time = work.timer.elapsed().user; // in microseconds
702 }
703 };
730 bool update_preconditioner = false,
731 optional<T> rho = nullopt,
732 optional<T> mu_eq = nullopt,
733 optional<T> mu_in = nullopt,
734 optional<T> manual_minimal_H_eigenvalue = nullopt)
735 {
737 box_constraints == true,
738 std::invalid_argument,
739 "wrong model setup: the QP object is designed without box "
740 "constraints, but the update does not include lower or upper box "
741 "inequalities.");
742 settings.update_preconditioner = update_preconditioner;
743 if (!work.is_initialized) {
744 init(H, g, A, b, C, l, u, update_preconditioner, rho, mu_eq, mu_in);
745 return;
746 }
747 // dense case
748 work.refactorize = false;
749 work.proximal_parameter_update = false;
750 if (settings.compute_timings) {
751 work.timer.stop();
752 work.timer.start();
753 }
754 PreconditionerStatus preconditioner_status;
755 if (update_preconditioner) {
757 } else {
758 preconditioner_status = proxsuite::proxqp::PreconditionerStatus::KEEP;
759 }
760 const bool matrix_update =
761 !(H == nullopt && g == nullopt && A == nullopt && b == nullopt &&
762 C == nullopt && u == nullopt && l == nullopt);
763 if (matrix_update) {
764 typedef optional<VecRef<T>> optional_VecRef;
766 g,
767 A,
768 b,
769 C,
770 l,
771 u,
772 optional_VecRef(nullopt),
773 optional_VecRef(nullopt),
774 model,
775 work,
776 box_constraints);
777 }
779 settings, results, work, rho, mu_eq, mu_in);
782 manual_minimal_H_eigenvalue, results, settings);
783 typedef optional<MatRef<T>> optional_MatRef;
784 typedef optional<VecRef<T>> optional_VecRef;
785 proxsuite::proxqp::dense::setup(/* avoid double assignation */
786 optional_MatRef(nullopt),
787 optional_VecRef(nullopt),
788 optional_MatRef(nullopt),
789 optional_VecRef(nullopt),
790 optional_MatRef(nullopt),
791 optional_VecRef(nullopt),
792 optional_VecRef(nullopt),
793 optional_VecRef(nullopt),
794 optional_VecRef(nullopt),
795 settings,
796 model,
797 work,
798 results,
799 box_constraints,
800 ruiz,
801 preconditioner_status,
802 hessian_type);
803
804 if (settings.compute_timings) {
805 results.info.setup_time = work.timer.elapsed().user; // in microseconds
806 }
807 };
838 optional<VecRef<T>> l_box,
839 optional<VecRef<T>> u_box,
840 bool update_preconditioner = false,
841 optional<T> rho = nullopt,
842 optional<T> mu_eq = nullopt,
843 optional<T> mu_in = nullopt,
844 optional<T> manual_minimal_H_eigenvalue = nullopt)
845 {
847 box_constraints == false && (l_box != nullopt || u_box != nullopt),
848 std::invalid_argument,
849 "wrong model setup: the QP object is designed without box "
850 "constraints, but the update includes lower or upper box inequalities.");
851 settings.update_preconditioner = update_preconditioner;
852 if (!work.is_initialized) {
853 init(H,
854 g,
855 A,
856 b,
857 C,
858 l,
859 u,
860 l_box,
861 u_box,
862 update_preconditioner,
863 rho,
864 mu_eq,
865 mu_in);
866 return;
867 }
868 // dense case
869 work.refactorize = false;
870 work.proximal_parameter_update = false;
871 if (settings.compute_timings) {
872 work.timer.stop();
873 work.timer.start();
874 }
875 PreconditionerStatus preconditioner_status;
876 if (update_preconditioner) {
878 } else {
879 preconditioner_status = proxsuite::proxqp::PreconditionerStatus::KEEP;
880 }
881 const bool matrix_update =
882 !(H == nullopt && g == nullopt && A == nullopt && b == nullopt &&
883 C == nullopt && u == nullopt && l == nullopt && u_box == nullopt &&
884 l_box == nullopt);
885 if (matrix_update) {
887 H, g, A, b, C, l, u, l_box, u_box, model, work, box_constraints);
888 }
890 settings, results, work, rho, mu_eq, mu_in);
893 manual_minimal_H_eigenvalue, results, settings);
894 typedef optional<MatRef<T>> optional_MatRef;
895 typedef optional<VecRef<T>> optional_VecRef;
896 proxsuite::proxqp::dense::setup(/* avoid double assignation */
897 optional_MatRef(nullopt),
898 optional_VecRef(nullopt),
899 optional_MatRef(nullopt),
900 optional_VecRef(nullopt),
901 optional_MatRef(nullopt),
902 optional_VecRef(nullopt),
903 optional_VecRef(nullopt),
904 optional_VecRef(nullopt),
905 optional_VecRef(nullopt),
906 settings,
907 model,
908 work,
909 results,
910 box_constraints,
911 ruiz,
912 preconditioner_status,
913 hessian_type);
914
915 if (settings.compute_timings) {
916 results.info.setup_time = work.timer.elapsed().user; // in microseconds
917 }
918 };
922 void solve()
923 {
924 qp_solve( //
925 settings,
926 model,
927 results,
928 work,
929 box_constraints,
930 dense_backend,
931 hessian_type,
932 ruiz);
933 };
943 {
945 qp_solve( //
946 settings,
947 model,
948 results,
949 work,
950 box_constraints,
951 dense_backend,
952 hessian_type,
953 ruiz);
954 };
958 void cleanup()
959 {
960 results.cleanup(settings);
961 work.cleanup(box_constraints);
962 }
963};
1000template<typename T>
1013 optional<T> eps_abs = nullopt,
1014 optional<T> eps_rel = nullopt,
1015 optional<T> rho = nullopt,
1016 optional<T> mu_eq = nullopt,
1017 optional<T> mu_in = nullopt,
1018 optional<bool> verbose = nullopt,
1019 bool compute_preconditioner = true,
1020 bool compute_timings = false,
1021 optional<isize> max_iter = nullopt,
1024 bool check_duality_gap = false,
1025 optional<T> eps_duality_gap_abs = nullopt,
1026 optional<T> eps_duality_gap_rel = nullopt,
1027 bool primal_infeasibility_solving = false,
1028 optional<T> manual_minimal_H_eigenvalue = nullopt)
1029{
1030 isize n(0);
1031 isize n_eq(0);
1032 isize n_in(0);
1033 if (H != nullopt) {
1034 n = H.value().rows();
1035 }
1036 if (A != nullopt) {
1037 n_eq = A.value().rows();
1038 }
1039 if (C != nullopt) {
1040 n_in = C.value().rows();
1041 }
1042
1043 QP<T> Qp(n, n_eq, n_in, false, DenseBackend::PrimalDualLDLT);
1044 Qp.settings.initial_guess = initial_guess;
1045 Qp.settings.check_duality_gap = check_duality_gap;
1046
1047 if (eps_abs != nullopt) {
1048 Qp.settings.eps_abs = eps_abs.value();
1049 }
1050 if (eps_rel != nullopt) {
1051 Qp.settings.eps_rel = eps_rel.value();
1052 }
1053 if (verbose != nullopt) {
1054 Qp.settings.verbose = verbose.value();
1055 }
1056 if (max_iter != nullopt) {
1057 Qp.settings.max_iter = max_iter.value();
1058 }
1059 if (eps_duality_gap_abs != nullopt) {
1060 Qp.settings.eps_duality_gap_abs = eps_duality_gap_abs.value();
1061 }
1062 if (eps_duality_gap_rel != nullopt) {
1063 Qp.settings.eps_duality_gap_rel = eps_duality_gap_rel.value();
1064 }
1065 Qp.settings.compute_timings = compute_timings;
1066 Qp.settings.primal_infeasibility_solving = primal_infeasibility_solving;
1067 if (manual_minimal_H_eigenvalue != nullopt) {
1068 Qp.init(H,
1069 g,
1070 A,
1071 b,
1072 C,
1073 l,
1074 u,
1075 compute_preconditioner,
1076 rho,
1077 mu_eq,
1078 mu_in,
1079 manual_minimal_H_eigenvalue.value());
1080 } else {
1081 Qp.init(
1082 H, g, A, b, C, l, u, compute_preconditioner, rho, mu_eq, mu_in, nullopt);
1083 }
1084 Qp.solve(x, y, z);
1085
1086 return Qp.results;
1087}
1130template<typename T>
1140 optional<VecRef<T>> l_box,
1141 optional<VecRef<T>> u_box,
1145 optional<T> eps_abs = nullopt,
1146 optional<T> eps_rel = nullopt,
1147 optional<T> rho = nullopt,
1148 optional<T> mu_eq = nullopt,
1149 optional<T> mu_in = nullopt,
1150 optional<bool> verbose = nullopt,
1151 bool compute_preconditioner = true,
1152 bool compute_timings = false,
1153 optional<isize> max_iter = nullopt,
1156 bool check_duality_gap = false,
1157 optional<T> eps_duality_gap_abs = nullopt,
1158 optional<T> eps_duality_gap_rel = nullopt,
1159 bool primal_infeasibility_solving = false,
1160 optional<T> manual_minimal_H_eigenvalue = nullopt)
1161{
1162 isize n(0);
1163 isize n_eq(0);
1164 isize n_in(0);
1165 if (H != nullopt) {
1166 n = H.value().rows();
1167 }
1168 if (A != nullopt) {
1169 n_eq = A.value().rows();
1170 }
1171 if (C != nullopt) {
1172 n_in = C.value().rows();
1173 }
1174
1175 QP<T> Qp(n, n_eq, n_in, true, DenseBackend::PrimalDualLDLT);
1176 Qp.settings.initial_guess = initial_guess;
1177 Qp.settings.check_duality_gap = check_duality_gap;
1178
1179 if (eps_abs != nullopt) {
1180 Qp.settings.eps_abs = eps_abs.value();
1181 }
1182 if (eps_rel != nullopt) {
1183 Qp.settings.eps_rel = eps_rel.value();
1184 }
1185 if (verbose != nullopt) {
1186 Qp.settings.verbose = verbose.value();
1187 }
1188 if (max_iter != nullopt) {
1189 Qp.settings.max_iter = max_iter.value();
1190 }
1191 if (eps_duality_gap_abs != nullopt) {
1192 Qp.settings.eps_duality_gap_abs = eps_duality_gap_abs.value();
1193 }
1194 if (eps_duality_gap_rel != nullopt) {
1195 Qp.settings.eps_duality_gap_rel = eps_duality_gap_rel.value();
1196 }
1197 Qp.settings.compute_timings = compute_timings;
1198 Qp.settings.primal_infeasibility_solving = primal_infeasibility_solving;
1199 if (manual_minimal_H_eigenvalue != nullopt) {
1200 Qp.init(H,
1201 g,
1202 A,
1203 b,
1204 C,
1205 l,
1206 u,
1207 l_box,
1208 u_box,
1209 compute_preconditioner,
1210 rho,
1211 mu_eq,
1212 mu_in,
1213 manual_minimal_H_eigenvalue.value());
1214 } else {
1215 Qp.init(H,
1216 g,
1217 A,
1218 b,
1219 C,
1220 l,
1221 u,
1222 l_box,
1223 u_box,
1224 compute_preconditioner,
1225 rho,
1226 mu_eq,
1227 mu_in,
1228 nullopt);
1229 }
1230 Qp.solve(x, y, z);
1231
1232 return Qp.results;
1233}
1234
1235template<typename T>
1236bool
1237operator==(const QP<T>& qp1, const QP<T>& qp2)
1238{
1239 bool value = qp1.model == qp2.model && qp1.settings == qp2.settings &&
1240 qp1.results == qp2.results &&
1242 return value;
1243}
1244
1245template<typename T>
1246bool
1247operator!=(const QP<T>& qp1, const QP<T>& qp2)
1248{
1249 return !(qp1 == qp2);
1250}
1251
1253template<typename T>
1255{
1260 std::vector<QP<T>> qp_vector;
1261 dense::isize m_size;
1262
1263 explicit BatchQP(size_t batch_size)
1264 {
1265 if (qp_vector.max_size() != batch_size) {
1266 qp_vector.clear();
1267 qp_vector.reserve(batch_size);
1268 }
1269 m_size = 0;
1270 }
1271
1275 QP<T>& init_qp_in_place(dense::isize dim,
1276 dense::isize n_eq,
1277 dense::isize n_in)
1278 {
1279 qp_vector.emplace_back(dim, n_eq, n_in);
1280 auto& qp = qp_vector.back();
1281 m_size++;
1282 return qp;
1283 };
1284
1288 void insert(const QP<T>& qp) { qp_vector.emplace_back(qp); };
1289
1293 QP<T>& get(isize i) { return qp_vector.at(size_t(i)); };
1294
1298 const QP<T>& get(isize i) const { return qp_vector.at(size_t(i)); };
1299
1303 QP<T>& operator[](isize i) { return get(i); };
1304
1308 const QP<T>& operator[](isize i) const { return get(i); };
1309
1310 dense::isize size() { return m_size; };
1311};
1312
1313} // namespace dense
1314} // namespace proxqp
1315} // namespace proxsuite
1316
1317#endif /* end of include guard PROXSUITE_PROXQP_DENSE_WRAPPER_HPP */
#define PROXSUITE_CHECK_ARGUMENT_SIZE(size, expected_size, message)
Definition macros.hpp:28
#define PROXSUITE_THROW_PRETTY(condition, exception, message)
Definition macros.hpp:18
void warm_start(optional< VecRef< T > > x_wm, optional< VecRef< T > > y_wm, optional< VecRef< T > > z_wm, Results< T > &results, Settings< T > &settings, Model< T > &model)
Definition helpers.hpp:717
void initial_guess(Workspace< T > &qpwork, Settings< T > &qpsettings, Model< T > &qpmodel, Results< T > &qpresults)
Definition helpers.hpp:342
Eigen::Ref< Mat< T, l > const > MatRef
Definition fwd.hpp:33
DenseBackend dense_backend_choice(DenseBackend _dense_backend, isize dim, isize n_eq, isize n_in, bool box_constraints)
This class defines the API of PROXQP solver with dense backend.
Definition wrapper.hpp:83
bool operator==(const Model< T > &model1, const Model< T > &model2)
Definition model.hpp:153
proxqp::Results< T > solve(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, optional< VecRef< T > > x=nullopt, optional< VecRef< T > > y=nullopt, optional< VecRef< T > > z=nullopt, optional< T > eps_abs=nullopt, optional< T > eps_rel=nullopt, optional< T > rho=nullopt, optional< T > mu_eq=nullopt, optional< T > mu_in=nullopt, optional< bool > verbose=nullopt, bool compute_preconditioner=true, bool compute_timings=false, optional< isize > max_iter=nullopt, proxsuite::proxqp::InitialGuessStatus initial_guess=proxsuite::proxqp::InitialGuessStatus::EQUALITY_CONSTRAINED_INITIAL_GUESS, bool check_duality_gap=false, optional< T > eps_duality_gap_abs=nullopt, optional< T > eps_duality_gap_rel=nullopt, bool primal_infeasibility_solving=false, optional< T > manual_minimal_H_eigenvalue=nullopt)
Definition wrapper.hpp:1002
void qp_solve(const Settings< T > &qpsettings, const Model< T > &qpmodel, Results< T > &qpresults, Workspace< T > &qpwork, const bool box_constraints, const DenseBackend &dense_backend, const HessianType &hessian_type, preconditioner::RuizEquilibration< T > &ruiz)
Definition solver.hpp:1090
void update_proximal_parameters(Settings< T > &settings, Results< T > &results, Workspace< T > &work, optional< T > rho_new, optional< T > mu_eq_new, optional< T > mu_in_new)
Definition helpers.hpp:680
void update_default_rho_with_minimal_Hessian_eigen_value(optional< T > manual_minimal_H_eigenvalue, Results< T > &results, Settings< T > &settings)
Definition helpers.hpp:176
void setup(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, optional< VecRef< T > > l_box, optional< VecRef< T > > u_box, Settings< T > &qpsettings, Model< T > &qpmodel, Workspace< T > &qpwork, Results< T > &qpresults, const bool box_constraints, preconditioner::RuizEquilibration< T > &ruiz, PreconditionerStatus preconditioner_status, const HessianType hessian_type)
Definition helpers.hpp:502
void update(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, optional< VecRef< T > > l_box, optional< VecRef< T > > u_box, Model< T > &model, Workspace< T > &work, const bool box_constraints)
Definition helpers.hpp:374
Eigen::Ref< Vec< T > const > VecRef
Definition fwd.hpp:26
bool operator!=(const Model< T > &model1, const Model< T > &model2)
Definition model.hpp:167
constexpr nullopt_t nullopt
Definition optional.hpp:42
This class stores all the results of PROXQP solvers with sparse and dense backends.
Definition results.hpp:68
This class defines the settings of PROXQP solvers with sparse and dense backends.
Definition settings.hpp:89
QP< T > & init_qp_in_place(dense::isize dim, dense::isize n_eq, dense::isize n_in)
Definition wrapper.hpp:1275
std::vector< QP< T > > qp_vector
Definition wrapper.hpp:1260
const QP< T > & operator[](isize i) const
Definition wrapper.hpp:1308
void insert(const QP< T > &qp)
Definition wrapper.hpp:1288
const QP< T > & get(isize i) const
Definition wrapper.hpp:1298
This class stores the model of the QP problem.
Definition model.hpp:24
HessianType which_hessian_type() const
Definition wrapper.hpp:336
DenseBackend which_dense_backend() const
Definition wrapper.hpp:335
QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints)
Definition wrapper.hpp:267
void update(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, optional< VecRef< T > > l_box, optional< VecRef< T > > u_box, bool update_preconditioner=false, optional< T > rho=nullopt, optional< T > mu_eq=nullopt, optional< T > mu_in=nullopt, optional< T > manual_minimal_H_eigenvalue=nullopt)
Definition wrapper.hpp:831
QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints, proxsuite::proxqp::HessianType _hessian_type, DenseBackend _dense_backend)
Definition wrapper.hpp:140
void solve(optional< VecRef< T > > x, optional< VecRef< T > > y, optional< VecRef< T > > z)
Definition wrapper.hpp:940
QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints, DenseBackend _dense_backend, proxsuite::proxqp::HessianType _hessian_type)
Definition wrapper.hpp:173
QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints, proxsuite::proxqp::HessianType _hessian_type)
Definition wrapper.hpp:205
void update(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, bool update_preconditioner=false, optional< T > rho=nullopt, optional< T > mu_eq=nullopt, optional< T > mu_in=nullopt, optional< T > manual_minimal_H_eigenvalue=nullopt)
Definition wrapper.hpp:723
QP(isize _dim, isize _n_eq, isize _n_in)
Definition wrapper.hpp:318
QP(isize _dim, isize _n_eq, isize _n_in, bool _box_constraints, DenseBackend _dense_backend)
Definition wrapper.hpp:237
bool is_box_constrained() const
Definition wrapper.hpp:334
void init(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, optional< VecRef< T > > l_box, optional< VecRef< T > > u_box, bool compute_preconditioner=true, optional< T > rho=nullopt, optional< T > mu_eq=nullopt, optional< T > mu_in=nullopt, optional< T > manual_minimal_H_eigenvalue=nullopt)
Definition wrapper.hpp:520
QP(isize _dim, isize _n_eq, isize _n_in, proxsuite::proxqp::HessianType _hessian_type)
Definition wrapper.hpp:293
preconditioner::RuizEquilibration< T > ruiz
Definition wrapper.hpp:129
void init(optional< MatRef< T > > H, optional< VecRef< T > > g, optional< MatRef< T > > A, optional< VecRef< T > > b, optional< MatRef< T > > C, optional< VecRef< T > > l, optional< VecRef< T > > u, bool compute_preconditioner=true, optional< T > rho=nullopt, optional< T > mu_eq=nullopt, optional< T > mu_in=nullopt, optional< T > manual_minimal_H_eigenvalue=nullopt)
Definition wrapper.hpp:354
This class defines the workspace of the dense solver.
Definition workspace.hpp:26