From 90fb612fc2f25bf092d7a7d45de329f913fe3bef Mon Sep 17 00:00:00 2001 From: salinasd Date: Fri, 30 Jan 2015 11:52:29 +0000 Subject: skbl constructor from list a simplices more efficient git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/trunk@448 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 3e3b30cd496da0a96df0379b11acf554aad6b5ae --- src/GudhUI/utils/Persistence_compute.h | 2 +- .../Skeleton_blocker_simple_traits.h | 2 + .../include/gudhi/Skeleton_blocker/internal/Trie.h | 199 ++ .../Skeleton_blockers_simplices_iterators.h | 87 +- .../include/gudhi/Skeleton_blocker_complex.h | 2450 ++++++++++---------- .../gudhi/Skeleton_blocker_simplifiable_complex.h | 11 +- src/Skeleton_blocker/test/CMakeLists.txt | 2 + .../test/TestSkeletonBlockerComplex.cpp | 119 +- src/common/include/gudhi/Utils.h | 2 +- 9 files changed, 1497 insertions(+), 1377 deletions(-) create mode 100644 src/Skeleton_blocker/include/gudhi/Skeleton_blocker/internal/Trie.h diff --git a/src/GudhUI/utils/Persistence_compute.h b/src/GudhUI/utils/Persistence_compute.h index 0ddc205a..50d340b8 100644 --- a/src/GudhUI/utils/Persistence_compute.h +++ b/src/GudhUI/utils/Persistence_compute.h @@ -87,7 +87,7 @@ public: Gudhi::persistent_cohomology::Persistent_cohomology< Gudhi::Simplex_tree<>,Gudhi::persistent_cohomology::Field_Zp > pcoh (st); pcoh.init_coefficients( params.p ); //initilizes the coefficient field for homology - pcoh.compute_persistent_cohomology( INFINITY ); //put params.min_persistence + pcoh.compute_persistent_cohomology( params.min_pers ); //put params.min_pers stream_<<"persistence: \n"; stream_<<"p dimension birth death: \n"; diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h index b8506b2c..52e454ea 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h @@ -77,6 +77,8 @@ struct Skeleton_blocker_simple_traits { : vertex(val) { } + operator int() const { return (int)vertex; } + boost_vertex_handle vertex; bool operator==(const Vertex_handle& other) const { diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/internal/Trie.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/internal/Trie.h new file mode 100644 index 00000000..03650f73 --- /dev/null +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/internal/Trie.h @@ -0,0 +1,199 @@ +/* + * Trie.h + * Created on: Jan 29, 2015 + * 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): David Salinas + * + * Copyright (C) 2014 INRIA Sophia Antipolis-Méditerranée (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 TRIE_H_ +#define TRIE_H_ + +#include +#include + +namespace Gudhi { + + +namespace skbl { + + +template +struct Trie{ + typedef typename SkeletonBlockerComplex::Vertex_handle Vertex_handle; + typedef typename SkeletonBlockerComplex::Simplex_handle Simplex_handle; + + Vertex_handle v; + std::vector > childs; + //std::vector > childs; -> use of deleted function +private: + const Trie* parent_; +public: + Trie():parent_(0){} + Trie(Vertex_handle v_):v(v_),parent_(0){} + + Trie(Vertex_handle v_,Trie* parent):v(v_),parent_(parent){} + + + bool operator==(const Trie& other) const{ + return (v == other.v) ; + } + + void add_child(Trie* child){ + if(child){ + std::shared_ptr ptr_to_add(child); + childs.push_back(ptr_to_add); + child->parent_ = this; + } + } + typedef typename Simplex_handle::Simplex_vertex_const_iterator Simplex_vertex_const_iterator; + + + Trie* make_trie(Simplex_vertex_const_iterator s_it,Simplex_vertex_const_iterator s_end){ + if(s_it == s_end) return 0; + else{ + Trie* res = new Trie(*s_it); + Trie* child = make_trie(++s_it,s_end); + res->add_child(child); + return res; + } + } +private: + //go down recursively in the tree while advancing the simplex iterator. + //when it reaches a leaf, it inserts the remaining that is not present + void add_simplex_helper(Simplex_vertex_const_iterator s_it,Simplex_vertex_const_iterator s_end){ + assert(*s_it == v); + ++s_it; + if(s_it==s_end) return ; + if(!is_leaf()){ + for(auto child : childs){ + if(child->v == *s_it) + return child->add_simplex_helper(s_it,s_end); + } + //s_it is not found and needs to be inserted + } + //not leaf -> remaining of s needs to be inserted + Trie* son_with_what_remains_of_s(make_trie(s_it,s_end)); + add_child(son_with_what_remains_of_s); + return; + } + + void maximal_faces_helper(std::vector& res) const{ + if(is_leaf()) res.push_back(simplex()); + else + for(auto child : childs) + child->maximal_faces_helper(res); + } + +public: + /** + * adds the simplex to the trie + */ + void add_simplex(const Simplex_handle& s){ + if(s.empty()) return; + assert(v==s.first_vertex()); + add_simplex_helper(s.begin(),s.end()); + } + + std::vector maximal_faces() const{ + std::vector res; + maximal_faces_helper(res); + return res; + } + + /** + * Goes to the root in the trie to consitute simplex + */ + void add_vertices_up_to_the_root(Simplex_handle& res) const{ + res.add_vertex(v); + if(parent_) + parent_->add_vertices_up_to_the_root(res); + } + + Simplex_handle simplex() const{ + Simplex_handle res; + add_vertices_up_to_the_root(res); + return res; + } + + bool is_leaf() const{ + return childs.empty(); + } + + bool is_root() const{ + return parent_==0; + } + + const Trie* parent() { + return parent_; + } + + void remove_leaf() { + assert(is_leaf); + if(!is_root()) + parent_->childs.erase(this); + } + + /** + * true iff the simplex corresponds to one node in the trie + */ + bool contains(const Simplex_handle& s) const{ + Trie const* current = this; + if(s.empty()) return true; + if(current->v != s.first_vertex()) return false; + auto s_pos = s.begin(); + ++s_pos; + while(s_pos != s.end() && current != 0){ + bool found = false; + for(const auto child : current->childs){ + if(child->v == *s_pos) { + ++s_pos; + current = child.get(); + found = true; + break; + } + } + if(!found) return false; + } + return current!=0; + } + + Trie* go_bottom_left(){ + if(is_leaf()) + return this; + else + return (*childs.begin())->go_bottom_left(); + } + + friend std::ostream& operator<<(std::ostream& stream, const Trie& trie){ + stream<< "T( "<< trie.v<< " "; + for(auto t : trie.childs) + stream << *t ; + stream<<")"; + return stream; + } +}; +} + +} + +#endif /* TRIE_H_ */ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_simplices_iterators.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_simplices_iterators.h index 0681d5da..1f5ca238 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_simplices_iterators.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_simplices_iterators.h @@ -32,7 +32,7 @@ #include "gudhi/Skeleton_blocker_link_complex.h" #include "gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h" - +#include "gudhi/Skeleton_blocker/internal/Trie.h" namespace Gudhi { @@ -68,91 +68,8 @@ class Simplex_around_vertex_iterator : // Link_vertex_handle == Complex_Vertex_handle but this renaming helps avoiding confusion -private: - - struct Trie{ - Vertex_handle v; - std::vector > childs; - private: - const Trie* parent_; - public: - - - //std::list > childs; -> use of deleted function - - Trie():parent_(0){} - Trie(Vertex_handle v_):v(v_),parent_(0){ - } - - Trie(Vertex_handle v_,Trie* parent):v(v_),parent_(parent){ - } - - - bool operator==(const Trie& other) const{ - return (v == other.v) ; - } - - void add_child(Trie* child){ - if(child){ - std::shared_ptr ptr_to_add(child); - childs.push_back(ptr_to_add); - child->parent_ = this; - } - } - - - friend std::ostream& operator<<(std::ostream& stream, const Trie& trie){ - stream<< "T( "<< trie.v<< " "; - for(auto t : trie.childs) - stream << *t ; - stream<<")"; - return stream; - } - - // goes to the root in the trie to consitute simplex - void add_vertices_up_to_the_root(Simplex_handle& res) const{ - res.add_vertex(v); - if(parent_) - parent_->add_vertices_up_to_the_root(res); - } - - Simplex_handle simplex() const{ - Simplex_handle res; - add_vertices_up_to_the_root(res); - return res; - } - - bool is_leaf() const{ - return childs.empty(); - } - - bool is_root() const{ - return parent_==0; - } - - const Trie* parent() { - return parent_; - } - - void remove_leaf() { - assert(is_leaf); - if(!is_root()) - parent_->childs.erase(this); - } - - private: - - - public: - - Trie* go_bottom_left(){ - if(is_leaf()) - return this; - else - return (*childs.begin())->go_bottom_left(); - } + typedef typename Gudhi::skbl::Trie Trie; - }; private: const Complex* complex; diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h index 7c0bc652..a013ba52 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h @@ -48,6 +48,8 @@ #include "gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h" #include "gudhi/Skeleton_blocker/internal/Top_faces.h" +#include "gudhi/Skeleton_blocker/internal/Trie.h" + #include "gudhi/Utils.h" namespace Gudhi { @@ -61,1275 +63,1199 @@ namespace skbl { */ template class Skeleton_blocker_complex { - template friend class Complex_vertex_iterator; - template friend class Complex_neighbors_vertices_iterator; - template friend class Complex_edge_iterator; - template friend class Complex_edge_around_vertex_iterator; - - template friend class Skeleton_blocker_link_complex; - template friend class Skeleton_blocker_link_superior; - template friend class Skeleton_blocker_sub_complex; - - public: - /** - * @brief The type of stored vertex node, specified by the template SkeletonBlockerDS - */ - typedef typename SkeletonBlockerDS::Graph_vertex Graph_vertex; - - /** - * @brief The type of stored edge node, specified by the template SkeletonBlockerDS - */ - typedef typename SkeletonBlockerDS::Graph_edge Graph_edge; - - typedef typename SkeletonBlockerDS::Root_vertex_handle Root_vertex_handle; - - /** - * @brief The type of an handle to a vertex of the complex. - */ - typedef typename SkeletonBlockerDS::Vertex_handle Vertex_handle; - typedef typename Root_vertex_handle::boost_vertex_handle boost_vertex_handle; - - /** - * @brief A ordered set of integers that represents a simplex. - */ - typedef Skeleton_blocker_simplex Simplex_handle; - typedef Skeleton_blocker_simplex Root_simplex_handle; - - /** - * @brief Handle to a blocker of the complex. - */ - typedef Simplex_handle* Blocker_handle; - - typedef typename Root_simplex_handle::Simplex_vertex_const_iterator Root_simplex_iterator; - typedef typename Simplex_handle::Simplex_vertex_const_iterator Simplex_handle_iterator; - - protected: - typedef typename boost::adjacency_list Graph; - // todo/remark : edges are not sorted, it heavily penalizes computation for SuperiorLink - // (eg Link with greater vertices) - // that burdens simplex iteration / complex initialization via list of simplices. - // to avoid that, one should modify the graph by storing two lists of adjacency for every - // vertex, the one with superior and the one with lower vertices, that way there is - // no more extra cost for computation of SuperiorLink - typedef typename boost::graph_traits::vertex_iterator boost_vertex_iterator; - typedef typename boost::graph_traits::edge_iterator boost_edge_iterator; - - protected: - typedef typename boost::graph_traits::adjacency_iterator boost_adjacency_iterator; - - public: - /** - * @brief Handle to an edge of the complex. - */ - typedef typename boost::graph_traits::edge_descriptor Edge_handle; - - protected: - typedef std::multimap BlockerMap; - typedef typename std::multimap::value_type BlockerPair; - typedef typename std::multimap::iterator BlockerMapIterator; - typedef typename std::multimap::const_iterator BlockerMapConstIterator; - - protected: - int num_vertices_; - int num_blockers_; - - typedef Skeleton_blocker_complex_visitor Visitor; - // typedef Visitor* Visitor_ptr; - Visitor* visitor; - - /** - * @details If 'x' is a Vertex_handle of a vertex in the complex then degree[x] = d is its degree. - * - * This quantity is updated when adding/removing edge. - * - * This is useful because the operation - * list.size() is done in linear time. - */ - std::vector degree_; - Graph skeleton; /** 1-skeleton of the simplicial complex. */ - - // todo remove!!! - - /** Each vertex can access to the blockers passing through it. */ - BlockerMap blocker_map_; - - public: - ///////////////////////////////////////////////////////////////////////////// - /** @name Constructors, Destructors - */ - //@{ - /** - *@brief constructs a simplicial complex with a given number of vertices and a visitor. - */ - explicit Skeleton_blocker_complex(int num_vertices_ = 0, Visitor* visitor_ = NULL) - : visitor(visitor_) { - clear(); - for (int i = 0; i < num_vertices_; ++i) { - add_vertex(); - } - } - - private: - /** - * this nested class is used for the next constructor that takes as an input a list of simplices. - * It stores a vector where the ith cell contains a set of i-simplices - * - */ - class Simplices_sets_from_list { - private: - typedef std::set Container_simplices; - typedef typename Container_simplices::iterator Simplices_iterator; - int dimension_; - std::vector simplices_; - - public: - explicit Simplices_sets_from_list(std::list& simplices) - : dimension_(-1) { - assert(!simplices.empty()); - - for (auto simplex = simplices.begin(); simplex != simplices.end(); - ++simplex) { - dimension_ = std::max(dimension_, static_cast(simplex->dimension())); - } - simplices_ = std::vector < Container_simplices > (dimension_ + 1); - - // compute k-simplices - for (auto simplex = simplices.begin(); simplex != simplices.end(); - ++simplex) { - simplices_[simplex->dimension()].insert(*simplex); - } - } - - Simplices_iterator begin(int k) { - assert(0 <= k && k <= dimension_); - return simplices_[k].begin(); - } - - Simplices_iterator end(int k) { - assert(0 <= k && k <= dimension_); - return simplices_[k].end(); - } - - Container_simplices& simplices(int k) { - return simplices_[k]; - } - - int dimension() { - return dimension_; - } - - bool contains(const Simplex_handle& simplex) const { - if (simplex.dimension() > dimension_) - return false; - else - return simplices_[simplex.dimension()].find(simplex) - != simplices_[simplex.dimension()].end(); - } - - void print() { - for (int i = 0; i < dimension_; ++i) { - std::cout << i << "-simplices" << std::endl; - auto l = simplices_[i]; - for (auto s : l) { - std::cout << s << std::endl; - } - } - } - }; - - void compute_next_expand(Simplices_sets_from_list& simplices, int dim, - std::list& next_expand) { - next_expand.clear(); - - // high_resolution_clock::time_point tbeginfor = high_resolution_clock::now(); - // auto durationlooplink = std::chrono::duration_cast( tbeginfor - tbeginfor ).count(); - - for (auto sigma = simplices.begin(dim); sigma != simplices.end(dim); - ++sigma) { - // high_resolution_clock::time_point tbeg = high_resolution_clock::now(); - Simplex_handle t(*sigma); - Skeleton_blocker_link_superior link(*this, t); - // xxx all time here, most likely because accessing superior edges needs passing through lower one - // currently - - // durationlooplink += std::chrono::duration_cast( high_resolution_clock::now() - tbeg ).count(); - - for (auto v : link.vertex_range()) { - Vertex_handle v_in_complex(*this->get_address(link.get_id(v))); - t.add_vertex(v_in_complex); - next_expand.push_back(t); - t.remove_vertex(v_in_complex); - } - } - // high_resolution_clock::time_point t2 = high_resolution_clock::now(); - // auto durationlooptotal = std::chrono::duration_cast( t2 - tbeginfor ).count(); - // DBGVALUE(durationlooptotal); - // DBGVALUE(durationlooplink); - } - - public: - /** - * @brief Constructor with a list of simplices. - * @details The list of simplices must be the list - * of simplices of a simplicial complex. - */ - Skeleton_blocker_complex(std::list& simplices, - Visitor* visitor_ = NULL) - : num_vertices_(0), - num_blockers_(0), - visitor(visitor_) { - Simplices_sets_from_list set_simplices(simplices); - - int dim = set_simplices.dimension(); - - // add 1-skeleton to the complex - for (auto v_it = set_simplices.begin(0); v_it != set_simplices.end(0); - ++v_it) - add_vertex(); - - for (auto e_it = set_simplices.begin(1); e_it != set_simplices.end(1); - ++e_it) { - Vertex_handle a = e_it->first_vertex(); - Vertex_handle b = e_it->last_vertex(); - assert(contains_vertex(a) && contains_vertex(b)); - add_edge(a, b); - } - - // then add blockers - for (int current_dim = 1; current_dim <= dim; ++current_dim) { - std::list expansion_simplices; - compute_next_expand(set_simplices, current_dim, expansion_simplices); - - for (const auto &simplex : expansion_simplices) { - if (!set_simplices.contains(simplex)) { - add_blocker(simplex); - } - } - } - } - - // We cannot use the default copy constructor since we need - // to make a copy of each of the blockers - Skeleton_blocker_complex(const Skeleton_blocker_complex& copy) { - visitor = NULL; - degree_ = copy.degree_; - skeleton = Graph(copy.skeleton); - num_vertices_ = copy.num_vertices_; - - num_blockers_ = 0; - // we copy the blockers - for (auto blocker : copy.const_blocker_range()) { - add_blocker(*blocker); - } - } - - /** - */ - Skeleton_blocker_complex& operator=(const Skeleton_blocker_complex& copy) { - clear(); - visitor = NULL; - degree_ = copy.degree_; - skeleton = Graph(copy.skeleton); - num_vertices_ = copy.num_vertices_; - - num_blockers_ = 0; - // we copy the blockers - for (auto blocker : copy.const_blocker_range()) - add_blocker(*blocker); - return *this; - } - - /** - * The destructor delete all blockers allocated. - */ - virtual ~Skeleton_blocker_complex() { - clear(); - } - - /** - * @details Clears the simplicial complex. After a call to this function, - * blockers are destroyed. The 1-skeleton and the set of blockers - * are both empty. - */ - virtual void clear() { - // xxx for now the responsabilty of freeing the visitor is for - // the user - visitor = NULL; - - degree_.clear(); - num_vertices_ = 0; - - remove_blockers(); - - skeleton.clear(); - } - - /** - *@brief allows to change the visitor. - */ - void set_visitor(Visitor* other_visitor) { - visitor = other_visitor; - } - - //@} - - ///////////////////////////////////////////////////////////////////////////// - /** @name Vertices operations - */ - //@{ - public: - /** - * @brief Return a local Vertex_handle of a vertex given a global one. - * @remark Assume that the vertex is present in the complex. - */ - Vertex_handle operator[](Root_vertex_handle global) const { - auto local(get_address(global)); - assert(local); - return *local; - } - - /** - * @brief Return the vertex node associated to local Vertex_handle. - * @remark Assume that the vertex is present in the complex. - */ - Graph_vertex& operator[](Vertex_handle address) { - assert( - 0 <= address.vertex && address.vertex < boost::num_vertices(skeleton)); - return skeleton[address.vertex]; - } - - /** - * @brief Return the vertex node associated to local Vertex_handle. - * @remark Assume that the vertex is present in the complex. - */ - const Graph_vertex& operator[](Vertex_handle address) const { - assert( - 0 <= address.vertex && address.vertex < boost::num_vertices(skeleton)); - return skeleton[address.vertex]; - } - - /** - * @brief Adds a vertex to the simplicial complex and returns its Vertex_handle. - */ - Vertex_handle add_vertex() { - Vertex_handle address(boost::add_vertex(skeleton)); - num_vertices_++; - (*this)[address].activate(); - // safe since we now that we are in the root complex and the field 'address' and 'id' - // are identical for every vertices - (*this)[address].set_id(Root_vertex_handle(address.vertex)); - degree_.push_back(0); - if (visitor) - visitor->on_add_vertex(address); - return address; - } - - /** - * @brief Remove a vertex from the simplicial complex - * @remark It just deactivates the vertex with a boolean flag but does not - * remove it from vertices from complexity issues. - */ - void remove_vertex(Vertex_handle address) { - assert(contains_vertex(address)); - // We remove b - boost::clear_vertex(address.vertex, skeleton); - (*this)[address].deactivate(); - num_vertices_--; - degree_[address.vertex] = -1; - if (visitor) - visitor->on_remove_vertex(address); - } - - /** - */ - bool contains_vertex(Vertex_handle u) const { - if (u.vertex < 0 || u.vertex >= boost::num_vertices(skeleton)) - return false; - return (*this)[u].is_active(); - } - - /** - */ - bool contains_vertex(Root_vertex_handle u) const { - boost::optional address = get_address(u); - return address && (*this)[*address].is_active(); - } - - /** - * @return true iff the simplicial complex contains all vertices - * of simplex sigma - */ - bool contains_vertices(const Simplex_handle & sigma) const { - for (auto vertex : sigma) - if (!contains_vertex(vertex)) - return false; - return true; - } - - /** - * @brief Given an Id return the address of the vertex having this Id in the complex. - * @remark For a simplicial complex, the address is the id but it may not be the case for a SubComplex. - */ - virtual boost::optional get_address( - Root_vertex_handle id) const { - boost::optional res; - if (id.vertex < boost::num_vertices(skeleton)) - res = Vertex_handle(id.vertex); // xxx - return res; - } - - /** - * return the id of a vertex of adress local present in the graph - */ - Root_vertex_handle get_id(Vertex_handle local) const { - assert(0 <= local.vertex && local.vertex < boost::num_vertices(skeleton)); - return (*this)[local].get_id(); - } - - /** - * @brief Convert an address of a vertex of a complex to the address in - * the current complex. - * @details - * If the current complex is a sub (or sup) complex of 'other', it converts - * the address of a vertex v expressed in 'other' to the address of the vertex - * v in the current one. - * @remark this methods uses Root_vertex_handle to identify the vertex and - * assumes the vertex is present in the current complex. - */ - Vertex_handle convert_handle_from_another_complex( - const Skeleton_blocker_complex& other, Vertex_handle vh_in_other) const { - auto vh_in_current_complex = get_address(other.get_id(vh_in_other)); - assert(vh_in_current_complex); - return *vh_in_current_complex; - } - - /** - * @brief return the graph degree of a vertex. - */ - int degree(Vertex_handle local) const { - assert(0 <= local.vertex && local.vertex < boost::num_vertices(skeleton)); - return degree_[local.vertex]; - } - - //@} - - ///////////////////////////////////////////////////////////////////////////// - /** @name Edges operations - */ - //@{ - public: - /** - * @brief return an edge handle if the two vertices forms - * an edge in the complex - */ - boost::optional operator[]( - const std::pair& ab) const { - boost::optional res; - std::pair edge_pair( - boost::edge(ab.first.vertex, ab.second.vertex, skeleton)); - if (edge_pair.second) - res = edge_pair.first; - return res; - } - - /** - * @brief returns the stored node associated to an edge - */ - Graph_edge& operator[](Edge_handle edge_handle) { - return skeleton[edge_handle]; - } - - /** - * @brief returns the stored node associated to an edge - */ - const Graph_edge& operator[](Edge_handle edge_handle) const { - return skeleton[edge_handle]; - } - - /** - * @brief returns the first vertex of an edge - * @details it assumes that the edge is present in the complex - */ - Vertex_handle first_vertex(Edge_handle edge_handle) const { - return static_cast(source(edge_handle, skeleton)); - } - - /** - * @brief returns the first vertex of an edge - * @details it assumes that the edge is present in the complex - */ - Vertex_handle second_vertex(Edge_handle edge_handle) const { - return static_cast(target(edge_handle, skeleton)); - } - - /** - * @brief returns the simplex made with the two vertices of the edge - * @details it assumes that the edge is present in the complex - - */ - Simplex_handle get_vertices(Edge_handle edge_handle) const { - auto edge((*this)[edge_handle]); - return Simplex_handle((*this)[edge.first()], (*this)[edge.second()]); - } - - /** - * @brief Adds an edge between vertices a and b and all its cofaces. - */ - Edge_handle add_edge(Vertex_handle a, Vertex_handle b) { - assert(contains_vertex(a) && contains_vertex(b)); - assert(a != b); - - auto edge_handle((*this)[std::make_pair(a, b)]); - // std::pair pair_descr_bool = (*this)[std::make_pair(a,b)]; - // Edge_handle edge_descr; - // bool edge_present = pair_descr_bool.second; - if (!edge_handle) { - edge_handle = boost::add_edge(a.vertex, b.vertex, skeleton).first; - (*this)[*edge_handle].setId(get_id(a), get_id(b)); - degree_[a.vertex]++; - degree_[b.vertex]++; - if (visitor) - visitor->on_add_edge(a, b); - } - return *edge_handle; - } - - /** - * @brief Adds all edges and their cofaces of a simplex to the simplicial complex. - */ - void add_edges(const Simplex_handle & sigma) { - Simplex_handle_iterator i, j; - for (i = sigma.begin(); i != sigma.end(); ++i) - for (j = i, j++; j != sigma.end(); ++j) - add_edge(*i, *j); - } - - /** - * @brief Removes an edge from the simplicial complex and all its cofaces. - * @details returns the former Edge_handle representing the edge - */ - virtual Edge_handle remove_edge(Vertex_handle a, Vertex_handle b) { - bool found; - Edge_handle edge; - tie(edge, found) = boost::edge(a.vertex, b.vertex, skeleton); - if (found) { - if (visitor) - visitor->on_remove_edge(a, b); - // if (heapCollapse.Contains(edge)) heapCollapse.Delete(edge); - boost::remove_edge(a.vertex, b.vertex, skeleton); - degree_[a.vertex]--; - degree_[b.vertex]--; - } - return edge; - } - - /** - * @brief Removes edge and its cofaces from the simplicial complex. - */ - void remove_edge(Edge_handle edge) { - assert(contains_vertex(first_vertex(edge))); - assert(contains_vertex(second_vertex(edge))); - remove_edge(first_vertex(edge), second_vertex(edge)); - } - - /** - * @brief The complex is reduced to its set of vertices. - * All the edges and blockers are removed. - */ - void keep_only_vertices() { - remove_blockers(); - - for (auto u : vertex_range()) { - while (this->degree(u) > 0) { - Vertex_handle v(*(adjacent_vertices(u.vertex, this->skeleton).first)); - this->remove_edge(u, v); - } - } - } - - /** - * @return true iff the simplicial complex contains an edge between - * vertices a and b - */ - bool contains_edge(Vertex_handle a, Vertex_handle b) const { - // if (a.vertex<0 || b.vertex <0) return false; - return boost::edge(a.vertex, b.vertex, skeleton).second; - } - - /** - * @return true iff the simplicial complex contains all vertices - * and all edges of simplex sigma - */ - bool contains_edges(const Simplex_handle & sigma) const { - for (auto i = sigma.begin(); i != sigma.end(); ++i) { - if (!contains_vertex(*i)) - return false; - for (auto j = i; ++j != sigma.end();) { - if (!contains_edge(*i, *j)) - return false; - } - } - return true; - } - //@} - - ///////////////////////////////////////////////////////////////////////////// - /** @name Blockers operations - */ - //@{ - /** - * @brief Adds the simplex to the set of blockers and - * returns a Blocker_handle toward it if was not present before and 0 otherwise. - */ - Blocker_handle add_blocker(const Simplex_handle& blocker) { - assert(blocker.dimension() > 1); - if (contains_blocker(blocker)) { - // std::cerr << "ATTEMPT TO ADD A BLOCKER ALREADY THERE ---> BLOCKER IGNORED" << endl; - return 0; - } else { - if (visitor) - visitor->on_add_blocker(blocker); - Blocker_handle blocker_pt = new Simplex_handle(blocker); - num_blockers_++; - auto vertex = blocker_pt->begin(); - while (vertex != blocker_pt->end()) { - blocker_map_.insert(BlockerPair(*vertex, blocker_pt)); - ++vertex; - } - return blocker_pt; - } - } - - protected: - /** - * @brief Adds the simplex to the set of blockers - */ - void add_blocker(Blocker_handle blocker) { - if (contains_blocker(*blocker)) { - // std::cerr << "ATTEMPT TO ADD A BLOCKER ALREADY THERE ---> BLOCKER IGNORED" << endl; - return; - } else { - if (visitor) - visitor->on_add_blocker(*blocker); - num_blockers_++; - auto vertex = blocker->begin(); - while (vertex != blocker->end()) { - blocker_map_.insert(BlockerPair(*vertex, blocker)); - ++vertex; - } - } - } - - protected: - /** - * Removes sigma from the blocker map of vertex v - */ - void remove_blocker(const Blocker_handle sigma, Vertex_handle v) { - Complex_blocker_around_vertex_iterator blocker; - for (blocker = blocker_range(v).begin(); blocker != blocker_range(v).end(); - ++blocker) { - if (*blocker == sigma) - break; - } - if (*blocker != sigma) { - std::cerr - << "bug ((*blocker).second == sigma) ie try to remove a blocker not present\n"; - assert(false); - } else { - blocker_map_.erase(blocker.current_position()); - } - } - - public: - /** - * @brief Removes the simplex from the set of blockers. - * @remark sigma has to belongs to the set of blockers - */ - void remove_blocker(const Blocker_handle sigma) { - for (auto vertex : *sigma) { - remove_blocker(sigma, vertex); - } - num_blockers_--; - } - - /** - * @brief Remove all blockers, in other words, it expand the simplicial - * complex to the smallest flag complex that contains it. - */ - void remove_blockers() { - // Desallocate the blockers - while (!blocker_map_.empty()) { - delete_blocker(blocker_map_.begin()->second); - } - num_blockers_ = 0; - blocker_map_.clear(); - } - - protected: - /** - * Removes the simplex sigma from the set of blockers. - * sigma has to belongs to the set of blockers - * - * @remark contrarily to delete_blockers does not call the destructor - */ - void remove_blocker(const Simplex_handle& sigma) { - assert(contains_blocker(sigma)); - for (auto vertex : sigma) - remove_blocker(sigma, vertex); - num_blockers_--; - } - - public: - /** - * Removes the simplex s from the set of blockers - * and desallocate s. - */ - void delete_blocker(Blocker_handle sigma) { - if (visitor) - visitor->on_delete_blocker(sigma); - remove_blocker(sigma); - delete sigma; - } - - /** - * @return true iff s is a blocker of the simplicial complex - */ - bool contains_blocker(const Blocker_handle s) const { - if (s->dimension() < 2) - return false; - - Vertex_handle a = s->first_vertex(); - - for (auto blocker : const_blocker_range(a)) { - if (s == *blocker) - return true; - } - return false; - } - - /** - * @return true iff s is a blocker of the simplicial complex - */ - bool contains_blocker(const Simplex_handle & s) const { - if (s.dimension() < 2) - return false; - - Vertex_handle a = s.first_vertex(); - - for (auto blocker : const_blocker_range(a)) { - if (s == *blocker) - return true; - } - return false; - } - - private: - /** - * @return true iff a blocker of the simplicial complex - * is a face of sigma. - */ - bool blocks(const Simplex_handle & sigma) const { - for (auto blocker : const_blocker_range()) { - if (sigma.contains(*blocker)) - return true; - } - return false; - } - - //@} - - protected: - /** - * @details Adds to simplex the neighbours of v e.g. \f$ n \leftarrow n \cup N(v) \f$. - * If keep_only_superior is true then only vertices that are greater than v are added. - */ - virtual void add_neighbours(Vertex_handle v, Simplex_handle & n, - bool keep_only_superior = false) const { - boost_adjacency_iterator ai, ai_end; - for (tie(ai, ai_end) = adjacent_vertices(v.vertex, skeleton); ai != ai_end; - ++ai) { - if (keep_only_superior) { - if (*ai > v.vertex) { - n.add_vertex(Vertex_handle(*ai)); - } - } else { - n.add_vertex(Vertex_handle(*ai)); - } - } - } - - /** - * @details Add to simplex res all vertices which are - * neighbours of alpha: ie \f$ res \leftarrow res \cup N(alpha) \f$. - * - * If 'keep_only_superior' is true then only vertices that are greater than alpha are added. - * todo revoir - * - */ - virtual void add_neighbours(const Simplex_handle &alpha, Simplex_handle & res, - bool keep_only_superior = false) const { - res.clear(); - auto alpha_vertex = alpha.begin(); - add_neighbours(*alpha_vertex, res, keep_only_superior); - for (alpha_vertex = (alpha.begin())++; alpha_vertex != alpha.end(); - ++alpha_vertex) - keep_neighbours(*alpha_vertex, res, keep_only_superior); - } - - /** - * @details Remove from simplex n all vertices which are - * not neighbours of v e.g. \f$ res \leftarrow res \cap N(v) \f$. - * If 'keep_only_superior' is true then only vertices that are greater than v are keeped. - */ - virtual void keep_neighbours(Vertex_handle v, Simplex_handle& res, - bool keep_only_superior = false) const { - Simplex_handle nv; - add_neighbours(v, nv, keep_only_superior); - res.intersection(nv); - } - - /** - * @details Remove from simplex all vertices which are - * neighbours of v eg \f$ res \leftarrow res \setminus N(v) \f$. - * If 'keep_only_superior' is true then only vertices that are greater than v are added. - */ - virtual void remove_neighbours(Vertex_handle v, Simplex_handle & res, - bool keep_only_superior = false) const { - Simplex_handle nv; - add_neighbours(v, nv, keep_only_superior); - res.difference(nv); - } - - public: - /** - * @brief Compute the local vertices of 's' in the current complex - * If one of them is not present in the complex then the return value is uninitialized. - * - * - */ - // xxx rename get_address et place un using dans sub_complex - boost::optional get_simplex_address( - const Root_simplex_handle& s) const { - boost::optional res; - - Simplex_handle s_address; - // Root_simplex_const_iterator i; - for (auto i = s.begin(); i != s.end(); ++i) { - boost::optional address = get_address(*i); - if (!address) - return res; - else - s_address.add_vertex(*address); - } - res = s_address; - return res; - } - - /** - * @brief returns a simplex with vertices which are the id of vertices of the - * argument. - */ - Root_simplex_handle get_id(const Simplex_handle& local_simplex) const { - Root_simplex_handle global_simplex; - for (auto x = local_simplex.begin(); x != local_simplex.end(); ++x) { - global_simplex.add_vertex(get_id(*x)); - } - return global_simplex; - } - - /** - * @brief returns true iff the simplex s belongs to the simplicial - * complex. - */ - virtual bool contains(const Simplex_handle & s) const { - if (s.dimension() == -1) { - return false; - } else if (s.dimension() == 0) { - return contains_vertex(s.first_vertex()); - } else { - return (contains_edges(s) && !blocks(s)); - } - } - - /* - * @brief returnrs true iff the complex is empty. - */ - bool empty() const { - return num_vertices() == 0; - } - - /* - * @brief returns the number of vertices in the complex. - */ - int num_vertices() const { - // remark boost::num_vertices(skeleton) counts deactivated vertices - return num_vertices_; - } - - /* - * @brief returns the number of edges in the complex. - * @details currently in O(n) - */ - // todo cache the value - int num_edges() const { - return boost::num_edges(skeleton); - } - - /* - * @brief returns the number of blockers in the complex. - */ - int num_blockers() const { - return num_blockers_; - } - - /* - * @brief returns true iff the graph of the 1-skeleton of the complex is complete. - */ - bool complete() const { - return (num_vertices() * (num_vertices() - 1)) / 2 == num_edges(); - } - - /** - * @brief returns the number of connected components in the graph of the 1-skeleton. - */ - int num_connected_components() const { - int num_vert_collapsed = skeleton.vertex_set().size() - num_vertices(); - std::vector component(skeleton.vertex_set().size()); - return boost::connected_components(this->skeleton, &component[0]) - - num_vert_collapsed; - } - - /** - * @brief %Test if the complex is a cone. - * @details Runs in O(n) where n is the number of vertices. - */ - bool is_cone() const { - if (num_vertices() == 0) - return false; - if (num_vertices() == 1) - return true; - for (auto vi : vertex_range()) { - // xxx todo faire une methode bool is_in_blocker(Vertex_handle) - if (blocker_map_.find(vi) == blocker_map_.end()) { - // no blocker passes through the vertex, we just need to - // check if the current vertex is linked to all others vertices of the complex - if (degree_[vi.vertex] == num_vertices() - 1) - return true; - } - } - return false; - } - - //@} - - ///////////////////////////////////////////////////////////////////////////// - /** @name Vertex iterators - */ - //@{ - typedef Complex_vertex_iterator CVI; // todo rename - - // - // @brief Range over the vertices of the simplicial complex. - // Methods .begin() and .end() return a Complex_vertex_iterator. - // - typedef boost::iterator_range< - Complex_vertex_iterator > Complex_vertex_range; - - /** - * @brief Returns a Complex_vertex_range over all vertices of the complex - */ - Complex_vertex_range vertex_range() const { - auto begin = Complex_vertex_iterator(this); - auto end = Complex_vertex_iterator(this, 0); - return Complex_vertex_range(begin, end); - } - - typedef boost::iterator_range< - Complex_neighbors_vertices_iterator > Complex_neighbors_vertices_range; - - /** - * @brief Returns a Complex_edge_range over all edges of the simplicial complex that passes trough v - */ - Complex_neighbors_vertices_range vertex_range(Vertex_handle v) const { - auto begin = Complex_neighbors_vertices_iterator( - this, v); - auto end = Complex_neighbors_vertices_iterator( - this, v, 0); - return Complex_neighbors_vertices_range(begin, end); - } - - //@} - - /** @name Edge iterators - */ - //@{ - - typedef boost::iterator_range< - Complex_edge_iterator>> Complex_edge_range; - - /** - * @brief Returns a Complex_edge_range over all edges of the simplicial complex - */ - Complex_edge_range edge_range() const { - auto begin = Complex_edge_iterator>(this); - auto end = Complex_edge_iterator>(this, 0); - return Complex_edge_range(begin, end); - } - - typedef boost::iterator_range >> - Complex_edge_around_vertex_range; - /** - * @brief Returns a Complex_edge_range over all edges of the simplicial complex that passes - * through 'v' - */ - Complex_edge_around_vertex_range edge_range(Vertex_handle v) const { - auto begin = Complex_edge_around_vertex_iterator>(this, v); - auto end = Complex_edge_around_vertex_iterator>(this, v, 0); - return Complex_edge_around_vertex_range(begin, end); - } - - //@} - - /** @name Triangles iterators - */ - //@{ - private: - typedef Skeleton_blocker_link_complex > Link; - typedef Skeleton_blocker_link_superior > Superior_link; - - public: - typedef Triangle_around_vertex_iterator Superior_triangle_around_vertex_iterator; - - typedef boost::iterator_range < Triangle_around_vertex_iterator > Complex_triangle_around_vertex_range; - - /** - * @brief Range over triangles around a vertex of the simplicial complex. - * Methods .begin() and .end() return a Triangle_around_vertex_iterator. - * - */ - Complex_triangle_around_vertex_range triangle_range(Vertex_handle v) const { - auto begin = Triangle_around_vertex_iterator(this, v); - auto end = Triangle_around_vertex_iterator(this, v, 0); - return Complex_triangle_around_vertex_range(begin, end); - } - - typedef boost::iterator_range > Complex_triangle_range; - - /** - * @brief Range over triangles of the simplicial complex. - * Methods .begin() and .end() return a Triangle_around_vertex_iterator. - * - */ - Complex_triangle_range triangle_range() const { - auto end = Triangle_iterator(this, 0); - if (empty()) { - return Complex_triangle_range(end, end); - } else { - auto begin = Triangle_iterator(this); - return Complex_triangle_range(begin, end); - } - } - - //@} - - /** @name Simplices iterators - */ - //@{ - typedef Simplex_around_vertex_iterator Complex_simplex_around_vertex_iterator; - - /** - * @brief Range over the simplices of the simplicial complex around a vertex. - * Methods .begin() and .end() return a Complex_simplex_around_vertex_iterator. - */ - typedef boost::iterator_range < Complex_simplex_around_vertex_iterator > Complex_simplex_around_vertex_range; - - /** - * @brief Returns a Complex_simplex_around_vertex_range over all the simplices around a vertex of the complex - */ - Complex_simplex_around_vertex_range simplex_range(Vertex_handle v) const { - assert(contains_vertex(v)); - return Complex_simplex_around_vertex_range( - Complex_simplex_around_vertex_iterator(this, v), - Complex_simplex_around_vertex_iterator(this, v, true)); - } - - // typedef Simplex_iterator Complex_simplex_iterator; - typedef Simplex_iterator Complex_simplex_iterator; - - typedef boost::iterator_range < Complex_simplex_iterator > Complex_simplex_range; - - /** - * @brief Returns a Complex_simplex_range over all the simplices of the complex - */ - Complex_simplex_range simplex_range() const { - Complex_simplex_iterator end(this, true); - if (empty()) { - return Complex_simplex_range(end, end); - } else { - Complex_simplex_iterator begin(this); - return Complex_simplex_range(begin, end); - } - } - - //@} - - /** @name Blockers iterators - */ - //@{ - private: - /** - * @brief Iterator over the blockers adjacent to a vertex - */ - typedef Blocker_iterator_around_vertex_internal< - typename std::multimap::iterator, - Blocker_handle> - Complex_blocker_around_vertex_iterator; - - /** - * @brief Iterator over (constant) blockers adjacent to a vertex - */ - typedef Blocker_iterator_around_vertex_internal< - typename std::multimap::const_iterator, - const Blocker_handle> - Const_complex_blocker_around_vertex_iterator; - - typedef boost::iterator_range Complex_blocker_around_vertex_range; - typedef boost::iterator_range Const_complex_blocker_around_vertex_range; - - public: - /** - * @brief Returns a range of the blockers of the complex passing through a vertex - */ - Complex_blocker_around_vertex_range blocker_range(Vertex_handle v) { - auto begin = Complex_blocker_around_vertex_iterator(blocker_map_.lower_bound(v)); - auto end = Complex_blocker_around_vertex_iterator(blocker_map_.upper_bound(v)); - return Complex_blocker_around_vertex_range(begin, end); - } - - /** - * @brief Returns a range of the blockers of the complex passing through a vertex - */ - Const_complex_blocker_around_vertex_range const_blocker_range(Vertex_handle v) const { - auto begin = Const_complex_blocker_around_vertex_iterator(blocker_map_.lower_bound(v)); - auto end = Const_complex_blocker_around_vertex_iterator(blocker_map_.upper_bound(v)); - return Const_complex_blocker_around_vertex_range(begin, end); - } - - private: - /** - * @brief Iterator over the blockers. - */ - typedef Blocker_iterator_internal< - typename std::multimap::iterator, - Blocker_handle> - Complex_blocker_iterator; - - /** - * @brief Iterator over the (constant) blockers. - */ - typedef Blocker_iterator_internal< - typename std::multimap::const_iterator, - const Blocker_handle> - Const_complex_blocker_iterator; - - typedef boost::iterator_range Complex_blocker_range; - typedef boost::iterator_range Const_complex_blocker_range; - - public: - /** - * @brief Returns a range of the blockers of the complex - */ - Complex_blocker_range blocker_range() { - auto begin = Complex_blocker_iterator(blocker_map_.begin(), blocker_map_.end() ); - auto end = Complex_blocker_iterator(blocker_map_.end() , blocker_map_.end() ); - return Complex_blocker_range(begin, end); - } - - /** - * @brief Returns a range of the blockers of the complex - */ - Const_complex_blocker_range const_blocker_range() const { - auto begin = Const_complex_blocker_iterator(blocker_map_.begin(), blocker_map_.end() ); - auto end = Const_complex_blocker_iterator(blocker_map_.end(), blocker_map_.end() ); - return Const_complex_blocker_range(begin, end); - } - - //@} - - ///////////////////////////////////////////////////////////////////////////// - /** @name Print and IO methods - */ - //@{ - public: - std::string to_string() const { - std::ostringstream stream; - stream << num_vertices() << " vertices:\n" << vertices_to_string() << std::endl; - stream << num_edges() << " edges:\n" << edges_to_string() << std::endl; - stream << num_blockers() << " blockers:\n" << blockers_to_string() << std::endl; - return stream.str(); - } - - std::string vertices_to_string() const { - std::ostringstream stream; - for (auto vertex : vertex_range()) { - stream << "(" << (*this)[vertex].get_id() << "),"; - } - stream << std::endl; - return stream.str(); - } - - std::string edges_to_string() const { - std::ostringstream stream; - for (auto edge : edge_range()) { - stream << "(" << (*this)[edge].first() << "," << (*this)[edge].second() << ")" << " id = " << (*this)[edge].index() << std::endl; - } - stream << std::endl; - return stream.str(); - } - - std::string blockers_to_string() const { - std::ostringstream stream; - for (auto bl : blocker_map_) { - stream << bl.first << " => " << bl.second << ":" << *bl.second << std::endl; - } - return stream.str(); - } - - //@} + template friend class Complex_vertex_iterator; + template friend class Complex_neighbors_vertices_iterator; + template friend class Complex_edge_iterator; + template friend class Complex_edge_around_vertex_iterator; + + template friend class Skeleton_blocker_link_complex; + template friend class Skeleton_blocker_link_superior; + template friend class Skeleton_blocker_sub_complex; + +public: + /** + * @brief The type of stored vertex node, specified by the template SkeletonBlockerDS + */ + typedef typename SkeletonBlockerDS::Graph_vertex Graph_vertex; + + /** + * @brief The type of stored edge node, specified by the template SkeletonBlockerDS + */ + typedef typename SkeletonBlockerDS::Graph_edge Graph_edge; + + typedef typename SkeletonBlockerDS::Root_vertex_handle Root_vertex_handle; + + /** + * @brief The type of an handle to a vertex of the complex. + */ + typedef typename SkeletonBlockerDS::Vertex_handle Vertex_handle; + typedef typename Root_vertex_handle::boost_vertex_handle boost_vertex_handle; + + /** + * @brief A ordered set of integers that represents a simplex. + */ + typedef Skeleton_blocker_simplex Simplex_handle; + typedef Skeleton_blocker_simplex Root_simplex_handle; + + /** + * @brief Handle to a blocker of the complex. + */ + typedef Simplex_handle* Blocker_handle; + + typedef typename Root_simplex_handle::Simplex_vertex_const_iterator Root_simplex_iterator; + typedef typename Simplex_handle::Simplex_vertex_const_iterator Simplex_handle_iterator; + +protected: + typedef typename boost::adjacency_list Graph; + // todo/remark : edges are not sorted, it heavily penalizes computation for SuperiorLink + // (eg Link with greater vertices) + // that burdens simplex iteration / complex initialization via list of simplices. + // to avoid that, one should modify the graph by storing two lists of adjacency for every + // vertex, the one with superior and the one with lower vertices, that way there is + // no more extra cost for computation of SuperiorLink + typedef typename boost::graph_traits::vertex_iterator boost_vertex_iterator; + typedef typename boost::graph_traits::edge_iterator boost_edge_iterator; + +protected: + typedef typename boost::graph_traits::adjacency_iterator boost_adjacency_iterator; + +public: + /** + * @brief Handle to an edge of the complex. + */ + typedef typename boost::graph_traits::edge_descriptor Edge_handle; + +protected: + typedef std::multimap BlockerMap; + typedef typename std::multimap::value_type BlockerPair; + typedef typename std::multimap::iterator BlockerMapIterator; + typedef typename std::multimap::const_iterator BlockerMapConstIterator; + +protected: + int num_vertices_; + int num_blockers_; + + typedef Skeleton_blocker_complex_visitor Visitor; + // typedef Visitor* Visitor_ptr; + Visitor* visitor; + + /** + * @details If 'x' is a Vertex_handle of a vertex in the complex then degree[x] = d is its degree. + * + * This quantity is updated when adding/removing edge. + * + * This is useful because the operation + * list.size() is done in linear time. + */ + std::vector degree_; + Graph skeleton; /** 1-skeleton of the simplicial complex. */ + + /** Each vertex can access to the blockers passing through it. */ + BlockerMap blocker_map_; + +public: + ///////////////////////////////////////////////////////////////////////////// + /** @name Constructors, Destructors + */ + //@{ + /** + *@brief constructs a simplicial complex with a given number of vertices and a visitor. + */ + explicit Skeleton_blocker_complex(int num_vertices_ = 0, Visitor* visitor_ = NULL) + : visitor(visitor_) { + clear(); + for (int i = 0; i < num_vertices_; ++i) { + add_vertex(); + } + } + +private: + typedef Trie> STrie; + template + std::vector compute_cofaces(SimpleHandleOutputIterator simplex_begin,SimpleHandleOutputIterator simplex_end){ + std::vector cofaces(num_vertices(),0); + for(auto i = 0u ; i < num_vertices(); ++i) + cofaces[i] = new STrie(Vertex_handle(i)); + for(auto s_it = simplex_begin; s_it != simplex_end; ++s_it){ + if(s_it->dimension()>=1) + cofaces[s_it->first_vertex()]->add_simplex(*s_it); + } + return cofaces; + } + + +public: + /** + * @brief Constructor with a list of simplices. + * @details is_flag_complex indicates if the complex is a flag complex or not (to know if blockers have to be computed or not). + */ + template + Skeleton_blocker_complex(SimpleHandleOutputIterator simplex_begin,SimpleHandleOutputIterator simplex_end,bool is_flag_complex = false,Visitor* visitor_ = NULL) + : num_vertices_(0), + num_blockers_(0), + visitor(visitor_){ + std::vector> edges; + + //first pass add vertices and store edges + for(auto s_it = simplex_begin; s_it != simplex_end; ++s_it){ + if(s_it->dimension()==0) add_vertex(); + if(s_it->dimension()==1) edges.emplace_back(s_it->first_vertex(),s_it->last_vertex()); + } + for(const auto& e : edges) + add_edge(e.first,e.second); + + if(!is_flag_complex){ + //need to compute blockers + std::vector cofaces(compute_cofaces(simplex_begin,simplex_end)); + std::vector max_faces; + + for(auto i = 0u ; i < num_vertices(); ++i){ + auto max_faces_i = cofaces[i]->maximal_faces(); + max_faces.insert(max_faces.end(),max_faces_i.begin(),max_faces_i.end()); + } + + // all time here + + // for every maximal face s, it picks its first vertex v and check for all nv \in N(v) + // if s union nv is in the complex, if not it is a blocker. + for(auto max_face : max_faces){ + if(max_face.dimension()>0){ + auto first_v = max_face.first_vertex(); + for(auto nv : vertex_range(first_v)){ + if(! max_face.contains(nv) && max_face.first_vertex()contains(max_face) + max_face.add_vertex(nv); + if(!cofaces[first_v]->contains(max_face)) + add_blocker(max_face); + max_face.remove_vertex(nv); + } + } + } + } + } + for(auto i = 0u ; i < num_vertices(); ++i) + delete(cofaces[i]); + } + } + + // We cannot use the default copy constructor since we need + // to make a copy of each of the blockers + Skeleton_blocker_complex(const Skeleton_blocker_complex& copy) { + visitor = NULL; + degree_ = copy.degree_; + skeleton = Graph(copy.skeleton); + num_vertices_ = copy.num_vertices_; + + num_blockers_ = 0; + // we copy the blockers + for (auto blocker : copy.const_blocker_range()) { + add_blocker(*blocker); + } + } + + /** + */ + Skeleton_blocker_complex& operator=(const Skeleton_blocker_complex& copy) { + clear(); + visitor = NULL; + degree_ = copy.degree_; + skeleton = Graph(copy.skeleton); + num_vertices_ = copy.num_vertices_; + + num_blockers_ = 0; + // we copy the blockers + for (auto blocker : copy.const_blocker_range()) + add_blocker(*blocker); + return *this; + } + + /** + * The destructor delete all blockers allocated. + */ + virtual ~Skeleton_blocker_complex() { + clear(); + } + + /** + * @details Clears the simplicial complex. After a call to this function, + * blockers are destroyed. The 1-skeleton and the set of blockers + * are both empty. + */ + virtual void clear() { + // xxx for now the responsabilty of freeing the visitor is for + // the user + visitor = NULL; + + degree_.clear(); + num_vertices_ = 0; + + remove_blockers(); + + skeleton.clear(); + } + + /** + *@brief allows to change the visitor. + */ + void set_visitor(Visitor* other_visitor) { + visitor = other_visitor; + } + + //@} + + ///////////////////////////////////////////////////////////////////////////// + /** @name Vertices operations + */ + //@{ +public: + /** + * @brief Return a local Vertex_handle of a vertex given a global one. + * @remark Assume that the vertex is present in the complex. + */ + Vertex_handle operator[](Root_vertex_handle global) const { + auto local(get_address(global)); + assert(local); + return *local; + } + + /** + * @brief Return the vertex node associated to local Vertex_handle. + * @remark Assume that the vertex is present in the complex. + */ + Graph_vertex& operator[](Vertex_handle address) { + assert( + 0 <= address.vertex && address.vertex < boost::num_vertices(skeleton)); + return skeleton[address.vertex]; + } + + /** + * @brief Return the vertex node associated to local Vertex_handle. + * @remark Assume that the vertex is present in the complex. + */ + const Graph_vertex& operator[](Vertex_handle address) const { + assert( + 0 <= address.vertex && address.vertex < boost::num_vertices(skeleton)); + return skeleton[address.vertex]; + } + + /** + * @brief Adds a vertex to the simplicial complex and returns its Vertex_handle. + */ + Vertex_handle add_vertex() { + Vertex_handle address(boost::add_vertex(skeleton)); + num_vertices_++; + (*this)[address].activate(); + // safe since we now that we are in the root complex and the field 'address' and 'id' + // are identical for every vertices + (*this)[address].set_id(Root_vertex_handle(address.vertex)); + degree_.push_back(0); + if (visitor) + visitor->on_add_vertex(address); + return address; + } + + /** + * @brief Remove a vertex from the simplicial complex + * @remark It just deactivates the vertex with a boolean flag but does not + * remove it from vertices from complexity issues. + */ + void remove_vertex(Vertex_handle address) { + assert(contains_vertex(address)); + // We remove b + boost::clear_vertex(address.vertex, skeleton); + (*this)[address].deactivate(); + num_vertices_--; + degree_[address.vertex] = -1; + if (visitor) + visitor->on_remove_vertex(address); + } + + /** + */ + bool contains_vertex(Vertex_handle u) const { + if (u.vertex < 0 || u.vertex >= boost::num_vertices(skeleton)) + return false; + return (*this)[u].is_active(); + } + + /** + */ + bool contains_vertex(Root_vertex_handle u) const { + boost::optional address = get_address(u); + return address && (*this)[*address].is_active(); + } + + /** + * @return true iff the simplicial complex contains all vertices + * of simplex sigma + */ + bool contains_vertices(const Simplex_handle & sigma) const { + for (auto vertex : sigma) + if (!contains_vertex(vertex)) + return false; + return true; + } + + /** + * @brief Given an Id return the address of the vertex having this Id in the complex. + * @remark For a simplicial complex, the address is the id but it may not be the case for a SubComplex. + */ + virtual boost::optional get_address( + Root_vertex_handle id) const { + boost::optional res; + if (id.vertex < boost::num_vertices(skeleton)) + res = Vertex_handle(id.vertex); // xxx + return res; + } + + /** + * return the id of a vertex of adress local present in the graph + */ + Root_vertex_handle get_id(Vertex_handle local) const { + assert(0 <= local.vertex && local.vertex < boost::num_vertices(skeleton)); + return (*this)[local].get_id(); + } + + /** + * @brief Convert an address of a vertex of a complex to the address in + * the current complex. + * @details + * If the current complex is a sub (or sup) complex of 'other', it converts + * the address of a vertex v expressed in 'other' to the address of the vertex + * v in the current one. + * @remark this methods uses Root_vertex_handle to identify the vertex and + * assumes the vertex is present in the current complex. + */ + Vertex_handle convert_handle_from_another_complex( + const Skeleton_blocker_complex& other, Vertex_handle vh_in_other) const { + auto vh_in_current_complex = get_address(other.get_id(vh_in_other)); + assert(vh_in_current_complex); + return *vh_in_current_complex; + } + + /** + * @brief return the graph degree of a vertex. + */ + int degree(Vertex_handle local) const { + assert(0 <= local.vertex && local.vertex < boost::num_vertices(skeleton)); + return degree_[local.vertex]; + } + + //@} + + ///////////////////////////////////////////////////////////////////////////// + /** @name Edges operations + */ + //@{ +public: + /** + * @brief return an edge handle if the two vertices forms + * an edge in the complex + */ + boost::optional operator[]( + const std::pair& ab) const { + boost::optional res; + std::pair edge_pair( + boost::edge(ab.first.vertex, ab.second.vertex, skeleton)); + if (edge_pair.second) + res = edge_pair.first; + return res; + } + + /** + * @brief returns the stored node associated to an edge + */ + Graph_edge& operator[](Edge_handle edge_handle) { + return skeleton[edge_handle]; + } + + /** + * @brief returns the stored node associated to an edge + */ + const Graph_edge& operator[](Edge_handle edge_handle) const { + return skeleton[edge_handle]; + } + + /** + * @brief returns the first vertex of an edge + * @details it assumes that the edge is present in the complex + */ + Vertex_handle first_vertex(Edge_handle edge_handle) const { + return static_cast(source(edge_handle, skeleton)); + } + + /** + * @brief returns the first vertex of an edge + * @details it assumes that the edge is present in the complex + */ + Vertex_handle second_vertex(Edge_handle edge_handle) const { + return static_cast(target(edge_handle, skeleton)); + } + + /** + * @brief returns the simplex made with the two vertices of the edge + * @details it assumes that the edge is present in the complex + + */ + Simplex_handle get_vertices(Edge_handle edge_handle) const { + auto edge((*this)[edge_handle]); + return Simplex_handle((*this)[edge.first()], (*this)[edge.second()]); + } + + /** + * @brief Adds an edge between vertices a and b and all its cofaces. + */ + Edge_handle add_edge(Vertex_handle a, Vertex_handle b) { + assert(contains_vertex(a) && contains_vertex(b)); + assert(a != b); + + auto edge_handle((*this)[std::make_pair(a, b)]); + // std::pair pair_descr_bool = (*this)[std::make_pair(a,b)]; + // Edge_handle edge_descr; + // bool edge_present = pair_descr_bool.second; + if (!edge_handle) { + edge_handle = boost::add_edge(a.vertex, b.vertex, skeleton).first; + (*this)[*edge_handle].setId(get_id(a), get_id(b)); + degree_[a.vertex]++; + degree_[b.vertex]++; + if (visitor) + visitor->on_add_edge(a, b); + } + return *edge_handle; + } + + /** + * @brief Adds all edges and their cofaces of a simplex to the simplicial complex. + */ + void add_edges(const Simplex_handle & sigma) { + Simplex_handle_iterator i, j; + for (i = sigma.begin(); i != sigma.end(); ++i) + for (j = i, j++; j != sigma.end(); ++j) + add_edge(*i, *j); + } + + /** + * @brief Removes an edge from the simplicial complex and all its cofaces. + * @details returns the former Edge_handle representing the edge + */ + virtual Edge_handle remove_edge(Vertex_handle a, Vertex_handle b) { + bool found; + Edge_handle edge; + tie(edge, found) = boost::edge(a.vertex, b.vertex, skeleton); + if (found) { + if (visitor) + visitor->on_remove_edge(a, b); + // if (heapCollapse.Contains(edge)) heapCollapse.Delete(edge); + boost::remove_edge(a.vertex, b.vertex, skeleton); + degree_[a.vertex]--; + degree_[b.vertex]--; + } + return edge; + } + + /** + * @brief Removes edge and its cofaces from the simplicial complex. + */ + void remove_edge(Edge_handle edge) { + assert(contains_vertex(first_vertex(edge))); + assert(contains_vertex(second_vertex(edge))); + remove_edge(first_vertex(edge), second_vertex(edge)); + } + + /** + * @brief The complex is reduced to its set of vertices. + * All the edges and blockers are removed. + */ + void keep_only_vertices() { + remove_blockers(); + + for (auto u : vertex_range()) { + while (this->degree(u) > 0) { + Vertex_handle v(*(adjacent_vertices(u.vertex, this->skeleton).first)); + this->remove_edge(u, v); + } + } + } + + /** + * @return true iff the simplicial complex contains an edge between + * vertices a and b + */ + bool contains_edge(Vertex_handle a, Vertex_handle b) const { + // if (a.vertex<0 || b.vertex <0) return false; + return boost::edge(a.vertex, b.vertex, skeleton).second; + } + + /** + * @return true iff the simplicial complex contains all vertices + * and all edges of simplex sigma + */ + bool contains_edges(const Simplex_handle & sigma) const { + for (auto i = sigma.begin(); i != sigma.end(); ++i) { + if (!contains_vertex(*i)) + return false; + for (auto j = i; ++j != sigma.end();) { + if (!contains_edge(*i, *j)) + return false; + } + } + return true; + } + //@} + + ///////////////////////////////////////////////////////////////////////////// + /** @name Blockers operations + */ + //@{ + /** + * @brief Adds the simplex to the set of blockers and + * returns a Blocker_handle toward it if was not present before and 0 otherwise. + */ + Blocker_handle add_blocker(const Simplex_handle& blocker) { + assert(blocker.dimension() > 1); + if (contains_blocker(blocker)) { + // std::cerr << "ATTEMPT TO ADD A BLOCKER ALREADY THERE ---> BLOCKER IGNORED" << endl; + return 0; + } else { + if (visitor) + visitor->on_add_blocker(blocker); + Blocker_handle blocker_pt = new Simplex_handle(blocker); + num_blockers_++; + auto vertex = blocker_pt->begin(); + while (vertex != blocker_pt->end()) { + blocker_map_.insert(BlockerPair(*vertex, blocker_pt)); + ++vertex; + } + return blocker_pt; + } + } + +protected: + /** + * @brief Adds the simplex to the set of blockers + */ + void add_blocker(Blocker_handle blocker) { + if (contains_blocker(*blocker)) { + // std::cerr << "ATTEMPT TO ADD A BLOCKER ALREADY THERE ---> BLOCKER IGNORED" << endl; + return; + } else { + if (visitor) + visitor->on_add_blocker(*blocker); + num_blockers_++; + auto vertex = blocker->begin(); + while (vertex != blocker->end()) { + blocker_map_.insert(BlockerPair(*vertex, blocker)); + ++vertex; + } + } + } + +protected: + /** + * Removes sigma from the blocker map of vertex v + */ + void remove_blocker(const Blocker_handle sigma, Vertex_handle v) { + Complex_blocker_around_vertex_iterator blocker; + for (blocker = blocker_range(v).begin(); blocker != blocker_range(v).end(); + ++blocker) { + if (*blocker == sigma) + break; + } + if (*blocker != sigma) { + std::cerr + << "bug ((*blocker).second == sigma) ie try to remove a blocker not present\n"; + assert(false); + } else { + blocker_map_.erase(blocker.current_position()); + } + } + +public: + /** + * @brief Removes the simplex from the set of blockers. + * @remark sigma has to belongs to the set of blockers + */ + void remove_blocker(const Blocker_handle sigma) { + for (auto vertex : *sigma) { + remove_blocker(sigma, vertex); + } + num_blockers_--; + } + + /** + * @brief Remove all blockers, in other words, it expand the simplicial + * complex to the smallest flag complex that contains it. + */ + void remove_blockers() { + // Desallocate the blockers + while (!blocker_map_.empty()) { + delete_blocker(blocker_map_.begin()->second); + } + num_blockers_ = 0; + blocker_map_.clear(); + } + +protected: + /** + * Removes the simplex sigma from the set of blockers. + * sigma has to belongs to the set of blockers + * + * @remark contrarily to delete_blockers does not call the destructor + */ + void remove_blocker(const Simplex_handle& sigma) { + assert(contains_blocker(sigma)); + for (auto vertex : sigma) + remove_blocker(sigma, vertex); + num_blockers_--; + } + +public: + /** + * Removes the simplex s from the set of blockers + * and desallocate s. + */ + void delete_blocker(Blocker_handle sigma) { + if (visitor) + visitor->on_delete_blocker(sigma); + remove_blocker(sigma); + delete sigma; + } + + /** + * @return true iff s is a blocker of the simplicial complex + */ + bool contains_blocker(const Blocker_handle s) const { + if (s->dimension() < 2) + return false; + + Vertex_handle a = s->first_vertex(); + + for (auto blocker : const_blocker_range(a)) { + if (s == *blocker) + return true; + } + return false; + } + + /** + * @return true iff s is a blocker of the simplicial complex + */ + bool contains_blocker(const Simplex_handle & s) const { + if (s.dimension() < 2) + return false; + + Vertex_handle a = s.first_vertex(); + + for (auto blocker : const_blocker_range(a)) { + if (s == *blocker) + return true; + } + return false; + } + +private: + /** + * @return true iff a blocker of the simplicial complex + * is a face of sigma. + */ + bool blocks(const Simplex_handle & sigma) const { + for (auto blocker : const_blocker_range()) { + if (sigma.contains(*blocker)) + return true; + } + return false; + } + + //@} + +protected: + /** + * @details Adds to simplex the neighbours of v e.g. \f$ n \leftarrow n \cup N(v) \f$. + * If keep_only_superior is true then only vertices that are greater than v are added. + */ + virtual void add_neighbours(Vertex_handle v, Simplex_handle & n, + bool keep_only_superior = false) const { + boost_adjacency_iterator ai, ai_end; + for (tie(ai, ai_end) = adjacent_vertices(v.vertex, skeleton); ai != ai_end; + ++ai) { + if (keep_only_superior) { + if (*ai > v.vertex) { + n.add_vertex(Vertex_handle(*ai)); + } + } else { + n.add_vertex(Vertex_handle(*ai)); + } + } + } + + /** + * @details Add to simplex res all vertices which are + * neighbours of alpha: ie \f$ res \leftarrow res \cup N(alpha) \f$. + * + * If 'keep_only_superior' is true then only vertices that are greater than alpha are added. + * todo revoir + * + */ + virtual void add_neighbours(const Simplex_handle &alpha, Simplex_handle & res, + bool keep_only_superior = false) const { + res.clear(); + auto alpha_vertex = alpha.begin(); + add_neighbours(*alpha_vertex, res, keep_only_superior); + for (alpha_vertex = (alpha.begin())++; alpha_vertex != alpha.end(); + ++alpha_vertex) + keep_neighbours(*alpha_vertex, res, keep_only_superior); + } + + /** + * @details Remove from simplex n all vertices which are + * not neighbours of v e.g. \f$ res \leftarrow res \cap N(v) \f$. + * If 'keep_only_superior' is true then only vertices that are greater than v are keeped. + */ + virtual void keep_neighbours(Vertex_handle v, Simplex_handle& res, + bool keep_only_superior = false) const { + Simplex_handle nv; + add_neighbours(v, nv, keep_only_superior); + res.intersection(nv); + } + + /** + * @details Remove from simplex all vertices which are + * neighbours of v eg \f$ res \leftarrow res \setminus N(v) \f$. + * If 'keep_only_superior' is true then only vertices that are greater than v are added. + */ + virtual void remove_neighbours(Vertex_handle v, Simplex_handle & res, + bool keep_only_superior = false) const { + Simplex_handle nv; + add_neighbours(v, nv, keep_only_superior); + res.difference(nv); + } + +public: + /** + * @brief Compute the local vertices of 's' in the current complex + * If one of them is not present in the complex then the return value is uninitialized. + * + * + */ + // xxx rename get_address et place un using dans sub_complex + boost::optional get_simplex_address( + const Root_simplex_handle& s) const { + boost::optional res; + + Simplex_handle s_address; + // Root_simplex_const_iterator i; + for (auto i = s.begin(); i != s.end(); ++i) { + boost::optional address = get_address(*i); + if (!address) + return res; + else + s_address.add_vertex(*address); + } + res = s_address; + return res; + } + + /** + * @brief returns a simplex with vertices which are the id of vertices of the + * argument. + */ + Root_simplex_handle get_id(const Simplex_handle& local_simplex) const { + Root_simplex_handle global_simplex; + for (auto x = local_simplex.begin(); x != local_simplex.end(); ++x) { + global_simplex.add_vertex(get_id(*x)); + } + return global_simplex; + } + + /** + * @brief returns true iff the simplex s belongs to the simplicial + * complex. + */ + virtual bool contains(const Simplex_handle & s) const { + if (s.dimension() == -1) { + return false; + } else if (s.dimension() == 0) { + return contains_vertex(s.first_vertex()); + } else { + return (contains_edges(s) && !blocks(s)); + } + } + + /* + * @brief returnrs true iff the complex is empty. + */ + bool empty() const { + return num_vertices() == 0; + } + + /* + * @brief returns the number of vertices in the complex. + */ + int num_vertices() const { + // remark boost::num_vertices(skeleton) counts deactivated vertices + return num_vertices_; + } + + /* + * @brief returns the number of edges in the complex. + * @details currently in O(n) + */ + // todo cache the value + int num_edges() const { + return boost::num_edges(skeleton); + } + + /* + * @brief returns the number of blockers in the complex. + */ + int num_blockers() const { + return num_blockers_; + } + + /* + * @brief returns true iff the graph of the 1-skeleton of the complex is complete. + */ + bool complete() const { + return (num_vertices() * (num_vertices() - 1)) / 2 == num_edges(); + } + + /** + * @brief returns the number of connected components in the graph of the 1-skeleton. + */ + int num_connected_components() const { + int num_vert_collapsed = skeleton.vertex_set().size() - num_vertices(); + std::vector component(skeleton.vertex_set().size()); + return boost::connected_components(this->skeleton, &component[0]) + - num_vert_collapsed; + } + + /** + * @brief %Test if the complex is a cone. + * @details Runs in O(n) where n is the number of vertices. + */ + bool is_cone() const { + if (num_vertices() == 0) + return false; + if (num_vertices() == 1) + return true; + for (auto vi : vertex_range()) { + // xxx todo faire une methode bool is_in_blocker(Vertex_handle) + if (blocker_map_.find(vi) == blocker_map_.end()) { + // no blocker passes through the vertex, we just need to + // check if the current vertex is linked to all others vertices of the complex + if (degree_[vi.vertex] == num_vertices() - 1) + return true; + } + } + return false; + } + + //@} + + ///////////////////////////////////////////////////////////////////////////// + /** @name Vertex iterators + */ + //@{ + typedef Complex_vertex_iterator CVI; // todo rename + + // + // @brief Range over the vertices of the simplicial complex. + // Methods .begin() and .end() return a Complex_vertex_iterator. + // + typedef boost::iterator_range< + Complex_vertex_iterator > Complex_vertex_range; + + /** + * @brief Returns a Complex_vertex_range over all vertices of the complex + */ + Complex_vertex_range vertex_range() const { + auto begin = Complex_vertex_iterator(this); + auto end = Complex_vertex_iterator(this, 0); + return Complex_vertex_range(begin, end); + } + + typedef boost::iterator_range< + Complex_neighbors_vertices_iterator > Complex_neighbors_vertices_range; + + /** + * @brief Returns a Complex_edge_range over all edges of the simplicial complex that passes trough v + */ + Complex_neighbors_vertices_range vertex_range(Vertex_handle v) const { + auto begin = Complex_neighbors_vertices_iterator( + this, v); + auto end = Complex_neighbors_vertices_iterator( + this, v, 0); + return Complex_neighbors_vertices_range(begin, end); + } + + //@} + + /** @name Edge iterators + */ + //@{ + + typedef boost::iterator_range< + Complex_edge_iterator>> Complex_edge_range; + + /** + * @brief Returns a Complex_edge_range over all edges of the simplicial complex + */ + Complex_edge_range edge_range() const { + auto begin = Complex_edge_iterator>(this); + auto end = Complex_edge_iterator>(this, 0); + return Complex_edge_range(begin, end); + } + + typedef boost::iterator_range >> + Complex_edge_around_vertex_range; + /** + * @brief Returns a Complex_edge_range over all edges of the simplicial complex that passes + * through 'v' + */ + Complex_edge_around_vertex_range edge_range(Vertex_handle v) const { + auto begin = Complex_edge_around_vertex_iterator>(this, v); + auto end = Complex_edge_around_vertex_iterator>(this, v, 0); + return Complex_edge_around_vertex_range(begin, end); + } + + //@} + + /** @name Triangles iterators + */ + //@{ +private: + typedef Skeleton_blocker_link_complex > Link; + typedef Skeleton_blocker_link_superior > Superior_link; + +public: + typedef Triangle_around_vertex_iterator Superior_triangle_around_vertex_iterator; + + typedef boost::iterator_range < Triangle_around_vertex_iterator > Complex_triangle_around_vertex_range; + + /** + * @brief Range over triangles around a vertex of the simplicial complex. + * Methods .begin() and .end() return a Triangle_around_vertex_iterator. + * + */ + Complex_triangle_around_vertex_range triangle_range(Vertex_handle v) const { + auto begin = Triangle_around_vertex_iterator(this, v); + auto end = Triangle_around_vertex_iterator(this, v, 0); + return Complex_triangle_around_vertex_range(begin, end); + } + + typedef boost::iterator_range > Complex_triangle_range; + + /** + * @brief Range over triangles of the simplicial complex. + * Methods .begin() and .end() return a Triangle_around_vertex_iterator. + * + */ + Complex_triangle_range triangle_range() const { + auto end = Triangle_iterator(this, 0); + if (empty()) { + return Complex_triangle_range(end, end); + } else { + auto begin = Triangle_iterator(this); + return Complex_triangle_range(begin, end); + } + } + + //@} + + /** @name Simplices iterators + */ + //@{ + typedef Simplex_around_vertex_iterator Complex_simplex_around_vertex_iterator; + + /** + * @brief Range over the simplices of the simplicial complex around a vertex. + * Methods .begin() and .end() return a Complex_simplex_around_vertex_iterator. + */ + typedef boost::iterator_range < Complex_simplex_around_vertex_iterator > Complex_simplex_around_vertex_range; + + /** + * @brief Returns a Complex_simplex_around_vertex_range over all the simplices around a vertex of the complex + */ + Complex_simplex_around_vertex_range simplex_range(Vertex_handle v) const { + assert(contains_vertex(v)); + return Complex_simplex_around_vertex_range( + Complex_simplex_around_vertex_iterator(this, v), + Complex_simplex_around_vertex_iterator(this, v, true)); + } + + // typedef Simplex_iterator Complex_simplex_iterator; + typedef Simplex_iterator Complex_simplex_iterator; + + typedef boost::iterator_range < Complex_simplex_iterator > Complex_simplex_range; + + /** + * @brief Returns a Complex_simplex_range over all the simplices of the complex + */ + Complex_simplex_range simplex_range() const { + Complex_simplex_iterator end(this, true); + if (empty()) { + return Complex_simplex_range(end, end); + } else { + Complex_simplex_iterator begin(this); + return Complex_simplex_range(begin, end); + } + } + + //@} + + /** @name Blockers iterators + */ + //@{ +private: + /** + * @brief Iterator over the blockers adjacent to a vertex + */ + typedef Blocker_iterator_around_vertex_internal< + typename std::multimap::iterator, + Blocker_handle> + Complex_blocker_around_vertex_iterator; + + /** + * @brief Iterator over (constant) blockers adjacent to a vertex + */ + typedef Blocker_iterator_around_vertex_internal< + typename std::multimap::const_iterator, + const Blocker_handle> + Const_complex_blocker_around_vertex_iterator; + + typedef boost::iterator_range Complex_blocker_around_vertex_range; + typedef boost::iterator_range Const_complex_blocker_around_vertex_range; + +public: + /** + * @brief Returns a range of the blockers of the complex passing through a vertex + */ + Complex_blocker_around_vertex_range blocker_range(Vertex_handle v) { + auto begin = Complex_blocker_around_vertex_iterator(blocker_map_.lower_bound(v)); + auto end = Complex_blocker_around_vertex_iterator(blocker_map_.upper_bound(v)); + return Complex_blocker_around_vertex_range(begin, end); + } + + /** + * @brief Returns a range of the blockers of the complex passing through a vertex + */ + Const_complex_blocker_around_vertex_range const_blocker_range(Vertex_handle v) const { + auto begin = Const_complex_blocker_around_vertex_iterator(blocker_map_.lower_bound(v)); + auto end = Const_complex_blocker_around_vertex_iterator(blocker_map_.upper_bound(v)); + return Const_complex_blocker_around_vertex_range(begin, end); + } + +private: + /** + * @brief Iterator over the blockers. + */ + typedef Blocker_iterator_internal< + typename std::multimap::iterator, + Blocker_handle> + Complex_blocker_iterator; + + /** + * @brief Iterator over the (constant) blockers. + */ + typedef Blocker_iterator_internal< + typename std::multimap::const_iterator, + const Blocker_handle> + Const_complex_blocker_iterator; + + typedef boost::iterator_range Complex_blocker_range; + typedef boost::iterator_range Const_complex_blocker_range; + +public: + /** + * @brief Returns a range of the blockers of the complex + */ + Complex_blocker_range blocker_range() { + auto begin = Complex_blocker_iterator(blocker_map_.begin(), blocker_map_.end() ); + auto end = Complex_blocker_iterator(blocker_map_.end() , blocker_map_.end() ); + return Complex_blocker_range(begin, end); + } + + /** + * @brief Returns a range of the blockers of the complex + */ + Const_complex_blocker_range const_blocker_range() const { + auto begin = Const_complex_blocker_iterator(blocker_map_.begin(), blocker_map_.end() ); + auto end = Const_complex_blocker_iterator(blocker_map_.end(), blocker_map_.end() ); + return Const_complex_blocker_range(begin, end); + } + + //@} + + ///////////////////////////////////////////////////////////////////////////// + /** @name Print and IO methods + */ + //@{ +public: + std::string to_string() const { + std::ostringstream stream; + stream << num_vertices() << " vertices:\n" << vertices_to_string() << std::endl; + stream << num_edges() << " edges:\n" << edges_to_string() << std::endl; + stream << num_blockers() << " blockers:\n" << blockers_to_string() << std::endl; + return stream.str(); + } + + std::string vertices_to_string() const { + std::ostringstream stream; + for (auto vertex : vertex_range()) { + stream << "(" << (*this)[vertex].get_id() << "),"; + } + stream << std::endl; + return stream.str(); + } + + std::string edges_to_string() const { + std::ostringstream stream; + for (auto edge : edge_range()) { + stream << "(" << (*this)[edge].first() << "," << (*this)[edge].second() << ")" << " id = " << (*this)[edge].index() << std::endl; + } + stream << std::endl; + return stream.str(); + } + + std::string blockers_to_string() const { + std::ostringstream stream; + for (auto bl : blocker_map_) { + stream << bl.first << " => " << bl.second << ":" << *bl.second << std::endl; + } + return stream.str(); + } + + //@} }; /** * build a simplicial complex from a collection * of top faces. + * return the total number of simplices */ template -unsigned make_complex_from_top_faces(Complex& complex, - SimplexHandleIterator begin, - SimplexHandleIterator end) { - typedef typename Complex::Simplex_handle Simplex_handle; - - int dimension = 0; - for (auto top_face = begin; top_face != end; ++top_face) - dimension = std::max(dimension, top_face->dimension()); - - std::vector < std::set - > simplices_per_dimension(dimension + 1); - - for (auto top_face = begin; top_face != end; ++top_face) { - register_faces(simplices_per_dimension, *top_face); - } - - // compute list of simplices - std::list simplices; - for (int dim = 0; dim <= dimension; ++dim) { - std::copy(simplices_per_dimension[dim].begin(), - simplices_per_dimension[dim].end(), - std::back_inserter(simplices)); - } - complex = Complex(simplices); - return simplices.size(); +unsigned make_complex_from_top_faces(Complex& complex,SimplexHandleIterator begin,SimplexHandleIterator end) { + typedef typename Complex::Simplex_handle Simplex_handle; + + std::vector simplices; + + for (auto top_face = begin; top_face != end; ++top_face) { + auto subfaces_topface = subfaces(*top_face); + simplices.insert(simplices.end(),subfaces_topface.begin(),subfaces_topface.end()); + } + complex = Complex(simplices.begin(),simplices.end()); + return simplices.size(); } } // namespace skbl diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h index 3bbfe738..db17fb28 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h @@ -78,11 +78,20 @@ class Skeleton_blocker_simplifiable_complex : public Skeleton_blocker_complex& simplices, Visitor* visitor_ = NULL) : Skeleton_blocker_complex(simplices, visitor_) { } + /** + * @brief Constructor with a list of simplices + * @details The list of simplices must be the list of simplices of a simplicial complex. + */ + template + explicit Skeleton_blocker_simplifiable_complex(SimpleHandleOutputIterator simplex_begin,SimpleHandleOutputIterator simplex_end,bool is_flag_complex = false,Visitor* visitor_ = NULL) : + Skeleton_blocker_complex(simplex_begin,simplex_end,is_flag_complex, visitor_) { } + + virtual ~Skeleton_blocker_simplifiable_complex() { } //@} diff --git a/src/Skeleton_blocker/test/CMakeLists.txt b/src/Skeleton_blocker/test/CMakeLists.txt index c341dcee..d70b5d7f 100644 --- a/src/Skeleton_blocker/test/CMakeLists.txt +++ b/src/Skeleton_blocker/test/CMakeLists.txt @@ -2,6 +2,8 @@ cmake_minimum_required(VERSION 2.6) project(GUDHIskbl) if(NOT MSVC) + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -g") + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} --coverage") set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} --coverage") set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} --coverage") diff --git a/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp b/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp index 865f6ff8..9392cc9e 100644 --- a/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp +++ b/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp @@ -1,24 +1,24 @@ - /* 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): David Salinas - * - * Copyright (C) 2014 INRIA Sophia Antipolis-Mediterranee (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 . - */ +/* 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): David Salinas + * + * Copyright (C) 2014 INRIA Sophia Antipolis-Mediterranee (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 . + */ #include #include #include @@ -30,6 +30,7 @@ #include "gudhi/Skeleton_blocker_link_complex.h" #include "gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h" #include "gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h" +#include "gudhi/Skeleton_blocker/internal/Trie.h" using namespace std; @@ -77,9 +78,9 @@ void build_complete(int n,Complex& complex){ for(int i=0;i=0;i--) -// for(int j=i-1;j>=0;j--) -// complex.add_edge(Vertex_handle(i),Vertex_handle(j)); + // for(int i=n-1;i>=0;i--) + // for(int j=i-1;j>=0;j--) + // complex.add_edge(Vertex_handle(i),Vertex_handle(j)); for(int i=0;i simplices(subfaces(simplex)); simplices.remove(simplex); - Complex complex(simplices); + Complex complex(simplices.begin(),simplices.end()); PRINT(complex.to_string()); @@ -746,6 +747,67 @@ bool test_constructor2(){ } +bool test_constructor3(){ + typedef Vertex_handle Vh; + typedef Simplex_handle Sh; + std::vector simplices; + auto subf(subfaces(Sh(Vh(0),Vh(1),Vh(2)))); + subf.pop_back(); //remove max face -> now a blocker 012 + simplices.insert(simplices.begin(),subf.begin(),subf.end()); + DBGCONT(simplices); + Complex complex(simplices.begin(),simplices.end()); + + DBGVALUE(complex.to_string()); + + return complex.num_vertices()==3 && complex.num_blockers()==1; +} + +bool test_constructor4(){ + typedef Vertex_handle Vh; + typedef Simplex_handle Sh; + std::vector simplices; + auto subf(subfaces(Sh(Vh(0),Vh(1),Vh(2),Vh(3)))); + // subf.pop_back(); //remove max face -> now a blocker 012 + simplices.insert(simplices.begin(),subf.begin(),subf.end()); + + simplices.push_back(Sh(Vh(4))); + simplices.push_back(Sh(Vh(4),Vh(1))); + simplices.push_back(Sh(Vh(4),Vh(0))); + + DBGCONT(simplices); + Complex complex(simplices.begin(),simplices.end()); + + DBGVALUE(complex.to_string()); + + return complex.num_vertices()==5 && complex.num_blockers()==1 && complex.num_edges()==8; +} + + + +bool test_constructor5(){ + typedef Vertex_handle Vh; + typedef Simplex_handle Sh; + std::vector simplices; + auto subf(subfaces(Sh(Vh(0),Vh(1),Vh(2)))); + simplices.insert(simplices.begin(),subf.begin(),subf.end()); + + simplices.push_back(Sh(Vh(3))); + simplices.push_back(Sh(Vh(3),Vh(1))); + simplices.push_back(Sh(Vh(3),Vh(2))); + simplices.push_back(Sh(Vh(4))); + simplices.push_back(Sh(Vh(4),Vh(1))); + simplices.push_back(Sh(Vh(4),Vh(0))); + simplices.push_back(Sh(Vh(5))); + simplices.push_back(Sh(Vh(5),Vh(2))); + simplices.push_back(Sh(Vh(5),Vh(0))); + + DBGCONT(simplices); + Complex complex(simplices.begin(),simplices.end()); + + DBGVALUE(complex.to_string()); + + return complex.num_vertices()==6 && complex.num_blockers()==3 && complex.num_edges()==9; +} int main (int argc, char *argv[]) @@ -778,6 +840,9 @@ int main (int argc, char *argv[]) tests_complex.add("test_constructor_list_simplices",test_constructor); tests_complex.add("test_constructor_list_simplices2",test_constructor2); + tests_complex.add("test_constructor_list_simplices3",test_constructor3); + tests_complex.add("test_constructor_list_simplices4",test_constructor4); + tests_complex.add("test_constructor_list_simplices5",test_constructor5); if(tests_complex.run()){ return EXIT_SUCCESS; diff --git a/src/common/include/gudhi/Utils.h b/src/common/include/gudhi/Utils.h index 7678685c..1626a0bf 100644 --- a/src/common/include/gudhi/Utils.h +++ b/src/common/include/gudhi/Utils.h @@ -25,7 +25,7 @@ #define PRINT(a) std::cerr << #a << ": " << (a) << " (DISP)"<