Merge pull request #509 from Hind-M/homology_coeff_field

Fix issue #502: check homology_coeff_field primality before computing persistence

7 files changed, 58 insertions, 17 deletions

diff --git a/src/Persistent_cohomology/example/CMakeLists.txt b/src/Persistent_cohomology/example/CMakeLists.txt
index c68c6524..3e7e9369 100644
--- a/src/Persistent_cohomology/example/CMakeLists.txt
+++ b/src/Persistent_cohomology/example/CMakeLists.txt
@@ -11,7 +11,7 @@ if (TBB_FOUND)
   target_link_libraries(persistence_from_simple_simplex_tree ${TBB_LIBRARIES})
 endif()
 add_test(NAME Persistent_cohomology_example_from_simple_simplex_tree COMMAND $<TARGET_FILE:persistence_from_simple_simplex_tree>
-    "1" "0")
+    "2" "0")
 
 if(TARGET Boost::program_options)
   add_executable(rips_persistence_step_by_step rips_persistence_step_by_step.cpp)
diff --git a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology/Field_Zp.h b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology/Field_Zp.h
index 0673625c..8ec89e41 100644
--- a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology/Field_Zp.h
+++ b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology/Field_Zp.h
@@ -13,6 +13,7 @@
 
 #include <utility>
 #include <vector>
+#include <stdexcept>
 
 namespace Gudhi {
 
@@ -34,14 +35,29 @@ class Field_Zp {
 
   void init(int charac) {
     assert(charac > 0);  // division by zero + non negative values
+
     Prime = charac;
+
+    // Check that the provided prime is less than the maximum allowed as int, calculation below, and 'plus_times_equal' function : 46337 ; i.e (max_prime-1)*max_prime <= INT_MAX
+    if(Prime > 46337)
+        throw std::invalid_argument("Maximum homology_coeff_field allowed value is 46337");
+
+    // Check for primality
+    if (Prime <= 1)
+        throw std::invalid_argument("homology_coeff_field must be a prime number");
+
     inverse_.clear();
     inverse_.reserve(charac);
     inverse_.push_back(0);
     for (int i = 1; i < Prime; ++i) {
       int inv = 1;
-      while (((inv * i) % Prime) != 1)
+      int mult = inv * i;
+      while ( (mult % Prime) != 1) {
         ++inv;
+        if(mult == Prime)
+            throw std::invalid_argument("homology_coeff_field must be a prime number");
+        mult = inv * i;
+      }
       inverse_.push_back(inv);
     }
   }
diff --git a/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp b/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp
index fe3f8517..041cb0fd 100644
--- a/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp
+++ b/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp
@@ -146,9 +146,14 @@ void test_rips_persistence_in_dimension(int dimension) {
   std::clog << "str_rips_persistence=" << str_rips_persistence << std::endl;
 }
 
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_0 )
+{
+  BOOST_CHECK_THROW(test_rips_persistence_in_dimension(0), std::invalid_argument);
+}
+
 BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_1 )
 {
-  test_rips_persistence_in_dimension(1);
+  BOOST_CHECK_THROW(test_rips_persistence_in_dimension(1), std::invalid_argument);
 }
 
 BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_2 )
@@ -161,15 +166,35 @@ BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_3 )
   test_rips_persistence_in_dimension(3);
 }
 
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_4 )
+{
+  BOOST_CHECK_THROW(test_rips_persistence_in_dimension(4), std::invalid_argument);
+}
+
 BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_5 )
 {
   test_rips_persistence_in_dimension(5);
 }
 
