Merge pull request #139 from fabinsch/topic/fix-cpp-example-update

add is_valid function for dense model and fix example + unittest

Author: jcarpent

bindings/python/src/expose-model.hpp

@@ -31,7 +31,10 @@ exposeDenseModel(pybind11::module_ m)
     .def_readonly("dim", &Model<T>::dim)
     .def_readonly("n_eq", &Model<T>::n_eq)
     .def_readonly("n_in", &Model<T>::n_in)
-    .def_readonly("n_total", &Model<T>::n_total);
+    .def_readonly("n_total", &Model<T>::n_total)
+    .def("is_valid",
+         &Model<T>::is_valid,
+         "Check if model is containing valid data.");
 }
 } // namespace python
 } // namespace dense

doc/2-PROXQP_API/2-ProxQP_api.md

@@ -305,7 +305,7 @@ In this table you have on the three columns from left to right: the name of the
 | default_rho                         | 1.E-6                              | Default rho parameter of result class (i.e., for each initial guess, except WARM_START_WITH_PREVIOUS_RESULT, after a new solve or update, the solver initializes rho to this value).
 | default_mu_eq                       | 1.E-3                              | Default mu_eq parameter of result class (i.e., for each initial guess, except WARM_START_WITH_PREVIOUS_RESULT, after a new solve or update, the solver initializes mu_eq to this value).
 | default_mu_in                       | 1.E-1                              | Default mu_in parameter of result class (i.e., for each initial guess, except WARM_START_WITH_PREVIOUS_RESULT, after a new solve or update, the solver initializes mu_in to this value).
-| compute_timings                     | True                               | If set to true, timings will be computed by the solver (setup time, solving time, and run time = setup time + solving time).
+| compute_timings                     | False                              | If set to true, timings will be computed by the solver (setup time, solving time, and run time = setup time + solving time).
 | max_iter                            | 10.000                             | Maximal number of authorized outer iterations.
 | max_iter_in                         | 1500                               | Maximal number of authorized inner iterations.
 | initial_guess                       | EQUALITY_CONSTRAINED_INITIAL_GUESS | Sets the initial guess option for initilizing x, y and z.

doc/3-ProxQP_solve.md

@@ -100,7 +100,7 @@ Different options are available for the solve function. In the table below you h
 | rho                                 | 1.E-6                                           | Proximal step size wrt primal variable.
 | VERBOSE                             | False                                           | If set to true, the solver prints information at each loop.
 | compute_preconditioner              | True                                            | If set to true, the preconditioner will be derived.
-| compute_timings                     | True                                            | If set to true, timings will be computed by the solver (setup time, solving time, and run time = setup time + solving time).
+| compute_timings                     | False                                           | If set to true, timings will be computed by the solver (setup time, solving time, and run time = setup time + solving time).
 | max_iter                            | 10.000                                          | Maximal number of authorized outer iterations.
 | initial_guess                       | EQUALITY_CONSTRAINED_INITIAL_GUESS              | Sets the initial guess option for initilizing x, y and z.
 

examples/cpp/solve_without_api.cpp

@@ -14,10 +14,12 @@ main()
   isize n = 10;
   isize n_eq(n / 4);
   isize n_in(n / 4);
-  T p = 0.15;            // level of sparsity
+  T p = 0.35;            // level of sparsity
   T conditioning = 10.0; // conditioning level for H
-  auto H = utils::rand::sparse_positive_definite_rand(n, conditioning, p);
+  auto H = utils::rand::sparse_positive_definite_rand(
+    n, conditioning, p); // upper triangular matrix
   Mat H_dense = Mat(H);
+  H_dense.template triangularView<Eigen::Lower>() = H_dense.transpose();
   Vec g = utils::rand::vector_rand<T>(n);
   auto A = utils::rand::sparse_matrix_rand<T>(n_eq, n, p);
   Mat A_dense = Mat(A);

examples/cpp/update_dense_qp.cpp

@@ -29,7 +29,7 @@ main()
   dense::Model<T> qp2 = utils::dense_strongly_convex_qp(
     dim, n_eq, n_in, sparsity_factor, strong_convexity_factor);
   // re update the model
