From 0423b7024dee787659b76fff4b4f659546a40aea Mon Sep 17 00:00:00 2001 From: vrouvrea Date: Thu, 29 Jun 2017 15:57:04 +0000 Subject: First working version of expansion insertion (different from rips expansion). git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/graph_expansion_with_blocker_oracle@2570 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 6a6bb4052c3111e783e9e138619f1e945c856708 --- src/Simplex_tree/example/CMakeLists.txt | 17 +++ src/Simplex_tree/example/block.cpp | 82 +++++++++++ src/Simplex_tree/example/simple_simplex_tree.cpp | 35 ++++- src/Simplex_tree/include/gudhi/Simplex_tree.h | 167 ++++++++++++++++------- 4 files changed, 251 insertions(+), 50 deletions(-) create mode 100644 src/Simplex_tree/example/block.cpp diff --git a/src/Simplex_tree/example/CMakeLists.txt b/src/Simplex_tree/example/CMakeLists.txt index 5dbbfcc0..d5f512b3 100644 --- a/src/Simplex_tree/example/CMakeLists.txt +++ b/src/Simplex_tree/example/CMakeLists.txt @@ -45,3 +45,20 @@ endif() add_test(NAME Simplex_tree_example_cech_complex_step_by_step COMMAND $ "${CMAKE_SOURCE_DIR}/data/points/alphacomplexdoc.off" "-r" "12." "-d" "3") install(TARGETS Simplex_tree_example_cech_complex_step_by_step DESTINATION bin) + + +# +# TO BE REMOVED !! +# + +#add_executable ( rips_step_by_step rips_step_by_step.cpp) +#target_link_libraries(rips_step_by_step ${Boost_PROGRAM_OPTIONS_LIBRARY}) +#if (TBB_FOUND) +# target_link_libraries(rips_step_by_step ${TBB_LIBRARIES}) +#endif() + + +add_executable ( block block.cpp ) +if (TBB_FOUND) + target_link_libraries(block ${TBB_LIBRARIES}) +endif() diff --git a/src/Simplex_tree/example/block.cpp b/src/Simplex_tree/example/block.cpp new file mode 100644 index 00000000..07ec3921 --- /dev/null +++ b/src/Simplex_tree/example/block.cpp @@ -0,0 +1,82 @@ +/* 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): Vincent Rouvreau + * + * Copyright (C) 2014 + * + * 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 . + */ + +#include +#include + +#include +#include // for pair +#include + +using Simplex_tree = Gudhi::Simplex_tree<>; +using Vertex_handle = Simplex_tree::Vertex_handle; +using Filtration_value = Simplex_tree::Filtration_value; +using typeVectorVertex = std::vector< Vertex_handle >; +using typePairSimplexBool = std::pair< Simplex_tree::Simplex_handle, bool >; + +int main(int argc, char * const argv[]) { + + // Construct the Simplex Tree + Simplex_tree simplexTree; + + simplexTree.insert_simplex({0, 1}); + simplexTree.insert_simplex({0, 2}); + simplexTree.insert_simplex({0, 3}); + simplexTree.insert_simplex({1, 2}); + simplexTree.insert_simplex({1, 3}); + simplexTree.insert_simplex({2, 3}); + simplexTree.insert_simplex({2, 4}); + simplexTree.insert_simplex({3, 6}); + simplexTree.insert_simplex({4, 5}); + simplexTree.insert_simplex({4, 6}); + simplexTree.insert_simplex({5, 6}); + simplexTree.insert_simplex({6}); + + std::cout << "********************************************************************\n"; + // Display the Simplex_tree - Can not be done in the middle of 2 inserts + std::cout << "* The complex contains " << simplexTree.num_simplices() << " simplices\n"; + std::cout << " - dimension " << simplexTree.dimension() << " - filtration " << simplexTree.filtration() << "\n"; + std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n"; + for (auto f_simplex : simplexTree.filtration_simplex_range()) { + std::cout << " " << "[" << simplexTree.filtration(f_simplex) << "] "; + for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) + std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + + simplexTree.expansion_with_blockers(3, [](){return true;}); + + simplexTree.initialize_filtration(); + std::cout << "********************************************************************\n"; + // Display the Simplex_tree - Can not be done in the middle of 2 inserts + std::cout << "* The complex contains " << simplexTree.num_simplices() << " simplices\n"; + std::cout << " - dimension " << simplexTree.dimension() << " - filtration " << simplexTree.filtration() << "\n"; + std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n"; + for (auto f_simplex : simplexTree.filtration_simplex_range()) { + std::cout << " " << "[" << simplexTree.filtration(f_simplex) << "] "; + for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) + std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + + return 0; +} diff --git a/src/Simplex_tree/example/simple_simplex_tree.cpp b/src/Simplex_tree/example/simple_simplex_tree.cpp index 60f9a35e..f27d7ab8 100644 --- a/src/Simplex_tree/example/simple_simplex_tree.cpp +++ b/src/Simplex_tree/example/simple_simplex_tree.cpp @@ -195,9 +195,8 @@ int main(int argc, char * const argv[]) { std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n"; for (auto f_simplex : simplexTree.filtration_simplex_range()) { std::cout << " " << "[" << simplexTree.filtration(f_simplex) << "] "; - for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) { - std::cout << static_cast(vertex) << " "; - } + for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) + std::cout << "(" << vertex << ")"; std::cout << std::endl; } // [0.1] 0 @@ -250,5 +249,35 @@ int main(int argc, char * const argv[]) { std::cout << "***+ YES IT IS!\n"; else std::cout << "***- NO IT ISN'T\n"; + + invSimplexVector = { 0, 1 }; + simplexFound = simplexTree.find({ 0, 1 }); + std::cout << "**************IS THE SIMPLEX {0,1} IN THE SIMPLEX TREE ?\n"; + if (simplexFound != simplexTree.null_simplex()) + std::cout << "***+ YES IT IS!\n"; + else + std::cout << "***- NO IT ISN'T\n"; + + std::cout << "**************COFACES OF {0,1} IN CODIMENSION 1 ARE\n"; + for (auto& simplex : simplexTree.cofaces_simplex_range(simplexTree.find({0,1}), 1)) { + for (auto vertex : simplexTree.simplex_vertex_range(simplex)) + std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + + std::cout << "**************STARS OF {0,1} ARE\n"; + for (auto& simplex : simplexTree.star_simplex_range(simplexTree.find({0,1}))) { + for (auto vertex : simplexTree.simplex_vertex_range(simplex)) + std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + + std::cout << "**************BOUNDARIES OF {0,1,2} ARE\n"; + for (auto& simplex : simplexTree.boundary_simplex_range(simplexTree.find({0,1,2}))) { + for (auto vertex : simplexTree.simplex_vertex_range(simplex)) + std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + return 0; } diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h index dbed47b8..ff31fd89 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h @@ -1007,36 +1007,11 @@ class Simplex_tree { * The Simplex_tree must contain no simplex of dimension bigger than * 1 when calling the method. */ void expansion(int max_dim) { - expansion_with_blockers(max_dim, - [](Simplex_handle origin_sh, - Simplex_handle dict1_sh, - Simplex_handle dict2_sh, - Siblings* siblings) { - // Default blocker is always insert with the maximal filtration value between - // origin, dict1 and dict2 - return std::make_pair(true, (std::max)({origin_sh->second.filtration(), - dict1_sh->second.filtration(), - dict2_sh->second.filtration()})); }); - } - - /** \brief Expands the Simplex_tree containing only its one skeleton - * until dimension max_dim. - * - * The expanded simplicial complex until dimension \f$d\f$ - * attached to a graph \f$G\f$ is the maximal simplicial complex of - * dimension at most \f$d\f$ admitting the graph \f$G\f$ as \f$1\f$-skeleton. - * The filtration value assigned to a simplex is the maximal filtration - * value of one of its edges. - * - * The Simplex_tree must contain no simplex of dimension bigger than - * 1 when calling the method. */ - template< typename Blocker > - void expansion_with_blockers(int max_dim, Blocker blocker_expansion_function) { dimension_ = max_dim; for (Dictionary_it root_it = root_.members_.begin(); root_it != root_.members_.end(); ++root_it) { if (has_children(root_it)) { - siblings_expansion_with_blockers(root_it->second.children(), max_dim - 1, blocker_expansion_function); + siblings_expansion(root_it->second.children(), max_dim - 1); } } dimension_ = max_dim - dimension_; @@ -1044,9 +1019,8 @@ class Simplex_tree { private: /** \brief Recursive expansion of the simplex tree.*/ - template< typename Blocker > - void siblings_expansion_with_blockers(Siblings* siblings, // must contain elements - int k, Blocker blocker_expansion_function) { + void siblings_expansion(Siblings * siblings, // must contain elements + int k) { if (dimension_ > k) { dimension_ = k; } @@ -1060,22 +1034,20 @@ class Simplex_tree { s_h != siblings->members().end(); ++s_h, ++next) { Simplex_handle root_sh = find_vertex(s_h->first); if (has_children(root_sh)) { - intersection_with_blockers( - inter, // output intersection - next, // begin - siblings->members().end(), // end - root_sh->second.children()->members().begin(), - root_sh->second.children()->members().end(), - s_h, siblings, - blocker_expansion_function); + intersection( + inter, // output intersection + next, // begin + siblings->members().end(), // end + root_sh->second.children()->members().begin(), + root_sh->second.children()->members().end(), + s_h->second.filtration()); if (inter.size() != 0) { Siblings * new_sib = new Siblings(siblings, // oncles s_h->first, // parent inter); // boost::container::ordered_unique_range_t - // As siblings_expansion_with_blockers is recusively called, inter must be cleared before inter.clear(); s_h->second.assign_children(new_sib); - siblings_expansion_with_blockers(new_sib, k - 1, blocker_expansion_function); + siblings_expansion(new_sib, k - 1); } else { // ensure the children property s_h->second.assign_children(siblings); @@ -1087,20 +1059,16 @@ class Simplex_tree { /** \brief Intersects Dictionary 1 [begin1;end1) with Dictionary 2 [begin2,end2) * and assigns the maximal possible Filtration_value to the Nodes. */ - template< typename Blocker > - static void intersection_with_blockers(std::vector >& intersection, + static void intersection(std::vector >& intersection, Dictionary_it begin1, Dictionary_it end1, Dictionary_it begin2, Dictionary_it end2, - Dictionary_it origin_sh, - Siblings* siblings, - Blocker blocker_expansion_function) { + Filtration_value filtration_) { if (begin1 == end1 || begin2 == end2) return; // ----->> while (true) { if (begin1->first == begin2->first) { - std::pair blocker_result = blocker_expansion_function(origin_sh, begin1, begin2, siblings); - if (blocker_result.first) - intersection.emplace_back(begin1->first, Node(nullptr, blocker_result.second)); + Filtration_value filt = (std::max)({begin1->second.filtration(), begin2->second.filtration(), filtration_}); + intersection.emplace_back(begin1->first, Node(nullptr, filt)); if (++begin1 == end1 || ++begin2 == end2) return; // ----->> } else if (begin1->first < begin2->first) { @@ -1113,6 +1081,111 @@ class Simplex_tree { } } + + + /*-------------------------------------------------------------------------------------------------------------------------*/ + /*-------------------------------------------------------------------------------------------------------------------------*/ + /*-------------------------------------------------------------------------------------------------------------------------*/ + + public: + /** \brief Expands the Simplex_tree containing only its one skeleton + * until dimension max_dim. + * + * The expanded simplicial complex until dimension \f$d\f$ + * attached to a graph \f$G\f$ is the maximal simplicial complex of + * dimension at most \f$d\f$ admitting the graph \f$G\f$ as \f$1\f$-skeleton. + * The filtration value assigned to a simplex is the maximal filtration + * value of one of its edges. + * + * The Simplex_tree must contain no simplex of dimension bigger than + * 1 when calling the method. */ + template< typename Blocker > + void expansion_with_blockers(int max_dim, Blocker blocker_expansion_function) { + dimension_ = max_dim; + // Loop must be from the end to the beginning, as higher dimension simplex are always on the left part of the tree + for (auto& simplex : boost::adaptors::reverse(root_.members())) { + if (has_children(&simplex)) { + std::cout << " *** root on " << static_cast(simplex.first) << std::endl; + siblings_expansion_with_blockers(simplex.second.children(), max_dim - 1, blocker_expansion_function); + } + } + dimension_ = max_dim - dimension_; + } + + private: + /** \brief Recursive expansion of the simplex tree.*/ + template< typename Blocker > + void siblings_expansion_with_blockers(Siblings* siblings, // must contain elements + int k, Blocker blocker_expansion_function) { + if (dimension_ > k) { + dimension_ = k; + } + if (k == 0) + return; + // No need to go deeper + if (siblings->members().size() < 2) + return; + // Reverse loop starting before the last one for 'next' to be the last one + for (auto simplex = siblings->members().rbegin() + 1; simplex != siblings->members().rend(); simplex++) { + auto next = siblings->members().rbegin(); + std::vector > intersection; + while(next != simplex) { + bool to_be_inserted = true; + std::cout << "to_be_inserted = " << to_be_inserted << " dim = " << k << " simplex = " << simplex->first << " - next = " << next->first << std::endl; + + for (auto& border : boundary_simplex_range(simplex)) { + to_be_inserted = to_be_inserted && find_child(border, next->first); + + for (auto vertex : simplex_vertex_range(border)) { + std::cout << "(" << vertex << ")"; + } + std::cout << " | "; + } + std::cout << std::endl; + if (to_be_inserted) { + std::cout << next->first << " to be inserted." << std::endl; + intersection.emplace_back(next->first, Node(nullptr, 0.0)); + } + + // loop until simplex is reached + next++; + } + if (intersection.size() != 0) { + // Reverse the order to insert + std::reverse(std::begin(intersection), std::end(intersection)); + Siblings * new_sib = new Siblings(siblings, // oncles + simplex->first, // parent + intersection); // boost::container::ordered_unique_range_t + // intersection must be cleared before the function to be called recursively + intersection.clear(); + simplex->second.assign_children(new_sib); + siblings_expansion_with_blockers(new_sib, k - 1, blocker_expansion_function); + } else { + // ensure the children property + simplex->second.assign_children(siblings); + intersection.clear(); + } + + } + + } + + /** \private Returns true if vh is a member of sh*/ + bool find_child(Simplex_handle sh, Vertex_handle vh) { + std::vector child = {vh}; + std::cout << "+" << vh; + for (auto vertex : simplex_vertex_range(sh)) { + std::cout << "+" << vertex; + child.push_back(vertex); + } + std::cout << " => " << (find(child) != null_simplex()) << "___ "; + return find(child) != null_simplex(); + } + + /*-------------------------------------------------------------------------------------------------------------------------*/ + /*-------------------------------------------------------------------------------------------------------------------------*/ + /*-------------------------------------------------------------------------------------------------------------------------*/ + public: /** \brief Write the hasse diagram of the simplicial complex in os. * -- cgit v1.2.3