Remove Rips_edge_list and review the interfaces

3 files changed, 117 insertions, 292 deletions

diff --git a/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h b/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
index e90d284d..e225f7db 100644
--- a/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
+++ b/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
@@ -12,7 +12,6 @@
 #ifndef FLAG_COMPLEX_SPARSE_MATRIX_H_
 #define FLAG_COMPLEX_SPARSE_MATRIX_H_
 
-#include <gudhi/Rips_edge_list.h>
 #include <gudhi/graph_simplicial_complex.h>
 
 #include <boost/functional/hash.hpp>
@@ -25,14 +24,14 @@
 #endif
 
 #include <iostream>
-#include <utility>
+#include <utility>  // for std::pair
 #include <vector>
-#include <queue>
 #include <unordered_map>
-#include <tuple>
-#include <list>
-#include <algorithm>
-
+#include <unordered_set>
+#include <set>
+#include <tuple>  // for std::tie
+#include <algorithm>  // for std::includes
+#include <iterator>  // for std::inserter
 
 namespace Gudhi {
 
@@ -48,28 +47,29 @@ namespace collapse {
  * A class to store the vertices v/s MaxSimplices Sparse Matrix and to perform collapse operations using the N^2()
  * Algorithm.
  *
- * \tparam Vertex_handle type must be a signed integer type. It admits a total order <.
+ * \tparam Vertex type must be a signed integer type. It admits a total order <.
  * \tparam Filtration type for the value of the filtration function. Must be comparable with <.
  */
-template<typename Vertex_handle, typename Filtration>
+template<typename Vertex, typename Filtration>
 class Flag_complex_sparse_matrix {
+ public:
+  using Vertex_handle = Vertex;
+  using Filtration_value = Filtration;
  private:
   using Edge = std::pair<Vertex_handle, Vertex_handle>;  // This is an ordered pair, An edge is stored with convention of the first
                                          // element being the smaller i.e {2,3} not {3,2}. However this is at the level
                                          // of row indices on actual vertex lables
-  using Filtered_edge = std::pair<Edge, Filtration>;
 
   using Row_index = std::size_t;
 
   using Map_vertex_to_index = std::unordered_map<Vertex_handle, Row_index>;
 
-  using Sparse_row_matrix = Eigen::SparseMatrix<Filtration, Eigen::RowMajor>;
+  using Sparse_row_matrix = Eigen::SparseMatrix<Filtration_value, Eigen::RowMajor>;
 
   using Row_indices_vector = std::vector<Row_index>;
 
  public:
-  using Filtered_sorted_edge_list = std::vector<std::tuple<Filtration, Vertex_handle, Vertex_handle>>;
-  using Filtration_value = Filtration;
+  using Filtered_edge = std::pair<Edge, Filtration_value>;
   using Proximity_graph = Gudhi::Proximity_graph<Flag_complex_sparse_matrix>;
 
  private:
@@ -113,7 +113,7 @@ class Flag_complex_sparse_matrix {
   //! Stores the Sparse matrix of Filtration values representing the Original Simplicial Complex.
   /*!
     \code
-    Sparse_row_matrix   = Eigen::SparseMatrix<Filtration, Eigen::RowMajor> ;
+    Sparse_row_matrix   = Eigen::SparseMatrix<Filtration_value, Eigen::RowMajor> ;
     \endcode
     ;
       */
@@ -206,7 +206,7 @@ class Flag_complex_sparse_matrix {
   }
 
   template<typename FilteredEdgeInsertion>
-  void set_edge_critical(Row_index indx, Filtration filt, FilteredEdgeInsertion filtered_edge_insert) {
+  void set_edge_critical(Row_index indx, Filtration_value filt, FilteredEdgeInsertion filtered_edge_insert) {
 #ifdef DEBUG_TRACES
     std::cout << "The curent index  with filtration value " << indx << ", " << filt << " is primary critical" <<
     std::endl;
@@ -294,7 +294,7 @@ class Flag_complex_sparse_matrix {
     return common;
   }
 
-  void insert_vertex(Vertex_handle vertex, Filtration filt_val) {
+  void insert_vertex(Vertex_handle vertex, Filtration_value filt_val) {
     auto rw = vertex_to_row_.find(vertex);
     if (rw == vertex_to_row_.end()) {
       // Initializing the diagonal element of the adjency matrix corresponding to rw_b.
@@ -306,7 +306,7 @@ class Flag_complex_sparse_matrix {
     }
   }
 
-  void insert_new_edges(Vertex_handle u, Vertex_handle v, Filtration filt_val)
+  void insert_new_edges(Vertex_handle u, Vertex_handle v, Filtration_value filt_val)
   {
     // The edge must not be added before, it should be a new edge.
     insert_vertex(u, filt_val);
@@ -340,23 +340,16 @@ class Flag_complex_sparse_matrix {
     <B>domination_indicator_</B> are initialised by init() function which is
     called at the begining of this. <br>
   */
-  template<typename DistanceMatrix>
-  Flag_complex_sparse_matrix(const DistanceMatrix& distance_matrix)
-  : rows(0) {
-    Vertex_handle num_vertices = std::distance(std::begin(distance_matrix), std::end(distance_matrix));
-
-    // This one is not part of the loop
-    vertices_.emplace(0);
-    // Only browse the lower part of the distance matrix
-    for (Vertex_handle line_index = 1; line_index < num_vertices; line_index++) {
-      for (Vertex_handle row_index = 0; row_index < line_index; row_index++) {
-#ifdef DEBUG_TRACES
-        std::cout << "Insert edge: fn[" << row_index << ", " << line_index << "] = "
-                  << distance_matrix[line_index][row_index] << std::endl;
-#endif  // DEBUG_TRACES
-        f_edge_vector_.push_back({{row_index, line_index}, distance_matrix[line_index][row_index]});
-      }
-      vertices_.emplace(line_index);
+  template<typename Filtered_edge_range>
+  Flag_complex_sparse_matrix(const Filtered_edge_range& filtered_edge_range)
+  : f_edge_vector_(filtered_edge_range.begin(), filtered_edge_range.end()),
+    rows(0) {
+    for (Filtered_edge filtered_edge : filtered_edge_range) {
+      Vertex_handle u;
+      Vertex_handle v;
+      std::tie(u,v) = std::get<0>(filtered_edge);
+      vertices_.emplace(u);
+      vertices_.emplace(v);
     }
   }
 
@@ -385,7 +378,6 @@ class Flag_complex_sparse_matrix {
     u_set_dominated_redges_.clear();
     critical_edge_indicator_.clear();
 
-    std::cout << "Sort it - " << f_edge_vector_.size() << std::endl;
     // Sort edges
     auto sort_by_filtration = [](const Filtered_edge& edge_a, const Filtered_edge& edge_b) -> bool
     {
@@ -397,9 +389,7 @@ class Flag_complex_sparse_matrix {
 #else
     std::stable_sort(f_edge_vector_.begin(), f_edge_vector_.end(), sort_by_filtration);
 #endif
-    std::cout << "Sorted" << std::endl;
 
-    std::cout << vertices_.size() << std::endl;
     // Initializing sparse_row_adjacency_matrix_, This is a row-major sparse matrix.
     sparse_row_adjacency_matrix_ = Sparse_row_matrix(vertices_.size(), vertices_.size());
 
@@ -408,7 +398,7 @@ class Flag_complex_sparse_matrix {
       Edge edge = std::get<0>(fec);
       Vertex_handle u = std::get<0>(edge);
       Vertex_handle v = std::get<1>(edge);
-      Filtration filt = std::get<1>(fec);
+      Filtration_value filt = std::get<1>(fec);
 
       // Inserts the edge in the sparse matrix to update the graph (G_i)
       insert_new_edges(u, v, filt);
diff --git a/src/Collapse/include/gudhi/Rips_edge_list.h b/src/Collapse/include/gudhi/Rips_edge_list.h
deleted file mode 100644
index b7c4dcff..00000000
--- a/src/Collapse/include/gudhi/Rips_edge_list.h
+++ /dev/null
@@ -1,184 +0,0 @@
-/*    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):       Clément Maria, Pawel Dlotko, Vincent Rouvreau Siddharth Pritam
- *
- *    Copyright (C) 2016  INRIA
- *
- *    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/>.
- */
-
-#ifndef RIPS_edge_list_H_
-#define RIPS_edge_list_H_
-
-#include <gudhi/Debug_utils.h>
-#include <gudhi/graph_simplicial_complex.h>
-#include <boost/graph/adjacency_list.hpp>
-
-#include <iostream>
-#include <vector>
-#include <map>
-#include <string>
-#include <limits>  // for numeric_limits
-#include <utility>  // for pair<>
-
-
-namespace Gudhi {
-
-namespace rips_edge_list {
-
-/**
- * \class Rips_complex
- * \brief Rips complex data structure.
- * 
- * \ingroup rips_complex
- * 
- * \details
- * The data structure is a one skeleton graph, or Rips graph, containing edges when the edge length is less or equal
- * to a given threshold. Edge length is computed from a user given point cloud with a given distance function, or a
- * distance matrix.
- * 
- * \tparam Filtration_value is the type used to store the filtration values of the simplicial complex.
- */
-template<typename Filtration_value>
-class Rips_edge_list {
- public:
-  /**
-   * \brief Type of the one skeleton graph stored inside the Rips complex structure.
-   */
-  // typedef typename boost::adjacency_list < boost::vecS, boost::vecS, boost::undirectedS
-  // , boost::property < vertex_filtration_t, Filtration_value >
-  // , boost::property < edge_filtration_t, Filtration_value >> OneSkeletonGraph;
-
- private:
-  typedef int Vertex_handle;
-
- public:
-  /** \brief Rips_complex constructor from a list of points.
-   *
-   * @param[in] points Range of points.
-   * @param[in] threshold Rips value.
-   * @param[in] distance distance function that returns a `Filtration_value` from 2 given points.
-   * 
-   * \tparam ForwardPointRange must be a range for which `std::begin` and `std::end` return input iterators on a
-   * point.
-   *
-   * \tparam Distance furnishes `operator()(const Point& p1, const Point& p2)`, where
-   * `Point` is a point from the `ForwardPointRange`, and that returns a `Filtration_value`.
-   */
-  template<typename ForwardPointRange, typename Distance >
-  Rips_edge_list(const ForwardPointRange& points, Filtration_value threshold, Distance distance) {
-    compute_proximity_graph(points, threshold, distance);
-  }
-
-  /** \brief Rips_complex constructor from a distance matrix.
-   *
-   * @param[in] distance_matrix Range of distances.
-   * @param[in] threshold Rips value.
-   * 
-   * \tparam DistanceMatrix must have a `size()` method and on which `distance_matrix[i][j]` returns
-   * the distance between points \f$i\f$ and \f$j\f$ as long as \f$ 0 \leqslant i < j \leqslant
-   * distance\_matrix.size().\f$
-   */
-  template<typename DistanceMatrix>
-  Rips_edge_list(const DistanceMatrix& distance_matrix, Filtration_value threshold) {
-    compute_proximity_graph(boost::irange((size_t)0, distance_matrix.size()), threshold,
-                            [&](size_t i, size_t j){return distance_matrix[j][i];});
-  }
-
-  /** \brief Initializes the egde list (one skeleton) complex from the Rips graph 
-   *
-   * \tparam EdgeListForRips must meet `EdgeListForRips` concept.
-   * 
-   * @param[in] edges EdgeListForRips to be created.
-   * @param[in] dim_max graph expansion for Rips until this given maximal dimension.
-   * @exception std::invalid_argument In debug mode, if `edges.num_vertices()` does not return 0.
-   *
-   */
-  template <typename EdgeListForRips>
-  void create_edges(EdgeListForRips& edge_list) {
-    GUDHI_CHECK(edges.num_vertices() == 0,
-                std::invalid_argument("Rips_complex::create_complex - edge list is not empty"));
-
-    // sort the tuple (filteration_valuem, (v1,v2){edge}) 
-    //By default the sort is done on the first element, so here it's filteration value.
-    std::sort(edges.begin(),edges.end());
-    for(size_t i = 0; i < edges.size(); i++ )
-      edge_list.emplace_back(edges.at(i));
-
-  }
-
- private:
-  /** \brief Computes the proximity graph of the points.
-   *
-   * If points contains n elements, the proximity graph is the graph with n vertices, and an edge [u,v] iff the
-   * distance function between points u and v is smaller than threshold.
-   *
-   * \tparam ForwardPointRange furnishes `.begin()` and `.end()`
-   * methods.
-   *
-   * \tparam Distance furnishes `operator()(const Point& p1, const Point& p2)`, where
-   * `Point` is a point from the `ForwardPointRange`, and that returns a `Filtration_value`.
-   */
-  template< typename ForwardPointRange, typename Distance >
-  void compute_proximity_graph(const ForwardPointRange& points, Filtration_value threshold,
-               Distance distance) {
-    edges.clear();
-    // Compute the proximity graph of the points.
-    // If points contains n elements, the proximity graph is the graph with n vertices, and an edge [u,v] iff the
-    // distance function between points u and v is smaller than threshold.
-    // --------------------------------------------------------------------------------------------
-    // Creates the vector of edges and its filtration values (returned by distance function)
-    Vertex_handle idx_u = 0;
-    for (auto it_u = std::begin(points); it_u != std::end(points); ++it_u, ++idx_u) {
-      Vertex_handle idx_v = idx_u + 1;
-      for (auto it_v = it_u + 1; it_v != std::end(points); ++it_v, ++idx_v) {
-        Filtration_value fil = distance(*it_u, *it_v);
-        if (fil <= threshold) {
-          edges.emplace_back(fil, idx_u, idx_v);
-        }
-      }
-    }
-
-    // --------------------------------------------------------------------------------------------
-    // Creates the proximity graph from edges and sets the property with the filtration value.
-    // Number of points is labeled from 0 to idx_u-1
-    // --------------------------------------------------------------------------------------------
-    // Do not use : rips_skeleton_graph_ = OneSkeletonGraph(...) -> deep copy of the graph (boost graph is not
-    // move-enabled)
-    // rips_skeleton_graph_.~OneSkeletonGraph();
-    // new(&rips_skeleton_graph_)OneSkeletonGraph(edges.begin(), edges.end(), edges_fil.begin(), idx_u);
-
-    // auto vertex_prop = boost::get(vertex_filtration_t(), rips_skeleton_graph_);
-
-    // using vertex_iterator = typename boost::graph_traits<OneSkeletonGraph>::vertex_iterator;
-    // vertex_iterator vi, vi_end;
-    // for (std::tie(vi, vi_end) = boost::vertices(rips_skeleton_graph_);
-    //      vi != vi_end; ++vi) {
-    //   boost::put(vertex_prop, *vi, 0.);
-    // }
-  }
-
- private:
-  // OneSkeletonGraph rips_skeleton_graph_;
-  std::vector< std::tuple < Filtration_value, Vertex_handle, Vertex_handle > > edges;
-  // std::vector< Filtration_value > edges_fil;
-};
-
-}  // namespace rips_complex
-
-}  // namespace Gudhi
-
-#endif  // RIPS_COMPLEX_H_
diff --git a/src/Collapse/test/collapse_unit_test.cpp b/src/Collapse/test/collapse_unit_test.cpp
index 3a07e088..1bec3810 100644
--- a/src/Collapse/test/collapse_unit_test.cpp
+++ b/src/Collapse/test/collapse_unit_test.cpp
@@ -8,67 +8,77 @@
  *      - YYYY/MM Author: Description of the modification
  */
 
-#include <iostream>
-#include <tuple>
-#include <vector>
-#include <array>
 
 #define BOOST_TEST_DYN_LINK
 #define BOOST_TEST_MODULE "collapse"
 #include <boost/test/unit_test.hpp>
 #include <boost/mpl/list.hpp>
 
-#include "gudhi/Flag_complex_sparse_matrix.h"
+#include <gudhi/Flag_complex_sparse_matrix.h>
+#include <gudhi/distance_functions.h>
+#include <gudhi/graph_simplicial_complex.h>
+
+#include <iostream>
+#include <tuple>
+#include <vector>
+#include <array>
+#include <cmath>
 
 using Filtration_value = float;
 using Vertex_handle = short;
-using Filtered_edge = std::tuple<Filtration_value, Vertex_handle, Vertex_handle>;
-using Filtered_sorted_edge_list = std::vector<Filtered_edge>;
 using Flag_complex_sparse_matrix = Gudhi::collapse::Flag_complex_sparse_matrix<Vertex_handle, Filtration_value>;
+using Filtered_edge = Flag_complex_sparse_matrix::Filtered_edge;
+using Filtered_edge_list = std::vector<Filtered_edge>;
 
-bool find_edge_in_list(const Filtered_edge& edge, const Filtered_sorted_edge_list& edge_list) {
+template<typename Filtered_edge_range>
+bool find_edge_in_list(const Filtered_edge& edge, const Filtered_edge_range& edge_list) {
   for (auto edge_from_list : edge_list) {
     if (edge_from_list == edge)
       return true;
   }
   return false;
 }
-/*
-void trace_and_check_collapse(const Filtered_sorted_edge_list& edges, const Filtered_sorted_edge_list& removed_edges) {
-  std::cout << "BEFORE COLLAPSE - Total number of edges: " << edges.size() << std::endl;
-  BOOST_CHECK(edges.size() > 0);
-  for (auto edge : edges) {
-    std::cout << "f[" << std::get<1>(edge) << ", " << std::get<2>(edge) << "] = " << std::get<0>(edge) << std::endl;
+
+template<typename Filtered_edge_range>
+void trace_and_check_collapse(const Filtered_edge_range& filtered_edges, const Filtered_edge_list& removed_edges) {
+  std::cout << "BEFORE COLLAPSE - Total number of edges: " << filtered_edges.size() << std::endl;
+  BOOST_CHECK(filtered_edges.size() > 0);
+  for (auto filtered_edge : filtered_edges) {
+    auto edge = std::get<0>(filtered_edge);
+    std::cout << "f[" << std::get<0>(edge) << ", " << std::get<1>(edge) << "] = " << std::get<1>(filtered_edge) << std::endl;
   }
 
   std::cout << "COLLAPSE - keep edges: " << std::endl;
-  Flag_complex_sparse_matrix flag_complex_sparse_matrix(edges);
-  Filtered_sorted_edge_list collapse_edges;
+  Flag_complex_sparse_matrix flag_complex_sparse_matrix(filtered_edges);
+  Filtered_edge_list collapse_edges;
   flag_complex_sparse_matrix.filtered_edge_collapse(
     [&collapse_edges](std::pair<Vertex_handle, Vertex_handle> edge, Filtration_value filtration) {
       std::cout << "f[" << std::get<0>(edge) << ", " << std::get<1>(edge) << "] = " << filtration << std::endl;
-        collapse_edges.push_back({filtration, std::get<0>(edge), std::get<1>(edge)});
+        collapse_edges.push_back({edge, filtration});
       });
   std::cout << "AFTER COLLAPSE - Total number of edges: " << collapse_edges.size() << std::endl;
-  BOOST_CHECK(collapse_edges.size() <= edges.size());
-  for (auto edge_from_collapse : collapse_edges) {
-    std::cout << "f[" << std::get<1>(edge_from_collapse) << ", " << std::get<2>(edge_from_collapse) << "] = "
-              << std::get<0>(edge_from_collapse) << std::endl;
+  BOOST_CHECK(collapse_edges.size() <= filtered_edges.size());
+  for (auto filtered_edge_from_collapse : collapse_edges) {
+    auto edge_from_collapse = std::get<0>(filtered_edge_from_collapse);
+    std::cout << "f[" << std::get<0>(edge_from_collapse) << ", " << std::get<1>(edge_from_collapse) << "] = "
+              << std::get<1>(filtered_edge_from_collapse) << std::endl;
     // Check each edge from collapse is in the input
-    BOOST_CHECK(find_edge_in_list(edge_from_collapse, edges));
+    BOOST_CHECK(find_edge_in_list(filtered_edge_from_collapse, filtered_edges));
   }
 
   std::cout << "CHECK COLLAPSE - Total number of removed edges: " << removed_edges.size() << std::endl;
-  for (auto removed_edge : removed_edges) {
-    std::cout << "f[" << std::get<1>(removed_edge) << ", " << std::get<2>(removed_edge) << "] = "
-              << std::get<0>(removed_edge) << std::endl;
+  for (auto removed_filtered_edge : removed_edges) {
+    auto removed_edge = std::get<0>(removed_filtered_edge);
+    std::cout << "f[" << std::get<0>(removed_edge) << ", " << std::get<1>(removed_edge) << "] = "
+              << std::get<1>(removed_filtered_edge) << std::endl;
     // Check each removed edge from collapse is in the input
-    BOOST_CHECK(!find_edge_in_list(removed_edge, collapse_edges));
+    BOOST_CHECK(!find_edge_in_list(removed_filtered_edge, collapse_edges));
   }
 
 }
 
 BOOST_AUTO_TEST_CASE(collapse) {
+  std::cout << "***** COLLAPSE *****" << std::endl;
   //  1   2
   //  o---o
   //  |   |
@@ -76,7 +86,7 @@ BOOST_AUTO_TEST_CASE(collapse) {
   //  |   |
   //  o---o
   //  0   3
-  Filtered_sorted_edge_list edges {{1., 0, 1}, {1., 1, 2}, {1., 2, 3}, {1., 3, 0}};
+  Filtered_edge_list edges {{{0, 1}, 1.}, {{1, 2}, 1.}, {{2, 3}, 1.}, {{3, 0}, 1.}};
   trace_and_check_collapse(edges, {});
   
   //  1   2
@@ -86,9 +96,9 @@ BOOST_AUTO_TEST_CASE(collapse) {
   //  |/ \|
   //  o---o
   //  0   3
-  edges.push_back({2., 0, 2});
-  edges.push_back({2., 1, 3});
-  trace_and_check_collapse(edges, {{2., 1, 3}});
+  edges.push_back({{0, 2}, 2.});
+  edges.push_back({{1, 3}, 2.});
+  trace_and_check_collapse(edges, {{{1, 3}, 2.}});
 
   //  1   2   4
   //  o---o---o
@@ -97,10 +107,10 @@ BOOST_AUTO_TEST_CASE(collapse) {
   //  |/ \|   |
   //  o---o---o
   //  0   3   5
-  edges.push_back({3., 2, 4});
-  edges.push_back({3., 4, 5});
-  edges.push_back({3., 5, 3});
-  trace_and_check_collapse(edges, {{2., 1, 3}});
+  edges.push_back({{2, 4}, 3.});
+  edges.push_back({{4, 5}, 3.});
+  edges.push_back({{5, 3}, 3.});
+  trace_and_check_collapse(edges, {{{1, 3}, 2.}});
 
   //  1   2   4
   //  o---o---o
@@ -109,9 +119,9 @@ BOOST_AUTO_TEST_CASE(collapse) {
   //  |/ \|/ \|
   //  o---o---o
   //  0   3   5
-  edges.push_back({4., 2, 5});
-  edges.push_back({4., 4, 3});
-  trace_and_check_collapse(edges, {{2., 1, 3}, {4., 4, 3}});
+  edges.push_back({{2, 5}, 4.});
+  edges.push_back({{4, 3}, 4.});
+  trace_and_check_collapse(edges, {{{1, 3}, 2.}, {{4, 3}, 4.}});
 
   //  1   2   4
   //  o---o---o
@@ -120,13 +130,27 @@ BOOST_AUTO_TEST_CASE(collapse) {
   //  |/ \|/ \|
   //  o---o---o
   //  0   3   5
-  edges.push_back({5., 1, 5});
-  edges.push_back({5., 0, 4});
-  trace_and_check_collapse(edges, {{2., 1, 3}, {4., 4, 3}, {5., 0, 4}});
-}*/
+  edges.push_back({{1, 5}, 5.});
+  edges.push_back({{0, 4}, 5.});
+  trace_and_check_collapse(edges, {{{1, 3}, 2.}, {{4, 3}, 4.}, {{0, 4}, 5.}});
+}
+
+BOOST_AUTO_TEST_CASE(collapse_from_array) {
+  std::cout << "***** COLLAPSE FROM ARRAY *****" << std::endl;
+  //  1   2
+  //  o---o
+  //  |\ /|
+  //  | x |
+  //  |/ \|
+  //  o---o
+  //  0   3
+  std::array<Filtered_edge, 6> f_edge_array = {{{{0, 1}, 1.}, {{1, 2}, 1.}, {{2, 3}, 1.}, {{3, 0}, 1.}, {{0, 2}, 2.}, {{1, 3}, 2.}}};
+  trace_and_check_collapse(f_edge_array, {{{1, 3}, 2.}});
+}
 
 
-BOOST_AUTO_TEST_CASE(collapse_from_distance_matrix) {
+BOOST_AUTO_TEST_CASE(collapse_from_proximity_graph) {
+  std::cout << "***** COLLAPSE FROM PROXIMITY GRAPH *****" << std::endl;
   //  1   2
   //  o---o
   //  |\ /|
@@ -134,39 +158,34 @@ BOOST_AUTO_TEST_CASE(collapse_from_distance_matrix) {
   //  |/ \|
   //  o---o
   //  0   3
-  // Lower diagonal distance matrix
-  std::array<std::array<double, 4>, 4> distance_matrix = {{{0., 0., 0., 0.},
-                                                           {1., 0., 0., 0.},
-                                                           {2., 1., 0., 0.},
-                                                           {1., 2., 1., 0.} }};
-  
-  std::cout << "COLLAPSE - keep edges: " << std::endl;
-  Flag_complex_sparse_matrix flag_complex_sparse_matrix(distance_matrix);
-  Filtered_sorted_edge_list collapse_edges;
+  std::vector<std::vector<Filtration_value>> point_cloud = {{0., 0.},
+                                                 {0., 1.},
+                                                 {1., 0.},
+                                                 {1., 1.} };
+
+  Filtration_value threshold = std::numeric_limits<Filtration_value>::infinity();
+  using Proximity_graph = Flag_complex_sparse_matrix::Proximity_graph;
+  Proximity_graph proximity_graph = Gudhi::compute_proximity_graph<Flag_complex_sparse_matrix>(point_cloud,
+                                                                                               threshold,
+                                                                                               Gudhi::Euclidean_distance());
+  Flag_complex_sparse_matrix flag_complex_sparse_matrix(proximity_graph);
+  Filtered_edge_list collapse_edges;
   flag_complex_sparse_matrix.filtered_edge_collapse(
     [&collapse_edges](std::pair<Vertex_handle, Vertex_handle> edge, Filtration_value filtration) {
       std::cout << "f[" << std::get<0>(edge) << ", " << std::get<1>(edge) << "] = " << filtration << std::endl;
-        collapse_edges.push_back({filtration, std::get<0>(edge), std::get<1>(edge)});
+        collapse_edges.push_back({edge, filtration});
       });
-  std::cout << "AFTER COLLAPSE - Total number of edges: " << collapse_edges.size() << std::endl;
   BOOST_CHECK(collapse_edges.size() == 5);
 
-  Filtered_sorted_edge_list edges {{1., 0, 1}, {1., 1, 2}, {1., 2, 3}, {1., 0, 3}, {2., 0, 2}, {2., 1, 3}};
-  
-  for (auto edge_from_collapse : collapse_edges) {
-    std::cout << "f[" << std::get<1>(edge_from_collapse) << ", " << std::get<2>(edge_from_collapse) << "] = "
-              << std::get<0>(edge_from_collapse) << std::endl;
-    // Check each edge from collapse is in the input
-    BOOST_CHECK(find_edge_in_list(edge_from_collapse, edges));
-  }
-
-  Filtered_sorted_edge_list removed_edges {{2., 1, 3}};
-  
-  std::cout << "CHECK COLLAPSE - Total number of removed edges: " << removed_edges.size() << std::endl;
-  for (auto removed_edge : removed_edges) {
-    std::cout << "f[" << std::get<1>(removed_edge) << ", " << std::get<2>(removed_edge) << "] = "
-              << std::get<0>(removed_edge) << std::endl;
-    // Check each removed edge from collapse is in the input
-    BOOST_CHECK(!find_edge_in_list(removed_edge, collapse_edges));
+  std::size_t filtration_is_edge_length_nb = 0;
+  std::size_t filtration_is_diagonal_length_nb = 0;
+  float epsilon = std::numeric_limits<Filtration_value>::epsilon();
+  for (auto filtered_edge : collapse_edges) {
+    if (std::get<1>(filtered_edge) == 1.)
+      filtration_is_edge_length_nb++;
+    if (std::fabs(std::get<1>(filtered_edge) - std::sqrt(2.)) <= epsilon)
+      filtration_is_diagonal_length_nb++;
   }
+  BOOST_CHECK(filtration_is_edge_length_nb == 4);
+  BOOST_CHECK(filtration_is_diagonal_length_nb == 1);
 }
\ No newline at end of file