diff options
author | vrouvrea <vrouvrea@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2016-10-20 10:04:05 +0000 |
---|---|---|
committer | vrouvrea <vrouvrea@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2016-10-20 10:04:05 +0000 |
commit | 8d656e33138ef8b3a7d86a7375c92646efc29511 (patch) | |
tree | 3711227c4c1b2a6e9f25dda1db8dafb8365063a0 /src/Persistent_cohomology | |
parent | 355dc2a0ae73f243fc0aa8d7c797509d8ba5e6b4 (diff) | |
parent | 30e538a98919004e36b3b382017884486919cb6e (diff) |
Merge last trunk modification
Fix doc issue
Still doc issue
git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/ST_cythonize@1739 636b058d-ea47-450e-bf9e-a15bfbe3eedb
Former-commit-id: 0a99345f06e93a3525691699a6fe1505979e8e8e
Diffstat (limited to 'src/Persistent_cohomology')
8 files changed, 233 insertions, 131 deletions
diff --git a/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h b/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h index c8081cac..433cfd3e 100644 --- a/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h +++ b/src/Persistent_cohomology/doc/Intro_persistent_cohomology.h @@ -46,7 +46,7 @@ namespace persistent_cohomology { composed of three elements: topological spaces, their homology groups and an evolution scheme. - <DT>Topological Spaces:</DT> + \section persistencetopolocalspaces Topological Spaces Topological spaces are represented by simplicial complexes. Let \f$V = \{1, \cdots ,|V|\}\f$ be a set of <EM>vertices</EM>. A <EM>simplex</EM> \f$\sigma\f$ is a subset of vertices @@ -84,7 +84,7 @@ namespace persistent_cohomology { <CODE>Filtration_value filtration (Simplex_handle)</CODE> that returns the value of the filtration on the simplex represented by the handle. - <DT>Homology:</DT> + \section persistencehomology Homology For a ring \f$\mathcal{R}\f$, the group of <EM>n-chains</EM>, denoted \f$\mathbf{C}_n(\mathbf{K},\mathcal{R})\f$, of \f$\mathbf{K}\f$ is the group of formal sums of @@ -111,7 +111,7 @@ namespace persistent_cohomology { We refer to \cite Munkres-elementsalgtop1984 for an introduction to homology theory and to \cite DBLP:books/daglib/0025666 for an introduction to persistent homology. - <DT>Indexing Scheme:</DT> + \section persistenceindexingscheme Indexing Scheme "Changing" a simplicial complex consists in applying a simplicial map. An <EM>indexing scheme</EM> is a directed graph together with a traversal order, such that two @@ -139,19 +139,63 @@ namespace persistent_cohomology { by increasing filtration values (breaking ties so as a simplex appears after its subsimplices of same filtration value) provides an indexing scheme. -\section Examples - We provide several example files: run these examples with -h for details on their use, and read the README file. +\section pcohexamples Examples -\li <CODE>rips_persistence.cpp</CODE> computes the Rips complex of a point cloud and its persistence diagram. +We provide several example files: run these examples with -h for details on their use, and read the README file. -\li <CODE>rips_multifield_persistence.cpp</CODE> computes the Rips complex of a point cloud and its persistence diagram -with a family of field coefficients. +\li <a href="_persistent_cohomology_2rips_persistence_8cpp-example.html"> +Persistent_cohomology/rips_persistence.cpp</a> computes the Rips complex of a point cloud and its persistence diagram. -\li <CODE>performance_rips_persistence.cpp</CODE> provides timings for the construction of the Rips complex on a set of -points sampling a Klein bottle in \f$\mathbb{R}^5\f$ with a simplex tree, its conversion to a +\li <a href="_persistent_cohomology_2rips_multifield_persistence_8cpp-example.html"> +Persistent_cohomology/rips_multifield_persistence.cpp</a> computes the Rips complex of a point cloud and its +persistence diagram with a family of field coefficients. + +\li <a href="_persistent_cohomology_2performance_rips_persistence_8cpp-example.html"> +Persistent_cohomology/performance_rips_persistence.cpp</a> provides timings for the construction of the Rips complex +on a set of points sampling a Klein bottle in \f$\mathbb{R}^5\f$ with a simplex tree, its conversion to a Hasse diagram and the computation of persistent homology and multi-field persistent homology for the different representations. +\li <a href="_persistent_cohomology_2alpha_complex_3d_persistence_8cpp-example.html"> +Persistent_cohomology/alpha_complex_3d_persistence.cpp</a> computes the persistent homology with +\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file. +\code $> ./alpha_complex_3d_persistence ../../data/points/tore3D_300.off 2 0.45 \endcode +\code Simplex_tree dim: 3 +2 0 0 inf +2 1 0.0682162 1.0001 +2 1 0.0934117 1.00003 +2 2 0.56444 1.03938 \endcode + +\li <a href="_persistent_cohomology_2alpha_complex_persistence_8cpp-example.html"> +Persistent_cohomology/alpha_complex_persistence.cpp</a> computes the persistent homology with +\f$\mathbb{Z}/p\mathbb{Z}\f$ coefficients of the alpha complex on points sampling from an OFF file. +\code $> ./alpha_complex_persistence -r 32 -p 2 -m 0.45 ../../data/points/tore3D_300.off \endcode +\code Alpha complex is of dimension 3 - 9273 simplices - 300 vertices. +Simplex_tree dim: 3 +2 0 0 inf +2 1 0.0682162 1.0001 +2 1 0.0934117 1.00003 +2 2 0.56444 1.03938 \endcode + +\li <a href="_persistent_cohomology_2periodic_alpha_complex_3d_persistence_8cpp-example.html"> +Persistent_cohomology/periodic_alpha_complex_3d_persistence.cpp</a> computes the persistent homology with +\f$\mathbb{Z}/2\mathbb{Z}\f$ coefficients of the periodic alpha complex on points sampling from an OFF file. +\code $> ./periodic_alpha_complex_3d_persistence ../../data/points/grid_10_10_10_in_0_1.off 3 1.0 \endcode +\code Periodic Delaunay computed. +Simplex_tree dim: 3 +3 0 0 inf +3 1 0.0025 inf +3 1 0.0025 inf +3 1 0.0025 inf +3 2 0.005 inf +3 2 0.005 inf +3 2 0.005 inf +3 3 0.0075 inf \endcode + +\li <a href="_persistent_cohomology_2plain_homology_8cpp-example.html"> +Persistent_cohomology/plain_homology.cpp</a> computes the plain homology of a simple simplicial complex without +filtration values. + \copyright GNU General Public License v3. */ diff --git a/src/Persistent_cohomology/example/CMakeLists.txt b/src/Persistent_cohomology/example/CMakeLists.txt index 5fb4ba12..758bd6b1 100644 --- a/src/Persistent_cohomology/example/CMakeLists.txt +++ b/src/Persistent_cohomology/example/CMakeLists.txt @@ -30,17 +30,13 @@ endif() add_test(plain_homology ${CMAKE_CURRENT_BINARY_DIR}/plain_homology) add_test(persistence_from_simple_simplex_tree ${CMAKE_CURRENT_BINARY_DIR}/persistence_from_simple_simplex_tree 1 0) -add_test(rips_persistence_3 ${CMAKE_CURRENT_BINARY_DIR}/rips_persistence ${CMAKE_SOURCE_DIR}/data/points/Kl.txt -r 0.2 -d 3 -p 3 -m 100) -add_test(rips_persistence_via_boundary_matrix_3 ${CMAKE_CURRENT_BINARY_DIR}/rips_persistence_via_boundary_matrix ${CMAKE_SOURCE_DIR}/data/points/tore3D_1307.txt -r 0.3 -d 3 -p 3 -m 100) -add_test(persistence_from_file_3_2_0 ${CMAKE_CURRENT_BINARY_DIR}/persistence_from_file ${CMAKE_SOURCE_DIR}/data/points/bunny_5000.st -p 2 -m 0) -add_test(persistence_from_file_3_3_100 ${CMAKE_CURRENT_BINARY_DIR}/persistence_from_file ${CMAKE_SOURCE_DIR}/data/points/bunny_5000.st -p 3 -m 100) +add_test(rips_persistence_3 ${CMAKE_CURRENT_BINARY_DIR}/rips_persistence ${CMAKE_SOURCE_DIR}/data/points/Kl.txt -r 0.16 -d 3 -p 3 -m 100) +add_test(rips_persistence_via_boundary_matrix_3 ${CMAKE_CURRENT_BINARY_DIR}/rips_persistence_via_boundary_matrix ${CMAKE_SOURCE_DIR}/data/points/Kl.txt -r 0.16 -d 3 -p 3 -m 100) +add_test(persistence_from_file_3_2_0 ${CMAKE_CURRENT_BINARY_DIR}/persistence_from_file ${CMAKE_SOURCE_DIR}/data/filtered_simplicial_complex/bunny_5000_complex.fsc -p 2 -m 0) +add_test(persistence_from_file_3_3_100 ${CMAKE_CURRENT_BINARY_DIR}/persistence_from_file ${CMAKE_SOURCE_DIR}/data/filtered_simplicial_complex/bunny_5000_complex.fsc -p 3 -m 100) if(GMP_FOUND) - message("GMP_LIBRARIES = ${GMP_LIBRARIES}") - if(GMPXX_FOUND) - message("GMPXX_LIBRARIES = ${GMPXX_LIBRARIES}") - add_executable(rips_multifield_persistence rips_multifield_persistence.cpp ) target_link_libraries(rips_multifield_persistence ${Boost_SYSTEM_LIBRARY} ${Boost_PROGRAM_OPTIONS_LIBRARY} ${GMPXX_LIBRARIES} ${GMP_LIBRARIES}) add_executable ( performance_rips_persistence performance_rips_persistence.cpp ) @@ -52,53 +48,37 @@ if(GMP_FOUND) add_test(rips_multifield_persistence_2_71 ${CMAKE_CURRENT_BINARY_DIR}/rips_multifield_persistence ${CMAKE_SOURCE_DIR}/data/points/Kl.txt -r 0.2 -d 3 -p 2 -q 71 -m 100) endif(GMPXX_FOUND) -else() - # message(WARNING "GMP not found.") endif(GMP_FOUND) - if(CGAL_FOUND) - add_executable(alpha_complex_3d_persistence alpha_complex_3d_persistence.cpp) - target_link_libraries(alpha_complex_3d_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) - - if (TBB_FOUND) - target_link_libraries(alpha_complex_3d_persistence ${TBB_LIBRARIES}) - endif(TBB_FOUND) - add_test(alpha_complex_3d_persistence_2_0_5 ${CMAKE_CURRENT_BINARY_DIR}/alpha_complex_3d_persistence ${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off 2 0.45) - - - if (NOT CGAL_VERSION VERSION_LESS 4.7.0) - message(STATUS "CGAL version: ${CGAL_VERSION}.") - - find_package(Eigen3 3.1.0) - if (EIGEN3_FOUND) - message(STATUS "Eigen3 version: ${EIGEN3_VERSION}.") - include( ${EIGEN3_USE_FILE} ) - - add_executable (alpha_complex_persistence alpha_complex_persistence.cpp) - target_link_libraries(alpha_complex_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY} ${Boost_PROGRAM_OPTIONS_LIBRARY}) - - add_executable(periodic_alpha_complex_3d_persistence periodic_alpha_complex_3d_persistence.cpp) - target_link_libraries(periodic_alpha_complex_3d_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) - - add_executable(custom_persistence_sort custom_persistence_sort.cpp) - target_link_libraries(custom_persistence_sort ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) - - if (TBB_FOUND) - target_link_libraries(alpha_complex_persistence ${TBB_LIBRARIES}) - target_link_libraries(periodic_alpha_complex_3d_persistence ${TBB_LIBRARIES}) - target_link_libraries(custom_persistence_sort ${TBB_LIBRARIES}) - endif(TBB_FOUND) - add_test(alpha_complex_persistence_2_0_45 ${CMAKE_CURRENT_BINARY_DIR}/alpha_complex_persistence ${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off -m 0.45 -p 2) - add_test(periodic_alpha_complex_3d_persistence_2_0 ${CMAKE_CURRENT_BINARY_DIR}/periodic_alpha_complex_3d_persistence ${CMAKE_SOURCE_DIR}/data/points/grid_10_10_10_in_0_1.off ${CMAKE_SOURCE_DIR}/data/points/iso_cuboid_3_in_0_1.txt 2 0) - add_test(custom_persistence_sort ${CMAKE_CURRENT_BINARY_DIR}/custom_persistence_sort) - - else() - message(WARNING "Eigen3 not found. Version 3.1.0 is required for Alpha shapes feature.") - endif(EIGEN3_FOUND) - else() - message(WARNING "CGAL version: ${CGAL_VERSION} is too old to compile Alpha shapes feature. Version 4.6.0 is required.") - endif () - else() - # message(WARNING "CGAL not found.") - endif(CGAL_FOUND) - +if(CGAL_FOUND) + add_executable(alpha_complex_3d_persistence alpha_complex_3d_persistence.cpp) + target_link_libraries(alpha_complex_3d_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) + + if (TBB_FOUND) + target_link_libraries(alpha_complex_3d_persistence ${TBB_LIBRARIES}) + endif(TBB_FOUND) + add_test(alpha_complex_3d_persistence_2_0_5 ${CMAKE_CURRENT_BINARY_DIR}/alpha_complex_3d_persistence ${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off 2 0.45) + + + if (NOT CGAL_VERSION VERSION_LESS 4.7.0) + if (EIGEN3_FOUND) + add_executable (alpha_complex_persistence alpha_complex_persistence.cpp) + target_link_libraries(alpha_complex_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY} ${Boost_PROGRAM_OPTIONS_LIBRARY}) + + add_executable(periodic_alpha_complex_3d_persistence periodic_alpha_complex_3d_persistence.cpp) + target_link_libraries(periodic_alpha_complex_3d_persistence ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) + + add_executable(custom_persistence_sort custom_persistence_sort.cpp) + target_link_libraries(custom_persistence_sort ${Boost_SYSTEM_LIBRARY} ${CGAL_LIBRARY}) + + if (TBB_FOUND) + target_link_libraries(alpha_complex_persistence ${TBB_LIBRARIES}) + target_link_libraries(periodic_alpha_complex_3d_persistence ${TBB_LIBRARIES}) + target_link_libraries(custom_persistence_sort ${TBB_LIBRARIES}) + endif(TBB_FOUND) + add_test(alpha_complex_persistence_2_0_45 ${CMAKE_CURRENT_BINARY_DIR}/alpha_complex_persistence ${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off -m 0.45 -p 2) + add_test(periodic_alpha_complex_3d_persistence_2_0 ${CMAKE_CURRENT_BINARY_DIR}/periodic_alpha_complex_3d_persistence ${CMAKE_SOURCE_DIR}/data/points/grid_10_10_10_in_0_1.off ${CMAKE_SOURCE_DIR}/data/points/iso_cuboid_3_in_0_1.txt 2 0) + add_test(custom_persistence_sort ${CMAKE_CURRENT_BINARY_DIR}/custom_persistence_sort) + endif(EIGEN3_FOUND) + endif (NOT CGAL_VERSION VERSION_LESS 4.7.0) +endif(CGAL_FOUND) diff --git a/src/Persistent_cohomology/example/alpha_complex_persistence.cpp b/src/Persistent_cohomology/example/alpha_complex_persistence.cpp index cb181936..2412569a 100644 --- a/src/Persistent_cohomology/example/alpha_complex_persistence.cpp +++ b/src/Persistent_cohomology/example/alpha_complex_persistence.cpp @@ -4,6 +4,8 @@ #include <gudhi/Alpha_complex.h> #include <gudhi/Persistent_cohomology.h> +// to construct a simplex_tree from alpha complex +#include <gudhi/Simplex_tree.h> #include <iostream> #include <string> @@ -30,35 +32,38 @@ int main(int argc, char **argv) { // Init of an alpha complex from an OFF file // ---------------------------------------------------------------------------- using Kernel = CGAL::Epick_d< CGAL::Dynamic_dimension_tag >; - Gudhi::alpha_complex::Alpha_complex<Kernel> alpha_complex_from_file(off_file_points, alpha_square_max_value); - - // ---------------------------------------------------------------------------- - // Display information about the alpha complex - // ---------------------------------------------------------------------------- - std::cout << "Alpha complex is of dimension " << alpha_complex_from_file.dimension() << - " - " << alpha_complex_from_file.num_simplices() << " simplices - " << - alpha_complex_from_file.num_vertices() << " vertices." << std::endl; - - // Sort the simplices in the order of the filtration - alpha_complex_from_file.initialize_filtration(); - - std::cout << "Simplex_tree dim: " << alpha_complex_from_file.dimension() << std::endl; - // Compute the persistence diagram of the complex - Gudhi::persistent_cohomology::Persistent_cohomology< Gudhi::alpha_complex::Alpha_complex<Kernel>, - Gudhi::persistent_cohomology::Field_Zp > pcoh(alpha_complex_from_file); - // initializes the coefficient field for homology - pcoh.init_coefficients(coeff_field_characteristic); - - pcoh.compute_persistent_cohomology(min_persistence); - - // Output the diagram in filediag - if (output_file_diag.empty()) { - pcoh.output_diagram(); - } else { - std::cout << "Result in file: " << output_file_diag << std::endl; - std::ofstream out(output_file_diag); - pcoh.output_diagram(out); - out.close(); + Gudhi::alpha_complex::Alpha_complex<Kernel> alpha_complex_from_file(off_file_points); + + Gudhi::Simplex_tree<> simplex; + if (alpha_complex_from_file.create_complex(simplex, alpha_square_max_value)) { + // ---------------------------------------------------------------------------- + // Display information about the alpha complex + // ---------------------------------------------------------------------------- + std::cout << "Simplicial complex is of dimension " << simplex.dimension() << + " - " << simplex.num_simplices() << " simplices - " << + simplex.num_vertices() << " vertices." << std::endl; + + // Sort the simplices in the order of the filtration + simplex.initialize_filtration(); + + std::cout << "Simplex_tree dim: " << simplex.dimension() << std::endl; + // Compute the persistence diagram of the complex + Gudhi::persistent_cohomology::Persistent_cohomology< Gudhi::Simplex_tree<>, + Gudhi::persistent_cohomology::Field_Zp > pcoh(simplex); + // initializes the coefficient field for homology + pcoh.init_coefficients(coeff_field_characteristic); + + pcoh.compute_persistent_cohomology(min_persistence); + + // Output the diagram in filediag + if (output_file_diag.empty()) { + pcoh.output_diagram(); + } else { + std::cout << "Result in file: " << output_file_diag << std::endl; + std::ofstream out(output_file_diag); + pcoh.output_diagram(out); + out.close(); + } } return 0; diff --git a/src/Persistent_cohomology/example/custom_persistence_sort.cpp b/src/Persistent_cohomology/example/custom_persistence_sort.cpp index 9af38611..64f2a4dc 100644 --- a/src/Persistent_cohomology/example/custom_persistence_sort.cpp +++ b/src/Persistent_cohomology/example/custom_persistence_sort.cpp @@ -27,6 +27,8 @@ #include <gudhi/Alpha_complex.h> #include <gudhi/Persistent_cohomology.h> +// to construct a simplex_tree from alpha complex +#include <gudhi/Simplex_tree.h> #include <iostream> #include <iterator> @@ -38,6 +40,9 @@ using Kernel = CGAL::Epick_d< CGAL::Dimension_tag<3> >; using Point = Kernel::Point_d; using Alpha_complex = Gudhi::alpha_complex::Alpha_complex<Kernel>; +using Simplex_tree = Gudhi::Simplex_tree<>; +using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology< Simplex_tree, + Gudhi::persistent_cohomology::Field_Zp >; std::vector<Point> random_points() { // Instanciate a random point generator @@ -60,7 +65,7 @@ std::vector<Point> random_points() { * Compare two intervals by dimension, then by length. */ struct cmp_intervals_by_dim_then_length { - explicit cmp_intervals_by_dim_then_length(Alpha_complex * sc) + explicit cmp_intervals_by_dim_then_length(Simplex_tree * sc) : sc_(sc) { } template<typename Persistent_interval> @@ -71,46 +76,62 @@ struct cmp_intervals_by_dim_then_length { else return (sc_->dimension(get < 0 > (p1)) > sc_->dimension(get < 0 > (p2))); } - Alpha_complex* sc_; + Simplex_tree* sc_; }; int main(int argc, char **argv) { std::vector<Point> points = random_points(); + std::cout << "Points size=" << points.size() << std::endl; // Alpha complex persistence computation from generated points - Alpha_complex alpha_complex_from_points(points, 0.6); - - using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology< Alpha_complex, - Gudhi::persistent_cohomology::Field_Zp >; - Persistent_cohomology pcoh(alpha_complex_from_points); - - // initializes the coefficient field for homology - Z/3Z - pcoh.init_coefficients(3); - pcoh.compute_persistent_cohomology(0.2); - - // Custom sort and output persistence - cmp_intervals_by_dim_then_length cmp(&alpha_complex_from_points); - auto persistent_pairs = pcoh.get_persistent_pairs(); - std::sort(std::begin(persistent_pairs), std::end(persistent_pairs), cmp); - for (auto pair : persistent_pairs) { - std::cout << alpha_complex_from_points.dimension(get<0>(pair)) << " " - << alpha_complex_from_points.filtration(get<0>(pair)) << " " - << alpha_complex_from_points.filtration(get<1>(pair)) << std::endl; + Alpha_complex alpha_complex_from_points(points); + std::cout << "alpha_complex_from_points" << std::endl; + + Simplex_tree simplex; + std::cout << "simplex" << std::endl; + if (alpha_complex_from_points.create_complex(simplex, 0.6)) { + std::cout << "simplex" << std::endl; + // ---------------------------------------------------------------------------- + // Display information about the alpha complex + // ---------------------------------------------------------------------------- + std::cout << "Simplicial complex is of dimension " << simplex.dimension() << + " - " << simplex.num_simplices() << " simplices - " << + simplex.num_vertices() << " vertices." << std::endl; + + // Sort the simplices in the order of the filtration + simplex.initialize_filtration(); + + std::cout << "Simplex_tree dim: " << simplex.dimension() << std::endl; + + Persistent_cohomology pcoh(simplex); + + // initializes the coefficient field for homology - Z/3Z + pcoh.init_coefficients(3); + pcoh.compute_persistent_cohomology(0.2); + + // Custom sort and output persistence + cmp_intervals_by_dim_then_length cmp(&simplex); + auto persistent_pairs = pcoh.get_persistent_pairs(); + std::sort(std::begin(persistent_pairs), std::end(persistent_pairs), cmp); + for (auto pair : persistent_pairs) { + std::cout << simplex.dimension(get<0>(pair)) << " " + << simplex.filtration(get<0>(pair)) << " " + << simplex.filtration(get<1>(pair)) << std::endl; + } + + // Persistent Betti numbers + std::cout << "The persistent Betti numbers in interval [0.40, 0.41] are : "; + for (int dim = 0; dim < simplex.dimension(); dim++) + std::cout << "b" << dim << " = " << pcoh.persistent_betti_number(dim, 0.40, 0.41) << " ; "; + std::cout << std::endl; + + // Betti numbers + std::vector<int> betti_numbers = pcoh.betti_numbers(); + std::cout << "The Betti numbers are : "; + for (std::size_t i = 0; i < betti_numbers.size(); i++) + std::cout << "b" << i << " = " << betti_numbers[i] << " ; "; + std::cout << std::endl; } - - // Persistent Betti numbers - std::cout << "The persistent Betti numbers in interval [0.40, 0.41] are : "; - for (int dim = 0; dim < alpha_complex_from_points.dimension(); dim++) - std::cout << "b" << dim << " = " << pcoh.persistent_betti_number(dim, 0.40, 0.41) << " ; "; - std::cout << std::endl; - - // Betti numbers - std::vector<int> betti_numbers = pcoh.betti_numbers(); - std::cout << "The Betti numbers are : "; - for (std::size_t i = 0; i < betti_numbers.size(); i++) - std::cout << "b" << i << " = " << betti_numbers[i] << " ; "; - std::cout << std::endl; - return 0; } diff --git a/src/Persistent_cohomology/example/plain_homology.cpp b/src/Persistent_cohomology/example/plain_homology.cpp index 2156ffd9..ae82c817 100644 --- a/src/Persistent_cohomology/example/plain_homology.cpp +++ b/src/Persistent_cohomology/example/plain_homology.cpp @@ -25,6 +25,7 @@ #include <iostream> #include <vector> +#include <cstdint> // for std::uint8_t using namespace Gudhi; @@ -39,6 +40,8 @@ struct MyOptions : Simplex_tree_options_full_featured { static const bool store_key = true; // I have few vertices typedef short Vertex_handle; + // Maximum number of simplices to compute persistence is 2^8 - 1 = 255. One is reserved for null_key + typedef std::uint8_t Simplex_key; }; typedef Simplex_tree<MyOptions> ST; diff --git a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h index 5b371f2b..108a53d7 100644 --- a/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h +++ b/src/Persistent_cohomology/include/gudhi/Persistent_cohomology.h @@ -41,6 +41,7 @@ #include <tuple> #include <algorithm> #include <string> +#include <stdexcept> // for std::out_of_range namespace Gudhi { @@ -89,6 +90,7 @@ class Persistent_cohomology { * * @param[in] cpx Complex for which the persistent homology is computed. * cpx is a model of FilteredComplex + * @exception std::out_of_range In case the number of simplices is more than Simplex_key type numeric limit. */ explicit Persistent_cohomology(Complex_ds& cpx) : cpx_(&cpx), @@ -106,6 +108,10 @@ class Persistent_cohomology { interval_length_policy(&cpx, 0), column_pool_(), // memory pools for the CAM cell_pool_() { + if (cpx_->num_simplices() > std::numeric_limits<Simplex_key>::max()) { + // num_simplices must be strictly lower than the limit, because a value is reserved for null_key. + throw std::out_of_range ("The number of simplices is more than Simplex_key type numeric limit."); + } Simplex_key idx_fil = 0; for (auto sh : cpx_->filtration_simplex_range()) { cpx_->assign_key(sh, idx_fil); diff --git a/src/Persistent_cohomology/test/betti_numbers_unit_test.cpp b/src/Persistent_cohomology/test/betti_numbers_unit_test.cpp index a4e00b45..40221005 100644 --- a/src/Persistent_cohomology/test/betti_numbers_unit_test.cpp +++ b/src/Persistent_cohomology/test/betti_numbers_unit_test.cpp @@ -12,7 +12,7 @@ #include "gudhi/Simplex_tree.h" #include "gudhi/Persistent_cohomology.h" -struct MyOptions : Gudhi::Simplex_tree_options_full_featured { +struct MiniSTOptions : Gudhi::Simplex_tree_options_full_featured { // Implicitly use 0 as filtration value for all simplices static const bool store_filtration = false; // The persistence algorithm needs this @@ -21,7 +21,7 @@ struct MyOptions : Gudhi::Simplex_tree_options_full_featured { typedef short Vertex_handle; }; -using Mini_simplex_tree = Gudhi::Simplex_tree<MyOptions>; +using Mini_simplex_tree = Gudhi::Simplex_tree<MiniSTOptions>; using Mini_st_persistence = Gudhi::persistent_cohomology::Persistent_cohomology<Mini_simplex_tree, Gudhi::persistent_cohomology::Field_Zp>; diff --git a/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp b/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp index 5c54ed5f..6efd749e 100644 --- a/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp +++ b/src/Persistent_cohomology/test/persistent_cohomology_unit_test.cpp @@ -4,6 +4,7 @@ #include <utility> // std::pair, std::make_pair #include <cmath> // float comparison #include <limits> +#include <cstdint> // for std::uint8_t #define BOOST_TEST_DYN_LINK #define BOOST_TEST_MODULE "persistent_cohomology" @@ -172,3 +173,45 @@ BOOST_AUTO_TEST_CASE( rips_persistent_cohomology_single_field_dim_5 ) // std::string str_rips_persistence = test_rips_persistence(6, 0); // TODO(VR): division by zero // std::string str_rips_persistence = test_rips_persistence(0, 0); + +/** SimplexTree minimal options to test the limits. + * + * Maximum number of simplices to compute persistence is <CODE>std::numeric_limits<std::uint8_t>::max()<\CODE> = 256.*/ +struct MiniSTOptions { + typedef linear_indexing_tag Indexing_tag; + typedef short Vertex_handle; + typedef double Filtration_value; + // Maximum number of simplices to compute persistence is 2^8 - 1 = 255. One is reserved for null_key + typedef std::uint8_t Simplex_key; + static const bool store_key = true; + static const bool store_filtration = false; + static const bool contiguous_vertices = false; +}; + +using Mini_simplex_tree = Gudhi::Simplex_tree<MiniSTOptions>; +using Mini_st_persistence = + Gudhi::persistent_cohomology::Persistent_cohomology<Mini_simplex_tree, Gudhi::persistent_cohomology::Field_Zp>; + +BOOST_AUTO_TEST_CASE( persistence_constructor_exception ) +{ + Mini_simplex_tree st; + + // To make number of simplices = 255 + const short simplex_0[] = {0, 1, 2, 3, 4, 5, 6, 7}; + st.insert_simplex_and_subfaces(simplex_0); + // FIXME: Remove this line + st.set_dimension(8); + + // Sort the simplices in the order of the filtration + st.initialize_filtration(); + + BOOST_CHECK(st.num_simplices() <= std::numeric_limits<MiniSTOptions::Simplex_key>::max()); + // Class for homology computation + BOOST_CHECK_NO_THROW(Mini_st_persistence pcoh(st)); + + st.insert_simplex({8}); + BOOST_CHECK(st.num_simplices() > std::numeric_limits<MiniSTOptions::Simplex_key>::max()); + // Class for homology computation + BOOST_CHECK_THROW(Mini_st_persistence pcoh2(st), std::out_of_range); + +} |