-  qp.update(qp2.H, qp2.g, qp2.A, qp2.b, qp2.C, qp2.u, qp2.l);
+  qp.update(qp2.H, qp2.g, qp2.A, qp2.b, qp2.C, qp2.l, qp2.u);
   // solve it
   qp.solve();
   // print an optimal solution x,y and z

include/proxsuite/proxqp/dense/helpers.hpp

@@ -264,6 +264,7 @@ update(optional<MatRef<T>> H,
   if (C != nullopt) {
     model.C = C.value();
   }
+  assert(model.is_valid());
 }
 /*!
  * Setups the QP solver model.
@@ -383,6 +384,7 @@ setup( //
     qpmodel.l = l.value();
   } // else qpmodel.l remains initialized to a matrix with zero elements or zero
     // shape
+  assert(qpmodel.is_valid());
 
   qpwork.H_scaled = qpmodel.H;
   qpwork.g_scaled = qpmodel.g;

include/proxsuite/proxqp/dense/model.hpp

@@ -84,6 +84,40 @@ struct Model
                                                    C_sparse, u, l };
     return res;
   }
+
+  bool is_valid()
+  {
+#define PROXSUITE_CHECK_SIZE(size, expected_size)                              \
+  if (size != 0) {                                                             \
+    if (!(size == expected_size))                                              \
+      return false;                                                            \
+  }
+
+    // check that all matrices and vectors of qpmodel have the correct size
+    // and that H and C have expected properties
+    PROXSUITE_CHECK_SIZE(g.size(), dim);
+    PROXSUITE_CHECK_SIZE(b.size(), n_eq);
+    PROXSUITE_CHECK_SIZE(u.size(), n_in);
+    PROXSUITE_CHECK_SIZE(l.size(), n_in);
+    if (H.size()) {
+      PROXSUITE_CHECK_SIZE(H.rows(), dim);
+      PROXSUITE_CHECK_SIZE(H.cols(), dim);
+      if (!H.isApprox(H.transpose(), 0.0))
+        return false;
+    }
+    if (A.size()) {
+      PROXSUITE_CHECK_SIZE(A.rows(), n_eq);
+      PROXSUITE_CHECK_SIZE(A.cols(), dim);
+    }
+    if (C.size()) {
+      PROXSUITE_CHECK_SIZE(C.rows(), n_in);
+      PROXSUITE_CHECK_SIZE(C.cols(), dim);
+      if (C.isZero())
+        return false;
+    }
+    return true;
+#undef PROXSUITE_CHECK_SIZE
+  }
 };
 } // namespace dense
 } // namespace proxqp

include/proxsuite/proxqp/dense/wrapper.hpp

@@ -412,7 +412,7 @@ solve(
   optional<T> mu_in = nullopt,
   optional<bool> verbose = nullopt,
   bool compute_preconditioner = true,
-  bool compute_timings = true,
+  bool compute_timings = false,
   optional<isize> max_iter = nullopt,
   proxsuite::proxqp::InitialGuessStatus initial_guess =
     proxsuite::proxqp::InitialGuessStatus::EQUALITY_CONSTRAINED_INITIAL_GUESS)
@@ -443,7 +443,7 @@ solve(
     Qp.settings.verbose = verbose.value();
   }
   if (max_iter != nullopt) {
-    Qp.settings.max_iter = verbose.value();
+    Qp.settings.max_iter = max_iter.value();
   }
   Qp.settings.compute_timings = compute_timings;
   Qp.init(H, g, A, b, C, l, u, compute_preconditioner, rho, mu_eq, mu_in);

include/proxsuite/proxqp/settings.hpp

@@ -182,7 +182,7 @@ struct Settings
                                           // once a problem
     bool update_preconditioner = true,
     bool compute_preconditioner = true,
-    bool compute_timings = true,
+    bool compute_timings = false,
     isize preconditioner_max_iter = 10,
     T preconditioner_accuracy = 1.e-3,
     T eps_primal_inf = 1.E-4,
@@ -226,78 +226,6 @@ struct Settings
     , sparse_backend(sparse_backend)
   {
   }
-  /*
- void set(
-           T alpha_bcl_ = 0.1,
-           T beta_bcl_ = 0.9,
-           T refactor_dual_feasibility_threshold_ = 1e-2,
-           T refactor_rho_threshold_ = 1e-7,
-           T mu_min_eq_ = 1e-9,
-           T mu_min_in_ = 1e-8,
-           T mu_max_eq_inv_ = 1e9,
-           T mu_max_in_inv_ = 1e8,
-           T mu_update_factor_ = 0.1,
-           T mu_update_inv_factor_ = 10,
-           T cold_reset_mu_eq_ = 1. / 1.1,
-           T cold_reset_mu_in_ = 1. / 1.1,
-           T cold_reset_mu_eq_inv_ = 1.1,
-           T cold_reset_mu_in_inv_ = 1.1,
-           T eps_abs_ = 1.e-8,
-           T eps_rel_ = 0,
-           isize max_iter_ = 10000,
-           isize max_iter_in_ = 1500,
-           isize safe_guard_ = 1.E4,
-           isize nb_iterative_refinement_ = 10,
-           T eps_refact_ = 1.e-6, // before eps_refact_=1.e-6
-           bool VERBOSE = false,
-           InitialGuessStatus initial_guess_ =
-             InitialGuessStatus::WARM_START_WITH_PREVIOUS_RESULT,
-           bool update_preconditioner_ = true,
-           bool compute_preconditioner_ = true,
-           bool compute_timings_ = true,
-           isize preconditioner_max_iter_ = 10,
-           T preconditioner_accuracy_ = 1.e-3,
-           T eps_primal_inf_ = 1.E-4,
-           T eps_dual_inf_ = 1.E-4,
-           bool bcl_update_ = true,
-           SparseBackend sparse_backend_ = SparseBackend::Automatic
- ){
-    alpha_bcl = alpha_bcl_;
-    beta_bcl = beta_bcl_ ;
-    refactor_dual_feasibility_threshold = refactor_dual_feasibility_threshold_;
-    refactor_rho_threshold = refactor_rho_threshold_;
-    mu_min_eq = mu_min_eq_;
-    mu_min_in = mu_min_in_;
-    mu_max_eq_inv = mu_max_eq_inv_;
-    mu_max_in_inv = mu_max_in_inv_;
-    mu_update_factor = mu_update_factor_;
-    mu_update_inv_factor = mu_update_inv_factor_;
-    cold_reset_mu_eq = cold_reset_mu_eq_;
-    cold_reset_mu_in = cold_reset_mu_in_;
-    cold_reset_mu_eq_inv = cold_reset_mu_eq_inv_;
-    cold_reset_mu_in_inv = cold_reset_mu_in_inv_;
-    eps_abs = eps_abs_;
-    eps_rel = eps_rel_;
-    max_iter = max_iter_;
-    max_iter_in = max_iter_in_;
-    safe_guard = safe_guard_;
-    nb_iterative_refinement = nb_iterative_refinement_;
-    eps_refact = eps_refact_;
-    verbose = VERBOSE;
-    initial_guess = initial_guess_;
-    update_preconditioner = update_preconditioner_;
-    compute_preconditioner = compute_preconditioner_;
-    compute_timings = compute_timings_;
-    preconditioner_max_iter = preconditioner_max_iter_;
-    preconditioner_accuracy = preconditioner_accuracy_;
-    eps_primal_inf = eps_primal_inf_;
-    eps_dual_inf = eps_dual_inf_;
-    bcl_update = bcl_update_;
-    sparse_backend = sparse_backend_;
-
-
- }
- */
 };
 
 } // namespace proxqp

include/proxsuite/proxqp/sparse/wrapper.hpp

@@ -705,7 +705,7 @@ solve(
   optional<T> mu_in = nullopt,
   optional<bool> verbose = nullopt,
   bool compute_preconditioner = true,
-  bool compute_timings = true,
+  bool compute_timings = false,
   optional<isize> max_iter = nullopt,
   proxsuite::proxqp::InitialGuessStatus initial_guess =
     proxsuite::proxqp::InitialGuessStatus::EQUALITY_CONSTRAINED_INITIAL_GUESS,
@@ -739,7 +739,7 @@ solve(
     Qp.settings.verbose = verbose.value();
   }
   if (max_iter != nullopt) {
-    Qp.settings.max_iter = verbose.value();
+    Qp.settings.max_iter = max_iter.value();
   }
   Qp.settings.compute_timings = compute_timings;
   Qp.settings.sparse_backend = sparse_backend;