From 39ffc7af64543f0b2ce17487771afbd010a97349 Mon Sep 17 00:00:00 2001 From: pdlotko Date: Fri, 27 Nov 2015 14:16:44 +0000 Subject: Adding changes to the cubcial complex class according to Marc and Vincent's comments. git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/bitmap@930 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: cb9134a5c18fae3e2e63bbaf8329dd168d08d3b3 --- .../example/Bitmap_cubical_complex.cpp | 147 +- src/Bitmap_cubical_complex/example/CMakeLists.txt | 2 +- .../example/Random_bitmap_cubical_complex.cpp | 174 +-- .../include/gudhi/Bitmap_cubical_complex.h | 1443 ++++++++++---------- .../include/gudhi/Bitmap_cubical_complex_base.h | 1352 ++++++++++-------- src/Bitmap_cubical_complex/include/gudhi/counter.h | 316 +++-- src/Bitmap_cubical_complex/test/Bitmap_test.cpp | 1261 ++++++++--------- 7 files changed, 2463 insertions(+), 2232 deletions(-) (limited to 'src/Bitmap_cubical_complex') diff --git a/src/Bitmap_cubical_complex/example/Bitmap_cubical_complex.cpp b/src/Bitmap_cubical_complex/example/Bitmap_cubical_complex.cpp index 37c16618..31da3609 100644 --- a/src/Bitmap_cubical_complex/example/Bitmap_cubical_complex.cpp +++ b/src/Bitmap_cubical_complex/example/Bitmap_cubical_complex.cpp @@ -1,71 +1,76 @@ -/* 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 Sophia-Saclay (France) - * - * 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 . - */ - -// for persistence algorithm -#include -#include -#include - -#include - -// standard stuff -#include -#include - -using namespace Gudhi; -using namespace Gudhi::persistent_cohomology; -using namespace std; - -int main(int argc, char** argv) { - cout << "This program computes persistent homology, by using Bitmap_cubical_complex class, of cubical complexes " - "provided in text files in Perseus style (the only numbed in the first line is a dimension D of a cubical " - "complex. In the lines I between 2 and D+1 there are numbers of top dimensional cells in the direction I. Let N " - "denote product of the numbers in the lines between 2 and D. In the lines D+2 to D+2+N there are filtrations of " - "top dimensional cells. We assume that the cells are in the lexicographical order. See CubicalOneSphere.txt or " - "CubicalTwoSphere.txt for example." << endl; - - int p = 2; - double min_persistence = 0; - - if (argc != 2) { - cout << "Wrong number of parameters. Please provide the name of a file with a Perseus style cubical complex at the " - "input. The program will now terminate.\n"; - return 1; - } - - Bitmap_cubical_complex b(argv[1]); - - - // Compute the persistence diagram of the complex - persistent_cohomology::Persistent_cohomology< Bitmap_cubical_complex, Field_Zp > pcoh(b); - pcoh.init_coefficients(p); // initializes the coefficient field for homology - pcoh.compute_persistent_cohomology(min_persistence); - - - stringstream ss; - ss << argv[1] << "_persistence"; - std::ofstream out((char*) ss.str().c_str()); - pcoh.output_diagram(out); - out.close(); - - return 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): Pawel Dlotko + * + * Copyright (C) 2015 INRIA Saclay (France) + * + * 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 . + */ + + +//for persistence algorithm +#include +#include +#include + + +using namespace Gudhi; +using namespace Gudhi::Cubical_complex; +using namespace Gudhi::persistent_cohomology; + +//standard stuff +#include +#include + +using namespace std; + +int main( int argc , char** argv ) +{ + cout << "This program computes persistent homology, by using bitmap_cubical_complex class,\ + of cubical complexes provided in text files in Perseus style (the only numbered in \ +the first line is a dimension D of a bitmap. In the lines I between 2 and D+1 there are\ + numbers of top dimensional cells in the direction I. Let N denote product \ +of the numbers in the lines between 2 and D. In the lines D+2 to D+2+N there are\ + filtrations of top dimensional cells. We assume that the cells are in the \ +lexicographical order. See CubicalOneSphere.txt or CubicalTwoSphere.txt for example." << endl; + + int p = 2; + double min_persistence = 0; + + if ( argc != 2 ) + { + cout << "Wrong number of parameters. Please provide the name of a file with a\ + Perseus style bitmap at the input. The program will now terminate.\n"; + return 1; + } + + Bitmap_cubical_complex b( argv[1] ); + + + // Compute the persistence diagram of the complex + persistent_cohomology::Persistent_cohomology< Bitmap_cubical_complex, Field_Zp > pcoh(b); + pcoh.init_coefficients( p ); //initilizes the coefficient field for homology + pcoh.compute_persistent_cohomology( min_persistence ); + + + stringstream ss; + ss << argv[1] << "_persistence"; + std::ofstream out((char*)ss.str().c_str()); + pcoh.output_diagram(out); + out.close(); + + return 0; +} diff --git a/src/Bitmap_cubical_complex/example/CMakeLists.txt b/src/Bitmap_cubical_complex/example/CMakeLists.txt index 05ef1319..dd252a79 100644 --- a/src/Bitmap_cubical_complex/example/CMakeLists.txt +++ b/src/Bitmap_cubical_complex/example/CMakeLists.txt @@ -1,5 +1,5 @@ cmake_minimum_required(VERSION 2.6) -project(GUDHISimplexTreeFromFile) +project(GUDHIBitmap) add_executable ( Bitmap_cubical_complex Bitmap_cubical_complex.cpp ) target_link_libraries(Bitmap_cubical_complex ${Boost_SYSTEM_LIBRARY}) diff --git a/src/Bitmap_cubical_complex/example/Random_bitmap_cubical_complex.cpp b/src/Bitmap_cubical_complex/example/Random_bitmap_cubical_complex.cpp index ac7557ce..60cfc113 100644 --- a/src/Bitmap_cubical_complex/example/Random_bitmap_cubical_complex.cpp +++ b/src/Bitmap_cubical_complex/example/Random_bitmap_cubical_complex.cpp @@ -1,81 +1,93 @@ -/* 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 Saclay (France) - * - * 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 . - */ - - -// for persistence algorithm -#include -#include -#include - -#include - -// standard stuff -#include -#include -#include -#include -#include - -using namespace Gudhi; -using namespace Gudhi::persistent_cohomology; -using namespace std; - -int main(int argc, char** argv) { - srand(time(0)); - - cout << "This program computes persistent homology, by using Bitmap_cubical_complex class, of cubical complexes. " - "The first parameter of the program is the dimension D of the cubical complex. The next D parameters are number " - "of top dimensional cubes in each dimension of the cubical complex. The program will create random cubical " - "complex of that sizes and compute persistent homology of it." << endl; - - int p = 2; - double min_persistence = 0; - - size_t dimensionOfBitmap = (size_t) atoi(argv[1]); - std::vector< unsigned > sizes; - size_t multipliers = 1; - for (size_t dim = 0; dim != dimensionOfBitmap; ++dim) { - unsigned sizeInThisDimension = (unsigned) atoi(argv[2 + dim]); - sizes.push_back(sizeInThisDimension); - multipliers *= sizeInThisDimension; - } - - std::vector< double > data; - for (size_t i = 0; i != multipliers; ++i) { - data.push_back(rand() / (double) RAND_MAX); - } - - Bitmap_cubical_complex b(sizes, data); - - // Compute the persistence diagram of the complex - persistent_cohomology::Persistent_cohomology< Bitmap_cubical_complex, Field_Zp > pcoh(b); - pcoh.init_coefficients(p); // initializes the coefficient field for homology - pcoh.compute_persistent_cohomology(min_persistence); - - stringstream ss; - ss << "randomComplex_persistence"; - std::ofstream out((char*) ss.str().c_str()); - pcoh.output_diagram(out); - out.close(); - - return 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): Pawel Dlotko + * + * Copyright (C) 2015 INRIA Saclay (France) + * + * 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 . + */ + + +//for persistence algorithm +#include +#include +#include + + +using namespace Gudhi; +using namespace Gudhi::Cubical_complex; +using namespace Gudhi::persistent_cohomology; + +//standard stuff +#include +#include +#include + +using namespace std; + +int main( int argc , char** argv ) +{ + srand( time(0) ); + + cout << "This program computes persistent homology, by using bitmap_cubical_complex class, of cubical complexes. \ +The first parameter of the program is the dimension D of the bitmap. \ +The next D parameters are number of top dimensional cubes in each dimension of the bitmap.\ +The program will create random cubical complex of that sizes and compute persistent homology of it." << endl; + + int p = 2; + double min_persistence = 0; + + if ( argc < 3 ) + { + cerr << "Wrong number of parameters, the program will now terminate\n"; + return 1; + } + + size_t dimensionOfBitmap = (size_t)atoi( argv[1] ); + std::vector< unsigned > sizes; + size_t multipliers = 1; + for ( size_t dim = 0 ; dim != dimensionOfBitmap ; ++dim ) + { + unsigned sizeInThisDimension = (unsigned)atoi( argv[2+dim] ); + sizes.push_back( sizeInThisDimension ); + multipliers *= sizeInThisDimension; + } + + std::vector< double > data; + for ( size_t i = 0 ; i != multipliers ; ++i ) + { + data.push_back( rand()/(double)RAND_MAX ); + } + + + + Bitmap_cubical_complex b( sizes , data ); + + + // Compute the persistence diagram of the complex + persistent_cohomology::Persistent_cohomology< Bitmap_cubical_complex, Field_Zp > pcoh(b); + pcoh.init_coefficients( p ); //initilizes the coefficient field for homology + pcoh.compute_persistent_cohomology( min_persistence ); + + + stringstream ss; + ss << "randomComplex_persistence"; + std::ofstream out((char*)ss.str().c_str()); + pcoh.output_diagram(out); + out.close(); + + return 0; +} diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h index 2f8cb0a3..f2c753d9 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h @@ -1,724 +1,719 @@ -/* 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 Sophia-Saclay (France) - * - * 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 -#include // for pair -#include // for sort -#include // for vector - -// global variable, was used just for debugging. -bool globalDbg = false; - -template -class Bitmap_cubical_complex : public Bitmap_cubical_complex_base { - public: - //******************************************************************************************************************// - // Typedefs and typenames - //******************************************************************************************************************// - friend class Simplex_handle; - typedef size_t Simplex_key; - typedef T Filtration_value; - - - //******************************************************************************************************************// - // Simplex handle class - //******************************************************************************************************************// - - /** - * Handle of a cell, required for compatibility with the function to compute persistence in Gudhi. Elements of this - * class are: the pointer to the bitmap B in which the considered cell is together with a position of this cell in B. - * Given this data, one can get all the information about the considered cell. - **/ - class Simplex_handle { - public: - Simplex_handle() { - if (globalDbg) { - std::cerr << "Simplex_handle()\n"; - } - this->b = 0; - this->position = 0; - } - - Simplex_handle(Bitmap_cubical_complex* b) { - if (globalDbg) { - std::cerr << "Simplex_handle(Bitmap_cubical_complex* b)\n"; - } - this->b = b; - this->position = 0; - } - - Simplex_handle(const Simplex_handle& org) : b(org.b) { - if (globalDbg) { - std::cerr << "Simplex_handle( const Simplex_handle& org )\n"; - } - this->position = org.position; - } - - Simplex_handle& operator=(const Simplex_handle& rhs) { - if (globalDbg) { - std::cerr << "Simplex_handle operator = \n"; - } - this->position = rhs.position; - this->b = rhs.b; - return *this; - } - - Simplex_handle(Bitmap_cubical_complex* b, Simplex_key position) { - if (globalDbg) { - std::cerr << "Simplex_handle(Bitmap_cubical_complex* b , Simplex_key position)\n"; - std::cerr << "Position : " << position << std::endl; - } - this->b = b; - this->position = position; - } - friend class Bitmap_cubical_complex; - private: - Bitmap_cubical_complex* b; - Simplex_key position; - // Assumption -- this field always keep the REAL position of simplex in the bitmap, no matter what keys have been. - // to deal with the keys, the class Bitmap_cubical_complex have extra vectors: keyAssociatedToSimplex and - // simplexAssociatedToKey that allow to move between actual cell and the key assigned to it. - }; - - - //******************************************************************************************************************// - // Constructors - //******************************************************************************************************************// - // Over here we need to definie various input types. I am proposing the following ones: - // Perseus style - // TODO(Pawel Dlotko): H5 files? - // TODO(Pawel Dlotko): binary files with little endiangs / big endians? - // TODO(Pawel Dlotko): constructor from a vector of elements of a type T. - - /** - * Constructor form a Perseus-style file. - **/ - Bitmap_cubical_complex(char* perseusStyleFile) : Bitmap_cubical_complex_base(perseusStyleFile) { - if (globalDbg) { - std::cerr << "Bitmap_cubical_complex( char* perseusStyleFile )\n"; - } - std::vector< size_t > keyAssociatedToSimplex(this->totalNumberOfCells + 1); - std::vector< size_t > simplexAssociatedToKey(this->totalNumberOfCells + 1); - - for (size_t i = 0; i != this->totalNumberOfCells; ++i) { - keyAssociatedToSimplex[i] = simplexAssociatedToKey[i] = i; - } - this->keyAssociatedToSimplex = keyAssociatedToSimplex; - this->simplexAssociatedToKey = simplexAssociatedToKey; - // 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->initializeElementsOrderedAccordingToFiltration(); - } - - /** - * 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(std::vector dimensions, std::vector topDimensionalCells) - : Bitmap_cubical_complex_base(dimensions, topDimensionalCells) { - std::vector< size_t > keyAssociatedToSimplex(this->totalNumberOfCells + 1); - std::vector< size_t > simplexAssociatedToKey(this->totalNumberOfCells + 1); - - for (size_t i = 0; i != this->totalNumberOfCells; ++i) { - keyAssociatedToSimplex[i] = simplexAssociatedToKey[i] = i; - } - this->keyAssociatedToSimplex = keyAssociatedToSimplex; - this->simplexAssociatedToKey = simplexAssociatedToKey; - // 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->initializeElementsOrderedAccordingToFiltration(); - } - - //******************************************************************************************************************// - // Other 'easy' functions - //******************************************************************************************************************// - - /** - * Returns number of all cubes in the complex. - **/ - size_t num_simplices()const { - return this->totalNumberOfCells; - } - - /** - * Returns a Simplex_handle to a cube that do not exist in this complex. - **/ - Simplex_handle null_simplex() { - return Simplex_handle(this, this->data.size()); - } - - /** - * Returns dimension of the complex. - **/ - size_t dimension() { - return this->sizes.size(); - } - - /** - * Return dimension of a cell pointed by the Simplex_handle. - **/ - size_t dimension(const Simplex_handle& sh) { - if (globalDbg) { - std::cerr << "int dimension(const Simplex_handle& sh)\n"; - } - if (sh.position != this->data.size()) return sh.b->get_dimension_of_a_cell(sh.position); - return std::numeric_limits::max(); - } - - /** - * Return the filtration of a cell pointed by the Simplex_handle. - **/ - T filtration(const Simplex_handle& sh) { - if (globalDbg) { - std::cerr << "T filtration(const Simplex_handle& sh)\n"; - } - // Returns the filtration value of a simplex. - if (sh.position != this->data.size()) return sh.b->data[ sh.position ]; - return INT_MAX; - } - - /** - * Return a key which is not a key of any cube in the considered data structure. - **/ - Simplex_key null_key() { - if (globalDbg) { - std::cerr << "Simplex_key null_key()\n"; - } - return this->data.size(); - } - - /** - * Return the key of a cube pointed by the Simplex_handle. - **/ - Simplex_key key(const Simplex_handle& sh) { - if (globalDbg) { - std::cerr << "Simplex_key key(const Simplex_handle& sh)\n"; - } - return sh.b->keyAssociatedToSimplex[ sh.position ]; - } - - /** - * 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"; - } - return Simplex_handle(this, this->simplexAssociatedToKey[ key ]); - } - - /** - * 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"; - } - this->keyAssociatedToSimplex[sh.position] = key; - this->simplexAssociatedToKey[key] = sh.position; - } - - /** - * Function called from a constructor. It is needed for Filtration_simplex_iterator to work. - **/ - void initializeElementsOrderedAccordingToFiltration(); - - - - //******************************************************************************************************************// - // Iterators - //******************************************************************************************************************// - - /** - * Boundary_simplex_iterator class allows iteration on boundary of each cube. - **/ - class Boundary_simplex_range; - - class Boundary_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > { - // Iterator on the simplices belonging to the boundary of a simplex. - // value_type must be 'Simplex_handle'. - public: - Boundary_simplex_iterator(Simplex_handle& sh) : sh(sh) { - if (globalDbg) { - std::cerr << "Boundary_simplex_iterator( Simplex_handle& sh )\n"; - } - this->position = 0; - this->boundaryElements = this->sh.b->get_boundary_of_a_cell(this->sh.position); - } - - Boundary_simplex_iterator operator++() { - if (globalDbg) { - std::cerr << "Boundary_simplex_iterator operator++()\n"; - } - ++this->position; - return *this; - } - - Boundary_simplex_iterator operator++(int) { - Boundary_simplex_iterator result = *this; - ++(*this); - return result; - } - - Boundary_simplex_iterator operator=(const Boundary_simplex_iterator& rhs) { - if (globalDbg) { - std::cerr << "Boundary_simplex_iterator operator =\n"; - } - this->sh = rhs.sh; - this->boundaryElements.clear(); - this->boundaryElementsinsert(this->boundaryElements.end(), - rhs.boundaryElements.begin(), rhs.boundaryElements.end()); - } - - bool operator==(const Boundary_simplex_iterator& rhs) { - if (globalDbg) { - std::cerr << "bool operator ==\n"; - } - if (this->position == rhs.position) { - if (this->boundaryElements.size() != rhs.boundaryElements.size())return false; - for (size_t i = 0; i != this->boundaryElements.size(); ++i) { - if (this->boundaryElements[i] != rhs.boundaryElements[i])return false; - } - return true; - } - return false; - } - - bool operator!=(const Boundary_simplex_iterator& rhs) { - if (globalDbg) { - std::cerr << "bool operator != \n"; - } - return !(*this == rhs); - } - - Simplex_handle operator*() { - if (globalDbg) { - std::cerr << "Simplex_handle operator*\n"; - } - return Simplex_handle(this->sh.b, this->boundaryElements[this->position]); - } - - friend class Boundary_simplex_range; - private: - Simplex_handle sh; - std::vector< size_t > boundaryElements; - size_t position; - }; - - /** - * Boundary_simplex_range class provides ranges for boundary iterators. - **/ - class Boundary_simplex_range { - // Range giving access to the simplices in the boundary of a simplex. - // .begin() and .end() return type Boundary_simplex_iterator. - public: - Boundary_simplex_range(const Simplex_handle& sh) : sh(sh) { } - - Boundary_simplex_iterator begin() { - if (globalDbg) { - std::cerr << "Boundary_simplex_iterator begin\n"; - } - Boundary_simplex_iterator it(this->sh); - return it; - } - - Boundary_simplex_iterator end() { - if (globalDbg) { - std::cerr << "Boundary_simplex_iterator end()\n"; - } - Boundary_simplex_iterator it(this->sh); - it.position = it.boundaryElements.size(); - return it; - } - - private: - Simplex_handle sh; - }; - - - /** - * 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< std::input_iterator_tag, Simplex_handle > { - // 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) { } - - 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; - } - - bool operator==(const Filtration_simplex_iterator& rhs) { - if (globalDbg) { - std::cerr << "bool operator == ( const Filtration_simplex_iterator& rhs )\n"; - } - if (this->position == rhs.position) { - return true; - } - return false; - } - - bool operator!=(const Filtration_simplex_iterator& rhs) { - 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 Simplex_handle(this->b, this->b->elementsOrderedAccordingToFiltration[ this->position ]); - } - - friend class Filtration_simplex_range; - private: - Bitmap_cubical_complex* b; - size_t position; - }; - - /** - * 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: - 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->elementsOrderedAccordingToFiltration.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) { - if (globalDbg) { - std::cerr << "Boundary_simplex_range boundary_simplex_range(Simplex_handle& sh)\n"; - } - // Returns a range giving access to all simplices of the boundary of a simplex, i.e. the set of - // codimension 1 subsimplices of the Simplex. - return Boundary_simplex_range(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(Pawel Dlotko): 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< size_t > bdry = this->get_boundary_of_a_cell(sh.position); - if (globalDbg) { - std::cerr << "std::pair endpoints( Simplex_handle sh )\n"; - std::cerr << "bdry.size() : " << bdry.size() << std::endl; - } - // 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 for which this method was called have less than two elements in the boundary."); - return std::make_pair(Simplex_handle(this, bdry[0]), Simplex_handle(this, bdry[1])); - } - - - /** - * Class needed for compatibility with Gudhi. Not useful for other purposes. - **/ - class Skeleton_simplex_range; - - class Skeleton_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > { - // 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, size_t d) : b(b), dimension(d) { - if (globalDbg) { - std::cerr << "Skeleton_simplex_iterator ( Bitmap_cubical_complex* b , 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), 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; - } - - bool operator==(const Skeleton_simplex_iterator& rhs) { - if (globalDbg) { - std::cerr << "bool operator ==\n"; - } - if (this->position == rhs.position) { - return true; - } - return false; - } - - bool operator!=(const Skeleton_simplex_iterator& rhs) { - 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 Simplex_handle(this->b, this->position); - } - - friend class Skeleton_simplex_range; - private: - Bitmap_cubical_complex* b; - size_t position; - int dimension; - }; - - /** - * 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: - Skeleton_simplex_range(Bitmap_cubical_complex* b, int 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; - int dimension; - }; - - /** - * Function needed for compatibility with Gudhi. Not useful for other purposes. - **/ - Skeleton_simplex_range skeleton_simplex_range(int dimension) { - if (globalDbg) { - std::cerr << "Skeleton_simplex_range skeleton_simplex_range( int dimension )\n"; - } - return Skeleton_simplex_range(this, dimension); - } - - - - //******************************************************************************************************************// - // functions used for debugging: - - /** - * Function used for debugging purposes. - **/ - void printKeyAssociatedToSimplex() { - for (size_t i = 0; i != this->data.size(); ++i) { - std::cerr << i << " -> " << this->simplexAssociatedToKey[i] << std::endl; - } - } - - /** - * Function used for debugging purposes. - **/ - size_t printRealPosition(const Simplex_handle& sh) { - return sh.position; - } - - private: - std::vector< size_t > keyAssociatedToSimplex; - std::vector< size_t > simplexAssociatedToKey; - // needed by Filtration_simplex_iterator. If this iterator is not used, this field is not initialized. - std::vector< size_t > elementsOrderedAccordingToFiltration; -}; - -template -bool compareElementsForElementsOrderedAccordingToFiltration(const std::pair< size_t, - std::pair< T, char > >& f, - const std::pair< size_t, - std::pair< T, char > >& s) { - if (globalDbg) { - std::cerr << "ompareElementsForElementsOrderedAccordingToFiltration\n"; - } - if (f.second.first < s.second.first) { - return true; - } else { - if (f.second.first > s.second.first) { - return false; - } else { - // in this case f.second.first == s.second.first, and we use dimension to compare: - if (f.second.second < s.second.second) { - return true; - } else { - if (f.second.second > s.second.second) { - return false; - } else { - // in this case, both the filtration value and the dimensions for those cells are the same. - // Since it may be nice to have a stable sorting procedure, in this case, we compare positions in the bitmap: - return ( f.first < s.first); - } - } - } - } -} - -template -void Bitmap_cubical_complex::initializeElementsOrderedAccordingToFiltration() { - if (globalDbg) { - std::cerr << "void Bitmap_cubical_complex::initializeElementsOrderedAccordingToFiltration() \n"; - } - // ( position , (filtration , dimension) ) - std::vector< std::pair< size_t, std::pair< T, char > > > dataOfElementsFromBitmap(this->data.size()); - for (size_t i = 0; i != this->data.size(); ++i) { - // TODO(Pawel Dlotko): This can be optimized by having a counter here. We do not need to re-compute the dimension - // for every cell from scratch - dataOfElementsFromBitmap[i] = std::make_pair(i, std::make_pair(this->data[i], this->get_dimension_of_a_cell(i))); - } - std::sort(dataOfElementsFromBitmap.begin(), dataOfElementsFromBitmap.end(), - compareElementsForElementsOrderedAccordingToFiltration); - - // Elements of bitmap ordered according to filtration then according to dimension then according to position in bitmap - std::vector< size_t > elements_of_bitmap_ordered(this->data.size()); - for (size_t i = 0; i != dataOfElementsFromBitmap.size(); ++i) { - elements_of_bitmap_ordered[i] = dataOfElementsFromBitmap[i].first; - } - this->elementsOrderedAccordingToFiltration = elements_of_bitmap_ordered; -} - - -//****************************************************************************************************************// -//****************************************************************************************************************// -//****************************************************************************************************************// -//****************************************************************************************************************// - -#endif // BITMAP_CUBICAL_COMPLEX_H_ + /* 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 Sophia-Saclay (France) + * + * 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 . + */ + + +#pragma once +#include +#include "Bitmap_cubical_complex_base.h" + + + +namespace Gudhi +{ + +namespace Cubical_complex +{ + +//global variable, was used just for debugging. +const bool globalDbg = false; + +template +class Bitmap_cubical_complex : public Bitmap_cubical_complex_base +{ +public: +//*********************************************// +//Typedefs and typenames +//*********************************************// + friend class Simplex_handle; + typedef size_t Simplex_key; + typedef T Filtration_value; + + +//*********************************************// +//Simplex handle class +//*********************************************// + /** + * Handle of a cell, required for compatibility with the function to compute persistence in Gudhi. + * Elements of this class are: the pointer to the bitmap B in which the considered cell is + * together with a position of this cell in B. Given this data, + * one can get all the information about the considered cell. + **/ + class Simplex_handle + { + public: + Simplex_handle() + { + if ( globalDbg ){cerr << "Simplex_handle()\n";} + this->b = 0; + this->position = 0; + } + + Simplex_handle(Bitmap_cubical_complex* b) + { + if ( globalDbg ) + { + cerr << "Simplex_handle(Bitmap_cubical_complex* b)\n"; + } + this->b = b; + this->position = 0; + } + + //Simplex_handle( const Simplex_handle& org ):b(org.b) + //{ + // if ( globalDbg ){cerr << "Simplex_handle( const Simplex_handle& org )\n";} + // this->position = org.position; + //} + + Simplex_handle operator = ( const Simplex_handle& rhs ) + { + if ( globalDbg ){cerr << "Simplex_handle operator = \n";} + this->position = rhs.position; + this->b = rhs.b; + return *this; + } + + Simplex_handle(Bitmap_cubical_complex* b , Simplex_key position) + { + if ( globalDbg ) + { + cerr << "Simplex_handle(Bitmap_cubical_complex* b , Simplex_key position)\n"; + cerr << "Position : " << position << endl; + } + this->b = b; + this->position = position; + } + friend class Bitmap_cubical_complex; + private: + Bitmap_cubical_complex* b; + Simplex_key position; + //Assumption -- field above always keep the REAL position of simplex in the bitmap, + //no matter what keys have been. + //to deal with the keys, the class Bitmap_cubical_complex have extra vectors: key_associated_to_simplex and + //simplex_associated_to_key that allow to move between actual cell and the key assigned to it. + }; + + +//*********************************************// +//Constructors +//*********************************************// + //Over here we need to definie various input types. I am proposing the following ones: + //Perseus style + //H5 files? TODO + //binary files with little endiangs / big endians? TODO + //constructor from a vector of elements of a type T. TODO + + /** + * Constructor form a Perseus-style file. + **/ + Bitmap_cubical_complex( const char* perseus_style_file ): + Bitmap_cubical_complex_base(perseus_style_file),key_associated_to_simplex(this->total_number_of_cells+1), + simplex_associated_to_key(this->total_number_of_cells+1) + { + if ( globalDbg ){cerr << "Bitmap_cubical_complex( const char* perseus_style_file )\n";} + for ( size_t i = 0 ; i != this->total_number_of_cells ; ++i ) + { + this->key_associated_to_simplex[i] = this->simplex_associated_to_key[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_elements_ordered_according_to_filtration(); + } + + + /** + * 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( std::vector& dimensions , std::vector& top_dimensional_cells ): + Bitmap_cubical_complex_base(dimensions,top_dimensional_cells), + key_associated_to_simplex(this->total_number_of_cells+1), + simplex_associated_to_key(this->total_number_of_cells+1) + { + for ( size_t i = 0 ; i != this->total_number_of_cells ; ++i ) + { + this->key_associated_to_simplex[i] = this->simplex_associated_to_key[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_elements_ordered_according_to_filtration(); + } + +//*********************************************// +//Other 'easy' functions +//*********************************************// + /** + * Returns number of all cubes in the complex. + **/ + 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. + **/ + Simplex_handle null_simplex() + { + return Simplex_handle(this,this->data.size()); + } + + /** + * Returns dimension of the complex. + **/ + size_t dimension() + { + return this->sizes.size(); + } + + /** + * Return dimension of a cell pointed by the Simplex_handle. + **/ + unsigned dimension(const Simplex_handle& sh) + { + if ( globalDbg ){cerr << "unsigned dimension(const Simplex_handle& sh)\n";} + if ( sh.position != this->data.size() ) return sh.b->get_dimension_of_a_cell( sh.position ); + return -1; + } + + /** + * Return the filtration of a cell pointed by the Simplex_handle. + **/ + T filtration(const Simplex_handle& sh) + { + if ( globalDbg ){cerr << "T filtration(const Simplex_handle& sh)\n";} + //Returns the filtration value of a simplex. + if ( sh.position != this->data.size() ) return sh.b->data[ sh.position ]; + return std::numeric_limits::max(); + } + + /** + * Return a key which is not a key of any cube in the considered data structure. + **/ + Simplex_key null_key() + { + if ( globalDbg ){cerr << "Simplex_key null_key()\n";} + return this->data.size(); + } + + /** + * Return the key of a cube pointed by the Simplex_handle. + **/ + Simplex_key key(const Simplex_handle& sh) + { + if ( globalDbg ){cerr << "Simplex_key key(const Simplex_handle& sh)\n";} + return sh.b->key_associated_to_simplex[ sh.position ]; + } + + /** + * Return the Simplex_handle given the key of the cube. + **/ + Simplex_handle simplex(Simplex_key key) + { + if ( globalDbg ){cerr << "Simplex_handle simplex(Simplex_key key)\n";} + return Simplex_handle( this , this->simplex_associated_to_key[ key ] ); + } + + /** + * Assign key to a cube pointed by the Simplex_handle + **/ + void assign_key(Simplex_handle& sh, Simplex_key key) + { + if ( globalDbg ){cerr << "void assign_key(Simplex_handle& sh, Simplex_key key)\n";} + this->key_associated_to_simplex[sh.position] = key; + this->simplex_associated_to_key[key] = sh.position; + } + + /** + * Function called from a constructor. It is needed for Filtration_simplex_iterator to work. + **/ + void initialize_elements_ordered_according_to_filtration(); + + + +//*********************************************// +//Iterators +//*********************************************// + + /** + * Boundary_simplex_iterator class allows iteration on boundary of each cube. + **/ + class Boundary_simplex_range; + class Boundary_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > + { + //Iterator on the simplices belonging to the boundary of a simplex. + //value_type must be 'Simplex_handle'. + public: + Boundary_simplex_iterator( Simplex_handle& sh ):sh(sh) + { + if ( globalDbg ){cerr << "Boundary_simplex_iterator( Simplex_handle& sh )\n";} + this->position = 0; + this->boundary_elements = this->sh.b->get_boundary_of_a_cell( this->sh.position ); + } + Boundary_simplex_iterator operator++() + { + if ( globalDbg ){cerr << "Boundary_simplex_iterator operator++()\n";} + ++this->position; + return *this; + } + Boundary_simplex_iterator operator++(int) + { + Boundary_simplex_iterator result = *this; + ++(*this); + return result; + } + Boundary_simplex_iterator operator =( const Boundary_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "Boundary_simplex_iterator operator =\n";} + this->sh = rhs.sh; + this->boundary_elements.clear(); + this->boundary_elementsinsert + (this->boundary_elements.end(), rhs.boundary_elements.begin(), rhs.boundary_elements.end()); + } + bool operator == ( const Boundary_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator ==\n";} + if ( this->position == rhs.position ) + { + if ( this->boundary_elements.size() != rhs.boundary_elements.size() )return false; + for ( size_t i = 0 ; i != this->boundary_elements.size() ; ++i ) + { + if ( this->boundary_elements[i] != rhs.boundary_elements[i] )return false; + } + return true; + } + return false; + } + + bool operator != ( const Boundary_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator != \n";} + return !(*this == rhs); + } + Simplex_handle operator*() + { + if ( globalDbg ){cerr << "Simplex_handle operator*\n";} + return Simplex_handle( this->sh.b , this->boundary_elements[this->position] ); + } + + friend class Boundary_simplex_range; + private: + Simplex_handle sh; + std::vector< size_t > boundary_elements; + size_t position; + }; + + + /** + * Boundary_simplex_range class provides ranges for boundary iterators. + **/ + class Boundary_simplex_range + { + //Range giving access to the simplices in the boundary of a simplex. + //.begin() and .end() return type Boundary_simplex_iterator. + public: + Boundary_simplex_range(const Simplex_handle& sh):sh(sh){}; + Boundary_simplex_iterator begin() + { + if ( globalDbg ){cerr << "Boundary_simplex_iterator begin\n";} + Boundary_simplex_iterator it( this->sh ); + return it; + } + Boundary_simplex_iterator end() + { + if ( globalDbg ){cerr << "Boundary_simplex_iterator end()\n";} + Boundary_simplex_iterator it( this->sh ); + it.position = it.boundary_elements.size(); + return it; + } + private: + Simplex_handle sh; + }; + + + /** + * 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< std::input_iterator_tag, Simplex_handle > + { + //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){}; + + Filtration_simplex_iterator operator++() + { + if ( globalDbg ){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 ){cerr << "Filtration_simplex_iterator operator =\n";} + this->b = rhs.b; + this->position = rhs.position; + } + bool operator == ( const Filtration_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator == ( const Filtration_simplex_iterator& rhs )\n";} + if ( this->position == rhs.position ) + { + return true; + } + return false; + } + + bool operator != ( const Filtration_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator != ( const Filtration_simplex_iterator& rhs )\n";} + return !(*this == rhs); + } + Simplex_handle operator*() + { + if ( globalDbg ){cerr << "Simplex_handle operator*()\n";} + return Simplex_handle( this->b , this->b->elements_ordered_according_to_filtration[ this->position ] ); + } + + friend class Filtration_simplex_range; + private: + Bitmap_cubical_complex* b; + size_t position; + }; + + + /** + * 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: + Filtration_simplex_range(Bitmap_cubical_complex* b):b(b){}; + Filtration_simplex_iterator begin() + { + if ( globalDbg ){cerr << "Filtration_simplex_iterator begin() \n";} + return Filtration_simplex_iterator( this->b ); + } + Filtration_simplex_iterator end() + { + if ( globalDbg ){cerr << "Filtration_simplex_iterator end()\n";} + Filtration_simplex_iterator it( this->b ); + it.position = this->b->elements_ordered_according_to_filtration.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) + { + if ( globalDbg ){cerr << "Boundary_simplex_range boundary_simplex_range(Simplex_handle& sh)\n";} + //Returns a range giving access to all simplices of the boundary of a simplex, + //i.e. the set of codimension 1 subsimplices of the Simplex. + return Boundary_simplex_range(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 ){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 -- 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< size_t > bdry = this->get_boundary_of_a_cell( sh.position ); + if ( globalDbg ) + { + cerr << "std::pair endpoints( Simplex_handle sh )\n"; + cerr << "bdry.size() : " << bdry.size() << endl; + } + //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( Simplex_handle(this,bdry[0]) , Simplex_handle(this,bdry[1]) ); + } + + + /** + * Class needed for compatibility with Gudhi. Not useful for other purposes. + **/ + class Skeleton_simplex_range; + class Skeleton_simplex_iterator : std::iterator< std::input_iterator_tag, Simplex_handle > + { + //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 , size_t d ):b(b),dimension(d) + { + if ( globalDbg ){cerr << "Skeleton_simplex_iterator ( Bitmap_cubical_complex* b , 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),dimension(0){}; + + Skeleton_simplex_iterator operator++() + { + if ( globalDbg ){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 ){cerr << "Skeleton_simplex_iterator operator =\n";} + this->b = rhs.b; + this->position = rhs.position; + } + bool operator == ( const Skeleton_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator ==\n";} + if ( this->position == rhs.position ) + { + return true; + } + return false; + } + + bool operator != ( const Skeleton_simplex_iterator& rhs ) + { + if ( globalDbg ){cerr << "bool operator != ( const Skeleton_simplex_iterator& rhs )\n";} + return !(*this == rhs); + } + Simplex_handle operator*() + { + if ( globalDbg ){cerr << "Simplex_handle operator*() \n";} + return Simplex_handle( this->b , this->position ); + } + + friend class Skeleton_simplex_range; + private: + Bitmap_cubical_complex* b; + size_t position; + unsigned dimension; + }; + /** + * 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: + Skeleton_simplex_range(Bitmap_cubical_complex* b , unsigned dimension):b(b),dimension(dimension){}; + Skeleton_simplex_iterator begin() + { + if ( globalDbg ){cerr << "Skeleton_simplex_iterator begin()\n";} + return Skeleton_simplex_iterator( this->b , this->dimension ); + } + Skeleton_simplex_iterator end() + { + if ( globalDbg ){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 ){cerr << "Skeleton_simplex_range skeleton_simplex_range( unsigned dimension )\n";} + return Skeleton_simplex_range( this , dimension ); + } + + + +//*********************************************// +//functions used for debugging: + /** + * Function used for debugging purposes. + **/ + //void printkey_associated_to_simplex() + //{ + // for ( size_t i = 0 ; i != this->data.size() ; ++i ) + // { + // cerr << i << " -> " << this->simplex_associated_to_key[i] << endl; + // } + //} + + /** + * Function used for debugging purposes. + **/ + size_t printRealPosition( const Simplex_handle& sh ) + { + return sh.position; + } + +private: + std::vector< size_t > key_associated_to_simplex; + std::vector< size_t > simplex_associated_to_key; + std::vector< size_t > elements_ordered_according_to_filtration; + //filed above is needed by Filtration_simplex_iterator. If this iterator is not used, this field is not initialized. +};//Bitmap_cubical_complex + +template +bool compare_elements_for_elements_ordered_according_to_filtration +( const std::pair< size_t , std::pair< T , char > >& f , const std::pair< size_t , std::pair< T , char > >& s ) +{ + if ( globalDbg ){cerr << "compare_elements_for_elements_ordered_according_to_filtration\n";} + if ( f.second.first < s.second.first ) + { + return true; + } + else + { + if ( f.second.first > s.second.first ) + { + return false; + } + else + { + //in this case f.second.first == s.second.first, and we use dimension to compare: + if ( f.second.second < s.second.second ) + { + return true; + } + else + { + if ( f.second.second > s.second.second ) + { + return false; + } + else + { + //in this case, both the filtration value and the dimensions for those cells are the same. + //Since it may be nice to have a stable sorting procedure, in this case, + //we compare positions in the bitmap: + return ( f.first < s.first ); + } + } + } + } +} + +template +void Bitmap_cubical_complex::initialize_elements_ordered_according_to_filtration() +{ + if ( globalDbg ) + { + cerr << "void Bitmap_cubical_complex::initialize_elements_ordered_according_to_filtration() \n"; + } + //( position , (filtration , dimension) ) + std::vector< std::pair< size_t , std::pair< T , char > > > data_of_elements_from_bitmap( this->data.size() ); + for ( size_t i = 0 ; i != this->data.size() ; ++i ) + { + //TODO -- this can be optimized by having a counter here. + //We do not need to re-compute the dimension for every cell from scratch + data_of_elements_from_bitmap[i] = + std::make_pair( i , std::make_pair( this->data[i] , this->get_dimension_of_a_cell(i) ) ); + } + std::sort( data_of_elements_from_bitmap.begin() , + data_of_elements_from_bitmap.end() , + compare_elements_for_elements_ordered_according_to_filtration ); + + std::vector< size_t > + elements_ordered_according_to_filtration_then_to_dimension_then_to_position + ( this->data.size() ); + for ( size_t i = 0 ; i != data_of_elements_from_bitmap.size() ; ++i ) + { + elements_ordered_according_to_filtration_then_to_dimension_then_to_position[i] + = data_of_elements_from_bitmap[i].first; + } + this->elements_ordered_according_to_filtration = + elements_ordered_according_to_filtration_then_to_dimension_then_to_position; +} + + +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// + + +} +} \ No newline at end of file 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 d9c91832..2c2bd481 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 @@ -1,593 +1,759 @@ -/* 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 Sophia-Saclay (France) - * - * 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_BASE_H_ -#define BITMAP_CUBICAL_COMPLEX_BASE_H_ - -#include - -#include -#include -#include -#include -#include -#include -#include -#include - -/** - * This is a class implementing a basic bitmap data structure to store cubical complexes. It implements only the most basic subroutines. - * The idea of the bitmap is the following. Our aim is to have a memory efficient data structure to store d-dimensional cubical complex C being a cubical decomposition - * of a rectangular region of a space. This is achieved by storing C as a vector of bits (this is where the name 'bitmap' came from). Each cell is represented by a single - * bit (in case of black and white bitmaps, or by a single element of a type T (here T is a filtration type of a bitmap, typically a double). All the informations needed for homology and - * persistent homology computations (like dimension of a cell, boundary and coboundary elements of a cell, are then obtained from the position of the element in C. - */ -template -class Bitmap_cubical_complex_base { - public: - /** - * There are a few constructors of a Bitmap_cubical_complex_base class. First one, that takes vector, creates an empty bitmap of a dimension equal to the number of elements in the - * input vector and size in the i-th dimension equal the number in the position i-of the input vector. - */ - Bitmap_cubical_complex_base(std::vector sizes_); - /** - * The second constructor takes as a input a Perseus style file. For more details, please consult the documentations of Perseus software as well as examples attached to this - * implementation. - **/ - Bitmap_cubical_complex_base(char* perseusStyleFile_); - /** - * The last constructor of a Bitmap_cubical_complex_base class accepts vector of dimensions (as the first one) together with vector of filtration values of top dimensional cells. - **/ - Bitmap_cubical_complex_base(std::vector dimensions_, std::vector topDimensionalCells_); - - /** - * The functions get_boundary_of_a_cell, get_coboundary_of_a_cell and get_cell_data are the basic functions that compute boundary / coboundary / dimension and the filtration - * value form a position of a cell in the structure of a bitmap. The input parameter of all of those function is a 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. - */ - inline std::vector< size_t > get_boundary_of_a_cell(size_t cell_); - /** - * The functions get_boundary_of_a_cell, get_coboundary_of_a_cell, get_dimension_of_a_cell and get_cell_data are the basic functions that compute boundary / coboundary / dimension and the filtration - * value form a position of a cell in the structure of a bitmap. The input parameter of all of those function is a 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. - **/ - inline std::vector< size_t > get_coboundary_of_a_cell(size_t cell_); - /** - * In the case of get_dimension_of_a_cell function, the output is a non-negative integer indicating the dimension of a cell. - **/ - inline unsigned get_dimension_of_a_cell(size_t cell_); - /** - * In the case of get_cell_data, the output parameter is a reference to the value of a cube in a given position. - **/ - inline T& get_cell_data(size_t cell_); - - - /** - * Typical input used to construct a baseBitmap class is a filtration given at the top dimensional cells. Then, there are a few ways one can pick the filtration of lower dimensional - * cells. The most typical one is by so called lower star filtration. This function is always called by any constructor which takes the top dimensional cells. If you use such a constructor, - * then there is no need to call this function. Call it only if you are putting the filtration of the cells by your own (for instance by using topDimensionalCellsIterator). - **/ - void impose_lower_star_filtration(); // assume that top dimensional cells are already set. - - /** - * Returns dimension of a complex. - **/ - inline unsigned dimension() { - return sizes.size(); - } - - /** - * Returns number of all cubes in the data structure. - **/ - inline unsigned size_of_bitmap() { - return this->data.size(); - } - - /** - * Writing to stream operator. - **/ - template - friend std::ostream& operator<<(std::ostream & os_, const Bitmap_cubical_complex_base& b_); - - // ITERATORS - - /** - * Iterator through all cells in the complex (in order they appear in the structure -- i.e. in lexicographical order). - **/ - typedef typename std::vector< T >::iterator all_cells_iterator; - - all_cells_iterator all_cells_begin()const { - return this->data.begin(); - } - - all_cells_iterator all_cells_end()const { - return this->data.end(); - } - - - typedef typename std::vector< T >::const_iterator all_cells_const_iterator; - - all_cells_const_iterator all_cells_const_begin()const { - return this->data.begin(); - } - - all_cells_const_iterator all_cells_const_end()const { - return this->data.end(); - } - - /** - * Iterator through top dimensional cells of the complex. The cells appear in order they are stored in the structure (i.e. in lexicographical order) - **/ - class Top_dimensional_cells_iterator : std::iterator< std::input_iterator_tag, double > { - public: - Top_dimensional_cells_iterator(Bitmap_cubical_complex_base& b_) : b(b_) { - for (size_t i = 0; i != b_.dimension(); ++i) { - this->counter.push_back(0); - } - } - - Top_dimensional_cells_iterator operator++() { - // first find first element of the counter that can be increased: - size_t dim = 0; - while ((dim != this->b.dimension()) && (this->counter[dim] == this->b.sizes[dim] - 1))++dim; - - if (dim != this->b.dimension()) { - ++this->counter[dim]; - for (size_t i = 0; i != dim; ++i) { - this->counter[i] = 0; - } - } else { - ++this->counter[0]; - } - return *this; - } - - Top_dimensional_cells_iterator operator++(int) { - Top_dimensional_cells_iterator result = *this; - ++(*this); - return result; - } - - Top_dimensional_cells_iterator operator=(const Top_dimensional_cells_iterator& rhs_) { - this->counter = rhs_.counter; - this->b = rhs_.b; - return *this; - } - - bool operator==(const Top_dimensional_cells_iterator& rhs_) { - if (&this->b != &rhs_.b)return false; - if (this->counter.size() != rhs_.counter.size())return false; - for (size_t i = 0; i != this->counter.size(); ++i) { - if (this->counter[i] != rhs_.counter[i])return false; - } - return true; - } - - bool operator!=(const Top_dimensional_cells_iterator& rhs_) { - return !(*this == rhs_); - } - - T& operator*() { - // given the counter, compute the index in the array and return this element. - unsigned index = 0; - for (size_t i = 0; i != this->counter.size(); ++i) { - index += (2 * this->counter[i] + 1) * this->b.multipliers[i]; - } - return this->b.data[index]; - } - - size_t computeIndexInBitmap() { - size_t index = 0; - for (size_t i = 0; i != this->counter.size(); ++i) { - index += (2 * this->counter[i] + 1) * this->b.multipliers[i]; - } - return index; - } - - void printCounter() { - for (size_t i = 0; i != this->counter.size(); ++i) { - std::cout << this->counter[i] << " "; - } - } - friend class Bitmap_cubical_complex_base; - protected: - std::vector< unsigned > counter; - Bitmap_cubical_complex_base& b; - }; - - Top_dimensional_cells_iterator top_dimensional_cells_begin() { - Top_dimensional_cells_iterator a(*this); - return a; - } - - Top_dimensional_cells_iterator top_dimensional_cells_end() { - Top_dimensional_cells_iterator a(*this); - for (size_t i = 0; i != this->dimension(); ++i) { - a.counter[i] = this->sizes[i] - 1; - } - a.counter[0]++; - return a; - } - - - //****************************************************************************************************************// - //****************************************************************************************************************// - //****************************************************************************************************************// - //****************************************************************************************************************// - - - //****************************************************************************************************************// - //****************************************************************************************************************// - //****************************************************************************************************************// - //****************************************************************************************************************// - - protected: - std::vector sizes; - std::vector multipliers; - std::vector data; - size_t totalNumberOfCells; - - void set_up_containers(std::vector sizes_) { - unsigned multiplier = 1; - for (size_t i = 0; i != sizes_.size(); ++i) { - this->sizes.push_back(sizes_[i]); - this->multipliers.push_back(multiplier); - // multiplier *= 2*(sizes[i]+1)+1; - multiplier *= 2 * sizes_[i] + 1; - } - // std::reverse( this->sizes.begin() , this->sizes.end() ); - std::vector data(multiplier); - std::fill(data.begin(), data.end(), INT_MAX); - this->totalNumberOfCells = multiplier; - this->data = data; - } - - size_t compute_position_in_bitmap(std::vector< int > counter_) { - size_t position = 0; - for (size_t i = 0; i != this->multipliers.size(); ++i) { - position += this->multipliers[i] * counter_[i]; - } - return position; - } - - std::vector compute_counter_for_given_cell(size_t cell_) { - std::vector counter; - for (size_t dim = this->sizes.size(); dim != 0; --dim) { - counter.push_back(cell_ / this->multipliers[dim - 1]); - cell_ = cell_ % this->multipliers[dim - 1]; - } - std::reverse(counter.begin(), counter.end()); - return counter; - } - - std::vector< size_t > generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions(std::vector< bool > directionsForPeriodicBCond_); -}; - -template -std::ostream& operator<<(std::ostream & out_, const Bitmap_cubical_complex_base& b_) { - // for ( typename bitmap::all_cells_const_iterator it = b.all_cells_const_begin() ; - // it != b.all_cells_const_end() ; ++it ) - for (typename Bitmap_cubical_complex_base::all_cells_const_iterator it = b_.all_cells_const_begin(); - it != b_.all_cells_const_end(); ++it) { - out_ << *it << " "; - } - return out_; -} - -template -Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(std::vector sizes_) { - this->set_up_containers(sizes_); -} - -template -Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(std::vector sizesInFollowingDirections_, - std::vector topDimensionalCells_) { - this->set_up_containers(sizesInFollowingDirections_); - - size_t numberOfTopDimensionalElements = 1; - for (size_t i = 0; i != sizesInFollowingDirections_.size(); ++i) { - numberOfTopDimensionalElements *= sizesInFollowingDirections_[i]; - } - if (numberOfTopDimensionalElements != topDimensionalCells_.size()) { - std::cerr << "Error in constructor Bitmap_cubical_complex_base( std::vector sizesInFollowingDirections_ , " - "std::vector topDimensionalCells_ ). Number of top dimensional elements that follow from " - "sizesInFollowingDirections vector is different than the size of topDimensionalCells vector." << std::endl; - throw("Error in constructor Bitmap_cubical_complex_base( std::vector sizesInFollowingDirections_ , " - "std::vector topDimensionalCells_ ). Number of top dimensional elements that follow from " - "sizesInFollowingDirections vector is different than the size of topDimensionalCells vector."); - } - - Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); - size_t index = 0; - for (it = this->top_dimensional_cells_begin(); it != this->top_dimensional_cells_end(); ++it) { - (*it) = topDimensionalCells_[index]; - ++index; - } - this->impose_lower_star_filtration(); -} - -template -Bitmap_cubical_complex_base::Bitmap_cubical_complex_base(char* perseusStyleFile_) { - bool dbg = false; - std::ifstream inFiltration, inIds; - inFiltration.open(perseusStyleFile_); - unsigned dimensionOfData; - inFiltration >> dimensionOfData; - - if (dbg) { - std::cerr << "dimensionOfData : " << dimensionOfData << std::endl; - } - - std::vector sizes; - for (size_t i = 0; i != dimensionOfData; ++i) { - int sizeInThisDimension; - inFiltration >> sizeInThisDimension; - sizeInThisDimension = abs(sizeInThisDimension); - sizes.push_back(sizeInThisDimension); - if (dbg) { - std::cerr << "sizeInThisDimension : " << sizeInThisDimension << std::endl; - } - } - this->set_up_containers(sizes); - - Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); - it = this->top_dimensional_cells_begin(); - - // TODO(Pawel Dlotko): Over here we also need to read id's of cell and put them to bitmapElement structure! - while (!inFiltration.eof()) { - double filtrationLevel; - inFiltration >> filtrationLevel; - if (dbg) { - std::cerr << "Cell of an index : " << it.computeIndexInBitmap() << " and dimension: " << - this->get_dimension_of_a_cell(it.computeIndexInBitmap()) << " get the value : " << - filtrationLevel << std::endl; - } - *it = filtrationLevel; - ++it; - } - inFiltration.close(); - this->impose_lower_star_filtration(); -} - -template -std::vector< size_t > Bitmap_cubical_complex_base::get_boundary_of_a_cell(size_t cell_) { - bool bdg = false; - // First of all, we need to take the list of coordinates in which the cell has nonzero length. - // We do it by using modified version to compute dimension of a cell: - std::vector< unsigned > dimensionsInWhichCellHasNonzeroLength; - unsigned dimension = 0; - size_t cell1 = cell_; - for (size_t i = this->multipliers.size(); i != 0; --i) { - unsigned position = cell1 / multipliers[i - 1]; - if (position % 2 == 1) { - dimensionsInWhichCellHasNonzeroLength.push_back(i - 1); - dimension++; - } - cell1 = cell1 % multipliers[i - 1]; - } - - if (bdg) { - std::cerr << "dimensionsInWhichCellHasNonzeroLength : \n"; - for (size_t i = 0; i != dimensionsInWhichCellHasNonzeroLength.size(); ++i) { - std::cerr << dimensionsInWhichCellHasNonzeroLength[i] << std::endl; - } - getchar(); - } - - std::vector< size_t > boundaryElements; - if (dimensionsInWhichCellHasNonzeroLength.size() == 0)return boundaryElements; - for (size_t i = 0; i != dimensionsInWhichCellHasNonzeroLength.size(); ++i) { - boundaryElements.push_back(cell_ - multipliers[ dimensionsInWhichCellHasNonzeroLength[i] ]); - boundaryElements.push_back(cell_ + multipliers[ dimensionsInWhichCellHasNonzeroLength[i] ]); - - if (bdg) std::cerr << "multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << - multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << std::endl; - if (bdg) std::cerr << "cell_ - multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << - cell_ - multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << std::endl; - if (bdg) std::cerr << "cell_ + multipliers[dimensionsInWhichCellHasNonzeroLength[i]] : " << - cell_ + multipliers[dimensionsInWhichCellHasNonzeroLength[i]] << std::endl; - } - return boundaryElements; -} - -template -std::vector< size_t > Bitmap_cubical_complex_base::get_coboundary_of_a_cell(size_t cell_) { - bool bdg = false; - // First of all, we need to take the list of coordinates in which the cell has nonzero length. - // We do it by using modified version to compute dimension of a cell: - std::vector< unsigned > dimensionsInWhichCellHasZeroLength; - unsigned dimension = 0; - size_t cell1 = cell_; - for (size_t i = this->multipliers.size(); i != 0; --i) { - unsigned position = cell1 / multipliers[i - 1]; - if (position % 2 == 0) { - dimensionsInWhichCellHasZeroLength.push_back(i - 1); - dimension++; - } - cell1 = cell1 % multipliers[i - 1]; - } - - std::vector counter = this->compute_counter_for_given_cell(cell_); - // reverse(counter.begin() , counter.end()); - - if (bdg) { - std::cerr << "dimensionsInWhichCellHasZeroLength : \n"; - for (size_t i = 0; i != dimensionsInWhichCellHasZeroLength.size(); ++i) { - std::cerr << dimensionsInWhichCellHasZeroLength[i] << std::endl; - } - std::cerr << "\n counter : " << std::endl; - for (size_t i = 0; i != counter.size(); ++i) { - std::cerr << counter[i] << std::endl; - } - getchar(); - } - - std::vector< size_t > coboundaryElements; - if (dimensionsInWhichCellHasZeroLength.size() == 0)return coboundaryElements; - for (size_t i = 0; i != dimensionsInWhichCellHasZeroLength.size(); ++i) { - if (bdg) { - std::cerr << "Dimension : " << i << std::endl; - if (counter[dimensionsInWhichCellHasZeroLength[i]] == 0) { - std::cerr << "In dimension : " << i << - " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; - } - if (counter[dimensionsInWhichCellHasZeroLength[i]] == 2 * this->sizes[dimensionsInWhichCellHasZeroLength[i]]) { - std::cerr << "In dimension : " << i << - " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; - } - } - - - if ((cell_ > multipliers[dimensionsInWhichCellHasZeroLength[i]]) && - (counter[dimensionsInWhichCellHasZeroLength[i]] != 0)) { - // if ( counter[dimensionsInWhichCellHasZeroLength[i]] != 0 ) - if (bdg)std::cerr << "Subtracting : " << cell_ - multipliers[dimensionsInWhichCellHasZeroLength[i]] << std::endl; - coboundaryElements.push_back(cell_ - multipliers[dimensionsInWhichCellHasZeroLength[i]]); - } - if ((cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]] < this->data.size()) && - (counter[dimensionsInWhichCellHasZeroLength[i]] != 2 * this->sizes[dimensionsInWhichCellHasZeroLength[i]])) { - // if ( counter[dimensionsInWhichCellHasZeroLength[i]] != 2*this->sizes[dimensionsInWhichCellHasZeroLength[i]] ) - coboundaryElements.push_back(cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]]); - if (bdg)std::cerr << "Adding : " << cell_ + multipliers[dimensionsInWhichCellHasZeroLength[i]] << std::endl; - } - } - return coboundaryElements; -} - -template -unsigned Bitmap_cubical_complex_base::get_dimension_of_a_cell(size_t cell_) { - bool dbg = false; - if (dbg)std::cerr << "\n\n\n Computing position o a cell of an index : " << cell_ << std::endl; - unsigned dimension = 0; - for (size_t i = this->multipliers.size(); i != 0; --i) { - unsigned position = cell_ / multipliers[i - 1]; - - if (dbg)std::cerr << "i-1 :" << i - 1 << std::endl; - if (dbg)std::cerr << "cell_ : " << cell_ << std::endl; - if (dbg)std::cerr << "position : " << position << std::endl; - if (dbg)std::cerr << "multipliers[" << i - 1 << "] = " << multipliers[i - 1] << std::endl; - if (dbg)getchar(); - - if (position % 2 == 1) { - if (dbg)std::cerr << "Nonzero length in this direction \n"; - dimension++; - } - cell_ = cell_ % multipliers[i - 1]; - } - return dimension; -} - -template -T& Bitmap_cubical_complex_base::get_cell_data(size_t cell_) { - return this->data[cell_]; -} - -template -void Bitmap_cubical_complex_base::impose_lower_star_filtration() { - bool dbg = false; - - // this vector will be used to check which elements have already been taken care of in imposing lower star filtration: - std::vector isThisCellConsidered(this->data.size(), false); - - std::vector indicesToConsider; - // we assume here that we already have a filtration on the top dimensional cells - // and we have to extend it to lower ones. - typename Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); - for (it = this->top_dimensional_cells_begin(); it != this->top_dimensional_cells_end(); ++it) { - indicesToConsider.push_back(it.computeIndexInBitmap()); - } - - while (indicesToConsider.size()) { - if (dbg) { - std::cerr << "indicesToConsider in this iteration \n"; - for (size_t i = 0; i != indicesToConsider.size(); ++i) { - std::cout << indicesToConsider[i] << " "; - } - getchar(); - } - std::vector newIndicesToConsider; - for (size_t i = 0; i != indicesToConsider.size(); ++i) { - std::vector bd = this->get_boundary_of_a_cell(indicesToConsider[i]); - for (size_t boundaryIt = 0; boundaryIt != bd.size(); ++boundaryIt) { - if (this->data[ bd[boundaryIt] ] > this->data[ indicesToConsider[i] ]) { - this->data[ bd[boundaryIt] ] = this->data[ indicesToConsider[i] ]; - } - if (isThisCellConsidered[ bd[boundaryIt] ] == false) { - newIndicesToConsider.push_back(bd[boundaryIt]); - isThisCellConsidered[ bd[boundaryIt] ] = true; - } - } - } - indicesToConsider.swap(newIndicesToConsider); - } -} - -template -std::vector< size_t > -Bitmap_cubical_complex_base::generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions(std::vector< bool > directionsForPeriodicBCond_) { - bool dbg = false; - if (this->sizes.size() != directionsForPeriodicBCond_.size()) - throw ("directionsForPeriodicBCond_ vector size is different from the size of the bitmap. " - "The program will now terminate \n"); - - std::vector sizes(this->sizes.size()); - for (size_t i = 0; i != this->sizes.size(); ++i)sizes[i] = 2 * this->sizes[i]; - - counter c(sizes); - - std::vector< size_t > result; - - for (size_t i = 0; i != this->data.size(); ++i) { - size_t position; - if (!c.isFinal()) { - position = i; - // result.push_back( i ); - } else { - std::vector< bool > finals = c.directionsOfFinals(); - bool jumpInPosition = false; - for (size_t dir = 0; dir != finals.size(); ++dir) { - if (finals[dir] == false)continue; - if (directionsForPeriodicBCond_[dir]) { - jumpInPosition = true; - } - } - if (jumpInPosition == true) { - // in this case this guy is final, so we need to find 'the opposite one' - position = compute_position_in_bitmap(c.findOpposite(directionsForPeriodicBCond_)); - } else { - position = i; - } - } - result.push_back(position); - if (dbg) { - std::cerr << " position : " << position << std::endl; - std::cerr << c << std::endl; - getchar(); - } - - c.increment(); - } - - return result; -} - -#endif // BITMAP_CUBICAL_COMPLEX_BASE_H_ +/* 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 Sophia-Saclay (France) + * + * 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 . + */ + +#pragma once + +#include +#include +#include +#include +#include +#include +#include +#include +#include "counter.h" + + +using namespace std; + +namespace Gudhi +{ + +namespace Cubical_complex +{ + + + +/** + * This is a class implementing a basic bitmap data structure to store cubical complexes. + * It implements only the most basic subroutines. + * The idea of the bitmap is the following. Our aim is to have a memory efficient + * data structure to store d-dimensional cubical complex + * C being a cubical decomposition + * of a rectangular region of a space. This is achieved by storing C as a + * vector of bits (this is where the name 'bitmap' came from). + * Each cell is represented by a single + * bit (in case of black and white bitmaps, or by a single element of a type T + * (here T is a filtration type of a bitmap, typically a double). + * All the informations needed for homology and + * persistent homology computations (like dimension of a cell, boundary and + * coboundary elements of a cell, are then obtained from the + * position of the element in C. + * The default filtration used in this implementation is the lower star filtration. + */ +template +class Bitmap_cubical_complex_base +{ +public: + /** + * There are a few constructors of a Bitmap_cubical_complex_base class. + * First one, that takes vector, creates an empty bitmap of a dimension equal + * the number of elements in the + * input vector and size in the i-th dimension equal the number in the position i-of the input vector. + */ + Bitmap_cubical_complex_base( std::vector& sizes ); + /** + * The second constructor takes as a input a Perseus style file. For more details, + * please consult the documentations of + * Perseus software as well as examples attached to this + * implementation. + **/ + Bitmap_cubical_complex_base( const char* perseus_style_file ); + /** + * The last constructor of a Bitmap_cubical_complex_base class accepts vector of dimensions (as the first one) + * together with vector of filtration values of top dimensional cells. + **/ + Bitmap_cubical_complex_base( std::vector& dimensions , const std::vector& top_dimensional_cells ); + + /** + * The functions get_boundary_of_a_cell, get_coboundary_of_a_cell, get_dimension_of_a_cell + * and get_cell_data are the basic + * functions that compute boundary / coboundary / dimension and the filtration + * value form a position of a cell in the structure of a bitmap. The input parameter of all of those function is a + * 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. + */ + inline std::vector< size_t > get_boundary_of_a_cell( size_t cell )const; + /** + * The functions get_coboundary_of_a_cell, get_coboundary_of_a_cell, + * get_dimension_of_a_cell and get_cell_data are the basic + * functions that compute boundary / coboundary / dimension and the filtration + * value form a position of a cell in the structure of a bitmap. + * The input parameter of all of those function is a 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. + **/ + inline std::vector< size_t > get_coboundary_of_a_cell( size_t cell )const; + /** + * In the case of get_dimension_of_a_cell function, the output is a non-negative integer + * indicating the dimension of a cell. + **/ + 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. + **/ + inline T& get_cell_data( size_t cell ); + + + /** + * Typical input used to construct a baseBitmap class is a filtration given at the top dimensional cells. + * Then, there are a few ways one can pick the filtration of lower dimensional + * cells. The most typical one is by so called lower star filtration. This function is always called by any + * constructor which takes the top dimensional cells. If you use such a constructor, + * then there is no need to call this function. Call it only if you are putting the filtration + * of the cells by your own (for instance by using Top_dimensional_cells_iterator). + **/ + void impose_lower_star_filtration();//assume that top dimensional cells are already set. + + /** + * Returns dimension of a complex. + **/ + inline unsigned dimension()const{ return sizes.size(); } + + /** + * Returns number of all cubes in the data structure. + **/ + inline unsigned size_of_bitmap()const + { + return this->data.size(); + } + + /** + * Writing to stream operator. + **/ + template + friend ostream& operator << ( ostream & os , const Bitmap_cubical_complex_base& b ); + + //ITERATORS + + /** + * Iterator through all cells in the complex (in order they appear in the structure -- i.e. + * in lexicographical order). + **/ + typedef typename std::vector< T >::iterator all_cells_iterator; + all_cells_iterator all_cells_begin()const + { + return this->data.begin(); + } + all_cells_iterator all_cells_end()const + { + return this->data.end(); + } + + + typedef typename std::vector< T >::const_iterator all_cells_const_iterator; + all_cells_const_iterator all_cells_const_begin()const + { + return this->data.begin(); + } + all_cells_const_iterator all_cells_const_end()const + { + return this->data.end(); + } + + /** + * Iterator through top dimensional cells of the complex. The cells appear in order they are stored + * in the structure (i.e. in lexicographical order) + **/ + class Top_dimensional_cells_iterator : std::iterator< std::input_iterator_tag, double > + { + public: + Top_dimensional_cells_iterator( Bitmap_cubical_complex_base& b ):b(b) + { + for ( size_t i = 0 ; i != b.dimension() ; ++i ) + { + this->counter.push_back(0); + } + } + Top_dimensional_cells_iterator operator++() + { + //first find first element of the counter that can be increased: + size_t dim = 0; + while ( ( dim != this->b.dimension() ) && ( this->counter[dim] == this->b.sizes[dim]-1 ) )++dim; + + if ( dim != this->b.dimension() ) + { + ++this->counter[dim]; + for ( size_t i = 0 ; i != dim ; ++i ) + { + this->counter[i] = 0; + } + } + else + { + ++this->counter[0]; + } + return *this; + } + Top_dimensional_cells_iterator operator++(int) + { + Top_dimensional_cells_iterator result = *this; + ++(*this); + return result; + } + Top_dimensional_cells_iterator operator =( const Top_dimensional_cells_iterator& rhs ) + { + this->counter = rhs.counter; + this->b = rhs.b; + return *this; + } + bool operator == ( const Top_dimensional_cells_iterator& rhs )const + { + if ( &this->b != &rhs.b )return false; + if ( this->counter.size() != rhs.counter.size() )return false; + for ( size_t i = 0 ; i != this->counter.size() ; ++i ) + { + if ( this->counter[i] != rhs.counter[i] )return false; + } + return true; + } + bool operator != ( const Top_dimensional_cells_iterator& rhs )const + { + return !(*this == rhs); + } + + T& operator*() + { + //given the counter, compute the index in the array and return this element. + unsigned index = 0; + for ( size_t i = 0 ; i != this->counter.size() ; ++i ) + { + index += (2*this->counter[i]+1)*this->b.multipliers[i]; + } + return this->b.data[index]; + } + + size_t compute_index_in_bitmap()const + { + size_t index = 0; + for ( size_t i = 0 ; i != this->counter.size() ; ++i ) + { + index += (2*this->counter[i]+1)*this->b.multipliers[i]; + } + return index; + } + + void print_counter()const + { + for ( size_t i = 0 ; i != this->counter.size() ; ++i ) + { + cout << this->counter[i] << " "; + } + } + friend class Bitmap_cubical_complex_base; + protected: + std::vector< unsigned > counter; + Bitmap_cubical_complex_base& b; + }; + Top_dimensional_cells_iterator top_dimensional_cells_begin() + { + Top_dimensional_cells_iterator a(*this); + return a; + } + Top_dimensional_cells_iterator top_dimensional_cells_end() + { + Top_dimensional_cells_iterator a(*this); + for ( size_t i = 0 ; i != this->dimension() ; ++i ) + { + a.counter[i] = this->sizes[i]-1; + } + a.counter[0]++; + return a; + } + + +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// + + +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// +//****************************************************************************************************************// + +protected: + std::vector sizes; + std::vector multipliers; + std::vector data; + size_t total_number_of_cells; + void set_up_containers( std::vector& sizes ) + { + unsigned multiplier = 1; + for ( size_t i = 0 ; i != sizes.size() ; ++i ) + { + this->sizes.push_back(sizes[i]); + this->multipliers.push_back(multiplier); + //multiplier *= 2*(sizes[i]+1)+1; + multiplier *= 2*sizes[i]+1; + } + //std::reverse( this->sizes.begin() , this->sizes.end() ); + std::vector data(multiplier); + std::fill( data.begin() , data.end() , std::numeric_limits::max() ); + this->total_number_of_cells = multiplier; + this->data = data; + } + + size_t compute_position_in_bitmap( std::vector< unsigned >& counter ) + { + size_t position = 0; + for ( size_t i = 0 ; i != this->multipliers.size() ; ++i ) + { + position += this->multipliers[i]*counter[i]; + } + return position; + } + + std::vector compute_counter_for_given_cell( size_t cell )const + { + std::vector counter; + for ( size_t dim = this->sizes.size() ; dim != 0 ; --dim ) + { + counter.push_back(cell/this->multipliers[dim-1]); + cell = cell%this->multipliers[dim-1]; + } + std::reverse( counter.begin() , counter.end() ); + return counter; + } + + std::vector< size_t > + generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions + ( std::vector< bool >& directions_for_periodic_b_cond ); +}; + + + + +template +ostream& operator << ( ostream & out , const Bitmap_cubical_complex_base& b ) +{ + for ( typename Bitmap_cubical_complex_base::all_cells_const_iterator + it = b.all_cells_const_begin() ; it != b.all_cells_const_end() ; ++it ) + { + out << *it << " "; + } + return out; +} + + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base +( std::vector& sizes ) +{ + this->set_up_containers( sizes ); +} + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base +( std::vector& sizes_in_following_directions , const std::vector& top_dimensional_cells ) +{ + this->set_up_containers( sizes_in_following_directions ); + + size_t number_of_top_dimensional_elements = 1; + for ( size_t i = 0 ; i != sizes_in_following_directions.size() ; ++i ) + { + number_of_top_dimensional_elements *= sizes_in_following_directions[i]; + } + if ( number_of_top_dimensional_elements != top_dimensional_cells.size() ) + { + cerr << + "Error in constructor\ + Bitmap_cubical_complex_base\ + ( std::vector sizes_in_following_directions , std::vector top_dimensional_cells ).\ + Number of top dimensional elements that follow from sizes_in_following_directions vector is different\ + than the size of top_dimensional_cells vector." << endl; + throw("Error in constructor Bitmap_cubical_complex_base( std::vector sizes_in_following_directions,\ + std::vector top_dimensional_cells )\ + . Number of top dimensional elements that follow from sizes_in_following_directions vector is different than the\ + size of top_dimensional_cells vector."); + } + + Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + size_t index = 0; + for ( it = this->top_dimensional_cells_begin() ; it != this->top_dimensional_cells_end() ; ++it ) + { + (*it) = top_dimensional_cells[index]; + ++index; + } + this->impose_lower_star_filtration(); +} + + +template +Bitmap_cubical_complex_base::Bitmap_cubical_complex_base( const char* perseus_style_file ) +{ + bool dbg = false; + ifstream inFiltration, inIds; + inFiltration.open( perseus_style_file ); + unsigned dimensionOfData; + inFiltration >> dimensionOfData; + + if (dbg){cerr << "dimensionOfData : " << dimensionOfData << endl;} + + std::vector sizes; + for ( size_t i = 0 ; i != dimensionOfData ; ++i ) + { + unsigned size_in_this_dimension; + inFiltration >> size_in_this_dimension; + size_in_this_dimension = abs(size_in_this_dimension); + sizes.push_back( size_in_this_dimension ); + if (dbg){cerr << "size_in_this_dimension : " << size_in_this_dimension << endl;} + } + this->set_up_containers( sizes ); + + Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + it = this->top_dimensional_cells_begin(); + + //TODO -- over here we also need to read id's of cell and put them to bitmapElement structure! + while ( !inFiltration.eof() ) + { + double filtrationLevel; + inFiltration >> filtrationLevel; + if ( dbg ) + { + cerr << "Cell of an index : " + << it.compute_index_in_bitmap() + << " and dimension: " + << this->get_dimension_of_a_cell(it.compute_index_in_bitmap()) + << " get the value : " << filtrationLevel << endl; + } + *it = filtrationLevel; + ++it; + } + inFiltration.close(); + this->impose_lower_star_filtration(); +} + + +template +std::vector< size_t > Bitmap_cubical_complex_base::get_boundary_of_a_cell( size_t cell )const +{ + bool bdg = false; + //first of all, we need to take the list of coordinates in which the cell has nonzero length. + //We do it by using modified version to compute dimension of a cell: + std::vector< unsigned > dimensions_in_which_cell_has_nonzero_length; + unsigned dimension = 0; + size_t cell1 = cell; + for ( size_t i = this->multipliers.size() ; i != 0 ; --i ) + { + unsigned position = cell1/multipliers[i-1]; + if ( position%2 == 1 ) + { + dimensions_in_which_cell_has_nonzero_length.push_back(i-1); + dimension++; + } + cell1 = cell1%multipliers[i-1]; + } + + if (bdg) + { + cerr << "dimensions_in_which_cell_has_nonzero_length : \n"; + for ( size_t i = 0 ; i != dimensions_in_which_cell_has_nonzero_length.size() ; ++i ) + { + cerr << dimensions_in_which_cell_has_nonzero_length[i] << endl; + } + getchar(); + } + + std::vector< size_t > boundary_elements; + if ( dimensions_in_which_cell_has_nonzero_length.size() == 0 )return boundary_elements; + for ( size_t i = 0 ; i != dimensions_in_which_cell_has_nonzero_length.size() ; ++i ) + { + boundary_elements.push_back( cell - multipliers[ dimensions_in_which_cell_has_nonzero_length[i] ] ); + boundary_elements.push_back( cell + multipliers[ dimensions_in_which_cell_has_nonzero_length[i] ] ); + + if (bdg) cerr << "multipliers[dimensions_in_which_cell_has_nonzero_length[i]] : " + << multipliers[dimensions_in_which_cell_has_nonzero_length[i]] << endl; + if (bdg) cerr << "cell - multipliers[dimensions_in_which_cell_has_nonzero_length[i]] : " + << cell - multipliers[dimensions_in_which_cell_has_nonzero_length[i]] << endl; + if (bdg) cerr << "cell + multipliers[dimensions_in_which_cell_has_nonzero_length[i]] : " + << cell + multipliers[dimensions_in_which_cell_has_nonzero_length[i]] << endl; + } + return boundary_elements; +} + + + + +template +std::vector< size_t > Bitmap_cubical_complex_base::get_coboundary_of_a_cell( size_t cell )const +{ + bool bdg = false; + //first of all, we need to take the list of coordinates in which the cell has nonzero length. + //We do it by using modified version to compute dimension of a cell: + std::vector< unsigned > dimensions_in_which_cell_has_zero_length; + unsigned dimension = 0; + size_t cell1 = cell; + for ( size_t i = this->multipliers.size() ; i != 0 ; --i ) + { + unsigned position = cell1/multipliers[i-1]; + if ( position%2 == 0 ) + { + dimensions_in_which_cell_has_zero_length.push_back(i-1); + dimension++; + } + cell1 = cell1%multipliers[i-1]; + } + + std::vector counter = this->compute_counter_for_given_cell( cell ); + //reverse(counter.begin() , counter.end()); + + if (bdg) + { + cerr << "dimensions_in_which_cell_has_zero_length : \n"; + for ( size_t i = 0 ; i != dimensions_in_which_cell_has_zero_length.size() ; ++i ) + { + cerr << dimensions_in_which_cell_has_zero_length[i] << endl; + } + cerr << "\n counter : " << endl; + for ( size_t i = 0 ; i != counter.size() ; ++i ) + { + cerr << counter[i] << endl; + } + getchar(); + } + + std::vector< size_t > coboundary_elements; + if ( dimensions_in_which_cell_has_zero_length.size() == 0 )return coboundary_elements; + for ( size_t i = 0 ; i != dimensions_in_which_cell_has_zero_length.size() ; ++i ) + { + if ( bdg ) + { + cerr << "Dimension : " << i << endl; + if (counter[dimensions_in_which_cell_has_zero_length[i]] == 0) + { + cerr << "In dimension : " << i + << " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; + } + if ( counter[dimensions_in_which_cell_has_zero_length[i]] + == + 2*this->sizes[dimensions_in_which_cell_has_zero_length[i]] ) + { + cerr << "In dimension : " << i + << " we cannot substract, since we will jump out of a Bitmap_cubical_complex_base \n"; + } + } + + + if ( (cell > multipliers[dimensions_in_which_cell_has_zero_length[i]]) + && (counter[dimensions_in_which_cell_has_zero_length[i]] != 0) ) + //if ( counter[dimensions_in_which_cell_has_zero_length[i]] != 0 ) + { + if ( bdg ) + { + cerr << "Subtracting : " << cell - multipliers[dimensions_in_which_cell_has_zero_length[i]] << endl; + } + coboundary_elements.push_back( cell - multipliers[dimensions_in_which_cell_has_zero_length[i]] ); + } + if ( + (cell + multipliers[dimensions_in_which_cell_has_zero_length[i]] < this->data.size()) && + (counter[dimensions_in_which_cell_has_zero_length[i]] + != + 2*this->sizes[dimensions_in_which_cell_has_zero_length[i]]) + ) + //if ( counter[dimensions_in_which_cell_has_zero_length[i]] != + //2*this->sizes[dimensions_in_which_cell_has_zero_length[i]] ) + { + coboundary_elements.push_back( cell + multipliers[dimensions_in_which_cell_has_zero_length[i]] ); + if ( bdg )cerr << "Adding : " << cell + multipliers[dimensions_in_which_cell_has_zero_length[i]] << endl; + } + } + return coboundary_elements; +} + + + + + + +template +unsigned Bitmap_cubical_complex_base::get_dimension_of_a_cell( size_t cell )const +{ + bool dbg = false; + if (dbg)cerr << "\n\n\n Computing position o a cell of an index : " << cell << endl; + unsigned dimension = 0; + for ( size_t i = this->multipliers.size() ; i != 0 ; --i ) + { + unsigned position = cell/multipliers[i-1]; + + if (dbg)cerr << "i-1 :" << i-1 << endl; + if (dbg)cerr << "cell : " << cell << endl; + if (dbg)cerr << "position : " << position << endl; + if (dbg)cerr << "multipliers["< +T& Bitmap_cubical_complex_base::get_cell_data( size_t cell ) +{ + return this->data[cell]; +} + + +template +void Bitmap_cubical_complex_base::impose_lower_star_filtration() +{ + bool dbg = false; + + //this vector will be used to check which elements have already been taken care of + //in imposing lower star filtration: + std::vector is_this_cell_considered( this->data.size() , false ); + + std::vector indices_to_consider; + //we assume here that we already have a filtration on the top dimensional cells and + //we have to extend it to lower ones. + typename Bitmap_cubical_complex_base::Top_dimensional_cells_iterator it(*this); + for ( it = this->top_dimensional_cells_begin() ; it != this->top_dimensional_cells_end() ; ++it ) + { + indices_to_consider.push_back( it.compute_index_in_bitmap() ); + } + + while ( indices_to_consider.size() ) + { + if ( dbg ) + { + cerr << "indices_to_consider in this iteration \n"; + for ( size_t i = 0 ; i != indices_to_consider.size() ; ++i ) + { + cout << indices_to_consider[i] << " "; + } + getchar(); + } + std::vector new_indices_to_consider; + for ( size_t i = 0 ; i != indices_to_consider.size() ; ++i ) + { + std::vector bd = this->get_boundary_of_a_cell( indices_to_consider[i] ); + for ( size_t boundaryIt = 0 ; boundaryIt != bd.size() ; ++boundaryIt ) + { + if ( this->data[ bd[boundaryIt] ] > this->data[ indices_to_consider[i] ] ) + { + this->data[ bd[boundaryIt] ] = this->data[ indices_to_consider[i] ]; + } + if ( is_this_cell_considered[ bd[boundaryIt] ] == false ) + { + new_indices_to_consider.push_back( bd[boundaryIt] ); + is_this_cell_considered[ bd[boundaryIt] ] = true; + } + } + } + indices_to_consider.swap(new_indices_to_consider); + } +} + + +template +bool compareFirstElementsOfTuples( const std::pair< std::pair< T , size_t > , char >& first , + const std::pair< std::pair< T , size_t > , char >& second ) +{ + if ( first.first.first < second.first.first ) + { + return true; + } + else + { + if ( first.first.first > second.first.first ) + { + return false; + } + //in this case first.first.first == second.first.first, so we need to compare dimensions + return first.second < second.second; + } +} + + + +template +std::vector< size_t > Bitmap_cubical_complex_base:: +generate_vector_of_shifts_for_bitmaps_with_periodic_boundary_conditions +( std::vector< bool >& directions_for_periodic_b_cond ) +{ + bool dbg = false; + if ( this->sizes.size() != directions_for_periodic_b_cond.size() ) + throw "directions_for_periodic_b_cond vector size is different from the size of the bitmap. Program terminate \n"; + + std::vector sizes( this->sizes.size() ); + for ( size_t i = 0 ; i != this->sizes.size() ; ++i )sizes[i] = 2*this->sizes[i]; + + counter c( sizes ); + + std::vector< size_t > result; + + for ( size_t i = 0 ; i != this->data.size() ; ++i ) + { + size_t position; + if ( !c.isFinal() ) + { + position = i; + //result.push_back( i ); + } + else + { + std::vector< bool > finals = c.directions_of_finals(); + bool jump_in_position = false; + for ( size_t dir = 0 ; dir != finals.size() ; ++dir ) + { + if ( finals[dir] == false )continue; + if ( directions_for_periodic_b_cond[dir] ) + { + jump_in_position = true; + } + } + if ( jump_in_position == true ) + { + //in this case this guy is final, so we need to find 'the opposite one' + position = compute_position_in_bitmap( c.find_opposite( directions_for_periodic_b_cond ) ); + } + else + { + position = i; + } + } + result.push_back( position ); + if ( dbg ) + { + cerr << " position : " << position << endl; + cerr << c << endl; + getchar(); + } + + c.increment(); + } + + return result; +} + +} + +} \ No newline at end of file diff --git a/src/Bitmap_cubical_complex/include/gudhi/counter.h b/src/Bitmap_cubical_complex/include/gudhi/counter.h index 9445a422..a5fda36f 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/counter.h +++ b/src/Bitmap_cubical_complex/include/gudhi/counter.h @@ -1,139 +1,177 @@ -/* 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 Sophia-Saclay (France) - * - * 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 COUNTER_H_ -#define COUNTER_H_ - -#include -#include - -/** - * This is an implementation of a simple counter. It is needed for the implementation of a bitmapCubicalComplex. - **/ - -class counter { - public: - /** - * Constructor of a counter class. It takes only the parameter which is the end value of the counter. - * The default beginning value is a vector of the same length as the end, filled-in with zeros. - **/ - counter(std::vector< int > endd) { - for (size_t i = 0; i != endd.size(); ++i) { - this->current.push_back(0); - this->begin.push_back(0); - this->end.push_back(endd[i]); - } - } - - /** - * Constructor of a counter class. It takes as the input beginn and endd vector. It assumes that begin vector is - * lexicographically below the end vector. - **/ - counter(std::vector< int > beginn, std::vector< int > endd) { - if (beginn.size() != endd.size()) - throw("In constructor of a counter, begin and end vectors do not have the same size. Program terminate"); - for (size_t i = 0; i != endd.size(); ++i) { - this->current.push_back(0); - this->begin.push_back(0); - this->end.push_back(endd[i]); - } - } - - /** - * Function to increment the counter. If the value returned by the function is true, then the incrementation process - * was successful. - * If the value of the function is false, that means, that the counter have reached its end-value. - **/ - bool increment() { - size_t i = 0; - while ((i != this->end.size()) && (this->current[i] == this->end[i])) { - ++i; - } - - if (i == this->end.size())return false; - ++this->current[i]; - for (size_t j = 0; j != i; ++j) { - this->current[j] = this->begin[j]; - } - return true; - } - - /** - * Function to check if we are at the end of counter. - **/ - bool isFinal() { - for (size_t i = 0; i != this->current.size(); ++i) { - if (this->current[i] == this->end[i])return true; - } - return false; - } - - /** - * Function required in the implementation of bitmapCubicalComplexWPeriodicBoundaryCondition. - * Its aim is to find an counter corresponding to the element the following boundary element is identified with - * when periodic boundary conditions are imposed. - **/ - std::vector< int > findOpposite(std::vector< bool > directionsForPeriodicBCond) { - std::vector< int > result; - for (size_t i = 0; i != this->current.size(); ++i) { - if ((this->current[i] == this->end[i]) && (directionsForPeriodicBCond[i] == true)) { - result.push_back(this->begin[i]); - } else { - result.push_back(this->current[i]); - } - } - return result; - } - - /** - * Function checking at which positions the current value of a counter is the final value of the counter. - **/ - std::vector< bool > directionsOfFinals() { - std::vector< bool > result; - for (size_t i = 0; i != this->current.size(); ++i) { - if (this->current[i] == this->end[i]) { - result.push_back(true); - } else { - result.push_back(false); - } - } - return result; - } - - /** - * Function to write counter to the stream. - **/ - friend std::ostream& operator<<(std::ostream& out, const counter& c) { - // std::cerr << "c.current.size() : " << c.current.size() << std::endl; - for (size_t i = 0; i != c.current.size(); ++i) { - out << c.current[i] << " "; - } - return out; - } - - private: - std::vector< int > begin; - std::vector< int > end; - std::vector< int > current; -}; - -#endif // COUNTER_H_ +/* 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 Sophia-Saclay (France) + * + * 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 . + */ + +#pragma once + +#include +#include + +using namespace std; + +namespace Gudhi +{ + +namespace Cubical_complex +{ + +/** +* This is an implementation of a counter being a vector of integers. +* The constructor of the class takes as an input two vectors W and V. +* It assumes that W < V coordinatewise. +* If the initial counter W is not specified, it is assumed to be vector of zeros. +* The class allows to iterate between W and V by using increment() function. +* The increment() function returns a bool value. +* The current counter reach the end counter V if the value returned by the increment function is FALSE. +* This class is needed for the implementation of a bitmapCubicalComplex. +**/ + +class counter +{ +public: + /** + * Constructor of a counter class. It takes only the parameter which is the end value of the counter. + * The default beginning value is a vector of the same length as the endd, filled-in with zeros. + **/ + counter(std::vector const& endd): begin(endd.size(),0), end(endd), current(endd.size(),0){} + //counter(std::vector< int >& endd) + //{ + // for ( size_t i = 0 ; i != endd.size() ; ++i ) + // { + // this->current.push_back(0); + // this->begin.push_back(0); + // this->end.push_back( endd[i] ); + // } + //} + + + /** + * Constructor of a counter class. It takes as the input beginn and end vector. + * It assumes that begin vector is lexicographically below the end vector. + **/ + counter(std::vector< unsigned >& beginn , std::vector< unsigned >& endd) + { + if ( beginn.size() != endd.size() ) + throw "In constructor of a counter, begin and end vectors do not have the same size. Program terminate"; + for ( size_t i = 0 ; i != endd.size() ; ++i ) + { + this->current.push_back(0); + this->begin.push_back(0); + this->end.push_back( endd[i] ); + } + } + + /** + * Function to increment the counter. If the value returned by the function is true, + * then the incrementation process was successful. + * If the value of the function is false, that means, that the counter have reached its end-value. + **/ + bool increment() + { + size_t i = 0; + while( (i != this->end.size()) && (this->current[i] == this->end[i]) ) + { + ++i; + } + + if ( i == this->end.size() )return false; + ++this->current[i]; + for ( size_t j = 0 ; j != i ; ++j ) + { + this->current[j] = this->begin[j]; + } + return true; + } + + /** + * Function to check if we are at the end of counter. + **/ + bool isFinal() + { + for ( size_t i = 0 ; i != this->current.size() ; ++i ) + { + if ( this->current[i] == this->end[i] )return true; + } + return false; + } + + /** + * Function required in the implementation of bitmapCubicalComplexWPeriodicBoundaryCondition. + * Its aim is to find an counter corresponding to the element the following + * boundary element is identified with when periodic boundary conditions are imposed. + **/ + std::vector< unsigned > find_opposite( std::vector< bool >& directionsForPeriodicBCond ) + { + std::vector< unsigned > result; + for ( size_t i = 0 ; i != this->current.size() ; ++i ) + { + if ( (this->current[i] == this->end[i]) && (directionsForPeriodicBCond[i] == true) ) + { + result.push_back( this->begin[i] ); + } + else + { + result.push_back( this->current[i] ); + } + } + return result; + } + + /** + * Function checking at which positions the current value of a counter is the final value of the counter. + **/ + std::vector< bool > directions_of_finals() + { + std::vector< bool > result; + for ( size_t i = 0 ; i != this->current.size() ; ++i ) + { + if ( this->current[i] == this->end[i] ) + { + result.push_back( true ); + } + else + { + result.push_back( false ); + } + } + return result; + } + + /** + * Function to write counter to the stream. + **/ + friend std::ostream& operator<<(std::ostream& out , const counter& c ) + { + //cerr << "c.current.size() : " << c.current.size() << endl; + for ( size_t i = 0 ; i != c.current.size() ; ++i ) + { + out << c.current[i] << " "; + } + return out; + } +private: + std::vector< unsigned > begin; + std::vector< unsigned > end; + std::vector< unsigned > current; +}; + +} +} \ No newline at end of file diff --git a/src/Bitmap_cubical_complex/test/Bitmap_test.cpp b/src/Bitmap_cubical_complex/test/Bitmap_test.cpp index 1c204bae..183b856a 100644 --- a/src/Bitmap_cubical_complex/test/Bitmap_test.cpp +++ b/src/Bitmap_cubical_complex/test/Bitmap_test.cpp @@ -1,623 +1,638 @@ -#include "gudhi/Bitmap_cubical_complex.h" - -#define BOOST_TEST_DYN_LINK -#define BOOST_TEST_MODULE "cubical_complex" -#include - -using namespace std; - -BOOST_AUTO_TEST_CASE(check_dimension) { - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - BOOST_CHECK(increasing.dimension() == 2); -} - -BOOST_AUTO_TEST_CASE(topDimensionalCellsIterator_test) { - std::vector< double > expectedFiltrationValues1; - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(100); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - expectedFiltrationValues1.push_back(0); - - std::vector< double > expectedFiltrationValues2; - expectedFiltrationValues2.push_back(1); - expectedFiltrationValues2.push_back(2); - expectedFiltrationValues2.push_back(3); - expectedFiltrationValues2.push_back(4); - expectedFiltrationValues2.push_back(5); - expectedFiltrationValues2.push_back(6); - expectedFiltrationValues2.push_back(7); - expectedFiltrationValues2.push_back(8); - expectedFiltrationValues2.push_back(9); - - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector< double > oneDimensionalCycle; - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(100); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - oneDimensionalCycle.push_back(0); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - Bitmap_cubical_complex hole(dimensions, oneDimensionalCycle); - - - int i = 0; - for (Bitmap_cubical_complex::Top_dimensional_cells_iterator it = increasing.top_dimensional_cells_begin(); it != increasing.top_dimensional_cells_end(); ++it) { - BOOST_CHECK(*it == expectedFiltrationValues2[i]); - ++i; - } - i = 0; - for (Bitmap_cubical_complex::Top_dimensional_cells_iterator it = hole.top_dimensional_cells_begin(); it != hole.top_dimensional_cells_end(); ++it) { - BOOST_CHECK(*it == expectedFiltrationValues1[i]); - ++i; - } -} - -BOOST_AUTO_TEST_CASE(compute_boundary_test_1) { - - std::vector boundary0; - std::vector boundary1; - boundary1.push_back(0); - boundary1.push_back(2); - std::vector boundary2; - std::vector boundary3; - boundary3.push_back(2); - boundary3.push_back(4); - std::vector boundary4; - std::vector boundary5; - boundary5.push_back(4); - boundary5.push_back(6); - std::vector boundary6; - std::vector boundary7; - boundary7.push_back(0); - boundary7.push_back(14); - std::vector boundary8; - boundary8.push_back(1); - boundary8.push_back(15); - boundary8.push_back(7); - boundary8.push_back(9); - std::vector boundary9; - boundary9.push_back(2); - boundary9.push_back(16); - std::vector boundary10; - boundary10.push_back(3); - boundary10.push_back(17); - boundary10.push_back(9); - boundary10.push_back(11); - std::vector boundary11; - boundary11.push_back(4); - boundary11.push_back(18); - std::vector boundary12; - boundary12.push_back(5); - boundary12.push_back(19); - boundary12.push_back(11); - boundary12.push_back(13); - std::vector boundary13; - boundary13.push_back(6); - boundary13.push_back(20); - std::vector boundary14; - std::vector boundary15; - boundary15.push_back(14); - boundary15.push_back(16); - std::vector boundary16; - std::vector boundary17; - boundary17.push_back(16); - boundary17.push_back(18); - std::vector boundary18; - std::vector boundary19; - boundary19.push_back(18); - boundary19.push_back(20); - std::vector boundary20; - std::vector boundary21; - boundary21.push_back(14); - boundary21.push_back(28); - std::vector boundary22; - boundary22.push_back(15); - boundary22.push_back(29); - boundary22.push_back(21); - boundary22.push_back(23); - std::vector boundary23; - boundary23.push_back(16); - boundary23.push_back(30); - std::vector boundary24; - boundary24.push_back(17); - boundary24.push_back(31); - boundary24.push_back(23); - boundary24.push_back(25); - std::vector boundary25; - boundary25.push_back(18); - boundary25.push_back(32); - std::vector boundary26; - boundary26.push_back(19); - boundary26.push_back(33); - boundary26.push_back(25); - boundary26.push_back(27); - std::vector boundary27; - boundary27.push_back(20); - boundary27.push_back(34); - std::vector boundary28; - std::vector boundary29; - boundary29.push_back(28); - boundary29.push_back(30); - std::vector boundary30; - std::vector boundary31; - boundary31.push_back(30); - boundary31.push_back(32); - std::vector boundary32; - std::vector boundary33; - boundary33.push_back(32); - boundary33.push_back(34); - std::vector boundary34; - std::vector boundary35; - boundary35.push_back(28); - boundary35.push_back(42); - std::vector boundary36; - boundary36.push_back(29); - boundary36.push_back(43); - boundary36.push_back(35); - boundary36.push_back(37); - std::vector boundary37; - boundary37.push_back(30); - boundary37.push_back(44); - std::vector boundary38; - boundary38.push_back(31); - boundary38.push_back(45); - boundary38.push_back(37); - boundary38.push_back(39); - std::vector boundary39; - boundary39.push_back(32); - boundary39.push_back(46); - std::vector boundary40; - boundary40.push_back(33); - boundary40.push_back(47); - boundary40.push_back(39); - boundary40.push_back(41); - std::vector boundary41; - boundary41.push_back(34); - boundary41.push_back(48); - std::vector boundary42; - std::vector boundary43; - boundary43.push_back(42); - boundary43.push_back(44); - std::vector boundary44; - std::vector boundary45; - boundary45.push_back(44); - boundary45.push_back(46); - std::vector boundary46; - std::vector boundary47; - boundary47.push_back(46); - boundary47.push_back(48); - std::vector boundary48; - std::vector< std::vector > boundaries; - boundaries.push_back(boundary0); - boundaries.push_back(boundary1); - boundaries.push_back(boundary2); - boundaries.push_back(boundary3); - boundaries.push_back(boundary4); - boundaries.push_back(boundary5); - boundaries.push_back(boundary6); - boundaries.push_back(boundary7); - boundaries.push_back(boundary8); - boundaries.push_back(boundary9); - boundaries.push_back(boundary10); - boundaries.push_back(boundary11); - boundaries.push_back(boundary12); - boundaries.push_back(boundary13); - boundaries.push_back(boundary14); - boundaries.push_back(boundary15); - boundaries.push_back(boundary16); - boundaries.push_back(boundary17); - boundaries.push_back(boundary18); - boundaries.push_back(boundary19); - boundaries.push_back(boundary20); - boundaries.push_back(boundary21); - boundaries.push_back(boundary22); - boundaries.push_back(boundary23); - boundaries.push_back(boundary24); - boundaries.push_back(boundary25); - boundaries.push_back(boundary26); - boundaries.push_back(boundary27); - boundaries.push_back(boundary28); - boundaries.push_back(boundary29); - boundaries.push_back(boundary30); - boundaries.push_back(boundary31); - boundaries.push_back(boundary32); - boundaries.push_back(boundary33); - boundaries.push_back(boundary34); - boundaries.push_back(boundary35); - boundaries.push_back(boundary36); - boundaries.push_back(boundary37); - boundaries.push_back(boundary38); - boundaries.push_back(boundary39); - boundaries.push_back(boundary40); - boundaries.push_back(boundary41); - boundaries.push_back(boundary42); - boundaries.push_back(boundary43); - boundaries.push_back(boundary44); - boundaries.push_back(boundary45); - boundaries.push_back(boundary46); - boundaries.push_back(boundary47); - boundaries.push_back(boundary48); - - - - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { - std::vector< size_t > bd = increasing.get_boundary_of_a_cell(i); - for (size_t j = 0; j != bd.size(); ++j) { - BOOST_CHECK(boundaries[i][j] == bd[j]); - } - } -} - -BOOST_AUTO_TEST_CASE(compute_boundary_test_2) { - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - - - std::vector coboundaryElements; - coboundaryElements.push_back(7); - coboundaryElements.push_back(1); - coboundaryElements.push_back(8); - coboundaryElements.push_back(9); - coboundaryElements.push_back(1); - coboundaryElements.push_back(3); - coboundaryElements.push_back(10); - coboundaryElements.push_back(11); - coboundaryElements.push_back(3); - coboundaryElements.push_back(5); - coboundaryElements.push_back(12); - coboundaryElements.push_back(13); - coboundaryElements.push_back(5); - coboundaryElements.push_back(8); - coboundaryElements.push_back(8); - coboundaryElements.push_back(10); - coboundaryElements.push_back(10); - coboundaryElements.push_back(12); - coboundaryElements.push_back(12); - coboundaryElements.push_back(7); - coboundaryElements.push_back(21); - coboundaryElements.push_back(15); - coboundaryElements.push_back(8); - coboundaryElements.push_back(22); - coboundaryElements.push_back(9); - coboundaryElements.push_back(23); - coboundaryElements.push_back(15); - coboundaryElements.push_back(17); - coboundaryElements.push_back(10); - coboundaryElements.push_back(24); - coboundaryElements.push_back(11); - coboundaryElements.push_back(25); - coboundaryElements.push_back(17); - coboundaryElements.push_back(19); - coboundaryElements.push_back(12); - coboundaryElements.push_back(26); - coboundaryElements.push_back(13); - coboundaryElements.push_back(27); - coboundaryElements.push_back(19); - coboundaryElements.push_back(22); - coboundaryElements.push_back(22); - coboundaryElements.push_back(24); - coboundaryElements.push_back(24); - coboundaryElements.push_back(26); - coboundaryElements.push_back(26); - coboundaryElements.push_back(21); - coboundaryElements.push_back(35); - coboundaryElements.push_back(29); - coboundaryElements.push_back(22); - coboundaryElements.push_back(36); - coboundaryElements.push_back(23); - coboundaryElements.push_back(37); - coboundaryElements.push_back(29); - coboundaryElements.push_back(31); - coboundaryElements.push_back(24); - coboundaryElements.push_back(38); - coboundaryElements.push_back(25); - coboundaryElements.push_back(39); - coboundaryElements.push_back(31); - coboundaryElements.push_back(33); - coboundaryElements.push_back(26); - coboundaryElements.push_back(40); - coboundaryElements.push_back(27); - coboundaryElements.push_back(41); - coboundaryElements.push_back(33); - coboundaryElements.push_back(36); - coboundaryElements.push_back(36); - coboundaryElements.push_back(38); - coboundaryElements.push_back(38); - coboundaryElements.push_back(40); - coboundaryElements.push_back(40); - coboundaryElements.push_back(35); - coboundaryElements.push_back(43); - coboundaryElements.push_back(36); - coboundaryElements.push_back(37); - coboundaryElements.push_back(43); - coboundaryElements.push_back(45); - coboundaryElements.push_back(38); - coboundaryElements.push_back(39); - coboundaryElements.push_back(45); - coboundaryElements.push_back(47); - coboundaryElements.push_back(40); - coboundaryElements.push_back(41); - coboundaryElements.push_back(47); - size_t number = 0; - for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { - std::vector< size_t > bd = increasing.get_coboundary_of_a_cell(i); - for (size_t j = 0; j != bd.size(); ++j) { - BOOST_CHECK(coboundaryElements[number] == bd[j]); - ++number; - } - - } -} - -BOOST_AUTO_TEST_CASE(compute_boundary_test_3) { - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - - std::vector dim; - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(2); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - dim.push_back(1); - dim.push_back(0); - - for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { - BOOST_CHECK(increasing.get_dimension_of_a_cell(i) == dim[i]); - } -} - -BOOST_AUTO_TEST_CASE(Filtration_simplex_iterator_test) { - std::vector< double > increasingFiltrationOfTopDimensionalCells; - increasingFiltrationOfTopDimensionalCells.push_back(1); - increasingFiltrationOfTopDimensionalCells.push_back(2); - increasingFiltrationOfTopDimensionalCells.push_back(3); - increasingFiltrationOfTopDimensionalCells.push_back(4); - increasingFiltrationOfTopDimensionalCells.push_back(5); - increasingFiltrationOfTopDimensionalCells.push_back(6); - increasingFiltrationOfTopDimensionalCells.push_back(7); - increasingFiltrationOfTopDimensionalCells.push_back(8); - increasingFiltrationOfTopDimensionalCells.push_back(9); - - std::vector dimensions; - dimensions.push_back(3); - dimensions.push_back(3); - - Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); - - std::vector< unsigned > dim; - dim.push_back(0); - dim.push_back(0); - dim.push_back(0); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - dim.push_back(0); - dim.push_back(1); - dim.push_back(1); - dim.push_back(2); - - std::vector fil; - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(1); - fil.push_back(2); - fil.push_back(2); - fil.push_back(2); - fil.push_back(2); - fil.push_back(2); - fil.push_back(2); - fil.push_back(3); - fil.push_back(3); - fil.push_back(3); - fil.push_back(3); - fil.push_back(3); - fil.push_back(3); - fil.push_back(4); - fil.push_back(4); - fil.push_back(4); - fil.push_back(4); - fil.push_back(4); - fil.push_back(4); - fil.push_back(5); - fil.push_back(5); - fil.push_back(5); - fil.push_back(5); - fil.push_back(6); - fil.push_back(6); - fil.push_back(6); - fil.push_back(6); - fil.push_back(7); - fil.push_back(7); - fil.push_back(7); - fil.push_back(7); - fil.push_back(7); - fil.push_back(7); - fil.push_back(8); - fil.push_back(8); - fil.push_back(8); - fil.push_back(8); - fil.push_back(9); - fil.push_back(9); - fil.push_back(9); - fil.push_back(9); - - - Bitmap_cubical_complex::Filtration_simplex_range range = increasing.filtration_simplex_range(); - size_t position = 0; - for (Bitmap_cubical_complex::Filtration_simplex_iterator it = range.begin(); it != range.end(); ++it) { - BOOST_CHECK(increasing.dimension(*it) == dim[position]); - BOOST_CHECK(increasing.filtration(*it) == fil[position]); - ++position; - } -} +#include +#include +#include +#include + +// standard stuff +#include +#include + +#define BOOST_TEST_DYN_LINK +#define BOOST_TEST_MODULE "cubical_complex" +#include + + +using namespace std; +using namespace Gudhi; +using namespace Gudhi::Cubical_complex; +using namespace Gudhi::persistent_cohomology; + + + +BOOST_AUTO_TEST_CASE(check_dimension) { + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + BOOST_CHECK(increasing.dimension() == 2); +} + +BOOST_AUTO_TEST_CASE(topDimensionalCellsIterator_test) { + std::vector< double > expectedFiltrationValues1; + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(100); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + expectedFiltrationValues1.push_back(0); + + std::vector< double > expectedFiltrationValues2; + expectedFiltrationValues2.push_back(1); + expectedFiltrationValues2.push_back(2); + expectedFiltrationValues2.push_back(3); + expectedFiltrationValues2.push_back(4); + expectedFiltrationValues2.push_back(5); + expectedFiltrationValues2.push_back(6); + expectedFiltrationValues2.push_back(7); + expectedFiltrationValues2.push_back(8); + expectedFiltrationValues2.push_back(9); + + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector< double > oneDimensionalCycle; + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(100); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + oneDimensionalCycle.push_back(0); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + Bitmap_cubical_complex hole(dimensions, oneDimensionalCycle); + + + int i = 0; + for (Bitmap_cubical_complex::Top_dimensional_cells_iterator + it = increasing.top_dimensional_cells_begin(); it != increasing.top_dimensional_cells_end(); ++it) { + BOOST_CHECK(*it == expectedFiltrationValues2[i]); + ++i; + } + i = 0; + for (Bitmap_cubical_complex::Top_dimensional_cells_iterator + it = hole.top_dimensional_cells_begin(); it != hole.top_dimensional_cells_end(); ++it) { + BOOST_CHECK(*it == expectedFiltrationValues1[i]); + ++i; + } +} + +BOOST_AUTO_TEST_CASE(compute_boundary_test_1) { + + std::vector boundary0; + std::vector boundary1; + boundary1.push_back(0); + boundary1.push_back(2); + std::vector boundary2; + std::vector boundary3; + boundary3.push_back(2); + boundary3.push_back(4); + std::vector boundary4; + std::vector boundary5; + boundary5.push_back(4); + boundary5.push_back(6); + std::vector boundary6; + std::vector boundary7; + boundary7.push_back(0); + boundary7.push_back(14); + std::vector boundary8; + boundary8.push_back(1); + boundary8.push_back(15); + boundary8.push_back(7); + boundary8.push_back(9); + std::vector boundary9; + boundary9.push_back(2); + boundary9.push_back(16); + std::vector boundary10; + boundary10.push_back(3); + boundary10.push_back(17); + boundary10.push_back(9); + boundary10.push_back(11); + std::vector boundary11; + boundary11.push_back(4); + boundary11.push_back(18); + std::vector boundary12; + boundary12.push_back(5); + boundary12.push_back(19); + boundary12.push_back(11); + boundary12.push_back(13); + std::vector boundary13; + boundary13.push_back(6); + boundary13.push_back(20); + std::vector boundary14; + std::vector boundary15; + boundary15.push_back(14); + boundary15.push_back(16); + std::vector boundary16; + std::vector boundary17; + boundary17.push_back(16); + boundary17.push_back(18); + std::vector boundary18; + std::vector boundary19; + boundary19.push_back(18); + boundary19.push_back(20); + std::vector boundary20; + std::vector boundary21; + boundary21.push_back(14); + boundary21.push_back(28); + std::vector boundary22; + boundary22.push_back(15); + boundary22.push_back(29); + boundary22.push_back(21); + boundary22.push_back(23); + std::vector boundary23; + boundary23.push_back(16); + boundary23.push_back(30); + std::vector boundary24; + boundary24.push_back(17); + boundary24.push_back(31); + boundary24.push_back(23); + boundary24.push_back(25); + std::vector boundary25; + boundary25.push_back(18); + boundary25.push_back(32); + std::vector boundary26; + boundary26.push_back(19); + boundary26.push_back(33); + boundary26.push_back(25); + boundary26.push_back(27); + std::vector boundary27; + boundary27.push_back(20); + boundary27.push_back(34); + std::vector boundary28; + std::vector boundary29; + boundary29.push_back(28); + boundary29.push_back(30); + std::vector boundary30; + std::vector boundary31; + boundary31.push_back(30); + boundary31.push_back(32); + std::vector boundary32; + std::vector boundary33; + boundary33.push_back(32); + boundary33.push_back(34); + std::vector boundary34; + std::vector boundary35; + boundary35.push_back(28); + boundary35.push_back(42); + std::vector boundary36; + boundary36.push_back(29); + boundary36.push_back(43); + boundary36.push_back(35); + boundary36.push_back(37); + std::vector boundary37; + boundary37.push_back(30); + boundary37.push_back(44); + std::vector boundary38; + boundary38.push_back(31); + boundary38.push_back(45); + boundary38.push_back(37); + boundary38.push_back(39); + std::vector boundary39; + boundary39.push_back(32); + boundary39.push_back(46); + std::vector boundary40; + boundary40.push_back(33); + boundary40.push_back(47); + boundary40.push_back(39); + boundary40.push_back(41); + std::vector boundary41; + boundary41.push_back(34); + boundary41.push_back(48); + std::vector boundary42; + std::vector boundary43; + boundary43.push_back(42); + boundary43.push_back(44); + std::vector boundary44; + std::vector boundary45; + boundary45.push_back(44); + boundary45.push_back(46); + std::vector boundary46; + std::vector boundary47; + boundary47.push_back(46); + boundary47.push_back(48); + std::vector boundary48; + std::vector< std::vector > boundaries; + boundaries.push_back(boundary0); + boundaries.push_back(boundary1); + boundaries.push_back(boundary2); + boundaries.push_back(boundary3); + boundaries.push_back(boundary4); + boundaries.push_back(boundary5); + boundaries.push_back(boundary6); + boundaries.push_back(boundary7); + boundaries.push_back(boundary8); + boundaries.push_back(boundary9); + boundaries.push_back(boundary10); + boundaries.push_back(boundary11); + boundaries.push_back(boundary12); + boundaries.push_back(boundary13); + boundaries.push_back(boundary14); + boundaries.push_back(boundary15); + boundaries.push_back(boundary16); + boundaries.push_back(boundary17); + boundaries.push_back(boundary18); + boundaries.push_back(boundary19); + boundaries.push_back(boundary20); + boundaries.push_back(boundary21); + boundaries.push_back(boundary22); + boundaries.push_back(boundary23); + boundaries.push_back(boundary24); + boundaries.push_back(boundary25); + boundaries.push_back(boundary26); + boundaries.push_back(boundary27); + boundaries.push_back(boundary28); + boundaries.push_back(boundary29); + boundaries.push_back(boundary30); + boundaries.push_back(boundary31); + boundaries.push_back(boundary32); + boundaries.push_back(boundary33); + boundaries.push_back(boundary34); + boundaries.push_back(boundary35); + boundaries.push_back(boundary36); + boundaries.push_back(boundary37); + boundaries.push_back(boundary38); + boundaries.push_back(boundary39); + boundaries.push_back(boundary40); + boundaries.push_back(boundary41); + boundaries.push_back(boundary42); + boundaries.push_back(boundary43); + boundaries.push_back(boundary44); + boundaries.push_back(boundary45); + boundaries.push_back(boundary46); + boundaries.push_back(boundary47); + boundaries.push_back(boundary48); + + + + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { + std::vector< size_t > bd = increasing.get_boundary_of_a_cell(i); + for (size_t j = 0; j != bd.size(); ++j) { + BOOST_CHECK(boundaries[i][j] == bd[j]); + } + } +} + +BOOST_AUTO_TEST_CASE(compute_boundary_test_2) { + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + + + std::vector coboundaryElements; + coboundaryElements.push_back(7); + coboundaryElements.push_back(1); + coboundaryElements.push_back(8); + coboundaryElements.push_back(9); + coboundaryElements.push_back(1); + coboundaryElements.push_back(3); + coboundaryElements.push_back(10); + coboundaryElements.push_back(11); + coboundaryElements.push_back(3); + coboundaryElements.push_back(5); + coboundaryElements.push_back(12); + coboundaryElements.push_back(13); + coboundaryElements.push_back(5); + coboundaryElements.push_back(8); + coboundaryElements.push_back(8); + coboundaryElements.push_back(10); + coboundaryElements.push_back(10); + coboundaryElements.push_back(12); + coboundaryElements.push_back(12); + coboundaryElements.push_back(7); + coboundaryElements.push_back(21); + coboundaryElements.push_back(15); + coboundaryElements.push_back(8); + coboundaryElements.push_back(22); + coboundaryElements.push_back(9); + coboundaryElements.push_back(23); + coboundaryElements.push_back(15); + coboundaryElements.push_back(17); + coboundaryElements.push_back(10); + coboundaryElements.push_back(24); + coboundaryElements.push_back(11); + coboundaryElements.push_back(25); + coboundaryElements.push_back(17); + coboundaryElements.push_back(19); + coboundaryElements.push_back(12); + coboundaryElements.push_back(26); + coboundaryElements.push_back(13); + coboundaryElements.push_back(27); + coboundaryElements.push_back(19); + coboundaryElements.push_back(22); + coboundaryElements.push_back(22); + coboundaryElements.push_back(24); + coboundaryElements.push_back(24); + coboundaryElements.push_back(26); + coboundaryElements.push_back(26); + coboundaryElements.push_back(21); + coboundaryElements.push_back(35); + coboundaryElements.push_back(29); + coboundaryElements.push_back(22); + coboundaryElements.push_back(36); + coboundaryElements.push_back(23); + coboundaryElements.push_back(37); + coboundaryElements.push_back(29); + coboundaryElements.push_back(31); + coboundaryElements.push_back(24); + coboundaryElements.push_back(38); + coboundaryElements.push_back(25); + coboundaryElements.push_back(39); + coboundaryElements.push_back(31); + coboundaryElements.push_back(33); + coboundaryElements.push_back(26); + coboundaryElements.push_back(40); + coboundaryElements.push_back(27); + coboundaryElements.push_back(41); + coboundaryElements.push_back(33); + coboundaryElements.push_back(36); + coboundaryElements.push_back(36); + coboundaryElements.push_back(38); + coboundaryElements.push_back(38); + coboundaryElements.push_back(40); + coboundaryElements.push_back(40); + coboundaryElements.push_back(35); + coboundaryElements.push_back(43); + coboundaryElements.push_back(36); + coboundaryElements.push_back(37); + coboundaryElements.push_back(43); + coboundaryElements.push_back(45); + coboundaryElements.push_back(38); + coboundaryElements.push_back(39); + coboundaryElements.push_back(45); + coboundaryElements.push_back(47); + coboundaryElements.push_back(40); + coboundaryElements.push_back(41); + coboundaryElements.push_back(47); + size_t number = 0; + for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { + std::vector< size_t > bd = increasing.get_coboundary_of_a_cell(i); + for (size_t j = 0; j != bd.size(); ++j) { + BOOST_CHECK(coboundaryElements[number] == bd[j]); + ++number; + } + + } +} + +BOOST_AUTO_TEST_CASE(compute_boundary_test_3) { + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + + std::vector dim; + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(2); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + dim.push_back(1); + dim.push_back(0); + + for (size_t i = 0; i != increasing.size_of_bitmap(); ++i) { + BOOST_CHECK(increasing.get_dimension_of_a_cell(i) == dim[i]); + } +} + +BOOST_AUTO_TEST_CASE(Filtration_simplex_iterator_test) { + std::vector< double > increasingFiltrationOfTopDimensionalCells; + increasingFiltrationOfTopDimensionalCells.push_back(1); + increasingFiltrationOfTopDimensionalCells.push_back(2); + increasingFiltrationOfTopDimensionalCells.push_back(3); + increasingFiltrationOfTopDimensionalCells.push_back(4); + increasingFiltrationOfTopDimensionalCells.push_back(5); + increasingFiltrationOfTopDimensionalCells.push_back(6); + increasingFiltrationOfTopDimensionalCells.push_back(7); + increasingFiltrationOfTopDimensionalCells.push_back(8); + increasingFiltrationOfTopDimensionalCells.push_back(9); + + std::vector dimensions; + dimensions.push_back(3); + dimensions.push_back(3); + + Bitmap_cubical_complex increasing(dimensions, increasingFiltrationOfTopDimensionalCells); + + std::vector< unsigned > dim; + dim.push_back(0); + dim.push_back(0); + dim.push_back(0); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + dim.push_back(0); + dim.push_back(1); + dim.push_back(1); + dim.push_back(2); + + std::vector fil; + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(1); + fil.push_back(2); + fil.push_back(2); + fil.push_back(2); + fil.push_back(2); + fil.push_back(2); + fil.push_back(2); + fil.push_back(3); + fil.push_back(3); + fil.push_back(3); + fil.push_back(3); + fil.push_back(3); + fil.push_back(3); + fil.push_back(4); + fil.push_back(4); + fil.push_back(4); + fil.push_back(4); + fil.push_back(4); + fil.push_back(4); + fil.push_back(5); + fil.push_back(5); + fil.push_back(5); + fil.push_back(5); + fil.push_back(6); + fil.push_back(6); + fil.push_back(6); + fil.push_back(6); + fil.push_back(7); + fil.push_back(7); + fil.push_back(7); + fil.push_back(7); + fil.push_back(7); + fil.push_back(7); + fil.push_back(8); + fil.push_back(8); + fil.push_back(8); + fil.push_back(8); + fil.push_back(9); + fil.push_back(9); + fil.push_back(9); + fil.push_back(9); + + + Bitmap_cubical_complex::Filtration_simplex_range range = increasing.filtration_simplex_range(); + size_t position = 0; + for (Bitmap_cubical_complex::Filtration_simplex_iterator it = range.begin(); it != range.end(); ++it) { + BOOST_CHECK(increasing.dimension(*it) == dim[position]); + BOOST_CHECK(increasing.filtration(*it) == fil[position]); + ++position; + } +} -- cgit v1.2.3