From 0a1d09caba097117e08c16d12bf682d611454c02 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Fri, 21 Oct 2022 16:54:09 +0200 Subject: Stores vertices on init_from_range and insert all of them in the Simplex_tree to avoid quadratic behavior --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 23 ++++++++++++++++++++--- 1 file changed, 20 insertions(+), 3 deletions(-) (limited to 'src/Alpha_complex/include/gudhi/Alpha_complex.h') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index aec8c1b1..d5d1a0f0 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -17,8 +17,7 @@ // to construct Alpha_complex from a OFF file of points #include -#include -#include // isnan, fmax +#include // isnan, fmax #include // for std::unique_ptr #include // for std::size_t @@ -45,6 +44,8 @@ #include // std::pair #include #include // for std::iota +#include // for std::sort +#include // for complex.insert_simplex_and_subfaces({...}, ...) // Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10 #if CGAL_VERSION_NR < 1041101000 @@ -146,6 +147,10 @@ class Alpha_complex { std::unique_ptr triangulation_; /** \brief Kernel for triangulation_ functions access.*/ A_kernel_d kernel_; + /** \brief Vertices to be inserted first by the create_complex method to avoid quadratic complexity. + * It isn't just [0, n) if some points have multiplicity (only one copy appears in the complex). + */ + std::vector vertices_; /** \brief Cache for geometric constructions: circumcenter and squared radius of a simplex.*/ std::vector cache_, old_cache_; @@ -279,6 +284,9 @@ class Alpha_complex { // structure to retrieve CGAL points from vertex handle - one vertex handle per point. // Needs to be constructed before as vertex handles arrives in no particular order. vertex_handle_to_iterator_.resize(point_cloud.size()); + // List of sorted unique vertices in the triangulation. We take advantage of the existing loop to construct it + // Vertices list avoids quadratic complexity with the Simplex_tree. We should not fill it up with Toplex_map e.g. + vertices_.reserve(triangulation_->number_of_vertices()); // Loop on triangulation vertices list for (CGAL_vertex_iterator vit = triangulation_->vertices_begin(); vit != triangulation_->vertices_end(); ++vit) { if (!triangulation_->is_infinite(*vit)) { @@ -286,8 +294,10 @@ class Alpha_complex { std::clog << "Vertex insertion - " << vit->data() << " -> " << vit->point() << std::endl; #endif // DEBUG_TRACES vertex_handle_to_iterator_[vit->data()] = vit; + vertices_.push_back(vit->data()); } } + std::sort(vertices_.begin(), vertices_.end()); // -------------------------------------------------------------------------------------------- } } @@ -384,12 +394,19 @@ class Alpha_complex { // -------------------------------------------------------------------------------------------- // Simplex_tree construction from loop on triangulation finite full cells list if (num_vertices() > 0) { + for (auto vertex : vertices_) { +#ifdef DEBUG_TRACES + std::clog << "SimplicialComplex insertion " << vertex << std::endl; +#endif // DEBUG_TRACES + complex.insert_simplex_and_subfaces({static_cast(vertex)}, std::numeric_limits::quiet_NaN()); + } + for (auto cit = triangulation_->finite_full_cells_begin(); cit != triangulation_->finite_full_cells_end(); ++cit) { Vector_vertex vertexVector; #ifdef DEBUG_TRACES - std::clog << "Simplex_tree insertion "; + std::clog << "SimplicialComplex insertion "; #endif // DEBUG_TRACES for (auto vit = cit->vertices_begin(); vit != cit->vertices_end(); ++vit) { if (*vit != nullptr) { -- cgit v1.2.3 From f4e90e943dd368b76f466c91075193124c191ef5 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Fri, 21 Oct 2022 17:16:01 +0200 Subject: Remove strange behavior with initializer list --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) (limited to 'src/Alpha_complex/include/gudhi/Alpha_complex.h') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index d5d1a0f0..7b837ae0 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -45,7 +45,6 @@ #include #include // for std::iota #include // for std::sort -#include // for complex.insert_simplex_and_subfaces({...}, ...) // Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10 #if CGAL_VERSION_NR < 1041101000 @@ -394,11 +393,13 @@ class Alpha_complex { // -------------------------------------------------------------------------------------------- // Simplex_tree construction from loop on triangulation finite full cells list if (num_vertices() > 0) { + std::vector one_vertex(1); for (auto vertex : vertices_) { #ifdef DEBUG_TRACES std::clog << "SimplicialComplex insertion " << vertex << std::endl; #endif // DEBUG_TRACES - complex.insert_simplex_and_subfaces({static_cast(vertex)}, std::numeric_limits::quiet_NaN()); + one_vertex[0] = vertex; + complex.insert_simplex_and_subfaces(one_vertex, std::numeric_limits::quiet_NaN()); } for (auto cit = triangulation_->finite_full_cells_begin(); -- cgit v1.2.3 From 11dec62e555e446fb34a1f5ba3d0c9aba1a7e956 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Thu, 3 Nov 2022 17:49:33 +0100 Subject: code review: Use std::ptrdiff_t as internal vertex handle type when not yet known --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) (limited to 'src/Alpha_complex/include/gudhi/Alpha_complex.h') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index 7b837ae0..a7372f19 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -101,13 +101,17 @@ template struct Is_Epeck_D> { static const bool val */ template, bool Weighted = false> class Alpha_complex { + private: + // Vertex_handle internal type (required by triangulation_ and vertices_). + using Internal_vertex_handle = std::ptrdiff_t; + public: /** \brief Geometric traits class that provides the geometric types and predicates needed by the triangulations.*/ using Geom_traits = std::conditional_t, Kernel>; // Add an int in TDS to save point index in the structure using TDS = CGAL::Triangulation_data_structure, + CGAL::Triangulation_vertex, CGAL::Triangulation_full_cell >; /** \brief A (Weighted or not) Delaunay triangulation of a set of points in \f$ \mathbb{R}^D\f$.*/ @@ -132,9 +136,6 @@ class Alpha_complex { // Vertex_iterator type from CGAL. using CGAL_vertex_iterator = typename Triangulation::Vertex_iterator; - // size_type type from CGAL. - using size_type = typename Triangulation::size_type; - // Structure to switch from simplex tree vertex handle to CGAL vertex iterator. using Vector_vertex_iterator = std::vector< CGAL_vertex_iterator >; @@ -149,7 +150,7 @@ class Alpha_complex { /** \brief Vertices to be inserted first by the create_complex method to avoid quadratic complexity. * It isn't just [0, n) if some points have multiplicity (only one copy appears in the complex). */ - std::vector vertices_; + std::vector vertices_; /** \brief Cache for geometric constructions: circumcenter and squared radius of a simplex.*/ std::vector cache_, old_cache_; @@ -261,11 +262,11 @@ class Alpha_complex { std::vector point_cloud(first, last); // Creates a vector {0, 1, ..., N-1} - std::vector indices(boost::counting_iterator(0), - boost::counting_iterator(point_cloud.size())); + std::vector indices(boost::counting_iterator(0), + boost::counting_iterator(point_cloud.size())); using Point_property_map = boost::iterator_property_map::iterator, - CGAL::Identity_property_map>; + CGAL::Identity_property_map>; using Search_traits_d = CGAL::Spatial_sort_traits_adapter_d; CGAL::spatial_sort(indices.begin(), indices.end(), Search_traits_d(std::begin(point_cloud))); -- cgit v1.2.3