/* 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): Pawel Dlotko * * Copyright (C) 2015 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 . */ #ifndef BITMAP_CUBICAL_COMPLEX_H_ #define BITMAP_CUBICAL_COMPLEX_H_ #include #include #ifdef GUDHI_USE_TBB #include #endif #include #include // for pair<> #include // for sort #include #include // for iota #include namespace Gudhi { namespace cubical_complex { // global variable, was used just for debugging. const bool globalDbg = false; template class is_before_in_filtration; /** * @brief Cubical complex represented as a bitmap. * @ingroup cubical_complex * @details This is a Bitmap_cubical_complex class. It joints a functionalities of Bitmap_cubical_complex_base and * Bitmap_cubical_complex_periodic_boundary_conditions_base classes into * Gudhi persistent homology engine. It is a template class that inherit from its template parameter. The template * parameter is supposed to be either Bitmap_cubical_complex_base or * Bitmap_cubical_complex_periodic_boundary_conditions_base class. **/ template class Bitmap_cubical_complex : public T { public: //*********************************************// // Typedefs and typenames //*********************************************// typedef std::size_t Simplex_key; typedef typename T::filtration_type Filtration_value; typedef Simplex_key Simplex_handle; //*********************************************// // Constructors //*********************************************// // Over here we need to define various input types. I am proposing the following ones: // Perseus style // TODO(PD) H5 files? // TODO(PD) binary files with little endiangs / big endians ? // TODO(PD) constructor from a vector of elements of a type T. ? /** * Constructor form a Perseus-style file. **/ Bitmap_cubical_complex(const char* perseus_style_file) : T(perseus_style_file), key_associated_to_simplex(this->total_number_of_cells + 1) { if (globalDbg) { std::cerr << "Bitmap_cubical_complex( const char* perseus_style_file )\n"; } for (std::size_t i = 0; i != this->total_number_of_cells; ++i) { this->key_associated_to_simplex[i] = i; } // we initialize this only once, in each constructor, when the bitmap is constructed. // If the user decide to change some elements of the bitmap, then this procedure need // to be called again. this->initialize_simplex_associated_to_key(); } /** * Constructor that requires vector of elements of type unsigned, which gives number of top dimensional cells * in the following directions and vector of element of a type T * with filtration on top dimensional cells. **/ Bitmap_cubical_complex(const std::vector& dimensions, const std::vector& top_dimensional_cells) : T(dimensions, top_dimensional_cells), key_associated_to_simplex(this->total_number_of_cells + 1) { for (std::size_t i = 0; i != this->total_number_of_cells; ++i) { this->key_associated_to_simplex[i] = i; } // we initialize this only once, in each constructor, when the bitmap is constructed. // If the user decide to change some elements of the bitmap, then this procedure need // to be called again. this->initialize_simplex_associated_to_key(); } /** * Constructor that requires vector of elements of type unsigned, which gives number of top dimensional cells * in the following directions and vector of element of a type Filtration_value * with filtration on top dimensional cells. The last parameter of the constructor is a vector of boolean of a length * equal to the dimension of cubical complex. * If the position i on this vector is true, then we impose periodic boundary conditions in this direction. **/ Bitmap_cubical_complex(const std::vector& dimensions, const std::vector& top_dimensional_cells, std::vector directions_in_which_periodic_b_cond_are_to_be_imposed) : T(dimensions, top_dimensional_cells, directions_in_which_periodic_b_cond_are_to_be_imposed), key_associated_to_simplex(this->total_number_of_cells + 1) { for (std::size_t i = 0; i != this->total_number_of_cells; ++i) { this->key_associated_to_simplex[i] = i; } // we initialize this only once, in each constructor, when the bitmap is constructed. // If the user decide to change some elements of the bitmap, then this procedure need // to be called again. this->initialize_simplex_associated_to_key(); } /** * Destructor of the Bitmap_cubical_complex class. **/ virtual ~Bitmap_cubical_complex() {} //*********************************************// // Other 'easy' functions //*********************************************// /** * Returns number of all cubes in the complex. **/ std::size_t num_simplices() const { return this->total_number_of_cells; } /** * Returns a Simplex_handle to a cube that do not exist in this complex. **/ static Simplex_handle null_simplex() { if (globalDbg) { std::cerr << "Simplex_handle null_simplex()\n"; } return std::numeric_limits::max(); } /** * Returns dimension of the complex. **/ inline std::size_t dimension() const { return this->sizes.size(); } /** * Return dimension of a cell pointed by the Simplex_handle. **/ inline unsigned dimension(Simplex_handle sh) const { if (globalDbg) { std::cerr << "unsigned dimension(const Simplex_handle& sh)\n"; } if (sh != null_simplex()) return this->get_dimension_of_a_cell(sh); return -1; } /** * Return the filtration of a cell pointed by the Simplex_handle. **/ Filtration_value filtration(Simplex_handle sh) { if (globalDbg) { std::cerr << "Filtration_value filtration(const Simplex_handle& sh)\n"; } // Returns the filtration value of a simplex. if (sh != null_simplex()) return this->data[sh]; return std::numeric_limits::infinity(); } /** * Return a key which is not a key of any cube in the considered data structure. **/ static Simplex_key null_key() { if (globalDbg) { std::cerr << "Simplex_key null_key()\n"; } return std::numeric_limits::max(); } /** * Return the key of a cube pointed by the Simplex_handle. **/ Simplex_key key(Simplex_handle sh) const { if (globalDbg) { std::cerr << "Simplex_key key(const Simplex_handle& sh)\n"; } if (sh != null_simplex()) { return this->key_associated_to_simplex[sh]; } return this->null_key(); } /** * Return the Simplex_handle given the key of the cube. **/ Simplex_handle simplex(Simplex_key key) { if (globalDbg) { std::cerr << "Simplex_handle simplex(Simplex_key key)\n"; } if (key != null_key()) { return this->simplex_associated_to_key[key]; } return null_simplex(); } /** * Assign key to a cube pointed by the Simplex_handle **/ void assign_key(Simplex_handle sh, Simplex_key key) { if (globalDbg) { std::cerr << "void assign_key(Simplex_handle& sh, Simplex_key key)\n"; } if (key == null_key()) return; this->key_associated_to_simplex[sh] = key; this->simplex_associated_to_key[key] = sh; } /** * Function called from a constructor. It is needed for Filtration_simplex_iterator to work. **/ void initialize_simplex_associated_to_key(); //*********************************************// // Iterators //*********************************************// /** * Boundary_simplex_range class provides ranges for boundary iterators. **/ typedef typename std::vector::iterator Boundary_simplex_iterator; typedef typename std::vector Boundary_simplex_range; /** * Filtration_simplex_iterator class provides an iterator though the whole structure in the order of filtration. * Secondary criteria for filtration are: * (1) Dimension of a cube (lower dimensional comes first). * (2) Position in the data structure (the ones that are earlies in the data structure comes first). **/ class Filtration_simplex_range; class Filtration_simplex_iterator : std::iterator { // Iterator over all simplices of the complex in the order of the indexing scheme. // 'value_type' must be 'Simplex_handle'. public: Filtration_simplex_iterator(Bitmap_cubical_complex* b) : b(b), position(0) {} Filtration_simplex_iterator() : b(NULL), position(0) {} Filtration_simplex_iterator operator++() { if (globalDbg) { std::cerr << "Filtration_simplex_iterator operator++\n"; } ++this->position; return (*this); } Filtration_simplex_iterator operator++(int) { Filtration_simplex_iterator result = *this; ++(*this); return result; } Filtration_simplex_iterator& operator=(const Filtration_simplex_iterator& rhs) { if (globalDbg) { std::cerr << "Filtration_simplex_iterator operator =\n"; } this->b = rhs.b; this->position = rhs.position; return (*this); } bool operator==(const Filtration_simplex_iterator& rhs) const { if (globalDbg) { std::cerr << "bool operator == ( const Filtration_simplex_iterator& rhs )\n"; } return (this->position == rhs.position); } bool operator!=(const Filtration_simplex_iterator& rhs) const { if (globalDbg) { std::cerr << "bool operator != ( const Filtration_simplex_iterator& rhs )\n"; } return !(*this == rhs); } Simplex_handle operator*() { if (globalDbg) { std::cerr << "Simplex_handle operator*()\n"; } return this->b->simplex_associated_to_key[this->position]; } friend class Filtration_simplex_range; private: Bitmap_cubical_complex* b; std::size_t position; }; /** * @brief Filtration_simplex_range provides the ranges for Filtration_simplex_iterator. **/ class Filtration_simplex_range { // Range over the simplices of the complex in the order of the filtration. // .begin() and .end() return type Filtration_simplex_iterator. public: typedef Filtration_simplex_iterator const_iterator; typedef Filtration_simplex_iterator iterator; Filtration_simplex_range(Bitmap_cubical_complex* b) : b(b) {} Filtration_simplex_iterator begin() { if (globalDbg) { std::cerr << "Filtration_simplex_iterator begin() \n"; } return Filtration_simplex_iterator(this->b); } Filtration_simplex_iterator end() { if (globalDbg) { std::cerr << "Filtration_simplex_iterator end()\n"; } Filtration_simplex_iterator it(this->b); it.position = this->b->simplex_associated_to_key.size(); return it; } private: Bitmap_cubical_complex* b; }; //*********************************************// // Methods to access iterators from the container: /** * boundary_simplex_range creates an object of a Boundary_simplex_range class * that provides ranges for the Boundary_simplex_iterator. **/ Boundary_simplex_range boundary_simplex_range(Simplex_handle sh) { return this->get_boundary_of_a_cell(sh); } /** * filtration_simplex_range creates an object of a Filtration_simplex_range class * that provides ranges for the Filtration_simplex_iterator. **/ Filtration_simplex_range filtration_simplex_range() { if (globalDbg) { std::cerr << "Filtration_simplex_range filtration_simplex_range()\n"; } // Returns a range over the simplices of the complex in the order of the filtration return Filtration_simplex_range(this); } //*********************************************// //*********************************************// // Elements which are in Gudhi now, but I (and in all the cases I asked also Marc) do not understand why they are // there. // TODO(PD) the file IndexingTag.h in the Gudhi library contains an empty structure, so // I understand that this is something that was planned (for simplicial maps?) // but was never finished. The only idea I have here is to use the same empty structure from // IndexingTag.h file, but only if the compiler needs it. If the compiler // do not need it, then I would rather not add here elements which I do not understand. // typedef Indexing_tag /** * Function needed for compatibility with Gudhi. Not useful for other purposes. **/ std::pair endpoints(Simplex_handle sh) { std::vector bdry = this->get_boundary_of_a_cell(sh); if (globalDbg) { std::cerr << "std::pair endpoints( Simplex_handle sh )\n"; std::cerr << "bdry.size() : " << bdry.size() << "\n"; } // this method returns two first elements from the boundary of sh. if (bdry.size() < 2) throw( "Error in endpoints in Bitmap_cubical_complex class. The cell have less than two elements in the " "boundary."); return std::make_pair(bdry[0], bdry[1]); } /** * Class needed for compatibility with Gudhi. Not useful for other purposes. **/ class Skeleton_simplex_range; class Skeleton_simplex_iterator : std::iterator { // Iterator over all simplices of the complex in the order of the indexing scheme. // 'value_type' must be 'Simplex_handle'. public: Skeleton_simplex_iterator(Bitmap_cubical_complex* b, std::size_t d) : b(b), dimension(d) { if (globalDbg) { std::cerr << "Skeleton_simplex_iterator ( Bitmap_cubical_complex* b , std::size_t d )\n"; } // find the position of the first simplex of a dimension d this->position = 0; while ((this->position != b->data.size()) && (this->b->get_dimension_of_a_cell(this->position) != this->dimension)) { ++this->position; } } Skeleton_simplex_iterator() : b(NULL), position(0), dimension(0) {} Skeleton_simplex_iterator operator++() { if (globalDbg) { std::cerr << "Skeleton_simplex_iterator operator++()\n"; } // increment the position as long as you did not get to the next element of the dimension dimension. ++this->position; while ((this->position != this->b->data.size()) && (this->b->get_dimension_of_a_cell(this->position) != this->dimension)) { ++this->position; } return (*this); } Skeleton_simplex_iterator operator++(int) { Skeleton_simplex_iterator result = *this; ++(*this); return result; } Skeleton_simplex_iterator& operator=(const Skeleton_simplex_iterator& rhs) { if (globalDbg) { std::cerr << "Skeleton_simplex_iterator operator =\n"; } this->b = rhs.b; this->position = rhs.position; this->dimension = rhs.dimension; return (*this); } bool operator==(const Skeleton_simplex_iterator& rhs) const { if (globalDbg) { std::cerr << "bool operator ==\n"; } return (this->position == rhs.position); } bool operator!=(const Skeleton_simplex_iterator& rhs) const { if (globalDbg) { std::cerr << "bool operator != ( const Skeleton_simplex_iterator& rhs )\n"; } return !(*this == rhs); } Simplex_handle operator*() { if (globalDbg) { std::cerr << "Simplex_handle operator*() \n"; } return this->position; } friend class Skeleton_simplex_range; private: Bitmap_cubical_complex* b; std::size_t position; unsigned dimension; }; /** * @brief Class needed for compatibility with Gudhi. Not useful for other purposes. **/ class Skeleton_simplex_range { // Range over the simplices of the complex in the order of the filtration. // .begin() and .end() return type Filtration_simplex_iterator. public: typedef Skeleton_simplex_iterator const_iterator; typedef Skeleton_simplex_iterator iterator; Skeleton_simplex_range(Bitmap_cubical_complex* b, unsigned dimension) : b(b), dimension(dimension) {} Skeleton_simplex_iterator begin() { if (globalDbg) { std::cerr << "Skeleton_simplex_iterator begin()\n"; } return Skeleton_simplex_iterator(this->b, this->dimension); } Skeleton_simplex_iterator end() { if (globalDbg) { std::cerr << "Skeleton_simplex_iterator end()\n"; } Skeleton_simplex_iterator it(this->b, this->dimension); it.position = this->b->data.size(); return it; } private: Bitmap_cubical_complex* b; unsigned dimension; }; /** * Function needed for compatibility with Gudhi. Not useful for other purposes. **/ Skeleton_simplex_range skeleton_simplex_range(unsigned dimension) { if (globalDbg) { std::cerr << "Skeleton_simplex_range skeleton_simplex_range( unsigned dimension )\n"; } return Skeleton_simplex_range(this, dimension); } friend class is_before_in_filtration; protected: std::vector key_associated_to_simplex; std::vector simplex_associated_to_key; }; // Bitmap_cubical_complex template void Bitmap_cubical_complex::initialize_simplex_associated_to_key() { if (globalDbg) { std::cerr << "void Bitmap_cubical_complex::initialize_elements_ordered_according_to_filtration() \n"; } this->simplex_associated_to_key = std::vector(this->data.size()); std::iota(std::begin(simplex_associated_to_key), std::end(simplex_associated_to_key), 0); #ifdef GUDHI_USE_TBB tbb::parallel_sort(simplex_associated_to_key.begin(), simplex_associated_to_key.end(), is_before_in_filtration(this)); #else std::sort(simplex_associated_to_key.begin(), simplex_associated_to_key.end(), is_before_in_filtration(this)); #endif // we still need to deal here with a key_associated_to_simplex: for (std::size_t i = 0; i != simplex_associated_to_key.size(); ++i) { this->key_associated_to_simplex[simplex_associated_to_key[i]] = i; } } template class is_before_in_filtration { public: explicit is_before_in_filtration(Bitmap_cubical_complex* CC) : CC_(CC) {} bool operator()(const typename Bitmap_cubical_complex::Simplex_handle& sh1, const typename Bitmap_cubical_complex::Simplex_handle& sh2) const { // Not using st_->filtration(sh1) because it uselessly tests for null_simplex. typedef typename T::filtration_type Filtration_value; Filtration_value fil1 = CC_->data[sh1]; Filtration_value fil2 = CC_->data[sh2]; if (fil1 != fil2) { return fil1 < fil2; } // in this case they are on the same filtration level, so the dimension decide. std::size_t dim1 = CC_->get_dimension_of_a_cell(sh1); std::size_t dim2 = CC_->get_dimension_of_a_cell(sh2); if (dim1 != dim2) { return dim1 < dim2; } // in this case both filtration and dimensions of the considered cubes are the same. To have stable sort, we simply // compare their positions in the bitmap: return sh1 < sh2; } protected: Bitmap_cubical_complex* CC_; }; } // namespace cubical_complex namespace Cubical_complex = cubical_complex; } // namespace Gudhi #endif // BITMAP_CUBICAL_COMPLEX_H_