Merge remote-tracking branch 'upstream/master' into cech_check_filt

3 files changed, 20 insertions, 22 deletions

diff --git a/src/Cech_complex/include/gudhi/Cech_complex_blocker.h b/src/Cech_complex/include/gudhi/Cech_complex_blocker.h
index 9b5c5add..3ea82826 100644
--- a/src/Cech_complex/include/gudhi/Cech_complex_blocker.h
+++ b/src/Cech_complex/include/gudhi/Cech_complex_blocker.h
@@ -54,6 +54,7 @@ class Cech_blocker {
 
   using Simplex_handle = typename SimplicialComplexForCech::Simplex_handle;
   using Filtration_value = typename SimplicialComplexForCech::Filtration_value;
+  using Simplex_key = typename SimplicialComplexForCech::Simplex_key;
 
   template<class PointIterator>
   Sphere get_sphere(PointIterator begin, PointIterator end) const {
@@ -70,7 +71,6 @@ class Cech_blocker {
    *  \return true if the simplex radius is greater than the Cech_complex max_radius*/
   bool operator()(Simplex_handle sh) {
     using Point_cloud =  std::vector<Point_d>;
-    CGAL::NT_converter<FT, Filtration_value> cast_to_fv;
     Filtration_value radius = 0;
     bool is_min_enclos_ball = false;
     Point_cloud points;
@@ -78,10 +78,10 @@ class Cech_blocker {
 
     // for each face of simplex sh, test outsider point is indeed inside enclosing ball, if yes, take it and exit loop, otherwise, new sphere is circumsphere of all vertices
     for (auto face_opposite_vertex : sc_ptr_->boundary_opposite_vertex_simplex_range(sh)) {
-        Sphere sph;
         auto k = sc_ptr_->key(face_opposite_vertex.first);
+        Simplex_key sph_key;
         if(k != sc_ptr_->null_key()) {
-            sph = cc_ptr_->get_cache().at(k);
+            sph_key = k;
         }
         else {
             for (auto vertex : sc_ptr_->simplex_vertex_range(face_opposite_vertex.first)) {
@@ -90,25 +90,22 @@ class Cech_blocker {
                 std::clog << "#(" << vertex << ")#";
 #endif  // DEBUG_TRACES
             }
-            sph = get_sphere(points.cbegin(), points.cend());
             // Put edge sphere in cache
-            sc_ptr_->assign_key(face_opposite_vertex.first, cc_ptr_->get_cache().size());
-            cc_ptr_->get_cache().push_back(sph);
+            sph_key = cc_ptr_->get_cache().size();
+            sc_ptr_->assign_key(face_opposite_vertex.first, sph_key);
+            cc_ptr_->get_cache().push_back(get_sphere(points.cbegin(), points.cend()));
             // Clear face points
             points.clear();
         }
         // Check if the minimal enclosing ball of current face contains the extra point/opposite vertex
+        Sphere const& sph = cc_ptr_->get_cache()[sph_key];
         if (kernel_.squared_distance_d_object()(sph.first, cc_ptr_->get_point(face_opposite_vertex.second)) <= sph.second) {
+            is_min_enclos_ball = true;
+            sc_ptr_->assign_key(sh, sph_key);
+            radius = sc_ptr_->filtration(face_opposite_vertex.first);
 #ifdef DEBUG_TRACES
             std::clog << "center: " << sph.first << ", radius: " <<  radius << std::endl;
 #endif  // DEBUG_TRACES
-            is_min_enclos_ball = true;
-#if CGAL_VERSION_NR >= 1050000000
-            if(exact_) CGAL::exact(sph.second);
-#endif
-            radius = std::sqrt(cast_to_fv(sph.second));
-            sc_ptr_->assign_key(sh, cc_ptr_->get_cache().size());
-            cc_ptr_->get_cache().push_back(sph);
             break;
         }
     }
@@ -121,6 +118,7 @@ class Cech_blocker {
 #if CGAL_VERSION_NR >= 1050000000
         if(exact_) CGAL::exact(sph.second);
 #endif
+        CGAL::NT_converter<FT, Filtration_value> cast_to_fv;
         radius = std::sqrt(cast_to_fv(sph.second));
 
         sc_ptr_->assign_key(sh, cc_ptr_->get_cache().size());
diff --git a/src/python/doc/installation.rst b/src/python/doc/installation.rst
index dd476054..4eefd415 100644
--- a/src/python/doc/installation.rst
+++ b/src/python/doc/installation.rst
@@ -175,7 +175,7 @@ A complete configuration would be :
     Scikit-learn version 1.0.1 
     POT version 0.8.0 
     HNSWlib found
-    PyKeOps version [pyKeOps]: 1.5 
+    PyKeOps version [pyKeOps]: 2.1
     EagerPy version 0.30.0 
     TensorFlow version 2.7.0 
     Sphinx version 4.3.0 
diff --git a/src/python/test/test_dtm.py b/src/python/test/test_dtm.py
index e46d616c..b276f041 100755
--- a/src/python/test/test_dtm.py
+++ b/src/python/test/test_dtm.py
@@ -91,11 +91,11 @@ def test_density():
 
 def test_dtm_overflow_warnings():
     pts = numpy.array([[10., 100000000000000000000000000000.], [1000., 100000000000000000000000000.]])
-
-    with warnings.catch_warnings(record=True) as w:
-        # TODO Test "keops" implementation as well when next version of pykeops (current is 1.5) is released (should fix the problem (cf. issue #543))
-        dtm = DistanceToMeasure(2, implementation="hnsw")
-        r = dtm.fit_transform(pts)
-        assert len(w) == 1
-        assert issubclass(w[0].category, RuntimeWarning)
-        assert "Overflow" in str(w[0].message)
+    impl_warn = ["keops", "hnsw"]
+    for impl in impl_warn:
+        with warnings.catch_warnings(record=True) as w:
+            dtm = DistanceToMeasure(2, implementation=impl)
+            r = dtm.fit_transform(pts)
+            assert len(w) == 1
+            assert issubclass(w[0].category, RuntimeWarning)
+            assert "Overflow" in str(w[0].message)