-// TODO(VR): not working from 6
-// std::string str_rips_persistence = test_rips_persistence(6, 0);
-// TODO(VR): division by zero
-// std::string str_rips_persistence = test_rips_persistence(0, 0);
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_6 )
+{
+  BOOST_CHECK_THROW(test_rips_persistence_in_dimension(6), std::invalid_argument);
+}
+
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_11 )
+{
+  test_rips_persistence_in_dimension(11);
+}
+
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_13 )
+{
+  test_rips_persistence_in_dimension(13);
+}
+
+BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_46349 )
+{
+  BOOST_CHECK_THROW(test_rips_persistence_in_dimension(46349), std::invalid_argument);
+}
 
 /** SimplexTree minimal options to test the limits.
  * 
diff --git a/src/python/gudhi/cubical_complex.pyx b/src/python/gudhi/cubical_complex.pyx
index 28fbe3af..97c69a2d 100644
--- a/src/python/gudhi/cubical_complex.pyx
+++ b/src/python/gudhi/cubical_complex.pyx
@@ -35,7 +35,7 @@ cdef extern from "Cubical_complex_interface.h" namespace "Gudhi":
 cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi":
     cdef cppclass Cubical_complex_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Cubical_complex::Cubical_complex_interface<>>":
         Cubical_complex_persistence_interface(Bitmap_cubical_complex_base_interface * st, bool persistence_dim_max) nogil
-        void compute_persistence(int homology_coeff_field, double min_persistence) nogil
+        void compute_persistence(int homology_coeff_field, double min_persistence) nogil except+
         vector[pair[int, pair[double, double]]] get_persistence() nogil
         vector[vector[int]] cofaces_of_cubical_persistence_pairs() nogil
         vector[int] betti_numbers() nogil
@@ -147,7 +147,7 @@ cdef class CubicalComplex:
         :func:`persistence` returns.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int.
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is
@@ -169,7 +169,7 @@ cdef class CubicalComplex:
         """This function computes and returns the persistence of the complex.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int.
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is
diff --git a/src/python/gudhi/periodic_cubical_complex.pyx b/src/python/gudhi/periodic_cubical_complex.pyx
index d353d2af..ef1d3080 100644
--- a/src/python/gudhi/periodic_cubical_complex.pyx
+++ b/src/python/gudhi/periodic_cubical_complex.pyx
@@ -32,7 +32,7 @@ cdef extern from "Cubical_complex_interface.h" namespace "Gudhi":
 cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi":
     cdef cppclass Periodic_cubical_complex_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Cubical_complex::Cubical_complex_interface<Gudhi::cubical_complex::Bitmap_cubical_complex_periodic_boundary_conditions_base<double>>>":
         Periodic_cubical_complex_persistence_interface(Periodic_cubical_complex_base_interface * st, bool persistence_dim_max) nogil
-        void compute_persistence(int homology_coeff_field, double min_persistence) nogil
+        void compute_persistence(int homology_coeff_field, double min_persistence) nogil except +
         vector[pair[int, pair[double, double]]] get_persistence() nogil
         vector[vector[int]] cofaces_of_cubical_persistence_pairs() nogil
         vector[int] betti_numbers() nogil
@@ -148,7 +148,7 @@ cdef class PeriodicCubicalComplex:
         :func:`persistence` returns.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int.
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is
@@ -170,7 +170,7 @@ cdef class PeriodicCubicalComplex:
         """This function computes and returns the persistence of the complex.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int.
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is
diff --git a/src/python/gudhi/simplex_tree.pxd b/src/python/gudhi/simplex_tree.pxd
index 3b8ea4f9..006a24ed 100644
--- a/src/python/gudhi/simplex_tree.pxd
+++ b/src/python/gudhi/simplex_tree.pxd
@@ -78,7 +78,7 @@ cdef extern from "Simplex_tree_interface.h" namespace "Gudhi":
 cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi":
     cdef cppclass Simplex_tree_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_full_featured>>":
         Simplex_tree_persistence_interface(Simplex_tree_interface_full_featured * st, bool persistence_dim_max) nogil
-        void compute_persistence(int homology_coeff_field, double min_persistence) nogil
+        void compute_persistence(int homology_coeff_field, double min_persistence) nogil except +
         vector[pair[int, pair[double, double]]] get_persistence() nogil
         vector[int] betti_numbers() nogil
         vector[int] persistent_betti_numbers(double from_value, double to_value) nogil
diff --git a/src/python/gudhi/simplex_tree.pyx b/src/python/gudhi/simplex_tree.pyx
index be08a3a1..9c51cb46 100644
--- a/src/python/gudhi/simplex_tree.pyx
+++ b/src/python/gudhi/simplex_tree.pyx
@@ -412,7 +412,7 @@ cdef class SimplexTree:
         """This function retrieves good values for extended persistence, and separate the diagrams into the Ordinary,
         Relative, Extended+ and Extended- subdiagrams.
 
-        :param homology_coeff_field: The homology coefficient field. Must be a prime number. Default value is 11.
+        :param homology_coeff_field: The homology coefficient field. Must be a prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int
         :param min_persistence: The minimum persistence value (i.e., the absolute value of the difference between the
             persistence diagram point coordinates) to take into account (strictly greater than min_persistence).
@@ -449,7 +449,7 @@ cdef class SimplexTree:
         """This function computes and returns the persistence of the simplicial complex.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number. Default value is 11.
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is
@@ -472,7 +472,7 @@ cdef class SimplexTree:
         when you do not want the list :func:`persistence` returns.
 
         :param homology_coeff_field: The homology coefficient field. Must be a
-            prime number. Default value is 11.
+            prime number. Default value is 11. Max is 46337.
         :type homology_coeff_field: int
         :param min_persistence: The minimum persistence value to take into
             account (strictly greater than min_persistence). Default value is