Merge last trunk modifications

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

Former-commit-id: 2a37c02ca51d0df28a981207db6ecc1bd0991ffa

8 files changed, 537 insertions, 59 deletions

diff --git a/src/Simplex_tree/doc/Intro_simplex_tree.h b/src/Simplex_tree/doc/Intro_simplex_tree.h
index f5b72ff6..769491d9 100644
--- a/src/Simplex_tree/doc/Intro_simplex_tree.h
+++ b/src/Simplex_tree/doc/Intro_simplex_tree.h
@@ -67,10 +67,13 @@ Information of the Simplex Tree:
 Number of vertices = 10   Number of simplices = 98 \endcode
  * 
  * \li <a href="_simplex_tree_2example_alpha_shapes_3_simplex_tree_from_off_file_8cpp-example.html">
- * Simplex_tree/example_alpha_shapes_3_simplex_tree_from_off_file.cpp</a> - Simplex tree is computed and displayed from a 3D alpha
- * complex (Requires CGAL, GMP and GMPXX to be installed)
- * 
+ * Simplex_tree/example_alpha_shapes_3_simplex_tree_from_off_file.cpp</a> - Simplex tree is computed and displayed
+ * from a 3D alpha complex (Requires CGAL, GMP and GMPXX to be installed).
  * 
+ * \li <a href="_simplex_tree_2graph_expansion_with_blocker_8cpp-example.html">
+ * Simplex_tree/graph_expansion_with_blocker.cpp</a> - Simple simplex tree construction from a one-skeleton graph with
+ * a simple blocker expansion method.
+ *
  * \subsection filteredcomplexeshassecomplex Hasse complex
  * The second one is the Hasse_complex. The Hasse complex is a data structure representing explicitly all co-dimension
  * 1 incidence relations in a complex. It is consequently faster when accessing the boundary of a simplex, but is less
diff --git a/src/Simplex_tree/example/CMakeLists.txt b/src/Simplex_tree/example/CMakeLists.txt
index e22cc92c..8bc4ad53 100644
--- a/src/Simplex_tree/example/CMakeLists.txt
+++ b/src/Simplex_tree/example/CMakeLists.txt
@@ -34,5 +34,12 @@ if(GMP_FOUND AND CGAL_FOUND)
       "${CMAKE_SOURCE_DIR}/data/points/bunny_5000.off")
 
   install(TARGETS Simplex_tree_example_alpha_shapes_3_from_off DESTINATION bin)
+endif()
 
+add_executable ( Simplex_tree_example_graph_expansion_with_blocker graph_expansion_with_blocker.cpp )
+if (TBB_FOUND)
+  target_link_libraries(Simplex_tree_example_graph_expansion_with_blocker ${TBB_LIBRARIES})
 endif()
