Merge pull request #272 from fabinsch/fix-is-valid

Check model is_valid up to eps

Author: jcarpent

include/proxsuite/proxqp/dense/helpers.hpp

@@ -112,9 +112,10 @@ estimate_minimal_eigen_value_of_symmetric_matrix(
   T power_iteration_accuracy,
   isize nb_power_iteration)
 {
-  PROXSUITE_THROW_PRETTY((!H.isApprox(H.transpose(), 0.0)),
-                         std::invalid_argument,
-                         "H is not symmetric.");
+  PROXSUITE_THROW_PRETTY(
+    (!H.isApprox(H.transpose(), std::numeric_limits<T>::epsilon())),
+    std::invalid_argument,
+    "H is not symmetric.");
   if (H.size()) {
     PROXSUITE_CHECK_ARGUMENT_SIZE(
       H.rows(),

include/proxsuite/proxqp/dense/model.hpp

@@ -120,9 +120,12 @@ struct Model
         H.rows(), dim, "H has not the expected number of rows.");
       PROXSUITE_CHECK_ARGUMENT_SIZE(
         H.cols(), dim, "H has not the expected number of cols.");
-      PROXSUITE_THROW_PRETTY((!H.isApprox(H.transpose(), 0.0)),
-                             std::invalid_argument,
-                             "H is not symmetric.");
+      PROXSUITE_THROW_PRETTY(
+        (!H.isApprox(
+          H.transpose(),
+          std::numeric_limits<typename decltype(H)::Scalar>::epsilon())),
+        std::invalid_argument,
+        "H is not symmetric.");
     }
     if (A.size()) {
       PROXSUITE_CHECK_ARGUMENT_SIZE(

include/proxsuite/proxqp/sparse/helpers.hpp

@@ -119,9 +119,10 @@ estimate_minimal_eigen_value_of_symmetric_matrix(SparseMat<T, I>& H,
                                                  T power_iteration_accuracy,
                                                  isize nb_power_iteration)
 {
-  PROXSUITE_THROW_PRETTY((!H.isApprox(H.transpose(), 0.0)),
-                         std::invalid_argument,
-                         "H is not symmetric.");
+  PROXSUITE_THROW_PRETTY(
+    (!H.isApprox(H.transpose(), std::numeric_limits<T>::epsilon())),
+    std::invalid_argument,
+    "H is not symmetric.");
   PROXSUITE_CHECK_ARGUMENT_SIZE(
     H.rows(),
     H.cols(),

test/src/dense_qp_wrapper.cpp

@@ -7564,3 +7564,27 @@ TEST_CASE(
                            minimal_eigenvalue) <= tol);
   }
 }
+
+DOCTEST_TEST_CASE("check that model.is_valid function for symmetric matrices "
+                  "works for epsilon precision")
+{
+  Eigen::Matrix<T, 3, 3> matrix = Eigen::Matrix<T, 3, 3>::Random();
+  Eigen::Matrix<T, 3, 3> symmetric_mat = matrix + matrix.transpose();
+
+  symmetric_mat(0, 1) =
+    symmetric_mat(1, 0) + std::numeric_limits<double>::epsilon();
+
+  // compare the two checks for symmetry with and without tolerance
+  bool is_symmetric_without_tolerance =
+    symmetric_mat.isApprox(symmetric_mat.transpose(), 0.0);
+  bool is_symmetric_with_tolerance = symmetric_mat.isApprox(
+    symmetric_mat.transpose(),
+    std::numeric_limits<typename decltype(symmetric_mat)::Scalar>::epsilon());
+  DOCTEST_CHECK(is_symmetric_without_tolerance == false);
+  DOCTEST_CHECK(is_symmetric_with_tolerance == true);
+
+  // initialize a model with a symmetric matrix as Hessian, this runs
+  // model.is_valid() that performs the check above
+  proxqp::dense::QP<T> qp(3, 0, 0);
+  qp.init(symmetric_mat, nullopt, nullopt, nullopt, nullopt, nullopt, nullopt);
+}