contiguous_vertices branch development merge.

git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/trunk@950 636b058d-ea47-450e-bf9e-a15bfbe3eedb

Former-commit-id: 1ec25e73a5d40695aa60175b378ec2ca64da480b

10 files changed, 171 insertions, 71 deletions

diff --git a/src/Persistent_cohomology/example/alpha_shapes_persistence.cpp b/src/Persistent_cohomology/example/alpha_shapes_persistence.cpp
index 6d5eebcf..92c0b065 100644
--- a/src/Persistent_cohomology/example/alpha_shapes_persistence.cpp
+++ b/src/Persistent_cohomology/example/alpha_shapes_persistence.cpp
@@ -66,7 +66,8 @@ typedef Alpha_shape_3::Edge Edge_3;
 typedef std::list<Alpha_shape_3::Vertex_handle> Vertex_list;
 
 // gudhi type definition
-typedef Simplex_tree<>::Vertex_handle Simplex_tree_vertex;
+typedef Simplex_tree<Simplex_tree_options_fast_persistence> ST;
+typedef ST::Vertex_handle Simplex_tree_vertex;
 typedef std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex > Alpha_shape_simplex_tree_map;
 typedef std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex> Alpha_shape_simplex_tree_pair;
 typedef std::vector< Simplex_tree_vertex > Simplex_tree_vector_vertex;