+add_test(NAME Simplex_tree_example_graph_expansion_with_blocker COMMAND $<TARGET_FILE:Simplex_tree_example_graph_expansion_with_blocker>)
+
+install(TARGETS Simplex_tree_example_graph_expansion_with_blocker DESTINATION bin)
diff --git a/src/Simplex_tree/example/graph_expansion_with_blocker.cpp b/src/Simplex_tree/example/graph_expansion_with_blocker.cpp
new file mode 100644
index 00000000..86bfb8cb
--- /dev/null
+++ b/src/Simplex_tree/example/graph_expansion_with_blocker.cpp
@@ -0,0 +1,79 @@
+/*    This file is part of the Gudhi Library. The Gudhi library
+ *    (Geometric Understanding in Higher Dimensions) is a generic C++
+ *    library for computational topology.
+ *
+ *    Author(s):       Vincent Rouvreau
+ *
+ *    Copyright (C) 2014
+ *
+ *    This program is free software: you can redistribute it and/or modify
+ *    it under the terms of the GNU General Public License as published by
+ *    the Free Software Foundation, either version 3 of the License, or
+ *    (at your option) any later version.
+ *
+ *    This program is distributed in the hope that it will be useful,
+ *    but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *    GNU General Public License for more details.
+ *
+ *    You should have received a copy of the GNU General Public License
+ *    along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <gudhi/Simplex_tree.h>
+
+#include <iostream>
+
+using Simplex_tree = Gudhi::Simplex_tree<>;
+using Simplex_handle = Simplex_tree::Simplex_handle;
+
+int main(int argc, char * const argv[]) {
+
+  // Construct the Simplex Tree with a 1-skeleton graph example
+  Simplex_tree simplexTree;
+
+  simplexTree.insert_simplex({0, 1}, 0.);
+  simplexTree.insert_simplex({0, 2}, 1.);
+  simplexTree.insert_simplex({0, 3}, 2.);
+  simplexTree.insert_simplex({1, 2}, 3.);
+  simplexTree.insert_simplex({1, 3}, 4.);
+  simplexTree.insert_simplex({2, 3}, 5.);
+  simplexTree.insert_simplex({2, 4}, 6.);
+  simplexTree.insert_simplex({3, 6}, 7.);
+  simplexTree.insert_simplex({4, 5}, 8.);
+  simplexTree.insert_simplex({4, 6}, 9.);
+  simplexTree.insert_simplex({5, 6}, 10.);
+  simplexTree.insert_simplex({6}, 10.);
+
+  simplexTree.expansion_with_blockers(3, [&](Simplex_handle sh){
+      bool result = false;
+      std::cout << "Blocker on [";
+      // User can loop on the vertices from the given simplex_handle i.e.
+      for (auto vertex : simplexTree.simplex_vertex_range(sh)) {
+        // We block the expansion, if the vertex '6' is in the given list of vertices
+        if (vertex == 6)
+          result = true;
+        std::cout << vertex << ", ";
+      }
+      std::cout << "] ( " << simplexTree.filtration(sh);
+      // User can re-assign a new filtration value directly in the blocker (default is the maximal value of boudaries)
+      simplexTree.assign_filtration(sh, simplexTree.filtration(sh) + 1.);
+
+      std::cout << " + 1. ) = " << result << std::endl;
+
+      return result;
+    });
+
+  std::cout << "********************************************************************\n";
+  std::cout << "* The complex contains " << simplexTree.num_simplices() << " simplices";
+  std::cout << "   - dimension " << simplexTree.dimension() << "\n";
+  std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
+  for (auto f_simplex : simplexTree.filtration_simplex_range()) {
+    std::cout << "   " << "[" << simplexTree.filtration(f_simplex) << "] ";
+    for (auto vertex : simplexTree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  return 0;
+}
diff --git a/src/Simplex_tree/example/simple_simplex_tree.cpp b/src/Simplex_tree/example/simple_simplex_tree.cpp
index 26d663fb..b6b65b88 100644
--- a/src/Simplex_tree/example/simple_simplex_tree.cpp
+++ b/src/Simplex_tree/example/simple_simplex_tree.cpp
@@ -185,18 +185,16 @@ int main(int argc, char * const argv[]) {
   }
 
   // ++ GENERAL VARIABLE SET
-  simplexTree.set_filtration(FOURTH_FILTRATION_VALUE);  // Max filtration value
 
   std::cout << "********************************************************************\n";
   // Display the Simplex_tree - Can not be done in the middle of 2 inserts
   std::cout << "* The complex contains " << simplexTree.num_simplices() << " simplices\n";
-  std::cout << "   - dimension " << simplexTree.dimension() << "   - filtration " << simplexTree.filtration() << "\n";
+  std::cout << "   - dimension " << simplexTree.dimension() << "\n";
   std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
   for (auto f_simplex : simplexTree.filtration_simplex_range()) {
     std::cout << "   " << "[" << simplexTree.filtration(f_simplex) << "] ";
-    for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) {
-      std::cout << static_cast<int>(vertex) << " ";
-    }
+    for (auto vertex : simplexTree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
     std::cout << std::endl;
   }
   //   [0.1] 0
@@ -249,5 +247,34 @@ int main(int argc, char * const argv[]) {
     std::cout << "***+ YES IT IS!\n";
   else
     std::cout << "***- NO IT ISN'T\n";
+
+  simplexFound = simplexTree.find({ 0, 1 });
+  std::cout << "**************IS THE SIMPLEX {0,1} IN THE SIMPLEX TREE ?\n";
+  if (simplexFound != simplexTree.null_simplex())
+    std::cout << "***+ YES IT IS!\n";
+  else
+    std::cout << "***- NO IT ISN'T\n";
+
+  std::cout << "**************COFACES OF {0,1} IN CODIMENSION 1 ARE\n";
+  for (auto& simplex : simplexTree.cofaces_simplex_range(simplexTree.find({0,1}), 1)) {
+    for (auto vertex : simplexTree.simplex_vertex_range(simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  std::cout << "**************STARS OF {0,1} ARE\n";
+  for (auto& simplex : simplexTree.star_simplex_range(simplexTree.find({0,1}))) {
+    for (auto vertex : simplexTree.simplex_vertex_range(simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  std::cout << "**************BOUNDARIES OF {0,1,2} ARE\n";
+  for (auto& simplex : simplexTree.boundary_simplex_range(simplexTree.find({0,1,2}))) {
+    for (auto vertex : simplexTree.simplex_vertex_range(simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
   return 0;
 }
diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h
index 54f5de13..b2d380ea 100644
--- a/src/Simplex_tree/include/gudhi/Simplex_tree.h
+++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h
@@ -289,7 +289,6 @@ class Simplex_tree {
   /** \brief Constructs an empty simplex tree. */
   Simplex_tree()
       : null_vertex_(-1),
-      threshold_(0),
       root_(nullptr, null_vertex_),
       filtration_vect_(),
       dimension_(-1) { }
