summaryrefslogtreecommitdiff
path: root/src/Collapse
diff options
context:
space:
mode:
authorROUVREAU Vincent <vincent.rouvreau@inria.fr>2020-04-09 21:46:42 +0200
committerROUVREAU Vincent <vincent.rouvreau@inria.fr>2020-04-09 21:46:42 +0200
commit9654c177078fc598c8a8424dd67d0742bf0defb9 (patch)
tree8e0d8f5ce711f51511f5ca4313ff4e24018a357f /src/Collapse
parent599e910811e1c4c743a61be65e089e798f578d4a (diff)
Use an output iterator for edge collapse return instead of storing it
Diffstat (limited to 'src/Collapse')
-rw-r--r--src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h20
-rw-r--r--src/Collapse/test/collapse_unit_test.cpp6
-rw-r--r--src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp46
-rw-r--r--src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp59
4 files changed, 37 insertions, 94 deletions
diff --git a/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h b/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
index d7014f2f..033c27a3 100644
--- a/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
+++ b/src/Collapse/include/gudhi/Flag_complex_sparse_matrix.h
@@ -55,7 +55,6 @@ using boolVector = std::vector<bool>;
using EdgeFiltVector = std::vector<EdgeFilt>;
typedef std::vector<std::tuple<double, Vertex, Vertex>> Filtered_sorted_edge_list;
-typedef std::unordered_map<Edge, bool, boost::hash<Edge>> u_edge_map;
typedef std::unordered_map<Edge, std::size_t, boost::hash<Edge>> u_edge_to_idx_map;
//! Class SparseMsMatrix
@@ -126,8 +125,6 @@ class Flag_complex_sparse_matrix {
// Vector of filtered edges, for edge-collapse, the indices of the edges are the row-indices.
EdgeFiltVector f_edge_vector;
- // List of non-dominated edges, the indices of the edges are the vertex lables!!.
- Filtered_sorted_edge_list critical_core_edges;
// Stores the indices from the sorted filtered edge vector.
// std::set<std::size_t> recurCriticalCoreIndcs;
@@ -214,7 +211,8 @@ class Flag_complex_sparse_matrix {
return edge_indices;
}
- void set_edge_critical(std::size_t indx, double filt) {
+ template<typename FilteredEdgeInsertion>
+ void set_edge_critical(std::size_t indx, double filt, FilteredEdgeInsertion filtered_edge_insert) {
#ifdef DEBUG_TRACES
std::cout << "The curent index with filtration value " << indx << ", " << filt << " is primary critical" <<
std::endl;
@@ -235,7 +233,7 @@ class Flag_complex_sparse_matrix {
std::cout << "The curent index became critical " << idx << std::endl;
#endif // DEBUG_TRACES
critical_edge_indicator.at(idx) = true;
- critical_core_edges.push_back({filt, u, v});
+ filtered_edge_insert({u, v}, filt);
std::set<std::size_t> inner_effected_indcs = three_clique_indices(idx);
for (auto inr_idx = inner_effected_indcs.rbegin(); inr_idx != inner_effected_indcs.rend(); inr_idx++) {
if (*inr_idx < idx) effectedIndcs.emplace(*inr_idx);
@@ -354,7 +352,8 @@ class Flag_complex_sparse_matrix {
}
// Performs edge collapse in a decreasing sequence of the filtration value.
- Filtered_sorted_edge_list filtered_edge_collapse() {
+ template<typename FilteredEdgeInsertion>
+ void filtered_edge_collapse(FilteredEdgeInsertion filtered_edge_insert) {
std::size_t endIdx = 0;
u_set_removed_redges.clear();
@@ -381,20 +380,15 @@ class Flag_complex_sparse_matrix {
if (not check_edge_domination(e)) {
critical_edge_indicator.at(endIdx) = true;
dominated_edge_indicator.at(endIdx) = false;
- critical_core_edges.push_back({filt, u, v});
+ filtered_edge_insert({u, v}, filt);
if (endIdx > 1)
- set_edge_critical(endIdx, filt);
+ set_edge_critical(endIdx, filt, filtered_edge_insert);
} else
dominated_edge_indicator.at(endIdx) = true;
endIdx++;
}
-#ifdef DEBUG_TRACES
- std::cout << "The total number of critical edges is: " << critical_core_edges.size() << std::endl;
- std::cout << "The total number of non-critical edges is: " << f_edge_vector.size() - critical_core_edges.size() << std::endl;
-#endif // DEBUG_TRACES
edgeCollapsed = true;
- return critical_core_edges;
}
void insert_vertex(const Vertex& vertex, double filt_val) {
diff --git a/src/Collapse/test/collapse_unit_test.cpp b/src/Collapse/test/collapse_unit_test.cpp
index 4857580c..29f33219 100644
--- a/src/Collapse/test/collapse_unit_test.cpp
+++ b/src/Collapse/test/collapse_unit_test.cpp
@@ -49,7 +49,11 @@ void trace_and_check_collapse(const Filtered_sorted_edge_list& edges, const Filt
}
Flag_complex_sparse_matrix flag_complex_sparse_matrix(edges);
- auto collapse_edges = flag_complex_sparse_matrix.filtered_edge_collapse();
+ Filtered_sorted_edge_list collapse_edges;
+ flag_complex_sparse_matrix.filtered_edge_collapse(
+ [&collapse_edges](std::pair<std::size_t, std::size_t> edge, double filtration) {
+ collapse_edges.push_back({std::get<0>(edge), std::get<1>(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) {
diff --git a/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp b/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp
index 7f5a9454..98a90892 100644
--- a/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp
+++ b/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp
@@ -1,5 +1,4 @@
#include <gudhi/Flag_complex_sparse_matrix.h>
-#include <gudhi/Rips_complex.h>
#include <gudhi/Simplex_tree.h>
#include <gudhi/Persistent_cohomology.h>
#include <gudhi/Rips_edge_list.h>
@@ -17,7 +16,6 @@ using Vector_of_points = std::vector<Point>;
using Simplex_tree = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>;
using Filtration_value = double;
-using Rips_complex = Gudhi::rips_complex::Rips_complex<Filtration_value>;
using Rips_edge_list = Gudhi::rips_edge_list::Rips_edge_list<Filtration_value>;
using Field_Zp = Gudhi::persistent_cohomology::Field_Zp;
using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<Simplex_tree, Field_Zp>;
@@ -62,32 +60,6 @@ void program_options(int argc, char* const argv[], double& min_persistence, doub
}
}
-class filt_edge_to_dist_matrix {
- public:
- template <class Distance_matrix, class Filtered_sorted_edge_list>
- filt_edge_to_dist_matrix(Distance_matrix& distance_mat, Filtered_sorted_edge_list& edge_filt,
- std::size_t number_of_points) {
- double inf = std::numeric_limits<double>::max();
- doubleVector distances;
- std::pair<std::size_t, std::size_t> e;
- for (std::size_t indx = 0; indx < number_of_points; indx++) {
- for (std::size_t j = 0; j <= indx; j++) {
- if (j == indx)
- distances.push_back(0);
- else
- distances.push_back(inf);
- }
- distance_mat.push_back(distances);
- distances.clear();
- }
-
- for (auto edIt = edge_filt.begin(); edIt != edge_filt.end(); edIt++) {
- e = std::minmax(std::get<1>(*edIt), std::get<2>(*edIt));
- distance_mat.at(std::get<1>(e)).at(std::get<0>(e)) = std::get<0>(*edIt);
- }
- }
-};
-
void program_options(int argc, char* argv[], std::string& csv_matrix_file, std::string& filediag,
Filtration_value& threshold, int& dim_max, int& p, Filtration_value& min_persistence);
@@ -133,16 +105,18 @@ int main(int argc, char* argv[]) {
std::cout << "Filtered edge collapse begins" << std::endl;
Flag_complex_sparse_matrix mat_filt_edge_coll(edge_t);
std::cout << "Matrix instansiated" << std::endl;
- Filtered_sorted_edge_list collapse_edges;
- collapse_edges = mat_filt_edge_coll.filtered_edge_collapse();
- filt_edge_to_dist_matrix(sparse_distances, collapse_edges, number_of_points);
- std::cout << "Total number of vertices after collapse in the sparse matrix are: " << mat_filt_edge_coll.num_vertices()
- << std::endl;
-
- Rips_complex rips_complex_after_collapse(sparse_distances, threshold);
Simplex_tree stree;
- rips_complex_after_collapse.create_complex(stree, dim_max);
+ mat_filt_edge_coll.filtered_edge_collapse(
+ [&stree](std::vector<std::size_t> edge, double filtration) {
+ // insert the 2 vertices with a 0. filtration value just like a Rips
+ stree.insert_simplex({edge[0]}, 0.);
+ stree.insert_simplex({edge[1]}, 0.);
+ // insert the edge
+ stree.insert_simplex(edge, filtration);
+ });
+
+ stree.expansion(dim_max);
std::cout << "The complex contains " << stree.num_simplices() << " simplices after collapse. \n";
std::cout << " and has dimension " << stree.dimension() << " \n";
diff --git a/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp b/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp
index a2840674..70d8d9c5 100644
--- a/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp
+++ b/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp
@@ -1,8 +1,6 @@
#include <gudhi/Flag_complex_sparse_matrix.h>
-#include <gudhi/Rips_complex.h>
#include <gudhi/Simplex_tree.h>
#include <gudhi/Persistent_cohomology.h>
-#include <gudhi/Rips_edge_list.h>
#include <gudhi/distance_functions.h>
#include <gudhi/reader_utils.h>
#include <gudhi/Points_off_io.h>
@@ -11,16 +9,18 @@
#include <boost/graph/adjacency_list.hpp>
#include <boost/program_options.hpp>
+#include<utility> // for std::pair
+#include<vector>
+
// Types definition
using Simplex_tree = Gudhi::Simplex_tree<>;
using Filtration_value = Simplex_tree::Filtration_value;
+using Vertex_handle = std::size_t; /*Simplex_tree::Vertex_handle;*/
using Point = std::vector<Filtration_value>;
using Vector_of_points = std::vector<Point>;
-using Rips_complex = Gudhi::rips_complex::Rips_complex<Filtration_value>;
-using Rips_edge_list = Gudhi::rips_edge_list::Rips_edge_list<Filtration_value>;
using Field_Zp = Gudhi::persistent_cohomology::Field_Zp;
using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<Simplex_tree, Field_Zp>;
using Distance_matrix = std::vector<std::vector<Filtration_value>>;
@@ -29,38 +29,10 @@ using Adjacency_list = boost::adjacency_list<boost::vecS, boost::vecS, boost::di
boost::property<Gudhi::vertex_filtration_t, double>,
boost::property<Gudhi::edge_filtration_t, double>>;
-
-class filt_edge_to_dist_matrix {
- public:
- template <class Distance_matrix, class Filtered_sorted_edge_list>
- filt_edge_to_dist_matrix(Distance_matrix& distance_mat, Filtered_sorted_edge_list& edge_filt,
- std::size_t number_of_points) {
- double inf = std::numeric_limits<double>::max();
- doubleVector distances;
- std::pair<std::size_t, std::size_t> e;
- for (std::size_t indx = 0; indx < number_of_points; indx++) {
- for (std::size_t j = 0; j <= indx; j++) {
- if (j == indx)
- distances.push_back(0);
- else
- distances.push_back(inf);
- }
- distance_mat.push_back(distances);
- distances.clear();
- }
-
- for (auto edIt = edge_filt.begin(); edIt != edge_filt.end(); edIt++) {
- e = std::minmax(std::get<1>(*edIt), std::get<2>(*edIt));
- distance_mat.at(std::get<1>(e)).at(std::get<0>(e)) = std::get<0>(*edIt);
- }
- }
-};
-
void program_options(int argc, char* argv[], std::string& off_file_points, std::string& filediag,
Filtration_value& threshold, int& dim_max, int& p, Filtration_value& min_persistence);
int main(int argc, char* argv[]) {
- typedef size_t Vertex_handle;
typedef std::vector<std::tuple<Filtration_value, Vertex_handle, Vertex_handle>> Filtered_sorted_edge_list;
auto the_begin = std::chrono::high_resolution_clock::now();
@@ -82,8 +54,6 @@ int main(int argc, char* argv[]) {
std::cout << min_persistence << ", " << threshold << ", " << dim_max
<< ", " << off_file_points << ", " << filediag << std::endl;
- Map map_empty;
-
Distance_matrix sparse_distances;
Gudhi::Points_off_reader<Point> off_reader(off_file_points);
@@ -120,18 +90,19 @@ int main(int argc, char* argv[]) {
std::cout << "Computing the one-skeleton for threshold: " << threshold << std::endl;
std::cout << "Matrix instansiated" << std::endl;
- Filtered_sorted_edge_list collapse_edges;
- collapse_edges = mat_filt_edge_coll.filtered_edge_collapse();
- filt_edge_to_dist_matrix(sparse_distances, collapse_edges, number_of_points);
- std::cout << "Total number of vertices after collapse in the sparse matrix are: " << mat_filt_edge_coll.num_vertices()
- << std::endl;
-
- // Rips_complex rips_complex_before_collapse(distances, threshold);
- Rips_complex rips_complex_after_collapse(sparse_distances, threshold);
Simplex_tree stree;
- rips_complex_after_collapse.create_complex(stree, dim_max);
-
+ mat_filt_edge_coll.filtered_edge_collapse(
+ [&stree](std::vector<std::size_t> edge, double filtration) {
+ // insert the 2 vertices with a 0. filtration value just like a Rips
+ stree.insert_simplex({edge[0]}, 0.);
+ stree.insert_simplex({edge[1]}, 0.);
+ // insert the edge
+ stree.insert_simplex(edge, filtration);
+ });
+
+ stree.expansion(dim_max);
+
std::cout << "The complex contains " << stree.num_simplices() << " simplices after collapse. \n";
std::cout << " and has dimension " << stree.dimension() << " \n";