@@ -184,7 +185,7 @@ int main(int argc, char * const argv[]) {
 
   // Loop on objects vector
   Vertex_list vertex_list;
-  Simplex_tree<> simplex_tree;
+  ST simplex_tree;
   Alpha_shape_simplex_tree_map map_cgal_simplex_tree;
   std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin();
   int dim_max = 0;
@@ -281,7 +282,7 @@ int main(int argc, char * const argv[]) {
 
   std::cout << "Simplex_tree dim: " << simplex_tree.dimension() << std::endl;
   // Compute the persistence diagram of the complex
-  Persistent_cohomology< Simplex_tree<>, Field_Zp > pcoh(simplex_tree);
+  Persistent_cohomology< ST, Field_Zp > pcoh(simplex_tree);
   // initializes the coefficient field for homology
   pcoh.init_coefficients(coeff_field_characteristic);
 
diff --git a/src/Persistent_cohomology/example/performance_rips_persistence.cpp b/src/Persistent_cohomology/example/performance_rips_persistence.cpp
index fc48d6b1..b4d282ac 100644
--- a/src/Persistent_cohomology/example/performance_rips_persistence.cpp
+++ b/src/Persistent_cohomology/example/performance_rips_persistence.cpp
@@ -86,7 +86,7 @@ int main(int argc, char * argv[]) {
   std::cout << "Compute Rips graph in " << elapsed_sec << " ms.\n";
 
   // Construct the Rips complex in a Simplex Tree
-  Simplex_tree<> st;
+  Simplex_tree<Simplex_tree_options_fast_persistence> st;
   start = std::chrono::system_clock::now();
 
   // insert the proximity graph in the simplex tree
diff --git a/src/Persistent_cohomology/example/persistence_from_file.cpp b/src/Persistent_cohomology/example/persistence_from_file.cpp
index 8eb8d0f3..67235467 100644
--- a/src/Persistent_cohomology/example/persistence_from_file.cpp
+++ b/src/Persistent_cohomology/example/persistence_from_file.cpp
@@ -54,7 +54,7 @@ int main(int argc, char * argv[]) {
       << std::endl;
   std::cout << "     - p=" << p << " - min_persistence=" << min_persistence << std::endl;
 
-  // Construct the Rips complex in a Simplex Tree
+  // Read the list of simplices from a file.
   Simplex_tree<> simplex_tree;
 
   std::ifstream simplex_tree_stream(simplex_tree_file);
diff --git a/src/Persistent_cohomology/example/rips_multifield_persistence.cpp b/src/Persistent_cohomology/example/rips_multifield_persistence.cpp
index 5277bf7a..c5cd775d 100644
--- a/src/Persistent_cohomology/example/rips_multifield_persistence.cpp
+++ b/src/Persistent_cohomology/example/rips_multifield_persistence.cpp
@@ -68,7 +68,8 @@ int main(int argc, char * argv[]) {
                                                , euclidean_distance<Point_t>);
 
   // Construct the Rips complex in a Simplex Tree
-  Simplex_tree<> st;
+  typedef Simplex_tree<Simplex_tree_options_fast_persistence> ST;
+  ST st;
   // insert the proximity graph in the simplex tree
   st.insert_graph(prox_graph);
   // expand the graph until dimension dim_max
@@ -78,7 +79,7 @@ int main(int argc, char * argv[]) {
   st.initialize_filtration();
 
   // Compute the persistence diagram of the complex
-  Persistent_cohomology< Simplex_tree<>, Multi_field > pcoh(st);
+  Persistent_cohomology<ST, Multi_field > pcoh(st);
   // initializes the coefficient field for homology
   pcoh.init_coefficients(min_p, max_p);
   // compute persistent homology, disgarding persistent features of life shorter than min_persistence
diff --git a/src/Persistent_cohomology/example/rips_persistence.cpp b/src/Persistent_cohomology/example/rips_persistence.cpp
index 9b1ef42f..2d926a0d 100644
--- a/src/Persistent_cohomology/example/rips_persistence.cpp
+++ b/src/Persistent_cohomology/example/rips_persistence.cpp
@@ -65,7 +65,8 @@ int main(int argc, char * argv[]) {
                                                , euclidean_distance<Point_t>);
 
   // Construct the Rips complex in a Simplex Tree
-  Simplex_tree<> st;
+  typedef Simplex_tree<Simplex_tree_options_fast_persistence> ST;
+  ST st;
   // insert the proximity graph in the simplex tree
   st.insert_graph(prox_graph);
   // expand the graph until dimension dim_max
@@ -78,7 +79,7 @@ int main(int argc, char * argv[]) {
   st.initialize_filtration();
 
   // Compute the persistence diagram of the complex
-  persistent_cohomology::Persistent_cohomology< Simplex_tree<>, Field_Zp > pcoh(st);
+  persistent_cohomology::Persistent_cohomology<ST, Field_Zp > pcoh(st);
   // initializes the coefficient field for homology
   pcoh.init_coefficients(p);
 
diff --git a/src/Simplex_tree/concept/SimplexTreeOptions.h b/src/Simplex_tree/concept/SimplexTreeOptions.h
index add3ebdd..d072cf34 100644
--- a/src/Simplex_tree/concept/SimplexTreeOptions.h
+++ b/src/Simplex_tree/concept/SimplexTreeOptions.h
@@ -37,5 +37,7 @@ struct SimplexTreeOptions {
   static const bool store_key;
   /// If true, each simplex has extra storage for one `Filtration_value`, and this value is propagated by operations like `Gudhi::Simplex_tree::expansion`. Without it, `Persistent_cohomology` degenerates to computing usual (non-persistent) cohomology.
   static const bool store_filtration;
+  /// If true, the list of vertices present in the complex must always be 0, ..., num_vertices-1, without any hole.
+  static constexpr bool contiguous_vertices;
 };
 
diff --git a/src/Simplex_tree/example/simplex_tree_from_alpha_shapes_3.cpp b/src/Simplex_tree/example/simplex_tree_from_alpha_shapes_3.cpp
index 45efe3ed..49d358ab 100644
--- a/src/Simplex_tree/example/simplex_tree_from_alpha_shapes_3.cpp
+++ b/src/Simplex_tree/example/simplex_tree_from_alpha_shapes_3.cpp
@@ -62,7 +62,8 @@ typedef Alpha_shape_3::Edge Edge;
 typedef std::list<Alpha_shape_3::Vertex_handle> Vertex_list;
 
 // gudhi type definition
-typedef Gudhi::Simplex_tree<>::Vertex_handle Simplex_tree_vertex;
+typedef Gudhi::Simplex_tree<> Simplex_tree;
+typedef Simplex_tree::Vertex_handle Simplex_tree_vertex;
 typedef std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex > Alpha_shape_simplex_tree_map;
 typedef std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex> Alpha_shape_simplex_tree_pair;
 typedef std::vector< Simplex_tree_vertex > Simplex_tree_vector_vertex;
@@ -161,7 +162,7 @@ int main(int argc, char * const argv[]) {
 
   // Loop on objects vector
   Vertex_list vertex_list;
-  Gudhi::Simplex_tree<> simplex_tree;
+  Simplex_tree simplex_tree;
   Alpha_shape_simplex_tree_map map_cgal_simplex_tree;
   std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin();
   for (auto object_iterator : the_objects) {
diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h
index ea61fa2e..408c0588 100644
--- a/src/Simplex_tree/include/gudhi/Simplex_tree.h
+++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h
@@ -78,15 +78,7 @@ namespace Gudhi {
  * @{
  */
 
-/// Model of SimplexTreeOptions.
-struct Simplex_tree_options_full_featured {
-  typedef linear_indexing_tag Indexing_tag;
-  typedef int Vertex_handle;
-  typedef double Filtration_value;
-  typedef int Simplex_key;
-  static const bool store_key = true;
-  static const bool store_filtration = true;
-};
+struct Simplex_tree_options_full_featured;
 
 /**
  * \brief Simplex Tree data structure for representing simplicial complexes.
@@ -568,9 +560,22 @@ class Simplex_tree {
   /** \brief Returns the Simplex_handle corresponding to the 0-simplex
    * representing the vertex with Vertex_handle v. */
   Simplex_handle find_vertex(Vertex_handle v) {
-    return root_.members_.begin() + v;
+    if (Options::contiguous_vertices) {
+      assert(contiguous_vertices());
+      return root_.members_.begin() + v;
+    } else {
+      return root_.members_.find(v);
+    }
+  }
+
+ public:
+  /** \private \brief Test if the vertices have contiguous numbering: 0, 1, etc.  */
+  bool contiguous_vertices() const {
+    if(root_.members_.empty()) return true;
+    if(root_.members_.begin()->first!=0) return false;
+    if(std::prev(root_.members_.end())->first!=root_.members_.size()-1) return false;
+    return true;
   }
-  //{ return root_.members_.find(v); }
 
  private:
   /** \brief Inserts a simplex represented by a vector of vertex.
@@ -1146,6 +1151,31 @@ std::istream& operator>>(std::istream & is, Simplex_tree<T...> & st) {
 
   return is;
 }
+
+/// Model of SimplexTreeOptions.
+struct Simplex_tree_options_full_featured {
+  typedef linear_indexing_tag Indexing_tag;
+  typedef int Vertex_handle;
+  typedef double Filtration_value;
+  typedef int Simplex_key;
+  static const bool store_key = true;
+  static const bool store_filtration = true;
+  static const bool contiguous_vertices = false;
+};
+
+/** Model of SimplexTreeOptions, faster than
+  `Simplex_tree_options_full_featured` but note the unsafe
+  `contiguous_vertices` option. */
+struct Simplex_tree_options_fast_persistence {
+  typedef linear_indexing_tag Indexing_tag;
+  typedef int Vertex_handle;
+  typedef float Filtration_value;
+  typedef int Simplex_key;
+  static const bool store_key = true;
+  static const bool store_filtration = true;
+  static const bool contiguous_vertices = true;
+};
+
 /** @} */  // end defgroup simplex_tree
 
 }  // namespace Gudhi
diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_iterators.h b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_iterators.h
index 372ef329..c5027f22 100644
--- a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_iterators.h
+++ b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_iterators.h
@@ -99,8 +99,7 @@ class Simplex_tree_boundary_simplex_iterator : public boost::iterator_facade<
 
 // any end() iterator
   explicit Simplex_tree_boundary_simplex_iterator(SimplexTree * st)
-      : last_(st->null_vertex()),
-        sib_(NULL) {
+      : sh_(st->null_simplex()) {
   }
 
   Simplex_tree_boundary_simplex_iterator(SimplexTree * st, Simplex_handle sh)
@@ -113,7 +112,7 @@ class Simplex_tree_boundary_simplex_iterator : public boost::iterator_facade<
     if (sib_ != NULL) {
       sh_ = sib_->find(next_);
     } else {
-      last_ = st->null_vertex();
+      sh_ = st->null_simplex();
     }  // vertex: == end()
   }
 
@@ -121,28 +120,42 @@ class Simplex_tree_boundary_simplex_iterator : public boost::iterator_facade<
   friend class boost::iterator_core_access;
 // valid when iterating along the SAME boundary.
   bool equal(Simplex_tree_boundary_simplex_iterator const& other) const {
-    return (sib_ == other.sib_ && last_ == other.last_);
+    return sh_ == other.sh_;
   }
 
   Simplex_handle const& dereference() const {
+    assert(sh_ != st_->null_simplex());
     return sh_;
   }
 
   void increment() {
     if (sib_ == NULL) {
-      last_ = st_->null_vertex();
+      sh_ = st_->null_simplex();
       return;
     }
 
     Siblings * for_sib = sib_;
-    for (auto rit = suffix_.rbegin(); rit != suffix_.rend(); ++rit) {
+    Siblings * new_sib = sib_->oncles();
+    auto rit = suffix_.rbegin();
+    if (SimplexTree::Options::contiguous_vertices
+	&& new_sib == nullptr
+	&& rit != suffix_.rend()) {
+      // We reached the root, use a short-cut to find a vertex. We could also
+      // optimize finding the second vertex of a segment, but people are
+      // expected to call endpoints().
+      assert(st_->contiguous_vertices());
+      sh_ = for_sib->members_.begin()+*rit;
+      for_sib = sh_->second.children();
+      ++rit;
+    }
+    for (; rit != suffix_.rend(); ++rit) {
       sh_ = for_sib->find(*rit);
       for_sib = sh_->second.children();
     }
     sh_ = for_sib->find(last_);  // sh_ points to the right simplex now
     suffix_.push_back(next_);
     next_ = sib_->parent();
-    sib_ = sib_->oncles();
+    sib_ = new_sib;
   }
 
   // Most of the storage should be moved to the range, iterators should be light.
diff --git a/src/Simplex_tree/test/simplex_tree_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
index fff00d77..874c3363 100644
--- a/src/Simplex_tree/test/simplex_tree_unit_test.cpp
+++ b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
@@ -4,10 +4,12 @@
 #include <utility> // std::pair, std::make_pair
 #include <cmath> // float comparison
 #include <limits>
+#include <functional> // greater
 
 #define BOOST_TEST_DYN_LINK
 #define BOOST_TEST_MODULE "simplex_tree"
 #include <boost/test/unit_test.hpp>
+#include <boost/mpl/list.hpp>
 
 //  ^
 // /!\ Nothing else from Simplex_tree shall be included to test includes are well defined.
@@ -15,26 +17,25 @@
 
 using namespace Gudhi;
 
-typedef Simplex_tree<> typeST;
-typedef std::pair<typeST::Simplex_handle, bool> typePairSimplexBool;
-typedef std::vector<Vertex_handle> typeVectorVertex;
-typedef std::pair<typeVectorVertex, Filtration_value> typeSimplex;
+typedef boost::mpl::list<Simplex_tree<>, Simplex_tree<Simplex_tree_options_fast_persistence>> list_of_tested_variants;
 
 const Vertex_handle DEFAULT_VERTEX_HANDLE = (const Vertex_handle) - 1;
 const Filtration_value DEFAULT_FILTRATION_VALUE = (const Filtration_value) 0.0;
 
+template<class typeST>
 void test_empty_simplex_tree(typeST& tst) {
   BOOST_CHECK(tst.null_vertex() == DEFAULT_VERTEX_HANDLE);
   BOOST_CHECK(tst.filtration() == DEFAULT_FILTRATION_VALUE);
   BOOST_CHECK(tst.num_vertices() == (size_t) 0);
   BOOST_CHECK(tst.num_simplices() == (size_t) 0);
-  typeST::Siblings* STRoot = tst.root();
-  BOOST_CHECK(STRoot != NULL);
-  BOOST_CHECK(STRoot->oncles() == NULL);
+  typename typeST::Siblings* STRoot = tst.root();
+  BOOST_CHECK(STRoot != nullptr);
+  BOOST_CHECK(STRoot->oncles() == nullptr);
   BOOST_CHECK(STRoot->parent() == DEFAULT_VERTEX_HANDLE);
   BOOST_CHECK(tst.dimension() == -1);
 }
 
+template<class typeST>
 void test_iterators_on_empty_simplex_tree(typeST& tst) {
   std::cout << "Iterator on vertices: " << std::endl;
   for (auto vertex : tst.complex_vertex_range()) {
@@ -56,8 +57,9 @@ void test_iterators_on_empty_simplex_tree(typeST& tst) {
   }
 }
 
-BOOST_AUTO_TEST_CASE(simplex_tree_when_empty) {
-  const Filtration_value DEFAULT_FILTRATION_VALUE = 0;
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_when_empty, typeST, list_of_tested_variants) {
+  typedef std::pair<typename typeST::Simplex_handle, bool> typePairSimplexBool;
+  typedef std::vector<Vertex_handle> typeVectorVertex;
 
   std::cout << "********************************************************************" << std::endl;
   std::cout << "TEST OF DEFAULT CONSTRUCTOR" << std::endl;
@@ -72,7 +74,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_when_empty) {
   BOOST_CHECK(simplexVectorEmpty.empty() == true);
   typePairSimplexBool returnEmptyValue = st.insert_simplex(simplexVectorEmpty,
                                                            DEFAULT_FILTRATION_VALUE);
-  BOOST_CHECK(returnEmptyValue.first == typeST::Simplex_handle(NULL));
+  BOOST_CHECK(returnEmptyValue.first == typename typeST::Simplex_handle(nullptr));
   BOOST_CHECK(returnEmptyValue.second == true);
 
   test_empty_simplex_tree(st);
@@ -84,7 +86,7 @@ bool AreAlmostTheSame(float a, float b) {
   return std::fabs(a - b) < std::numeric_limits<float>::epsilon();
 }
 
-BOOST_AUTO_TEST_CASE(simplex_tree_from_file) {
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_from_file, typeST, list_of_tested_variants) {
   // TEST OF INSERTION
   std::cout << "********************************************************************" << std::endl;
   std::cout << "TEST OF SIMPLEX TREE FROM A FILE" << std::endl;
@@ -101,7 +103,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_from_file) {
   // Check
   BOOST_CHECK(st.num_simplices() == 143353);
   BOOST_CHECK(st.dimension() == 3);
-  BOOST_CHECK(st.filtration() == 0.4);
+  BOOST_CHECK(AreAlmostTheSame(st.filtration(), 0.4));
 
   int previous_size = 0;
   for (auto f_simplex : st.filtration_simplex_range()) {
@@ -119,6 +121,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_from_file) {
   simplex_tree_stream.close();
 }
 
+template<class typeST, class typeSimplex>
 void test_simplex_tree_contains(typeST& simplexTree, typeSimplex& simplex, int pos) {
   auto f_simplex = simplexTree.filtration_simplex_range().begin() + pos;
 
@@ -135,16 +138,18 @@ void test_simplex_tree_contains(typeST& simplexTree, typeSimplex& simplex, int p
   }
 }
 
+template<class typeST, class typePairSimplexBool>
 void test_simplex_tree_insert_returns_true(const typePairSimplexBool& returnValue) {
   BOOST_CHECK(returnValue.second == true);
-  typeST::Simplex_handle shReturned = returnValue.first; // Simplex_handle = boost::container::flat_map< Vertex_handle, Node >::iterator
-  BOOST_CHECK(shReturned != typeST::Simplex_handle(NULL));
+  typename typeST::Simplex_handle shReturned = returnValue.first; // Simplex_handle = boost::container::flat_map< Vertex_handle, Node >::iterator
+  BOOST_CHECK(shReturned != typename typeST::Simplex_handle(nullptr));
 }
 
 // Global variables
 Filtration_value max_fil = DEFAULT_FILTRATION_VALUE;
 int dim_max = -1;
 
+template<class typeST, class Filtration_value>
 void set_and_test_simplex_tree_dim_fil(typeST& simplexTree, int vectorSize, const Filtration_value& fil) {
   if (vectorSize > dim_max + 1) {
     dim_max = vectorSize - 1;
@@ -173,11 +178,17 @@ void set_and_test_simplex_tree_dim_fil(typeST& simplexTree, int vectorSize, cons
   BOOST_CHECK(simplexTree.num_simplices() == num_simp);
 }
 
-BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_insertion, typeST, list_of_tested_variants) {
+  typedef std::pair<typename typeST::Simplex_handle, bool> typePairSimplexBool;
+  typedef std::vector<Vertex_handle> typeVectorVertex;
+  typedef std::pair<typeVectorVertex, Filtration_value> typeSimplex;
   const Filtration_value FIRST_FILTRATION_VALUE = 0.1;
   const Filtration_value SECOND_FILTRATION_VALUE = 0.2;
   const Filtration_value THIRD_FILTRATION_VALUE = 0.3;
   const Filtration_value FOURTH_FILTRATION_VALUE = 0.4;
+  // reset since we run the test several times
+  dim_max = -1;
+  max_fil = DEFAULT_FILTRATION_VALUE;
 
   // TEST OF INSERTION
   std::cout << "********************************************************************" << std::endl;
@@ -191,7 +202,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex firstSimplex = std::make_pair(firstSimplexVector, Filtration_value(FIRST_FILTRATION_VALUE));
   typePairSimplexBool returnValue = st.insert_simplex(firstSimplex.first, firstSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, firstSimplexVector.size(), firstSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 1);
 
@@ -202,7 +213,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex secondSimplex = std::make_pair(secondSimplexVector, Filtration_value(FIRST_FILTRATION_VALUE));
   returnValue = st.insert_simplex(secondSimplex.first, secondSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, secondSimplexVector.size(), secondSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 2);
 
@@ -213,7 +224,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex thirdSimplex = std::make_pair(thirdSimplexVector, Filtration_value(SECOND_FILTRATION_VALUE));
   returnValue = st.insert_simplex(thirdSimplex.first, thirdSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, thirdSimplexVector.size(), thirdSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 2); // Not incremented !!
 
@@ -224,7 +235,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex fourthSimplex = std::make_pair(fourthSimplexVector, Filtration_value(FIRST_FILTRATION_VALUE));
   returnValue = st.insert_simplex(fourthSimplex.first, fourthSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, fourthSimplexVector.size(), fourthSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 3);
 
@@ -235,7 +246,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex fifthSimplex = std::make_pair(fifthSimplexVector, Filtration_value(SECOND_FILTRATION_VALUE));
   returnValue = st.insert_simplex(fifthSimplex.first, fifthSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, fifthSimplexVector.size(), fifthSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 3); // Not incremented !!
 
@@ -246,7 +257,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex sixthSimplex = std::make_pair(sixthSimplexVector, Filtration_value(SECOND_FILTRATION_VALUE));
   returnValue = st.insert_simplex(sixthSimplex.first, sixthSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, sixthSimplexVector.size(), sixthSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 3); // Not incremented !!
 
@@ -257,7 +268,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex seventhSimplex = std::make_pair(seventhSimplexVector, Filtration_value(THIRD_FILTRATION_VALUE));
   returnValue = st.insert_simplex(seventhSimplex.first, seventhSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, seventhSimplexVector.size(), seventhSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 3); // Not incremented !!
 
@@ -268,7 +279,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex eighthSimplex = std::make_pair(eighthSimplexVector, Filtration_value(FIRST_FILTRATION_VALUE));
   returnValue = st.insert_simplex(eighthSimplex.first, eighthSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, eighthSimplexVector.size(), eighthSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 4);
 
@@ -279,7 +290,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   typeSimplex ninethSimplex = std::make_pair(ninethSimplexVector, Filtration_value(SECOND_FILTRATION_VALUE));
   returnValue = st.insert_simplex(ninethSimplex.first, ninethSimplex.second);
 
-  test_simplex_tree_insert_returns_true(returnValue);
+  test_simplex_tree_insert_returns_true<typeST>(returnValue);
   set_and_test_simplex_tree_dim_fil(st, ninethSimplexVector.size(), ninethSimplex.second);
   BOOST_CHECK(st.num_vertices() == (size_t) 4); // Not incremented !!
 
@@ -292,8 +303,8 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
   returnValue = st.insert_simplex(tenthSimplex.first, tenthSimplex.second);
 
   BOOST_CHECK(returnValue.second == false);
-  typeST::Simplex_handle shReturned = returnValue.first; // Simplex_handle = boost::container::flat_map< Vertex_handle, Node >::iterator
-  BOOST_CHECK(shReturned == typeST::Simplex_handle(NULL));
+  typename typeST::Simplex_handle shReturned = returnValue.first; // Simplex_handle = boost::container::flat_map< Vertex_handle, Node >::iterator
+  BOOST_CHECK(shReturned == typename typeST::Simplex_handle(nullptr));
   BOOST_CHECK(st.num_vertices() == (size_t) 4); // Not incremented !!
   BOOST_CHECK(st.dimension() == dim_max);
   BOOST_CHECK(AreAlmostTheSame(st.filtration(), max_fil));
@@ -307,7 +318,7 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
 
   BOOST_CHECK(returnValue.second == false);
   shReturned = returnValue.first; // Simplex_handle = boost::container::flat_map< Vertex_handle, Node >::iterator
-  BOOST_CHECK(shReturned == typeST::Simplex_handle(NULL));
+  BOOST_CHECK(shReturned == typename typeST::Simplex_handle(nullptr));
   BOOST_CHECK(st.num_vertices() == (size_t) 4); // Not incremented !!
   BOOST_CHECK(st.dimension() == dim_max);
   BOOST_CHECK(AreAlmostTheSame(st.filtration(), max_fil));
@@ -362,9 +373,10 @@ BOOST_AUTO_TEST_CASE(simplex_tree_insertion) {
 
 }
 
-bool sort_in_decr_order (Vertex_handle i,Vertex_handle j) { return (i>j); }
-
-BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
+BOOST_AUTO_TEST_CASE_TEMPLATE(NSimplexAndSubfaces_tree_insertion, typeST, list_of_tested_variants) {
+  typedef std::pair<typename typeST::Simplex_handle, bool> typePairSimplexBool;
+  typedef std::vector<Vertex_handle> typeVectorVertex;
+  typedef std::pair<typeVectorVertex, Filtration_value> typeSimplex;
   std::cout << "********************************************************************" << std::endl;
   std::cout << "TEST OF RECURSIVE INSERTION" << std::endl;
   typeST st;
@@ -382,7 +394,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Check it is well inserted
   BOOST_CHECK(true == returnValue.second);
   position = 0;
-  std::sort(SimplexVector1.begin(), SimplexVector1.end(), sort_in_decr_order);
+  std::sort(SimplexVector1.begin(), SimplexVector1.end(), std::greater<Vertex_handle>());
   for (auto vertex : st.simplex_vertex_range(returnValue.first)) {
     // Check returned Simplex_handle
     std::cout << "vertex = " << vertex << " | vector[" << position << "] = " << SimplexVector1[position] << std::endl;
@@ -401,7 +413,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Check it is well inserted
   BOOST_CHECK(true == returnValue.second);
   position = 0;
-  std::sort(SimplexVector2.begin(), SimplexVector2.end(), sort_in_decr_order);
+  std::sort(SimplexVector2.begin(), SimplexVector2.end(), std::greater<Vertex_handle>());
   for (auto vertex : st.simplex_vertex_range(returnValue.first)) {
     // Check returned Simplex_handle
     std::cout << "vertex = " << vertex << " | vector[" << position << "] = " << SimplexVector2[position] << std::endl;
@@ -420,7 +432,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Check it is well inserted
   BOOST_CHECK(true == returnValue.second);
   position = 0;
-  std::sort(SimplexVector3.begin(), SimplexVector3.end(), sort_in_decr_order);
+  std::sort(SimplexVector3.begin(), SimplexVector3.end(), std::greater<Vertex_handle>());
   for (auto vertex : st.simplex_vertex_range(returnValue.first)) {
     // Check returned Simplex_handle
     std::cout << "vertex = " << vertex << " | vector[" << position << "] = " << SimplexVector3[position] << std::endl;
@@ -450,7 +462,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Check it is well inserted
   BOOST_CHECK(true == returnValue.second);
   position = 0;
-  std::sort(SimplexVector5.begin(), SimplexVector5.end(), sort_in_decr_order);
+  std::sort(SimplexVector5.begin(), SimplexVector5.end(), std::greater<Vertex_handle>());
   for (auto vertex : st.simplex_vertex_range(returnValue.first)) {
     // Check returned Simplex_handle
     std::cout << "vertex = " << vertex << " | vector[" << position << "] = " << SimplexVector5[position] << std::endl;
@@ -469,7 +481,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Check it is well inserted
   BOOST_CHECK(true == returnValue.second);
   position = 0;
-  std::sort(SimplexVector6.begin(), SimplexVector6.end(), sort_in_decr_order);
+  std::sort(SimplexVector6.begin(), SimplexVector6.end(), std::greater<Vertex_handle>());
   for (auto vertex : st.simplex_vertex_range(returnValue.first)) {
     // Check returned Simplex_handle
     std::cout << "vertex = " << vertex << " | vector[" << position << "] = " << SimplexVector6[position] << std::endl;
@@ -509,7 +521,7 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   // Find in the simplex_tree
   // ------------------------------------------------------------------------------------------------------------------
   typeVectorVertex simpleSimplexVector{1};
-  Simplex_tree<>::Simplex_handle simplexFound = st.find(simpleSimplexVector);
+  typename typeST::Simplex_handle simplexFound = st.find(simpleSimplexVector);
   std::cout << "**************IS THE SIMPLEX {1} IN THE SIMPLEX TREE ?\n";
   if (simplexFound != st.null_simplex())
     std::cout << "***+ YES IT IS!\n";
@@ -570,14 +582,15 @@ BOOST_AUTO_TEST_CASE(NSimplexAndSubfaces_tree_insertion) {
   }
 }
 
-void test_cofaces(typeST& st, std::vector<Vertex_handle> expected, int dim, std::vector<typeST::Simplex_handle> res) {
-  typeST::Cofaces_simplex_range cofaces;
+template<class typeST, class Vertex_handle>
+void test_cofaces(typeST& st, const std::vector<Vertex_handle>& expected, int dim, const std::vector<typename typeST::Simplex_handle>& res) {
+  typename typeST::Cofaces_simplex_range cofaces;
   if (dim == 0)
     cofaces = st.star_simplex_range(st.find(expected));
   else
     cofaces = st.cofaces_simplex_range(st.find(expected), dim);
   for (auto simplex = cofaces.begin(); simplex != cofaces.end(); ++simplex) {
-    typeST::Simplex_vertex_range rg = st.simplex_vertex_range(*simplex);
+    typename typeST::Simplex_vertex_range rg = st.simplex_vertex_range(*simplex);
     for (auto vertex = rg.begin(); vertex != rg.end(); ++vertex) {
       std::cout << "(" << *vertex << ")";
     }
@@ -586,7 +599,8 @@ void test_cofaces(typeST& st, std::vector<Vertex_handle> expected, int dim, std:
   }
 }
 
-BOOST_AUTO_TEST_CASE(coface_on_simplex_tree) {
+BOOST_AUTO_TEST_CASE_TEMPLATE(coface_on_simplex_tree, typeST, list_of_tested_variants) {
+  typedef std::vector<Vertex_handle> typeVectorVertex;
   std::cout << "********************************************************************" << std::endl;
   std::cout << "TEST COFACE ALGORITHM" << std::endl;
   typeST st;
@@ -616,7 +630,7 @@ BOOST_AUTO_TEST_CASE(coface_on_simplex_tree) {
   st.set_dimension(3);
 
   std::vector<Vertex_handle> simplex_result;
-  std::vector<typeST::Simplex_handle> result;
+  std::vector<typename typeST::Simplex_handle> result;
   std::cout << "First test - Star of (3):" << std::endl;
 
   simplex_result = {3};
@@ -684,7 +698,8 @@ BOOST_AUTO_TEST_CASE(coface_on_simplex_tree) {
 
 }
 
-BOOST_AUTO_TEST_CASE(copy_move_on_simplex_tree) {
+BOOST_AUTO_TEST_CASE_TEMPLATE(copy_move_on_simplex_tree, typeST, list_of_tested_variants) {
+  typedef std::vector<Vertex_handle> typeVectorVertex;
   std::cout << "********************************************************************" << std::endl;
   std::cout << "TEST COPY MOVE CONSTRUCTORS" << std::endl;
   typeST st;
@@ -744,3 +759,39 @@ BOOST_AUTO_TEST_CASE(copy_move_on_simplex_tree) {
   
   std::cout << "Printing st once again- address = " << &st << std::endl;
 }
+
+template<class typeST>
+void test_simplex_is_vertex(typeST& st, typename typeST::Simplex_handle sh, typename typeST::Vertex_handle v) {
+  BOOST_CHECK(st.dimension(sh) == 0);
+  auto&& r = st.simplex_vertex_range(sh);
+  auto i = std::begin(r);
+  BOOST_CHECK(*i == v);
+  BOOST_CHECK(++i == std::end(r));
+}
+
+BOOST_AUTO_TEST_CASE(non_contiguous) {
+  typedef Simplex_tree<> typeST;
+  typedef typeST::Vertex_handle Vertex_handle;
+  typedef typeST::Simplex_handle Simplex_handle;
+  std::cout << "********************************************************************" << std::endl;
+  std::cout << "TEST NON-CONTIGUOUS VERTICES" << std::endl;
+  typeST st;
+  Vertex_handle e[] = {3,-7};
+  std::cout << "Insert" << std::endl;
+  st.insert_simplex_and_subfaces(e);
+  BOOST_CHECK(st.num_vertices() == 2);
+  BOOST_CHECK(st.num_simplices() == 3);
+  std::cout << "Find" << std::endl;
+  Simplex_handle sh = st.find(e);
+  BOOST_CHECK(sh != st.null_simplex());
+  std::cout << "Endpoints" << std::endl;
+  auto p = st.endpoints(sh);
+  test_simplex_is_vertex(st, p.first, 3);
+  test_simplex_is_vertex(st, p.second, -7);
+  std::cout << "Boundary" << std::endl;
+  auto&& b = st.boundary_simplex_range(sh);
+  auto i = std::begin(b);
+  test_simplex_is_vertex(st, *i, -7);
+  test_simplex_is_vertex(st, *++i, 3);
+  BOOST_CHECK(++i == std::end(b));
+}