@@ -297,7 +296,6 @@ class Simplex_tree {
   /** \brief User-defined copy constructor reproduces the whole tree structure. */
   Simplex_tree(const Simplex_tree& simplex_source)
       : null_vertex_(simplex_source.null_vertex_),
-      threshold_(simplex_source.threshold_),
       root_(nullptr, null_vertex_ , simplex_source.root_.members_),
       filtration_vect_(),
       dimension_(simplex_source.dimension_) {
@@ -323,12 +321,10 @@ class Simplex_tree {
   /** \brief User-defined move constructor moves the whole tree structure. */
   Simplex_tree(Simplex_tree && old)
       : null_vertex_(std::move(old.null_vertex_)),
-      threshold_(std::move(old.threshold_)),
       root_(std::move(old.root_)),
       filtration_vect_(std::move(old.filtration_vect_)),
       dimension_(std::move(old.dimension_)) {
     old.dimension_ = -1;
-    old.threshold_ = 0;
     old.root_ = Siblings(nullptr, null_vertex_);
   }
 
@@ -356,7 +352,6 @@ class Simplex_tree {
   /** \brief Checks if two simplex trees are equal. */
   bool operator==(Simplex_tree& st2) {
     if ((null_vertex_ != st2.null_vertex_) ||
-        (threshold_ != st2.threshold_) ||
         (dimension_ != st2.dimension_))
       return false;
     return rec_equal(&root_, &st2.root_);
@@ -407,14 +402,14 @@ class Simplex_tree {
 
   /** \brief Returns the filtration value of a simplex.
    *
-   * Called on the null_simplex, returns INFINITY.
+   * Called on the null_simplex, it returns infinity.
    * If SimplexTreeOptions::store_filtration is false, returns 0.
    */
   static Filtration_value filtration(Simplex_handle sh) {
     if (sh != null_simplex()) {
       return sh->second.filtration();
     } else {
-      return INFINITY;
+      return std::numeric_limits<Filtration_value>::infinity();
     }
   }
 
@@ -427,11 +422,6 @@ class Simplex_tree {
     sh->second.assign_filtration(fv);
   }
 
-  /** \brief Returns an upper bound of the filtration values of the simplices. */
-  Filtration_value filtration() const {
-    return threshold_;
-  }
-
   /** \brief Returns a Simplex_handle different from all Simplex_handles
    * associated to the simplices in the simplicial complex.
    *
@@ -770,11 +760,6 @@ class Simplex_tree {
     return &root_;
   }
 
-  /** Set an upper bound for the filtration values. */
-  void set_filtration(Filtration_value fil) {
-      threshold_ = fil;
-  }
-
   /** Set a dimension for the simplicial complex. */
   void set_dimension(int dimension) {
     dimension_ = dimension;
@@ -1095,6 +1080,118 @@ class Simplex_tree {
   }
 
  public:
+  /** \brief Expands a simplex tree containing only a graph. Simplices corresponding to cliques in the graph are added
+   * incrementally, faces before cofaces, unless the simplex has dimension larger than `max_dim` or `block_simplex`
+   * returns true for this simplex.
+   *
+   * @param[in] max_dim Expansion maximal dimension value.
+   * @param[in] block_simplex Blocker oracle. Its concept is <CODE>bool block_simplex(Simplex_handle sh)</CODE>
+   *
+   * The function identifies a candidate simplex whose faces are all already in the complex, inserts
+   * it with a filtration value corresponding to the maximum of the filtration values of the faces, then calls
+   * `block_simplex` on a `Simplex_handle` for this new simplex. If `block_simplex` returns true, the simplex is
+   * removed, otherwise it is kept. Note that the evaluation of `block_simplex` is a good time to update the
+   * filtration value of the simplex if you want a customized value. The algorithm then proceeds with the next
+   * candidate.
+   *
+   * @warning several candidates of the same dimension may be inserted simultaneously before calling `block_simplex`,
+   * so if you examine the complex in `block_simplex`, you may hit a few simplices of the same dimension that have not
+   * been vetted by `block_simplex` yet, or have already been rejected but not yet removed.
+   */
+  template< typename Blocker >
+  void expansion_with_blockers(int max_dim, Blocker block_simplex) {
+    // Loop must be from the end to the beginning, as higher dimension simplex are always on the left part of the tree
+    for (auto& simplex : boost::adaptors::reverse(root_.members())) {
+      if (has_children(&simplex)) {
+        siblings_expansion_with_blockers(simplex.second.children(), max_dim, max_dim - 1, block_simplex);
+      }
+    }
+  }
+
+ private:
+  /** \brief Recursive expansion with blockers of the simplex tree.*/
+  template< typename Blocker >
+  void siblings_expansion_with_blockers(Siblings* siblings, int max_dim, int k, Blocker block_simplex) {
+    if (dimension_ < max_dim - k) {
+      dimension_ = max_dim - k;
+    }
+    if (k == 0)
+      return;
+    // No need to go deeper
+    if (siblings->members().size() < 2)
+      return;
+    // Reverse loop starting before the last one for 'next' to be the last one
+    for (auto simplex = siblings->members().rbegin() + 1; simplex != siblings->members().rend(); simplex++) {
+      std::vector<std::pair<Vertex_handle, Node> > intersection;
+      for(auto next = siblings->members().rbegin(); next != simplex; next++) {
+        bool to_be_inserted = true;
+        Filtration_value filt = simplex->second.filtration();
+        // If all the boundaries are present, 'next' needs to be inserted
+        for (Simplex_handle border : boundary_simplex_range(simplex)) {
+          Simplex_handle border_child = find_child(border, next->first);
+          if (border_child == null_simplex()) {
+            to_be_inserted=false;
+            break;
+          }
+          filt = std::max(filt, filtration(border_child));
+        }
+        if (to_be_inserted) {
+          intersection.emplace_back(next->first, Node(nullptr, filt));
+        }
+      }
+      if (intersection.size() != 0) {
+        // Reverse the order to insert
+        Siblings * new_sib = new Siblings(siblings,  // oncles
+                                          simplex->first,  // parent
+                                          boost::adaptors::reverse(intersection));  // boost::container::ordered_unique_range_t
+        std::vector<Simplex_handle> blocked_new_sib_list;
+        // As all intersections are inserted, we can call the blocker function on all new_sib members
+        for (auto new_sib_member = new_sib->members().begin();
+             new_sib_member != new_sib->members().end();
+             new_sib_member++) {
+           bool blocker_result = block_simplex(new_sib_member);
+           // new_sib member has been blocked by the blocker function
+           // add it to the list to be removed - do not perform it while looping on it
+           if (blocker_result)
+             blocked_new_sib_list.push_back(new_sib_member);
+        }
+        bool removed = false;
+        for (auto& blocked_new_sib_member : blocked_new_sib_list){
+          removed = removed || remove_maximal_simplex(blocked_new_sib_member);
+        }
+        if (removed) {
+          // ensure the children property
+          simplex->second.assign_children(siblings);
+        } else {
+          // ensure recursive call
+          simplex->second.assign_children(new_sib);
+          siblings_expansion_with_blockers(new_sib, max_dim, k - 1, block_simplex);
+        }
+      } else {
+        // ensure the children property
+        simplex->second.assign_children(siblings);
+      }
+    }
+  }
+
+  /* \private Returns the Simplex_handle composed of the vertex list (from the Simplex_handle), plus the given
+   * Vertex_handle if the Vertex_handle is found in the Simplex_handle children list.
+   * Returns null_simplex() if it does not exist
+  */
+  Simplex_handle find_child(Simplex_handle sh, Vertex_handle vh) const {
+    if (!has_children(sh))
+      return null_simplex();
+
+    Simplex_handle child = sh->second.children()->find(vh);
+    // Specific case of boost::flat_map does not find, returns boost::flat_map::end()
+    // in simplex tree we want a null_simplex()
+    if (child == sh->second.children()->members().end())
+      return null_simplex();
+
+    return child;
+  }
+
+ public:
   /** \brief Write the hasse diagram of the simplicial complex in os.
    *
    * Each row in the file correspond to a simplex. A line is written:
@@ -1240,14 +1337,16 @@ class Simplex_tree {
  public:
   /** \brief Remove a maximal simplex.
    * @param[in] sh Simplex handle on the maximal simplex to remove.
+   * @return a boolean value that is an implementation detail, and that the user is supposed to ignore
    * \pre Please check the simplex has no coface before removing it.
    * \exception std::invalid_argument In debug mode, if sh has children.
    * \post Be aware that removing is shifting data in a flat_map (`initialize_filtration()` to be done).
    * \post Note that the dimension of the simplicial complex may be lower after calling `remove_maximal_simplex()`
    * than it was before. However, `upper_bond_dimension()` will return the old value, which remains a valid upper
    * bound. If you care, you can call `dimension()` to recompute the exact dimension.
+   * \internal @return true if the leaf's branch has no other leaves (branch's children has been re-assigned), false otherwise.
    */
-  void remove_maximal_simplex(Simplex_handle sh) {
+  bool remove_maximal_simplex(Simplex_handle sh) {
     // Guarantee the simplex has no children
     GUDHI_CHECK(!has_children(sh),
                 std::invalid_argument("Simplex_tree::remove_maximal_simplex - argument has children"));
@@ -1265,13 +1364,13 @@ class Simplex_tree {
       delete child;
       // dimension may need to be lowered
       dimension_to_be_lowered_ = true;
+      return true;
     }
+    return false;
   }
 
  private:
   Vertex_handle null_vertex_;
-  /** \brief Upper bound on the filtration values of the simplices.*/
-  Filtration_value threshold_;
   /** \brief Total number of simplices in the complex, without the empty simplex.*/
   /** \brief Set of simplex tree Nodes representing the vertices.*/
   Siblings root_;
@@ -1300,17 +1399,19 @@ std::istream& operator>>(std::istream & is, Simplex_tree<T...> & st) {
   typedef Simplex_tree<T...> ST;
   std::vector<typename ST::Vertex_handle> simplex;
   typename ST::Filtration_value fil;
-  typename ST::Filtration_value max_fil = 0;
+  int max_dim = -1;
   while (read_simplex(is, simplex, fil)) {
     // read all simplices in the file as a list of vertices
-    if (max_fil < fil) {
-      max_fil = fil;
+    // Warning : simplex_size needs to be casted in int - Can be 0
+    int dim = static_cast<int> (simplex.size() - 1);
+    if (max_dim < dim) {
+      max_dim = dim;
     }
     // insert every simplex in the simplex tree
     st.insert_simplex(simplex, fil);
     simplex.clear();
   }
-  st.set_filtration(max_fil);
+  st.set_dimension(max_dim);
 
   return is;
 }
diff --git a/src/Simplex_tree/test/simplex_tree_graph_expansion_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_graph_expansion_unit_test.cpp
new file mode 100644
index 00000000..19ce3321
--- /dev/null
+++ b/src/Simplex_tree/test/simplex_tree_graph_expansion_unit_test.cpp
@@ -0,0 +1,235 @@
+#include <iostream>
+#include <fstream>
+#include <string>
+#include <algorithm>
+#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.
+#include "gudhi/Simplex_tree.h"
+
+using namespace Gudhi;
+
+typedef boost::mpl::list<Simplex_tree<>, Simplex_tree<Simplex_tree_options_fast_persistence>> list_of_tested_variants;
+
+
+bool AreAlmostTheSame(float a, float b) {
+  return std::fabs(a - b) < std::numeric_limits<float>::epsilon();
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_expansion_with_blockers_3, typeST, list_of_tested_variants) {
+  using Simplex_handle = typename typeST::Simplex_handle;
+  // Construct the Simplex Tree with a 1-skeleton graph example
+  typeST simplex_tree;
+
+  simplex_tree.insert_simplex({0, 1}, 0.);
+  simplex_tree.insert_simplex({0, 2}, 1.);
+  simplex_tree.insert_simplex({0, 3}, 2.);
+  simplex_tree.insert_simplex({1, 2}, 3.);
+  simplex_tree.insert_simplex({1, 3}, 4.);
+  simplex_tree.insert_simplex({2, 3}, 5.);
+  simplex_tree.insert_simplex({2, 4}, 6.);
+  simplex_tree.insert_simplex({3, 6}, 7.);
+  simplex_tree.insert_simplex({4, 5}, 8.);
+  simplex_tree.insert_simplex({4, 6}, 9.);
+  simplex_tree.insert_simplex({5, 6}, 10.);
+  simplex_tree.insert_simplex({6}, 10.);
+
+  simplex_tree.expansion_with_blockers(3, [&](Simplex_handle sh){
+      bool result = false;
+      std::cout << "Blocker on [";
+      // User can loop on the vertices from the given simplex_handle i.e.
+      for (auto vertex : simplex_tree.simplex_vertex_range(sh)) {
+        // We block the expansion, if the vertex '6' is in the given list of vertices
+        if (vertex == 6)
+          result = true;
+        std::cout << vertex << ", ";
+      }
+      std::cout << "] ( " << simplex_tree.filtration(sh);
+      // User can re-assign a new filtration value directly in the blocker (default is the maximal value of boudaries)
+      simplex_tree.assign_filtration(sh, simplex_tree.filtration(sh) + 1.);
+
+      std::cout << " + 1. ) = " << result << std::endl;
+
+      return result;
+    });
+
+  std::cout << "********************************************************************\n";
+  std::cout << "simplex_tree_expansion_with_blockers_3\n";
+  std::cout << "********************************************************************\n";
+  std::cout << "* The complex contains " << simplex_tree.num_simplices() << " simplices";
+  std::cout << " - dimension " << simplex_tree.dimension() << "\n";
+  std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
+  for (auto f_simplex : simplex_tree.filtration_simplex_range()) {
+    std::cout << "   " << "[" << simplex_tree.filtration(f_simplex) << "] ";
+    for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  BOOST_CHECK(simplex_tree.num_simplices() == 23);
+  BOOST_CHECK(simplex_tree.dimension() == 3);
+  // {4, 5, 6} shall be blocked
+  BOOST_CHECK(simplex_tree.find({4, 5, 6}) == simplex_tree.null_simplex());
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2})), 4.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,3})), 5.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,2,3})), 6.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({1,2,3})), 6.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2,3})), 7.));
+
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_expansion_with_blockers_2, typeST, list_of_tested_variants) {
+  using Simplex_handle = typename typeST::Simplex_handle;
+  // Construct the Simplex Tree with a 1-skeleton graph example
+  typeST simplex_tree;
+
+  simplex_tree.insert_simplex({0, 1}, 0.);
+  simplex_tree.insert_simplex({0, 2}, 1.);
+  simplex_tree.insert_simplex({0, 3}, 2.);
+  simplex_tree.insert_simplex({1, 2}, 3.);
+  simplex_tree.insert_simplex({1, 3}, 4.);
+  simplex_tree.insert_simplex({2, 3}, 5.);
+  simplex_tree.insert_simplex({2, 4}, 6.);
+  simplex_tree.insert_simplex({3, 6}, 7.);
+  simplex_tree.insert_simplex({4, 5}, 8.);
+  simplex_tree.insert_simplex({4, 6}, 9.);
+  simplex_tree.insert_simplex({5, 6}, 10.);
+  simplex_tree.insert_simplex({6}, 10.);
+
+  simplex_tree.expansion_with_blockers(2, [&](Simplex_handle sh){
+      bool result = false;
+      std::cout << "Blocker on [";
+      // User can loop on the vertices from the given simplex_handle i.e.
+      for (auto vertex : simplex_tree.simplex_vertex_range(sh)) {
+        // We block the expansion, if the vertex '6' is in the given list of vertices
+        if (vertex == 6)
+          result = true;
+        std::cout << vertex << ", ";
+      }
+      std::cout << "] ( " << simplex_tree.filtration(sh);
+      // User can re-assign a new filtration value directly in the blocker (default is the maximal value of boudaries)
+      simplex_tree.assign_filtration(sh, simplex_tree.filtration(sh) + 1.);
+
+      std::cout << " + 1. ) = " << result << std::endl;
+
+      return result;
+    });
+
+  std::cout << "********************************************************************\n";
+  std::cout << "simplex_tree_expansion_with_blockers_2\n";
+  std::cout << "********************************************************************\n";
+  std::cout << "* The complex contains " << simplex_tree.num_simplices() << " simplices";
+  std::cout << " - dimension " << simplex_tree.dimension() << "\n";
+  std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
+  for (auto f_simplex : simplex_tree.filtration_simplex_range()) {
+    std::cout << "   " << "[" << simplex_tree.filtration(f_simplex) << "] ";
+    for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  BOOST_CHECK(simplex_tree.num_simplices() == 22);
+  BOOST_CHECK(simplex_tree.dimension() == 2);
+  // {4, 5, 6} shall be blocked
+  BOOST_CHECK(simplex_tree.find({4, 5, 6}) == simplex_tree.null_simplex());
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2})), 4.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,3})), 5.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,2,3})), 6.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({1,2,3})), 6.));
+  BOOST_CHECK(simplex_tree.find({0,1,2,3}) == simplex_tree.null_simplex());
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_expansion, typeST, list_of_tested_variants) {
+  // Construct the Simplex Tree with a 1-skeleton graph example
+  typeST simplex_tree;
+
+  simplex_tree.insert_simplex({0, 1}, 0.);
+  simplex_tree.insert_simplex({0, 2}, 1.);
+  simplex_tree.insert_simplex({0, 3}, 2.);
+  simplex_tree.insert_simplex({1, 2}, 3.);
+  simplex_tree.insert_simplex({1, 3}, 4.);
+  simplex_tree.insert_simplex({2, 3}, 5.);
+  simplex_tree.insert_simplex({2, 4}, 6.);
+  simplex_tree.insert_simplex({3, 6}, 7.);
+  simplex_tree.insert_simplex({4, 5}, 8.);
+  simplex_tree.insert_simplex({4, 6}, 9.);
+  simplex_tree.insert_simplex({5, 6}, 10.);
+  simplex_tree.insert_simplex({6}, 10.);
+
+  simplex_tree.expansion(3);
+  std::cout << "********************************************************************\n";
+  std::cout << "simplex_tree_expansion_3\n";
+  std::cout << "********************************************************************\n";
+  std::cout << "* The complex contains " << simplex_tree.num_simplices() << " simplices";
+  std::cout << " - dimension " << simplex_tree.dimension() << "\n";
+  std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
+  for (auto f_simplex : simplex_tree.filtration_simplex_range()) {
+    std::cout << "   " << "[" << simplex_tree.filtration(f_simplex) << "] ";
+    for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  BOOST_CHECK(simplex_tree.num_simplices() == 24);
+  BOOST_CHECK(simplex_tree.dimension() == 3);
+
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({4,5,6})), 10.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2})), 3.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,3})), 4.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,2,3})), 5.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({1,2,3})), 5.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2,3})), 5.));
+
+}
+
+BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_expansion_2, typeST, list_of_tested_variants) {
+  // Construct the Simplex Tree with a 1-skeleton graph example
+  typeST simplex_tree;
+
+  simplex_tree.insert_simplex({0, 1}, 0.);
+  simplex_tree.insert_simplex({0, 2}, 1.);
+  simplex_tree.insert_simplex({0, 3}, 2.);
+  simplex_tree.insert_simplex({1, 2}, 3.);
+  simplex_tree.insert_simplex({1, 3}, 4.);
+  simplex_tree.insert_simplex({2, 3}, 5.);
+  simplex_tree.insert_simplex({2, 4}, 6.);
+  simplex_tree.insert_simplex({3, 6}, 7.);
+  simplex_tree.insert_simplex({4, 5}, 8.);
+  simplex_tree.insert_simplex({4, 6}, 9.);
+  simplex_tree.insert_simplex({5, 6}, 10.);
+  simplex_tree.insert_simplex({6}, 10.);
+
+  simplex_tree.expansion(2);
+
+  std::cout << "********************************************************************\n";
+  std::cout << "simplex_tree_expansion_2\n";
+  std::cout << "********************************************************************\n";
+  std::cout << "* The complex contains " << simplex_tree.num_simplices() << " simplices";
+  std::cout << " - dimension " << simplex_tree.dimension() << "\n";
+  std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n";
+  for (auto f_simplex : simplex_tree.filtration_simplex_range()) {
+    std::cout << "   " << "[" << simplex_tree.filtration(f_simplex) << "] ";
+    for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex))
+      std::cout << "(" << vertex << ")";
+    std::cout << std::endl;
+  }
+
+  BOOST_CHECK(simplex_tree.num_simplices() == 23);
+  BOOST_CHECK(simplex_tree.dimension() == 2);
+
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({4,5,6})), 10.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,2})), 3.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,1,3})), 4.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({0,2,3})), 5.));
+  BOOST_CHECK(AreAlmostTheSame(simplex_tree.filtration(simplex_tree.find({1,2,3})), 5.));
+  BOOST_CHECK(simplex_tree.find({0,1,2,3}) == simplex_tree.null_simplex());
+}
diff --git a/src/Simplex_tree/test/simplex_tree_iostream_operator_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_iostream_operator_unit_test.cpp
index fed764d9..ecb9f025 100644
--- a/src/Simplex_tree/test/simplex_tree_iostream_operator_unit_test.cpp
+++ b/src/Simplex_tree/test/simplex_tree_iostream_operator_unit_test.cpp
@@ -33,13 +33,11 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(iostream_operator, Stree_type, list_of_tested_vari
   st.insert_simplex_and_subfaces({3, 4, 5}, 3.0);
   st.insert_simplex_and_subfaces({3, 0}, 2.0);
   st.insert_simplex_and_subfaces({2, 1, 0}, 3.0);
-  // FIXME
-  st.set_filtration(4.0);
 
   st.initialize_filtration();
   // Display the Simplex_tree
-  std::cout << "The ORIGINAL complex contains " << st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << st.dimension() << "   - filtration " << st.filtration() << std::endl;
+  std::cout << "The ORIGINAL complex contains " << st.num_simplices() << " simplices - dimension = "
+            << st.dimension() << std::endl;
   std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
   for (auto f_simplex : st.filtration_simplex_range()) {
     std::cout << "   " << "[" << st.filtration(f_simplex) << "] ";
@@ -67,8 +65,8 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(iostream_operator, Stree_type, list_of_tested_vari
   simplex_tree_istream >> read_st;
 
   // Display the Simplex_tree
-  std::cout << "The READ complex contains " << read_st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << read_st.dimension() << "   - filtration " << read_st.filtration() << std::endl;
+  std::cout << "The READ complex contains " << read_st.num_simplices() << " simplices - dimension = "
+            << read_st.dimension() << std::endl;
   std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
   for (auto f_simplex : read_st.filtration_simplex_range()) {
     std::cout << "   " << "[" << read_st.filtration(f_simplex) << "] ";
@@ -95,9 +93,8 @@ BOOST_AUTO_TEST_CASE(mini_iostream_operator) {
 
   st.initialize_filtration();
   // Display the Simplex_tree
-  std::cout << "The ORIGINAL complex contains " << st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << st.dimension() << "   - filtration " << st.filtration() << std::endl;
-  std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
+  std::cout << "The ORIGINAL complex contains " << st.num_simplices() << " simplices - dimension = "
+            << st.dimension() << std::endl;
   for (auto f_simplex : st.filtration_simplex_range()) {
     std::cout << "   " << "[" << st.filtration(f_simplex) << "] ";
     for (auto vertex : st.simplex_vertex_range(f_simplex)) {
@@ -124,8 +121,8 @@ BOOST_AUTO_TEST_CASE(mini_iostream_operator) {
   simplex_tree_istream >> read_st;
 
   // Display the Simplex_tree
-  std::cout << "The READ complex contains " << read_st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << read_st.dimension() << "   - filtration " << read_st.filtration() << std::endl;
+  std::cout << "The READ complex contains " << read_st.num_simplices() << " simplices - dimension = "
+            << read_st.dimension() << std::endl;
   std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
   for (auto f_simplex : read_st.filtration_simplex_range()) {
     std::cout << "   " << "[" << read_st.filtration(f_simplex) << "] ";
diff --git a/src/Simplex_tree/test/simplex_tree_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
index f67ff010..580d610a 100644
--- a/src/Simplex_tree/test/simplex_tree_unit_test.cpp
+++ b/src/Simplex_tree/test/simplex_tree_unit_test.cpp
@@ -26,7 +26,6 @@ void test_empty_simplex_tree(typeST& tst) {
   typedef typename typeST::Vertex_handle Vertex_handle;
   const Vertex_handle DEFAULT_VERTEX_VALUE = Vertex_handle(- 1);
   BOOST_CHECK(tst.null_vertex() == DEFAULT_VERTEX_VALUE);
-  BOOST_CHECK(tst.filtration() == 0.0);
   BOOST_CHECK(tst.num_vertices() == (size_t) 0);
   BOOST_CHECK(tst.num_simplices() == (size_t) 0);
   typename typeST::Siblings* STRoot = tst.root();
@@ -86,6 +85,40 @@ bool AreAlmostTheSame(float a, float b) {
   return std::fabs(a - b) < std::numeric_limits<float>::epsilon();
 }
 
+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;
+  typeST st;
+
+  std::string inputFile("simplex_tree_for_unit_test.txt");
+  std::ifstream simplex_tree_stream(inputFile.c_str());
+  simplex_tree_stream >> st;
+
+  // Display the Simplex_tree
+  std::cout << "The complex contains " << st.num_simplices() << " simplices" << std::endl;
+  std::cout << "   - dimension " << st.dimension() << std::endl;
+
+  // Check
+  BOOST_CHECK(st.num_simplices() == 143353);
+  BOOST_CHECK(st.dimension() == 3);
+
+  int previous_size = 0;
+  for (auto f_simplex : st.filtration_simplex_range()) {
+    // Size of simplex
+    int size = 0;
+    for (auto vertex : st.simplex_vertex_range(f_simplex)) {
+      // Remove warning
+      (void) vertex;
+      size++;
+    }
+    BOOST_CHECK(AreAlmostTheSame(st.filtration(f_simplex), (0.1 * size))); // Specific test: filtration = 0.1 * simplex_size
+    BOOST_CHECK(previous_size <= size); // Check list is sorted (because of sorted filtrations in simplex_tree.txt)
+    previous_size = size;
+  }
+  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;
@@ -112,19 +145,18 @@ void test_simplex_tree_insert_returns_true(const typePairSimplexBool& returnValu
 }
 
 // Global variables
-double max_fil = 0.0;
+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 (fil > max_fil) {
-    max_fil = fil;
-    simplexTree.set_filtration(max_fil);
-    std::cout << "   set_and_test_simplex_tree_dim_fil - max_fil=" << max_fil
+  if (vectorSize > dim_max + 1) {
+    dim_max = vectorSize - 1;
+    simplexTree.set_dimension(dim_max);
+    std::cout << "   set_and_test_simplex_tree_dim_fil - dim_max=" << dim_max
         << std::endl;
   }
 
-  BOOST_CHECK(simplexTree.dimension() >= vectorSize - 1);
-  BOOST_CHECK(AreAlmostTheSame(simplexTree.filtration(), max_fil));
+  BOOST_CHECK(simplexTree.dimension() == dim_max);
 
   // Another way to count simplices:
   size_t num_simp = 0;
@@ -147,7 +179,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_insertion, typeST, list_of_tested_var
   const Filtration_value THIRD_FILTRATION_VALUE = 0.3;
   const Filtration_value FOURTH_FILTRATION_VALUE = 0.4;
   // reset since we run the test several times
-  max_fil = 0.0;
+  dim_max = -1;
 
   // TEST OF INSERTION
   std::cout << "********************************************************************" << std::endl;
@@ -267,7 +299,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_insertion, typeST, list_of_tested_var
   BOOST_CHECK(shReturned == typename typeST::Simplex_handle(nullptr));
   std::cout << "st.num_vertices()=" << st.num_vertices() << std::endl;
   BOOST_CHECK(st.num_vertices() == (size_t) 4); // Not incremented !!
-  BOOST_CHECK(AreAlmostTheSame(st.filtration(), max_fil));
+  BOOST_CHECK(st.dimension() == dim_max);
 
   // ++ ELEVENTH
   std::cout << "   - INSERT (2,1,0) (already inserted)" << std::endl;
@@ -281,9 +313,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_insertion, typeST, list_of_tested_var
   shReturned = returnValue.first;
   BOOST_CHECK(shReturned == typename typeST::Simplex_handle(nullptr));
   BOOST_CHECK(st.num_vertices() == (size_t) 4); // Not incremented !!
-  std::cout << "   - INSERT (2,1,0) (already inserted)" << std::endl;
-  BOOST_CHECK(st.dimension() == 2);
-  BOOST_CHECK(AreAlmostTheSame(st.filtration(), max_fil));
+  BOOST_CHECK(st.dimension() == dim_max);
 
   /* Inserted simplex:        */
   /*    1                     */
@@ -323,7 +353,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_insertion, typeST, list_of_tested_var
 
   // Display the Simplex_tree - Can not be done in the middle of 2 inserts
   std::cout << "The complex contains " << st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << st.dimension() << "   - filtration " << st.filtration() << std::endl;
+  std::cout << "   - dimension " << st.dimension() << std::endl;
   std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
   for (auto f_simplex : st.filtration_simplex_range()) {
     std::cout << "   " << "[" << st.filtration(f_simplex) << "] ";
@@ -533,7 +563,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(NSimplexAndSubfaces_tree_insertion, typeST, list_o
 
   // Display the Simplex_tree - Can not be done in the middle of 2 inserts
   std::cout << "The complex contains " << st.num_simplices() << " simplices" << std::endl;
-  std::cout << "   - dimension " << st.dimension() << "   - filtration " << st.filtration() << std::endl;
+  std::cout << "   - dimension " << st.dimension() << std::endl;
   std::cout << std::endl << std::endl << "Iterator on Simplices in the filtration, with [filtration value]:" << std::endl;
   for (auto f_simplex : st.filtration_simplex_range()) {
     std::cout << "   " << "[" << st.filtration(f_simplex) << "] ";
@@ -708,7 +738,6 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(copy_move_on_simplex_tree, typeST, list_of_tested_
   typeST st_empty;
   // Check st has been emptied by the move
   BOOST_CHECK(st == st_empty);
-  BOOST_CHECK(st.filtration() == 0);
   BOOST_CHECK(st.dimension() == -1);
   BOOST_CHECK(st.num_simplices() == 0);
   BOOST_CHECK(st.num_vertices() == (size_t)0);