diff options
Diffstat (limited to 'src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h')
| -rw-r--r-- | src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h | 113 |
1 files changed, 97 insertions, 16 deletions
diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h index 0442ac34..8e8d3a52 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h @@ -32,7 +32,8 @@ #include <algorithm> #include <iterator> #include <limits> -#include <utility> // for pair<> +#include <utility> +#include <stdexcept> namespace Gudhi { @@ -100,6 +101,8 @@ class Bitmap_cubical_complex_base { * non-negative integer, indicating a position of a cube in the data structure. * In the case of functions that compute (co)boundary, the output is a vector if non-negative integers pointing to * the positions of (co)boundary element of the input cell. + * The boundary elements are guaranteed to be returned so that the + * incidence coefficients of boundary elements are alternating. */ virtual inline std::vector< size_t > get_boundary_of_a_cell(size_t cell)const; /** @@ -112,13 +115,79 @@ class Bitmap_cubical_complex_base { * In the case of functions that compute (co)boundary, the output is a vector if * non-negative integers pointing to the * positions of (co)boundary element of the input cell. + * Note that unlike in the case of boundary, over here the elements are + * not guaranteed to be returned with alternating incidence numbers. + * **/ virtual inline std::vector< size_t > get_coboundary_of_a_cell(size_t cell)const; + /** + * This procedure compute incidence numbers between cubes. For a cube \f$A\f$ of + * dimension n and a cube \f$B \subset A\f$ of dimension n-1, an incidence + * between \f$A\f$ and \f$B\f$ is the integer with which \f$B\f$ appears in the boundary of \f$A\f$. + * Note that first parameter is a cube of dimension n, + * and the second parameter is an adjusted cube in dimension n-1. + * Given \f$A = [b_1,e_1] \times \ldots \ [b_{j-1},e_{j-1}] \times [b_{j},e_{j}] \times [b_{j+1},e_{j+1}] \times \ldots \times [b_{n},e_{n}] \f$ + * such that \f$ b_{j} \neq e_{j} \f$ + * and \f$B = [b_1,e_1] \times \ldots \ [b_{j-1},e_{j-1}] \times [a,a] \times [b_{j+1},e_{j+1}] \times \ldots \times [b_{n},e_{n}] \f$ + * where \f$ a = b_{j}\f$ or \f$ a = e_{j}\f$, the incidence between \f$A\f$ and \f$B\f$ + * computed by this procedure is given by formula: + * \f$ c\ (-1)^{\sum_{i=1}^{j-1} dim [b_{i},e_{i}]} \f$ + * Where \f$ dim [b_{i},e_{i}] = 0 \f$ if \f$ b_{i}=e_{i} \f$ and 1 in other case. + * c is -1 if \f$ a = b_{j}\f$ and 1 if \f$ a = e_{j}\f$. + * @exception std::logic_error In case when the cube \f$B\f$ is not n-1 + * dimensional face of a cube \f$A\f$. + **/ + virtual int compute_incidence_between_cells( size_t coface , size_t face )const + { + + //first get the counters for coface and face: + std::vector<unsigned> coface_counter = this->compute_counter_for_given_cell( coface ); + std::vector<unsigned> face_counter = this->compute_counter_for_given_cell( face ); + + //coface_counter and face_counter should agree at all positions except from one: + int number_of_position_in_which_counters_do_not_agree = -1; + size_t number_of_full_faces_that_comes_before = 0; + for ( size_t i = 0 ; i != coface_counter.size() ; ++i ) + { + if ( (coface_counter[i]%2 == 1)&&(number_of_position_in_which_counters_do_not_agree==-1) ) + { + ++number_of_full_faces_that_comes_before; + } + if ( coface_counter[i] != face_counter[i] ) + { + if ( number_of_position_in_which_counters_do_not_agree != -1 ) + { + std::cout << "Cells given to compute_incidence_between_cells procedure do not form a pair of coface-face.\n"; + throw std::logic_error("Cells given to compute_incidence_between_cells procedure do not form a pair of coface-face."); + } + number_of_position_in_which_counters_do_not_agree = i; + } + } + + int incidence = 1; + if ( number_of_full_faces_that_comes_before%2 )incidence = -1; + //if the face cell is on the right from coface cell: + if ( coface_counter[number_of_position_in_which_counters_do_not_agree]+1 == + face_counter[number_of_position_in_which_counters_do_not_agree] + ) + { + incidence *= -1; + } + + return incidence; + } + + /** * In the case of get_dimension_of_a_cell function, the output is a non-negative integer * indicating the dimension of a cell. + * Note that unlike in the case of boundary, over here the elements are + * not guaranteed to be returned with alternating incidence numbers. + * To compute incidence between cells use compute_incidence_between_cells + * procedure **/ inline unsigned get_dimension_of_a_cell(size_t cell)const; + /** * In the case of get_cell_data, the output parameter is a reference to the value of a cube in a given position. * This allows reading and changing the value of filtration. Note that if the value of a filtration is changed, the @@ -289,7 +358,7 @@ class Bitmap_cubical_complex_base { * boundary_simplex_range creates an object of a Boundary_simplex_range class * that provides ranges for the Boundary_simplex_iterator. **/ - Boundary_range boundary_range(size_t sh) { + Boundary_range boundary_range(size_t sh) { return this->get_boundary_of_a_cell(sh); } @@ -668,21 +737,32 @@ Bitmap_cubical_complex_base<T>::Bitmap_cubical_complex_base(const char* perseus_ } template <typename T> -std::vector< size_t > Bitmap_cubical_complex_base<T>::get_boundary_of_a_cell(size_t cell)const { +std::vector< size_t > Bitmap_cubical_complex_base<T>::get_boundary_of_a_cell(size_t cell)const { std::vector< size_t > boundary_elements; // Speed traded of for memory. Check if it is better in practice. boundary_elements.reserve(this->dimension()*2); - size_t cell1 = cell; + size_t sum_of_dimensions = 0; + size_t cell1 = cell; for (size_t i = this->multipliers.size(); i != 0; --i) { unsigned position = cell1 / this->multipliers[i - 1]; if (position % 2 == 1) { - boundary_elements.push_back(cell - this->multipliers[ i - 1 ]); - boundary_elements.push_back(cell + this->multipliers[ i - 1 ]); + if ( sum_of_dimensions%2 ) + { + boundary_elements.push_back(cell + this->multipliers[ i - 1 ]); + boundary_elements.push_back(cell - this->multipliers[ i - 1 ]); + } + else + { + boundary_elements.push_back(cell - this->multipliers[ i - 1 ]); + boundary_elements.push_back(cell + this->multipliers[ i - 1 ]); + } + ++sum_of_dimensions; } cell1 = cell1 % this->multipliers[i - 1]; - } + } + return boundary_elements; } @@ -690,23 +770,24 @@ template <typename T> std::vector< size_t > Bitmap_cubical_complex_base<T>::get_coboundary_of_a_cell(size_t cell)const { std::vector<unsigned> counter = this->compute_counter_for_given_cell(cell); std::vector< size_t > coboundary_elements; - size_t cell1 = cell; + size_t cell1 = cell; for (size_t i = this->multipliers.size(); i != 0; --i) { unsigned position = cell1 / this->multipliers[i - 1]; if (position % 2 == 0) { - if ((cell > this->multipliers[i - 1]) && (counter[i - 1] != 0)) { - coboundary_elements.push_back(cell - this->multipliers[i - 1]); - } - if ( - (cell + this->multipliers[i - 1] < this->data.size()) && (counter[i - 1] != 2 * this->sizes[i - 1])) { - coboundary_elements.push_back(cell + this->multipliers[i - 1]); - } + if ((cell > this->multipliers[i - 1]) && (counter[i - 1] != 0)) { + coboundary_elements.push_back(cell - this->multipliers[i - 1]); + } + if ( + (cell + this->multipliers[i - 1] < this->data.size()) && (counter[i - 1] != 2 * this->sizes[i - 1])) { + coboundary_elements.push_back(cell + this->multipliers[i - 1]); + } } cell1 = cell1 % this->multipliers[i - 1]; - } + } return coboundary_elements; } + template <typename T> unsigned Bitmap_cubical_complex_base<T>::get_dimension_of_a_cell(size_t cell)const { bool dbg = false; |