From c45a4897ee4bee39801efc7ce176f33b25a75fb0 Mon Sep 17 00:00:00 2001 From: vrouvrea Date: Tue, 20 Jan 2015 10:38:41 +0000 Subject: cpplint fixes git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/TDA_dev_1.1.0@416 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 3644577a2a7a8798f5f53d33c49c4862e3ec6787 --- .../example/Skeleton_blocker_iteration.cpp | 3 +- .../include/gudhi/Skeleton_blocker.h | 18 +- .../Skeleton_blocker_complex_visitor.h | 202 +- .../Skeleton_blocker_link_superior.h | 124 +- .../Skeleton_blocker/Skeleton_blocker_off_io.h | 110 +- .../Skeleton_blocker_simple_geometric_traits.h | 143 +- .../Skeleton_blocker_simple_traits.h | 302 +-- .../Skeleton_blocker/Skeleton_blocker_simplex.h | 734 +++--- .../Skeleton_blocker_sub_complex.h | 520 ++-- .../iterators/Skeleton_blockers_edges_iterators.h | 2 +- .../include/gudhi/Skeleton_blocker_complex.h | 2650 ++++++++++---------- .../gudhi/Skeleton_blocker_geometric_complex.h | 223 +- .../include/gudhi/Skeleton_blocker_link_complex.h | 550 ++-- .../gudhi/Skeleton_blocker_simplifiable_complex.h | 1090 ++++---- src/Skeleton_blocker/test/TestSimplifiable.cpp | 45 +- .../test/TestSkeletonBlockerComplex.cpp | 10 +- 16 files changed, 3316 insertions(+), 3410 deletions(-) (limited to 'src/Skeleton_blocker') diff --git a/src/Skeleton_blocker/example/Skeleton_blocker_iteration.cpp b/src/Skeleton_blocker/example/Skeleton_blocker_iteration.cpp index 832ae0ec..221c7881 100644 --- a/src/Skeleton_blocker/example/Skeleton_blocker_iteration.cpp +++ b/src/Skeleton_blocker/example/Skeleton_blocker_iteration.cpp @@ -21,6 +21,7 @@ */ #include + #include #include #include @@ -29,8 +30,6 @@ #include "gudhi/Skeleton_blocker.h" -//#include "gudhi/Skeleton_blocker_complex.h" -//#include "gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h" using namespace std; using namespace Gudhi; diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker.h index 8f7e1590..1d1a9bba 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker.h @@ -19,8 +19,8 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ -#ifndef GUDHI_SKELETON_BLOCKER_H -#define GUDHI_SKELETON_BLOCKER_H +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_H_ #include "gudhi/Skeleton_blocker_complex.h" #include "gudhi/Skeleton_blocker_geometric_complex.h" @@ -31,11 +31,11 @@ #include "gudhi/Skeleton_blocker/Skeleton_blocker_simple_geometric_traits.h" -#include "gudhi/Utils.h" //xxx +#include "gudhi/Utils.h" // xxx -namespace Gudhi{ -namespace skbl{ +namespace Gudhi { +namespace skbl { /** \defgroup skbl Skeleton-Blocker \author David Salinas @@ -189,12 +189,12 @@ their collaboration to write the two initial papers \copyright GNU General Public License v3. \verbatim Contact: David Salinas, david.salinas@inria.fr \endverbatim */ -/** @} */ // end defgroup -} -} +/** @} */ // end defgroup +} // namespace skbl +} // namespace Gudhi -#endif +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h index 3afd4cc7..829ab1e8 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h @@ -1,31 +1,30 @@ - /* 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 . - */ -#ifndef GUDHI_COMPLEXVISITOR_H_ -#define GUDHI_COMPLEXVISITOR_H_ - +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_COMPLEX_VISITOR_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_COMPLEX_VISITOR_H_ #include "gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h" -namespace Gudhi{ +namespace Gudhi { namespace skbl { // todo rajouter les const @@ -34,94 +33,103 @@ namespace skbl { *@class Skeleton_blocker_complex_visitor *@brief Interface for a visitor of a simplicial complex. */ -template +template class Skeleton_blocker_complex_visitor { -public: - virtual ~Skeleton_blocker_complex_visitor(){}; - - virtual void on_add_vertex(Vertex_handle) = 0; - virtual void on_remove_vertex(Vertex_handle) = 0; - - virtual void on_add_edge(Vertex_handle a,Vertex_handle b) = 0; - virtual void on_remove_edge(Vertex_handle a,Vertex_handle b) = 0; - - /** - * @brief Called when an edge changes, during the contraction of - * an edge - */ - virtual void on_changed_edge(Vertex_handle a,Vertex_handle b) = 0; - - /** - * @brief Called when performing an edge contraction when - * an edge bx is replaced by an edge ax (not already present). - * Precisely, this methods is called this way in contract_edge : - * add_edge(a,x) - * on_swaped_edge(a,b,x) - * remove_edge(b,x) - */ - virtual void on_swaped_edge(Vertex_handle a,Vertex_handle b,Vertex_handle x)=0; - virtual void on_add_blocker(const Skeleton_blocker_simplex&) = 0; - virtual void on_delete_blocker(const Skeleton_blocker_simplex*) = 0; + public: + virtual ~Skeleton_blocker_complex_visitor() {} + + virtual void on_add_vertex(Vertex_handle) = 0; + virtual void on_remove_vertex(Vertex_handle) = 0; + + virtual void on_add_edge(Vertex_handle a, Vertex_handle b) = 0; + virtual void on_remove_edge(Vertex_handle a, Vertex_handle b) = 0; + + /** + * @brief Called when an edge changes, during the contraction of + * an edge + */ + virtual void on_changed_edge(Vertex_handle a, Vertex_handle b) = 0; + + /** + * @brief Called when performing an edge contraction when + * an edge bx is replaced by an edge ax (not already present). + * Precisely, this methods is called this way in contract_edge : + * add_edge(a,x) + * on_swaped_edge(a,b,x) + * remove_edge(b,x) + */ + virtual void on_swaped_edge(Vertex_handle a, Vertex_handle b, + Vertex_handle x)=0; + virtual void on_add_blocker( + const Skeleton_blocker_simplex&) = 0; + virtual void on_delete_blocker( + const Skeleton_blocker_simplex*) = 0; }; - /** *@class Dummy_complex_visitor *@brief A dummy visitor of a simplicial complex that does nothing * */ -template -class Dummy_complex_visitor : public Skeleton_blocker_complex_visitor{ -public: - void on_add_vertex(Vertex_handle) {} - void on_remove_vertex(Vertex_handle){} - void on_add_edge(Vertex_handle a,Vertex_handle b){} - void on_remove_edge(Vertex_handle a,Vertex_handle b){} - void on_changed_edge(Vertex_handle a,Vertex_handle b){} - void on_swaped_edge(Vertex_handle a,Vertex_handle b,Vertex_handle x){} - void on_add_blocker(const Skeleton_blocker_simplex&){} - void on_delete_blocker(const Skeleton_blocker_simplex*){} - +template +class Dummy_complex_visitor : public Skeleton_blocker_complex_visitor< + Vertex_handle> { + public: + void on_add_vertex(Vertex_handle) { + } + void on_remove_vertex(Vertex_handle) { + } + void on_add_edge(Vertex_handle a, Vertex_handle b) { + } + void on_remove_edge(Vertex_handle a, Vertex_handle b) { + } + void on_changed_edge(Vertex_handle a, Vertex_handle b) { + } + void on_swaped_edge(Vertex_handle a, Vertex_handle b, Vertex_handle x) { + } + void on_add_blocker(const Skeleton_blocker_simplex&) { + } + void on_delete_blocker(const Skeleton_blocker_simplex*) { + } }; - - /** *@class Print_complex_visitor *@brief A visitor that prints the visit information. * */ -template -class Print_complex_visitor : public Skeleton_blocker_complex_visitor{ -public: - void on_add_vertex(Vertex_handle v) { - std::cerr << "on_add_vertex:"<& b){ - std::cerr << "on_add_blocker:"<* b){ - std::cerr << "on_delete_blocker:"< +class Print_complex_visitor : public Skeleton_blocker_complex_visitor< + Vertex_handle> { + public: + void on_add_vertex(Vertex_handle v) { + std::cerr << "on_add_vertex:" << v << std::endl; + } + void on_remove_vertex(Vertex_handle v) { + std::cerr << "on_remove_vertex:" << v << std::endl; + } + void on_add_edge(Vertex_handle a, Vertex_handle b) { + std::cerr << "on_add_edge:" << a << "," << b << std::endl; + } + void on_remove_edge(Vertex_handle a, Vertex_handle b) { + std::cerr << "on_remove_edge:" << a << "," << b << std::endl; + } + void on_changed_edge(Vertex_handle a, Vertex_handle b) { + std::cerr << "on_changed_edge:" << a << "," << b << std::endl; + } + void on_swaped_edge(Vertex_handle a, Vertex_handle b, Vertex_handle x) { + std::cerr << "on_swaped_edge:" << a << "," << b << "," << x << std::endl; + } + void on_add_blocker(const Skeleton_blocker_simplex& b) { + std::cerr << "on_add_blocker:" << b << std::endl; + } + void on_delete_blocker(const Skeleton_blocker_simplex* b) { + std::cerr << "on_delete_blocker:" << b << std::endl; + } }; -} +} // namespace skbl -} // namespace GUDHI +} // namespace Gudhi -#endif /* GUDHI_COMPLEXVISITOR_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_COMPLEX_VISITOR_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h index 4e1c5178..17d58956 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_link_superior.h @@ -1,81 +1,77 @@ - /* 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 . - */ -#ifndef GUDHI_SKELETON_BLOCKER_LINK_SUPERIOR_H_ -#define GUDHI_SKELETON_BLOCKER_LINK_SUPERIOR_H_ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_LINK_SUPERIOR_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_LINK_SUPERIOR_H_ #include "gudhi/Skeleton_blocker_link_complex.h" -namespace Gudhi{ +namespace Gudhi { namespace skbl { template class Skeleton_blocker_sub_complex; - /** * \brief Class representing the link of a simplicial complex encoded by a skeleton/blockers pair. * It computes only vertices greater than the simplex used to build the link. */ template -class Skeleton_blocker_link_superior : public Skeleton_blocker_link_complex -{ - typedef typename ComplexType::Edge_handle Edge_handle; - - typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; - -public: - - typedef typename ComplexType::Vertex_handle Vertex_handle; - typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; - typedef typename ComplexType::Simplex_handle Simplex_handle; - typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; - typedef typename ComplexType::BlockerMap BlockerMap; - typedef typename ComplexType::BlockerPair BlockerPair; - typedef typename ComplexType::BlockerMapIterator BlockerMapIterator; - typedef typename ComplexType::BlockerMapConstIterator BlockerMapConstIterator; - typedef typename ComplexType::Simplex_handle::Simplex_vertex_const_iterator AddressSimplexConstIterator; - typedef typename ComplexType::Root_simplex_handle::Simplex_vertex_const_iterator IdSimplexConstIterator; - - - Skeleton_blocker_link_superior() - :Skeleton_blocker_link_complex(true) - { - } - - Skeleton_blocker_link_superior(const ComplexType & parent_complex, Simplex_handle& alpha_parent_adress) - :Skeleton_blocker_link_complex(parent_complex,alpha_parent_adress,true) - { - } - - Skeleton_blocker_link_superior(const ComplexType & parent_complex, Vertex_handle a_parent_adress) - :Skeleton_blocker_link_complex(parent_complex,a_parent_adress,true) - { - } - +class Skeleton_blocker_link_superior : public Skeleton_blocker_link_complex< + ComplexType> { + typedef typename ComplexType::Edge_handle Edge_handle; + + typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; + + public: + typedef typename ComplexType::Vertex_handle Vertex_handle; + typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; + typedef typename ComplexType::Simplex_handle Simplex_handle; + typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; + typedef typename ComplexType::BlockerMap BlockerMap; + typedef typename ComplexType::BlockerPair BlockerPair; + typedef typename ComplexType::BlockerMapIterator BlockerMapIterator; + typedef typename ComplexType::BlockerMapConstIterator BlockerMapConstIterator; + typedef typename ComplexType::Simplex_handle::Simplex_vertex_const_iterator AddressSimplexConstIterator; + typedef typename ComplexType::Root_simplex_handle::Simplex_vertex_const_iterator IdSimplexConstIterator; + + Skeleton_blocker_link_superior() + : Skeleton_blocker_link_complex(true) { + } + + Skeleton_blocker_link_superior(const ComplexType & parent_complex, + Simplex_handle& alpha_parent_adress) + : Skeleton_blocker_link_complex(parent_complex, + alpha_parent_adress, true) { + } + + Skeleton_blocker_link_superior(const ComplexType & parent_complex, + Vertex_handle a_parent_adress) + : Skeleton_blocker_link_complex(parent_complex, + a_parent_adress, true) { + } }; -} - -} // namespace GUDHI +} // namespace skbl +} // namespace Gudhi -#endif /* SKELETON_BLOCKER_LINK_SUPERIOR_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_LINK_SUPERIOR_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_off_io.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_off_io.h index fd44fba5..bc832d5d 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_off_io.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_off_io.h @@ -19,8 +19,11 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ -#ifndef SKELETON_BLOCKER_OFF_IO_H_ -#define SKELETON_BLOCKER_OFF_IO_H_ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_OFF_IO_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_OFF_IO_H_ + +#include +#include #include "gudhi/Off_reader.h" @@ -32,71 +35,70 @@ namespace skbl { *@brief Off reader visitor that can be passed to Off_reader to read a Skeleton_blocker_complex. */ template -class Skeleton_blocker_off_visitor_reader{ - Complex& complex_; - typedef typename Complex::Vertex_handle Vertex_handle; - typedef typename Complex::Point Point; +class Skeleton_blocker_off_visitor_reader { + Complex& complex_; + typedef typename Complex::Vertex_handle Vertex_handle; + typedef typename Complex::Point Point; - const bool load_only_points_; + const bool load_only_points_; -public: - Skeleton_blocker_off_visitor_reader(Complex& complex,bool load_only_points = false): - complex_(complex), - load_only_points_(load_only_points){} + public: + explicit Skeleton_blocker_off_visitor_reader(Complex& complex, bool load_only_points = false): + complex_(complex), + load_only_points_(load_only_points) {} - void init(int dim,int num_vertices,int num_faces,int num_edges){ - //todo do an assert to check that this number are correctly read - //todo reserve size for vector points - } + void init(int dim, int num_vertices, int num_faces, int num_edges) { + // todo do an assert to check that this number are correctly read + // todo reserve size for vector points + } - void point(const std::vector& point){ - complex_.add_vertex(point); - } + void point(const std::vector& point) { + complex_.add_vertex(point); + } - void maximal_face(const std::vector& face){ - if (!load_only_points_){ - for(size_t i = 0 ; i < face.size();++i) - for(size_t j = i+1 ; j < face.size();++j){ - complex_.add_edge(Vertex_handle(face[i]),Vertex_handle(face[j])); - } - } - } + void maximal_face(const std::vector& face) { + if (!load_only_points_) { + for (size_t i = 0; i < face.size(); ++i) + for (size_t j = i+1; j < face.size(); ++j) { + complex_.add_edge(Vertex_handle(face[i]), Vertex_handle(face[j])); + } + } + } - void done(){ - } + void done() {} }; /** *@brief Class that allows to load a Skeleton_blocker_complex from an off file. */ template -class Skeleton_blocker_off_reader{ -public: - /** - * name_file : file to read - * read_complex : complex that will receive the file content - * read_only_points : specify true if only the points must be read - */ - Skeleton_blocker_off_reader(const std::string & name_file,Complex& read_complex,bool read_only_points = false):valid_(false){ - std::ifstream stream(name_file); - if(stream.is_open()){ - Skeleton_blocker_off_visitor_reader off_visitor(read_complex,read_only_points); - Off_reader off_reader(stream); - valid_ = off_reader.read(off_visitor); - } - } - - /** - * return true iff reading did not meet problems. - */ - bool is_valid() const{ - return valid_; - } - -private: - bool valid_; +class Skeleton_blocker_off_reader { + public: + /** + * name_file : file to read + * read_complex : complex that will receive the file content + * read_only_points : specify true if only the points must be read + */ + Skeleton_blocker_off_reader(const std::string & name_file, Complex& read_complex, bool read_only_points = false):valid_(false) { + std::ifstream stream(name_file); + if (stream.is_open()) { + Skeleton_blocker_off_visitor_reader off_visitor(read_complex, read_only_points); + Off_reader off_reader(stream); + valid_ = off_reader.read(off_visitor); + } + } + + /** + * return true iff reading did not meet problems. + */ + bool is_valid() const { + return valid_; + } + + private: + bool valid_; }; } // namespace skbl @@ -105,4 +107,4 @@ private: } // namespace Gudhi -#endif /* SKELETON_BLOCKER_OFF_IO_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_OFF_IO_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_geometric_traits.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_geometric_traits.h index a24bea25..d3a5b9d8 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_geometric_traits.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simple_geometric_traits.h @@ -1,35 +1,35 @@ - /* 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 . - */ -#ifndef GUDHI_SKELETON_BLOCKERS_SIMPLE_GEOMETRIC_TRAITS_H_ -#define GUDHI_SKELETON_BLOCKERS_SIMPLE_GEOMETRIC_TRAITS_H_ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_GEOMETRIC_TRAITS_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_GEOMETRIC_TRAITS_H_ #include #include -#include "gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h" -namespace Gudhi{ +#include "gudhi/Skeleton_blocker/Skeleton_blocker_simple_traits.h" -namespace skbl{ +namespace Gudhi { +namespace skbl { /** * @extends SkeletonBlockerGeometricDS @@ -38,48 +38,57 @@ namespace skbl{ * can be passed as a template argument to Skeleton_blocker_geometric_complex */ template -struct Skeleton_blocker_simple_geometric_traits : public skbl::Skeleton_blocker_simple_traits { -public: - - typedef GeometryTrait GT; - typedef typename GT::Point Point; - typedef typename Skeleton_blocker_simple_traits::Root_vertex_handle Root_vertex_handle; - typedef typename Skeleton_blocker_simple_traits::Graph_vertex Simple_vertex; - - /** - * @brief Vertex with a point attached. - */ - class Simple_geometric_vertex : public Simple_vertex{ - template friend class Skeleton_blocker_geometric_complex; - private: - Point point_; - - Point& point(){ return point_; } - const Point& point() const { return point_; } - }; - - - class Simple_geometric_edge : public Skeleton_blocker_simple_traits::Graph_edge{ - int index_; - public: - Simple_geometric_edge():Skeleton_blocker_simple_traits::Graph_edge(),index_(-1){} - /** - * @brief Allows to modify the index of the edge. - * The indices of the edge are used to store heap information - * in the edge contraction algorithm. - */ - int& index(){return index_;} - int index() const {return index_;} - }; - - - typedef Simple_geometric_vertex Graph_vertex; - typedef Skeleton_blocker_simple_traits::Graph_edge Graph_edge; +struct Skeleton_blocker_simple_geometric_traits : + public skbl::Skeleton_blocker_simple_traits { + public: + typedef GeometryTrait GT; + typedef typename GT::Point Point; + typedef typename Skeleton_blocker_simple_traits::Root_vertex_handle Root_vertex_handle; + typedef typename Skeleton_blocker_simple_traits::Graph_vertex Simple_vertex; + + /** + * @brief Vertex with a point attached. + */ + class Simple_geometric_vertex : public Simple_vertex { + template friend class Skeleton_blocker_geometric_complex; + private: + Point point_; + + Point& point() { + return point_; + } + const Point& point() const { + return point_; + } + }; + + class Simple_geometric_edge : + public Skeleton_blocker_simple_traits::Graph_edge { + int index_; + public: + Simple_geometric_edge() + : Skeleton_blocker_simple_traits::Graph_edge(), + index_(-1) { + } + /** + * @brief Allows to modify the index of the edge. + * The indices of the edge are used to store heap information + * in the edge contraction algorithm. + */ + int& index() { + return index_; + } + int index() const { + return index_; + } + }; + + typedef Simple_geometric_vertex Graph_vertex; + typedef Skeleton_blocker_simple_traits::Graph_edge Graph_edge; }; +} // namespace skbl -} - -} // namespace GUDHI +} // namespace Gudhi -#endif /* GUDHI_SKELETON_BLOCKERS_SIMPLE_GEOMETRIC_TRAITS_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_GEOMETRIC_TRAITS_H_ 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 e17ea9b5..b8506b2c 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 @@ -1,32 +1,32 @@ - /* 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 . - */ -#ifndef GUDHI_SKELETON_BLOCKERS_SIMPLE_TRAITS_H_ -#define GUDHI_SKELETON_BLOCKERS_SIMPLE_TRAITS_H_ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_TRAITS_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_TRAITS_H_ #include #include #include "Skeleton_blocker_simplex.h" -namespace Gudhi{ +namespace Gudhi { namespace skbl { @@ -36,124 +36,144 @@ namespace skbl { * @brief Simple traits that is a model of SkeletonBlockerDS and * can be passed as a template argument to Skeleton_blocker_complex */ -struct Skeleton_blocker_simple_traits{ - /** - * @brief Global and local handle similar to boost subgraphs. - * Vertices are stored in a vector. - * For the root simplicial complex, the local and global descriptors are the same. - * For a subcomplex L and one of its vertices 'v', the local descriptor of 'v' is its position in - * the vertex vector of the subcomplex L whereas its global descriptor is the position of 'v' - * in the vertex vector of the root simplicial complex. - */ - struct Root_vertex_handle{ - typedef int boost_vertex_handle; - explicit Root_vertex_handle(boost_vertex_handle val = -1 ):vertex(val){} - boost_vertex_handle vertex; - - bool operator!=( const Root_vertex_handle& other) const{ - return ! (this->vertex == other.vertex); - } - - bool operator==( const Root_vertex_handle& other) const{ - return this->vertex == other.vertex; - } - - bool operator<( const Root_vertex_handle& other) const{ - return this->vertex < other.vertex; - } - - friend std::ostream& operator << (std::ostream& o, const Root_vertex_handle & v) - { - o << v.vertex; - return o; - } - }; - - struct Vertex_handle{ - typedef int boost_vertex_handle; - Vertex_handle(boost_vertex_handle val=-1):vertex(val){} - - boost_vertex_handle vertex; - - bool operator==( const Vertex_handle& other) const{ - return this->vertex == other.vertex; - } - - bool operator!=( const Vertex_handle& other) const{ - return this->vertex != other.vertex; - } - - bool operator<(const Vertex_handle& other) const{ - return this->vertex < other.vertex; - } - - friend std::ostream& operator << (std::ostream& o, const Vertex_handle & v) - { - o << v.vertex; - return o; - } - }; - - - - class Graph_vertex { - bool is_active_; - Root_vertex_handle id_; - public: - virtual ~Graph_vertex(){} - - void activate(){is_active_=true;} - void deactivate(){is_active_=false;} - bool is_active() const{return is_active_;} - void set_id(Root_vertex_handle i){id_=i;} - Root_vertex_handle get_id() const{return id_;} - - virtual std::string to_string() const { - std::ostringstream res; - res<< id_; - return res.str(); - } - - friend std::ostream& operator << (std::ostream& o, const Graph_vertex & v){ - o << v.to_string(); - return o; - } - }; - - class Graph_edge { - Root_vertex_handle a_; - Root_vertex_handle b_; - int index_; - public: - - Graph_edge():a_(-1),b_(-1),index_(-1){} - - int& index(){return index_;} - int index() const {return index_;} - - void setId(Root_vertex_handle a,Root_vertex_handle b){ - a_ = a; - b_ = b; - } - - Root_vertex_handle first() const { - return a_; - } - - Root_vertex_handle second() const { - return b_; - } - - friend std::ostream& operator << (std::ostream& o, const Graph_edge & v){ - o << "("<boost subgraphs. + * Vertices are stored in a vector. + * For the root simplicial complex, the local and global descriptors are the same. + * For a subcomplex L and one of its vertices 'v', the local descriptor of 'v' is its position in + * the vertex vector of the subcomplex L whereas its global descriptor is the position of 'v' + * in the vertex vector of the root simplicial complex. + */ + struct Root_vertex_handle { + typedef int boost_vertex_handle; + explicit Root_vertex_handle(boost_vertex_handle val = -1) + : vertex(val) { + } + boost_vertex_handle vertex; + + bool operator!=(const Root_vertex_handle& other) const { + return !(this->vertex == other.vertex); + } + + bool operator==(const Root_vertex_handle& other) const { + return this->vertex == other.vertex; + } + + bool operator<(const Root_vertex_handle& other) const { + return this->vertex < other.vertex; + } + + friend std::ostream& operator <<(std::ostream& o, + const Root_vertex_handle & v) { + o << v.vertex; + return o; + } + }; + + struct Vertex_handle { + typedef int boost_vertex_handle; + explicit Vertex_handle(boost_vertex_handle val = -1) + : vertex(val) { + } + + boost_vertex_handle vertex; + + bool operator==(const Vertex_handle& other) const { + return this->vertex == other.vertex; + } + + bool operator!=(const Vertex_handle& other) const { + return this->vertex != other.vertex; + } + + bool operator<(const Vertex_handle& other) const { + return this->vertex < other.vertex; + } + + friend std::ostream& operator <<(std::ostream& o, const Vertex_handle & v) { + o << v.vertex; + return o; + } + }; + + class Graph_vertex { + bool is_active_; + Root_vertex_handle id_; + + public: + virtual ~Graph_vertex() { + } + + void activate() { + is_active_ = true; + } + void deactivate() { + is_active_ = false; + } + bool is_active() const { + return is_active_; + } + void set_id(Root_vertex_handle i) { + id_ = i; + } + Root_vertex_handle get_id() const { + return id_; + } + + virtual std::string to_string() const { + std::ostringstream res; + res << id_; + return res.str(); + } + + friend std::ostream& operator <<(std::ostream& o, const Graph_vertex & v) { + o << v.to_string(); + return o; + } + }; + + class Graph_edge { + Root_vertex_handle a_; + Root_vertex_handle b_; + int index_; + + public: + Graph_edge() + : a_(-1), + b_(-1), + index_(-1) { + } + + int& index() { + return index_; + } + int index() const { + return index_; + } + + void setId(Root_vertex_handle a, Root_vertex_handle b) { + a_ = a; + b_ = b; + } + + Root_vertex_handle first() const { + return a_; + } + + Root_vertex_handle second() const { + return b_; + } + + friend std::ostream& operator <<(std::ostream& o, const Graph_edge & v) { + o << "(" << v.a_ << "," << v.b_ << " - id = " << v.index(); + return o; + } + }; }; -} +} // namespace skbl -} // namespace GUDHI +} // namespace Gudhi -#endif /* GUDHI_SKELETON_BLOCKERS_SIMPLE_TRAITS_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLE_TRAITS_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h index fb3d9470..be997b6b 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h @@ -1,35 +1,37 @@ - /* 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 . - */ - -#ifndef GUDHI_SKELETON_BLOCKER_SIMPLEX_H -#define GUDHI_SKELETON_BLOCKER_SIMPLEX_H - -#include -#include -#include -#include +/* 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 . + */ + +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLEX_H_ + +#include +#include +#include +#include #include +#include +#include -namespace Gudhi{ +namespace Gudhi { namespace skbl { @@ -43,355 +45,337 @@ namespace skbl { * * */ -template +template class Skeleton_blocker_simplex { - -private : - std::set simplex_set; - -public: - typedef typename T::boost_vertex_handle boost_vertex_handle; - - typedef T Vertex_handle; - - - typedef typename std::set::const_iterator Simplex_vertex_const_iterator; - typedef typename std::set::iterator Simplex_vertex_iterator; - - /** @name Constructors / Destructors / Initialization - */ - //@{ - - -// Skeleton_blocker_simplex():simplex_set() {} - - /** - * Clear the simplex - */ - void clear() { - simplex_set.clear(); - } - - - Skeleton_blocker_simplex(std::initializer_list& list) { - for_each(list.begin(),list.end(),add_vertex); - } - - - template - Skeleton_blocker_simplex(Args... args){ - add_vertices(args...); - } - - - template - void add_vertices(T v,Args... args){ - add_vertex(v); - add_vertices(args...); - } - - void add_vertices(T v){ - add_vertex(v); - } - - void add_vertices(){ - } - - /** - * Initialize a simplex with a string such as {0,1,2} - */ - Skeleton_blocker_simplex(std::string token){ - clear(); - if ((token[0] == '{') && (token[token.size()-1] == '}' ) ) - { - token.erase(0,1); - token.erase(token.size()-1,1); - while(token.size()!=0 ){ - int coma_position=token.find_first_of(','); - //cout << "coma_position:"< v; - v.reserve(std::min(simplex_set.size(), a.simplex_set.size())); - - set_intersection(simplex_set.begin(),simplex_set.end(), - a.simplex_set.begin(),a.simplex_set.end(), - std::back_inserter(v)); - clear(); - for (auto i:v) - simplex_set.insert(i); - } - - /** - * Substracts a from the simplex: - * \f$ (*this) \leftarrow (*this) \setminus a \f$ - */ - void difference(const Skeleton_blocker_simplex & a){ - std::vector v; - v.reserve(simplex_set.size()); - - set_difference(simplex_set.begin(),simplex_set.end(), - a.simplex_set.begin(),a.simplex_set.end(), - std::back_inserter(v)); - - clear(); - for (auto i:v) - simplex_set.insert(i); - } - - /** - * Add vertices of a to the simplex: - * \f$ (*this) \leftarrow (*this) \cup a \f$ - */ - void union_vertices(const Skeleton_blocker_simplex & a){ - std::vector v; - v.reserve(simplex_set.size() + a.simplex_set.size()); - - set_union(simplex_set.begin(),simplex_set.end(), - a.simplex_set.begin(),a.simplex_set.end(), - std::back_inserter(v)); - clear(); - simplex_set.insert(v.begin(),v.end()); - } - - typename std::set::const_iterator begin() const{ - return simplex_set.cbegin(); - } - - typename std::set::const_iterator end() const{ - return simplex_set.cend(); - } - - typename std::set::iterator begin(){ - return simplex_set.begin(); - } - - typename std::set::iterator end(){ - return simplex_set.end(); - } - - - - //@} - - /** @name Queries - */ - //@{ - - /** - * Returns the dimension of the simplex. - */ - int dimension() const{ - return (simplex_set.size() - 1); - } - - bool empty() const{ - return simplex_set.empty(); - } - - /** - * Returns the first vertex of the (oriented) simplex. - * - * Be careful : assumes the simplex is non-empty. - */ - T first_vertex() const{ - assert(!empty()); - return *(begin()); - } - - /** - * Returns the last vertex of the (oriented) simplex. - * - * Be careful : assumes the simplex is non-empty. - */ - T last_vertex() const - { - assert(!empty()); - return *(simplex_set.rbegin()); - } - /** - * @return true iff the simplex contains the simplex a, i.e. iff \f$ a \subset (*this) \f$. - */ - bool contains(const Skeleton_blocker_simplex & a) const{ - return includes(simplex_set.cbegin(),simplex_set.cend(),a.simplex_set.cbegin(),a.simplex_set.cend()); - } - - /** - * @return true iff the simplex contains the difference \f$ a \setminus b \f$. - */ - bool contains_difference(const Skeleton_blocker_simplex& a, const Skeleton_blocker_simplex& b) const{ - auto first1 = begin(); - auto last1 = end(); - - auto first2 = a.begin(); - auto last2 = a.end(); - - while (first2!=last2) { - // we ignore vertices of b - if(b.contains(*first2)){ - ++first2; - } - else{ - if ( (first1==last1) || (*first2<*first1) ) return false; - if (!(*first1<*first2)) ++first2; - ++first1; - } - } - return true; - } - - /** - * @return true iff the simplex contains the difference \f$ a \setminus \{ x \} \f$. - */ - bool contains_difference(const Skeleton_blocker_simplex& a, T x) const{ - auto first1 = begin(); - auto last1 = end(); - - auto first2 = a.begin(); - auto last2 = a.end(); - - while (first2!=last2) { - // we ignore vertices of b - if(x == *first2){ - ++first2; - } - else{ - if ( (first1==last1) || (*first2<*first1) ) return false; - if (!(*first1<*first2)) ++first2; - ++first1; - } - } - return true; - } - - /** - * @return true iff the simplex contains the difference \f$ a \setminus \{ x \} \f$. - */ - bool contains_difference(const Skeleton_blocker_simplex& a, T x, T y) const{ - auto first1 = begin(); - auto last1 = end(); - - auto first2 = a.begin(); - auto last2 = a.end(); - - while (first2!=last2) { - // we ignore vertices of b - if(x == *first2 || y == *first2){ - ++first2; - } - else{ - if ( (first1==last1) || (*first2<*first1) ) return false; - if (!(*first1<*first2)) ++first2; - ++first1; - } - } - return true; - } - - - /** - * @return true iff the simplex contains the vertex v, i.e. iff \f$ v \in (*this) \f$. - */ - bool contains(T v) const{ - return (simplex_set.find(v) != simplex_set.end()); - } - - /** - * @return \f$ (*this) \cap a = \emptyset \f$. - */ - bool disjoint(const Skeleton_blocker_simplex& a) const{ - std::vector v; - v.reserve(std::min(simplex_set.size(), a.simplex_set.size())); - - set_intersection(simplex_set.cbegin(),simplex_set.cend(), - a.simplex_set.cbegin(),a.simplex_set.cend(), - std::back_inserter(v)); - - return (v.size()==0); - } - - - bool operator==(const Skeleton_blocker_simplex& other) const{ - return (this->simplex_set == other.simplex_set); - } - - bool operator!=(const Skeleton_blocker_simplex& other) const{ - return (this->simplex_set != other.simplex_set); - } - - bool operator<(const Skeleton_blocker_simplex& other) const{ - return (std::lexicographical_compare(this->simplex_set.begin(),this->simplex_set.end(), - other.begin(),other.end())); - } - - //@} - - - - - - /** - * Display a simplex - */ - friend std::ostream& operator << (std::ostream& o, const Skeleton_blocker_simplex & sigma) - { - bool first = true; - o << "{"; - for(auto i : sigma) - { - if(first) first = false ; - else o << ","; - o << i; - } - o << "}"; - return o; - } - - + private: + std::set simplex_set; + + public: + typedef typename T::boost_vertex_handle boost_vertex_handle; + + typedef T Vertex_handle; + + typedef typename std::set::const_iterator Simplex_vertex_const_iterator; + typedef typename std::set::iterator Simplex_vertex_iterator; + + /** @name Constructors / Destructors / Initialization + */ + //@{ + + // Skeleton_blocker_simplex():simplex_set() {} + /** + * Clear the simplex + */ + void clear() { + simplex_set.clear(); + } + + Skeleton_blocker_simplex(std::initializer_list& list) { + for_each(list.begin(), list.end(), add_vertex); + } + + template + explicit Skeleton_blocker_simplex(Args ... args) { + add_vertices(args...); + } + + template + void add_vertices(T v, Args ... args) { + add_vertex(v); + add_vertices(args...); + } + + void add_vertices(T v) { + add_vertex(v); + } + + void add_vertices() { + } + + /** + * Initialize a simplex with a string such as {0,1,2} + */ + explicit Skeleton_blocker_simplex(std::string token) { + clear(); + if ((token[0] == '{') && (token[token.size() - 1] == '}')) { + token.erase(0, 1); + token.erase(token.size() - 1, 1); + while (token.size() != 0) { + int coma_position = token.find_first_of(','); + // cout << "coma_position:"< v; + v.reserve(std::min(simplex_set.size(), a.simplex_set.size())); + + set_intersection(simplex_set.begin(), simplex_set.end(), + a.simplex_set.begin(), a.simplex_set.end(), + std::back_inserter(v)); + clear(); + for (auto i : v) + simplex_set.insert(i); + } + + /** + * Substracts a from the simplex: + * \f$ (*this) \leftarrow (*this) \setminus a \f$ + */ + void difference(const Skeleton_blocker_simplex & a) { + std::vector v; + v.reserve(simplex_set.size()); + + set_difference(simplex_set.begin(), simplex_set.end(), + a.simplex_set.begin(), a.simplex_set.end(), + std::back_inserter(v)); + + clear(); + for (auto i : v) + simplex_set.insert(i); + } + + /** + * Add vertices of a to the simplex: + * \f$ (*this) \leftarrow (*this) \cup a \f$ + */ + void union_vertices(const Skeleton_blocker_simplex & a) { + std::vector v; + v.reserve(simplex_set.size() + a.simplex_set.size()); + + set_union(simplex_set.begin(), simplex_set.end(), a.simplex_set.begin(), + a.simplex_set.end(), std::back_inserter(v)); + clear(); + simplex_set.insert(v.begin(), v.end()); + } + + typename std::set::const_iterator begin() const { + return simplex_set.cbegin(); + } + + typename std::set::const_iterator end() const { + return simplex_set.cend(); + } + + typename std::set::iterator begin() { + return simplex_set.begin(); + } + + typename std::set::iterator end() { + return simplex_set.end(); + } + + //@} + + /** @name Queries + */ + //@{ + /** + * Returns the dimension of the simplex. + */ + int dimension() const { + return (simplex_set.size() - 1); + } + + bool empty() const { + return simplex_set.empty(); + } + + /** + * Returns the first vertex of the (oriented) simplex. + * + * Be careful : assumes the simplex is non-empty. + */ + T first_vertex() const { + assert(!empty()); + return *(begin()); + } + + /** + * Returns the last vertex of the (oriented) simplex. + * + * Be careful : assumes the simplex is non-empty. + */ + T last_vertex() const { + assert(!empty()); + return *(simplex_set.rbegin()); + } + /** + * @return true iff the simplex contains the simplex a, i.e. iff \f$ a \subset (*this) \f$. + */ + bool contains(const Skeleton_blocker_simplex & a) const { + return includes(simplex_set.cbegin(), simplex_set.cend(), + a.simplex_set.cbegin(), a.simplex_set.cend()); + } + + /** + * @return true iff the simplex contains the difference \f$ a \setminus b \f$. + */ + bool contains_difference(const Skeleton_blocker_simplex& a, + const Skeleton_blocker_simplex& b) const { + auto first1 = begin(); + auto last1 = end(); + + auto first2 = a.begin(); + auto last2 = a.end(); + + while (first2 != last2) { + // we ignore vertices of b + if (b.contains(*first2)) { + ++first2; + } else { + if ((first1 == last1) || (*first2 < *first1)) + return false; + if (!(*first1 < *first2)) + ++first2; + ++first1; + } + } + return true; + } + + /** + * @return true iff the simplex contains the difference \f$ a \setminus \{ x \} \f$. + */ + bool contains_difference(const Skeleton_blocker_simplex& a, T x) const { + auto first1 = begin(); + auto last1 = end(); + + auto first2 = a.begin(); + auto last2 = a.end(); + + while (first2 != last2) { + // we ignore vertices of b + if (x == *first2) { + ++first2; + } else { + if ((first1 == last1) || (*first2 < *first1)) + return false; + if (!(*first1 < *first2)) + ++first2; + ++first1; + } + } + return true; + } + + /** + * @return true iff the simplex contains the difference \f$ a \setminus \{ x \} \f$. + */ + bool contains_difference(const Skeleton_blocker_simplex& a, T x, T y) const { + auto first1 = begin(); + auto last1 = end(); + + auto first2 = a.begin(); + auto last2 = a.end(); + + while (first2 != last2) { + // we ignore vertices of b + if (x == *first2 || y == *first2) { + ++first2; + } else { + if ((first1 == last1) || (*first2 < *first1)) + return false; + if (!(*first1 < *first2)) + ++first2; + ++first1; + } + } + return true; + } + + /** + * @return true iff the simplex contains the vertex v, i.e. iff \f$ v \in (*this) \f$. + */ + bool contains(T v) const { + return (simplex_set.find(v) != simplex_set.end()); + } + + /** + * @return \f$ (*this) \cap a = \emptyset \f$. + */ + bool disjoint(const Skeleton_blocker_simplex& a) const { + std::vector v; + v.reserve(std::min(simplex_set.size(), a.simplex_set.size())); + + set_intersection(simplex_set.cbegin(), simplex_set.cend(), + a.simplex_set.cbegin(), a.simplex_set.cend(), + std::back_inserter(v)); + + return (v.size() == 0); + } + + bool operator==(const Skeleton_blocker_simplex& other) const { + return (this->simplex_set == other.simplex_set); + } + + bool operator!=(const Skeleton_blocker_simplex& other) const { + return (this->simplex_set != other.simplex_set); + } + + bool operator<(const Skeleton_blocker_simplex& other) const { + return (std::lexicographical_compare(this->simplex_set.begin(), + this->simplex_set.end(), other.begin(), + other.end())); + } + + //@} + + /** + * Display a simplex + */ + friend std::ostream& operator <<(std::ostream& o, + const Skeleton_blocker_simplex & sigma) { + bool first = true; + o << "{"; + for (auto i : sigma) { + if (first) + first = false; + else + o << ","; + o << i; + } + o << "}"; + return o; + } }; -} - -} // namespace GUDHI +} // namespace skbl -#endif +} // namespace Gudhi +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SIMPLEX_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h index 9b4253d1..e906df75 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h @@ -1,33 +1,36 @@ - /* 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 . - */ - -#ifndef GUDHI_SKELETON_BLOCKER_SUB_COMPLEX_H -#define GUDHI_SKELETON_BLOCKER_SUB_COMPLEX_H +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ + +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_ + +#include +#include #include "gudhi/Skeleton_blocker_complex.h" #include "gudhi/Skeleton_blocker/Skeleton_blocker_simplex.h" #include "gudhi/Utils.h" -namespace Gudhi{ +namespace Gudhi { namespace skbl { @@ -52,266 +55,251 @@ namespace skbl { * */ template -class Skeleton_blocker_sub_complex : public ComplexType -{ - -protected: - - template friend class Skeleton_blocker_link_complex; - - typedef typename ComplexType::Graph Graph; - typedef typename ComplexType::Edge_handle Edge_handle; - - typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; - - -public: - using ComplexType::add_vertex; - using ComplexType::add_edge; - using ComplexType::add_blocker; - - typedef typename ComplexType::Vertex_handle Vertex_handle; - typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; - typedef typename ComplexType::Simplex_handle Simplex_handle; - typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; - - -protected: - - - ///** - //* @brief Returns true iff the simplex formed by all vertices contained in 'addresses_sigma_in_link' - //* but 'vertex_to_be_ignored' is in 'link' - //*/ - /* - template friend bool - proper_face_in_union( - Skeleton_blocker_sub_complex & link, - std::vector > & addresses_sigma_in_link, - int vertex_to_be_ignored);*/ - - /** - * @brief Determines whether all proper faces of simplex 'sigma' belong to 'link1' \cup 'link2' - * where 'link1' and 'link2' are subcomplexes of the same complex of type ComplexType - */ - // template friend bool - // proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2){ - - //template friend bool - //proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2); - - typedef std::map IdAddressMap; - typedef typename IdAddressMap::value_type AddressPair; - typedef typename IdAddressMap::iterator IdAddressMapIterator; - typedef typename IdAddressMap::const_iterator IdAddressMapConstIterator; - std::map adresses; - - -public: - /** - * Add a vertex 'global' of K to L. When added to L, this vertex will receive - * another number, addresses(global), its local adress. - * return the adress where the vertex lay on L. - * The vertex corresponding to 'global' must not be already present - * in the complex. - */ - Vertex_handle add_vertex(Root_vertex_handle global){ - assert(!this->contains_vertex(global)); - Vertex_handle address(boost::add_vertex(this->skeleton)); - this->num_vertices_++; - (*this)[address].activate(); - (*this)[address].set_id(global); - adresses.insert(AddressPair(global, address)); - this->degree_.push_back(0); - return address; - } - - - /** - * Add an edge (v1_root,v2_root) to the sub-complex. - * It assumes that both vertices corresponding to v1_root and v2_root are present - * in the sub-complex. - */ - void add_edge(Root_vertex_handle v1_root, Root_vertex_handle v2_root){ - auto v1_sub(this->get_address(v1_root)); - auto v2_sub(this->get_address(v2_root)); - assert(v1_sub && v2_sub); - this->ComplexType::add_edge(*v1_sub, *v2_sub); - } - - /** - * Add a blocker to the sub-complex. - * It assumes that all vertices of blocker_root are present - * in the sub-complex. - */ - void add_blocker(const Root_simplex_handle& blocker_root){ - auto blocker_sub = this->get_address(blocker_root); - assert(blocker_sub); - this->add_blocker(new Simplex_handle(*blocker_sub)); - } - - - -public: - - /** - * Constructs the restricted complex of 'parent_complex' to - * vertices of 'simplex'. - */ - void make_restricted_complex(const ComplexType & parent_complex, const Simplex_handle& simplex){ - this->clear(); - // add vertices to the sub complex - for (auto x : simplex){ - assert(parent_complex.contains_vertex(x)); - auto x_local = this->add_vertex(parent_complex[x].get_id()); - (*this)[x_local] = parent_complex[x]; - } - - // add edges to the sub complex - for (auto x : simplex){ - // x_neigh is the neighbor of x intersected with vertices_simplex - Simplex_handle x_neigh; - parent_complex.add_neighbours(x, x_neigh, true); - x_neigh.intersection(simplex); - for (auto y : x_neigh){ - this->add_edge(parent_complex[x].get_id(), parent_complex[y].get_id()); - } - } - - // add blockers to the sub complex - for (auto blocker : parent_complex.const_blocker_range()){ - // check if it is the first time we encounter the blocker - if (simplex.contains(*blocker)){ - Root_simplex_handle blocker_root(parent_complex.get_id(*(blocker))); - Simplex_handle blocker_restr(*(this->get_simplex_address(blocker_root))); - this->add_blocker(new Simplex_handle(blocker_restr)); - } - } - } - - - - void clear(){ - adresses.clear(); - ComplexType::clear(); - } - - /** - * Compute the local vertex in L corresponding to the vertex global in K. - * runs in O(log n) if n = num_vertices() - */ - boost::optional get_address(Root_vertex_handle global) const{ - boost::optional res; - IdAddressMapConstIterator it = adresses.find(global); - if (it == adresses.end()) res.reset(); - else res = (*it).second; - return res; - } - - // /** - // * Allocates a simplex in L corresponding to the simplex s in K - // * with its local adresses and returns an AddressSimplex. - // */ - // boost::optional get_address(const Root_simplex_handle & s) const; - - -//private: - /** - * same as get_address except that it will return a simplex in any case. - * The vertices that were not found are not added. - */ - // @remark should be private but problem with VS - std::vector > get_addresses(const Root_simplex_handle & s) const{ - std::vector > res; - for (auto i : s) - { - res.push_back(get_address(i)); - } - return res; - } - +class Skeleton_blocker_sub_complex : public ComplexType { + protected: + template friend class Skeleton_blocker_link_complex; + + typedef typename ComplexType::Graph Graph; + typedef typename ComplexType::Edge_handle Edge_handle; + + typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; + + public: + using ComplexType::add_vertex; + using ComplexType::add_edge; + using ComplexType::add_blocker; + + typedef typename ComplexType::Vertex_handle Vertex_handle; + typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; + typedef typename ComplexType::Simplex_handle Simplex_handle; + typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; + + protected: + ///** + //* @brief Returns true iff the simplex formed by all vertices contained in 'addresses_sigma_in_link' + //* but 'vertex_to_be_ignored' is in 'link' + //*/ + /* + template friend bool + proper_face_in_union( + Skeleton_blocker_sub_complex & link, + std::vector > & addresses_sigma_in_link, + int vertex_to_be_ignored);*/ + + /** + * @brief Determines whether all proper faces of simplex 'sigma' belong to 'link1' \cup 'link2' + * where 'link1' and 'link2' are subcomplexes of the same complex of type ComplexType + */ + // template friend bool + // proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2){ + // template friend bool + // proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2); + typedef std::map IdAddressMap; + typedef typename IdAddressMap::value_type AddressPair; + typedef typename IdAddressMap::iterator IdAddressMapIterator; + typedef typename IdAddressMap::const_iterator IdAddressMapConstIterator; + std::map adresses; + + public: + /** + * Add a vertex 'global' of K to L. When added to L, this vertex will receive + * another number, addresses(global), its local adress. + * return the adress where the vertex lay on L. + * The vertex corresponding to 'global' must not be already present + * in the complex. + */ + Vertex_handle add_vertex(Root_vertex_handle global) { + assert(!this->contains_vertex(global)); + Vertex_handle address(boost::add_vertex(this->skeleton)); + this->num_vertices_++; + (*this)[address].activate(); + (*this)[address].set_id(global); + adresses.insert(AddressPair(global, address)); + this->degree_.push_back(0); + return address; + } + + /** + * Add an edge (v1_root,v2_root) to the sub-complex. + * It assumes that both vertices corresponding to v1_root and v2_root are present + * in the sub-complex. + */ + void add_edge(Root_vertex_handle v1_root, Root_vertex_handle v2_root) { + auto v1_sub(this->get_address(v1_root)); + auto v2_sub(this->get_address(v2_root)); + assert(v1_sub && v2_sub); + this->ComplexType::add_edge(*v1_sub, *v2_sub); + } + + /** + * Add a blocker to the sub-complex. + * It assumes that all vertices of blocker_root are present + * in the sub-complex. + */ + void add_blocker(const Root_simplex_handle& blocker_root) { + auto blocker_sub = this->get_address(blocker_root); + assert(blocker_sub); + this->add_blocker(new Simplex_handle(*blocker_sub)); + } + + public: + /** + * Constructs the restricted complex of 'parent_complex' to + * vertices of 'simplex'. + */ + void make_restricted_complex(const ComplexType & parent_complex, + const Simplex_handle& simplex) { + this->clear(); + // add vertices to the sub complex + for (auto x : simplex) { + assert(parent_complex.contains_vertex(x)); + auto x_local = this->add_vertex(parent_complex[x].get_id()); + (*this)[x_local] = parent_complex[x]; + } + + // add edges to the sub complex + for (auto x : simplex) { + // x_neigh is the neighbor of x intersected with vertices_simplex + Simplex_handle x_neigh; + parent_complex.add_neighbours(x, x_neigh, true); + x_neigh.intersection(simplex); + for (auto y : x_neigh) { + this->add_edge(parent_complex[x].get_id(), parent_complex[y].get_id()); + } + } + + // add blockers to the sub complex + for (auto blocker : parent_complex.const_blocker_range()) { + // check if it is the first time we encounter the blocker + if (simplex.contains(*blocker)) { + Root_simplex_handle blocker_root(parent_complex.get_id(*(blocker))); + Simplex_handle blocker_restr( + *(this->get_simplex_address(blocker_root))); + this->add_blocker(new Simplex_handle(blocker_restr)); + } + } + } + + void clear() { + adresses.clear(); + ComplexType::clear(); + } + + /** + * Compute the local vertex in L corresponding to the vertex global in K. + * runs in O(log n) if n = num_vertices() + */ + boost::optional get_address(Root_vertex_handle global) const { + boost::optional < Vertex_handle > res; + IdAddressMapConstIterator it = adresses.find(global); + if (it == adresses.end()) + res.reset(); + else + res = (*it).second; + return res; + } + + // /** + // * Allocates a simplex in L corresponding to the simplex s in K + // * with its local adresses and returns an AddressSimplex. + // */ + // boost::optional get_address(const Root_simplex_handle & s) const; + +// private: + /** + * same as get_address except that it will return a simplex in any case. + * The vertices that were not found are not added. + */ + // @remark should be private but problem with VS + std::vector > get_addresses( + const Root_simplex_handle & s) const { + std::vector < boost::optional > res; + for (auto i : s) { + res.push_back(get_address(i)); + } + return res; + } }; - /** * @remark remarque perte de temps a creer un nouveau simplexe a chaque fois * alors qu'on pourrait utiliser a la place de 'addresses_sigma_in_link' - * un simplex avec des valeurs spéciales ComplexDS::null_vertex par exemple + * un simplex avec des valeurs sp�ciales ComplexDS::null_vertex par exemple * pour indiquer qu'un vertex n'appartient pas au complex */ template bool proper_face_in_union( - Skeleton_blocker_sub_complex & link, - std::vector > & addresses_sigma_in_link, - int vertex_to_be_ignored) -{ - // we test that all vertices of 'addresses_sigma_in_link' but 'vertex_to_be_ignored' - // are in link1 if it is the case we construct the corresponding simplex - bool vertices_sigma_are_in_link = true; - typename ComplexType::Simplex_handle sigma_in_link; - for (int i = 0; i < addresses_sigma_in_link.size(); ++i){ - if (i != vertex_to_be_ignored){ - if (!addresses_sigma_in_link[i]){ - vertices_sigma_are_in_link = false; - break; - } - else sigma_in_link.add_vertex(*addresses_sigma_in_link[i]); - } - } - // If one of vertices of the simplex is not in the complex then it returns false - // Otherwise, it tests if the simplex is in the complex - return vertices_sigma_are_in_link && link.contains(sigma_in_link); + Skeleton_blocker_sub_complex & link, + std::vector > & addresses_sigma_in_link, + int vertex_to_be_ignored) { + // we test that all vertices of 'addresses_sigma_in_link' but 'vertex_to_be_ignored' + // are in link1 if it is the case we construct the corresponding simplex + bool vertices_sigma_are_in_link = true; + typename ComplexType::Simplex_handle sigma_in_link; + for (int i = 0; i < addresses_sigma_in_link.size(); ++i) { + if (i != vertex_to_be_ignored) { + if (!addresses_sigma_in_link[i]) { + vertices_sigma_are_in_link = false; + break; + } else { + sigma_in_link.add_vertex(*addresses_sigma_in_link[i]); + } + } + } + // If one of vertices of the simplex is not in the complex then it returns false + // Otherwise, it tests if the simplex is in the complex + return vertices_sigma_are_in_link && link.contains(sigma_in_link); } /* - template - bool - proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2) - { - typedef typename ComplexType::Vertex_handle Vertex_handle; - std::vector > addresses_sigma_in_link1 = link1.get_addresses(sigma); - std::vector > addresses_sigma_in_link2 = link2.get_addresses(sigma); - - for (int current_index = 0; current_index < addresses_sigma_in_link1.size(); ++current_index) - { - - if (!proper_face_in_union(link1, addresses_sigma_in_link1, current_index) - && !proper_face_in_union(link2, addresses_sigma_in_link2, current_index)){ - return false; - } - } - return true; - }*/ - - + template + bool + proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2) + { + typedef typename ComplexType::Vertex_handle Vertex_handle; + std::vector > addresses_sigma_in_link1 = link1.get_addresses(sigma); + std::vector > addresses_sigma_in_link2 = link2.get_addresses(sigma); + + for (int current_index = 0; current_index < addresses_sigma_in_link1.size(); ++current_index) + { + + if (!proper_face_in_union(link1, addresses_sigma_in_link1, current_index) + && !proper_face_in_union(link2, addresses_sigma_in_link2, current_index)){ + return false; + } + } + return true; + }*/ // Remark: this function should be friend in order to leave get_adresses private -// however doing so seemes currently not possible due to a visual studio bug c2668 +// however doing so seemes currently not possible due to a visual studio bug c2668 // "the compiler does not support partial ordering of template functions as specified in the C++ Standard" // http://www.serkey.com/error-c2668-ambiguous-call-to-overloaded-function-bb45ft.html template -bool -proper_faces_in_union(Skeleton_blocker_simplex & sigma, Skeleton_blocker_sub_complex & link1, Skeleton_blocker_sub_complex & link2) -{ - typedef typename ComplexType::Vertex_handle Vertex_handle; - std::vector > addresses_sigma_in_link1 = link1.get_addresses(sigma); - std::vector > addresses_sigma_in_link2 = link2.get_addresses(sigma); - - for (int current_index = 0; current_index < addresses_sigma_in_link1.size(); ++current_index) - { - - if (!proper_face_in_union(link1, addresses_sigma_in_link1, current_index) - && !proper_face_in_union(link2, addresses_sigma_in_link2, current_index)){ - return false; - } - } - return true; +bool proper_faces_in_union( + Skeleton_blocker_simplex & sigma, + Skeleton_blocker_sub_complex & link1, + Skeleton_blocker_sub_complex & link2) { + typedef typename ComplexType::Vertex_handle Vertex_handle; + std::vector < boost::optional > addresses_sigma_in_link1 = + link1.get_addresses(sigma); + std::vector < boost::optional > addresses_sigma_in_link2 = + link2.get_addresses(sigma); + + for (int current_index = 0; current_index < addresses_sigma_in_link1.size(); + ++current_index) { + if (!proper_face_in_union(link1, addresses_sigma_in_link1, current_index) + && !proper_face_in_union(link2, addresses_sigma_in_link2, + current_index)) { + return false; + } + } + return true; } -} // namespace skbl - -} // namespace GUDHI +} // namespace skbl +} // namespace Gudhi -#endif +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SKELETON_BLOCKER_SUB_COMPLEX_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_edges_iterators.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_edges_iterators.h index 918e3473..feae3d57 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_edges_iterators.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker/iterators/Skeleton_blockers_edges_iterators.h @@ -85,7 +85,7 @@ public: } Edge_handle dereference() const{ - return *(*complex)[std::make_pair(v,*current_)]; + return *(*complex)[std::make_pair(v,static_cast(*current_))]; } private: diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h index 15425384..7c0bc652 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_complex.h @@ -20,23 +20,25 @@ * along with this program. If not, see . */ -#ifndef GUDHI_SKELETON_BLOCKER_COMPLEX_H -#define GUDHI_SKELETON_BLOCKER_COMPLEX_H +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_COMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_COMPLEX_H_ +#include +#include +#include +#include +#include +#include +#include +#include #include #include #include #include -#include #include -#include -#include -#include -#include -#include -#include -#include +#include +#include #include "gudhi/Skeleton_blocker/iterators/Skeleton_blockers_iterators.h" #include "gudhi/Skeleton_blocker_link_complex.h" @@ -48,1385 +50,1291 @@ #include "gudhi/Skeleton_blocker/internal/Top_faces.h" #include "gudhi/Utils.h" - namespace Gudhi { namespace skbl { - /** *@class Skeleton_blocker_complex *@brief Abstract Simplicial Complex represented with a skeleton/blockers pair. *@ingroup 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 - < boost::setS, //edges - boost::vecS, // vertices - boost::undirectedS, - Graph_vertex, - Graph_edge - > 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. - */ - Skeleton_blocker_complex(int num_vertices_ = 0,Visitor* visitor_=NULL):visitor(visitor_){ - clear(); - for (int i=0; i Container_simplices; - typedef typename Container_simplices::iterator Simplices_iterator; - int dimension_; - std::vector simplices_; - - public: - 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_,(int)simplex->dimension()); - } - simplices_ = std::vector(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"<& 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 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 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_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< " << bl.second << ":"<<*bl.second <<"\n"; - } - return stream.str(); - } - - //@} - +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(); + } + + //@} }; - - /** * build a simplicial complex from a collection * of top faces. */ -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(); +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(); } -} // namespace skbl - -} // namespace GUDHI - +} // namespace skbl +} // namespace Gudhi -#endif +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_COMPLEX_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_geometric_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_geometric_complex.h index 4dbabd60..e2918901 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_geometric_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_geometric_complex.h @@ -1,35 +1,32 @@ - /* 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 . - */ -#ifndef SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ -#define SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ - +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ #include "gudhi/Utils.h" #include "gudhi/Skeleton_blocker_simplifiable_complex.h" #include "gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h" - - -namespace Gudhi{ +namespace Gudhi { namespace skbl { @@ -39,95 +36,91 @@ namespace skbl { * @ingroup skbl */ template -class Skeleton_blocker_geometric_complex : public Skeleton_blocker_simplifiable_complex -{ -public: - typedef typename SkeletonBlockerGeometricDS::GT GT; - - typedef Skeleton_blocker_simplifiable_complex SimplifiableSkeletonblocker; - - typedef typename SimplifiableSkeletonblocker::Vertex_handle Vertex_handle; - typedef typename SimplifiableSkeletonblocker::Root_vertex_handle Root_vertex_handle; - typedef typename SimplifiableSkeletonblocker::Edge_handle Edge_handle; - typedef typename SimplifiableSkeletonblocker::Simplex_handle Simplex_handle; - - - typedef typename SimplifiableSkeletonblocker::Graph_vertex Graph_vertex; - - typedef typename SkeletonBlockerGeometricDS::Point Point; - - - /** - * @brief Add a vertex to the complex with its associated point. - */ - Vertex_handle add_vertex(const Point& point){ - Vertex_handle ad = SimplifiableSkeletonblocker::add_vertex(); - (*this)[ad].point() = point; - return ad; - } - - - /** - * @brief Returns the Point associated to the vertex v. - */ - const Point& point(Vertex_handle v) const{ - assert(this->contains_vertex(v)); - return (*this)[v].point(); - } - - /** - * @brief Returns the Point associated to the vertex v. - */ - Point& point(Vertex_handle v) { - assert(this->contains_vertex(v)); - return (*this)[v].point(); - } - - const Point& point(Root_vertex_handle global_v) const{ - Vertex_handle local_v ( (*this)[global_v]) ; - assert(this->contains_vertex(local_v)); - return (*this)[local_v].point(); - } - - Point& point(Root_vertex_handle global_v) { - Vertex_handle local_v ( (*this)[global_v]) ; - assert(this->contains_vertex(local_v)); - return (*this)[local_v].point(); - } - - - typedef Skeleton_blocker_link_complex Geometric_link; - - /** - * Constructs the link of 'simplex' with points coordinates. - */ - Geometric_link link(const Simplex_handle& simplex) const{ - Geometric_link link(*this,simplex); - //we now add the point info - add_points_to_link(link); - return link; - } - - /** - * Constructs the link of 'simplex' with points coordinates. - */ - Geometric_link link(Edge_handle edge) const{ - Geometric_link link(*this,edge); - //we now add the point info - add_points_to_link(link); - return link; - } -private: - void add_points_to_link(Geometric_link& link) const{ - for(Vertex_handle v : link.vertex_range()){ - Root_vertex_handle v_root(link.get_id(v)); - link.point(v) = (*this).point(v_root); - } - } +class Skeleton_blocker_geometric_complex : + public Skeleton_blocker_simplifiable_complex { + public: + typedef typename SkeletonBlockerGeometricDS::GT GT; + + typedef Skeleton_blocker_simplifiable_complex SimplifiableSkeletonblocker; + + typedef typename SimplifiableSkeletonblocker::Vertex_handle Vertex_handle; + typedef typename SimplifiableSkeletonblocker::Root_vertex_handle Root_vertex_handle; + typedef typename SimplifiableSkeletonblocker::Edge_handle Edge_handle; + typedef typename SimplifiableSkeletonblocker::Simplex_handle Simplex_handle; + + typedef typename SimplifiableSkeletonblocker::Graph_vertex Graph_vertex; + + typedef typename SkeletonBlockerGeometricDS::Point Point; + + /** + * @brief Add a vertex to the complex with its associated point. + */ + Vertex_handle add_vertex(const Point& point) { + Vertex_handle ad = SimplifiableSkeletonblocker::add_vertex(); + (*this)[ad].point() = point; + return ad; + } + + /** + * @brief Returns the Point associated to the vertex v. + */ + const Point& point(Vertex_handle v) const { + assert(this->contains_vertex(v)); + return (*this)[v].point(); + } + + /** + * @brief Returns the Point associated to the vertex v. + */ + Point& point(Vertex_handle v) { + assert(this->contains_vertex(v)); + return (*this)[v].point(); + } + + const Point& point(Root_vertex_handle global_v) const { + Vertex_handle local_v((*this)[global_v]); + assert(this->contains_vertex(local_v)); + return (*this)[local_v].point(); + } + + Point& point(Root_vertex_handle global_v) { + Vertex_handle local_v((*this)[global_v]); + assert(this->contains_vertex(local_v)); + return (*this)[local_v].point(); + } + + typedef Skeleton_blocker_link_complex Geometric_link; + + /** + * Constructs the link of 'simplex' with points coordinates. + */ + Geometric_link link(const Simplex_handle& simplex) const { + Geometric_link link(*this, simplex); + // we now add the point info + add_points_to_link(link); + return link; + } + + /** + * Constructs the link of 'simplex' with points coordinates. + */ + Geometric_link link(Edge_handle edge) const { + Geometric_link link(*this, edge); + // we now add the point info + add_points_to_link(link); + return link; + } + private: + void add_points_to_link(Geometric_link& link) const { + for (Vertex_handle v : link.vertex_range()) { + Root_vertex_handle v_root(link.get_id(v)); + link.point(v) = (*this).point(v_root); + } + } }; -} +} // namespace skbl -} // namespace GUDHI +} // namespace Gudhi -#endif /* SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_GEOMETRIC_COMPLEX_H_ diff --git a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_link_complex.h b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_link_complex.h index f3ebfa60..725ecce5 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_link_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_link_complex.h @@ -1,37 +1,36 @@ - /* 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 . - */ -#ifndef GUDHI_SKELETON_BLOCKERS_LINK_COMPLEX_H_ -#define GUDHI_SKELETON_BLOCKERS_LINK_COMPLEX_H_ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_LINK_COMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_LINK_COMPLEX_H_ #include "gudhi/Utils.h" #include "gudhi/Skeleton_blocker_complex.h" -namespace Gudhi{ +namespace Gudhi { namespace skbl { template class Skeleton_blocker_sub_complex; - /** * \brief Class representing the link of a simplicial complex encoded by a skeleton/blockers pair. * It inherits from Skeleton_blocker_sub_complex because such complex is a sub complex of a @@ -39,252 +38,263 @@ template class Skeleton_blocker_sub_complex; * \ingroup skbl */ template -class Skeleton_blocker_link_complex : public Skeleton_blocker_sub_complex -{ - template friend class Skeleton_blocker_link_superior; - typedef typename ComplexType::Edge_handle Edge_handle; - - typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; - - -private: - - bool only_superior_vertices_; - -public: - typedef typename ComplexType::Vertex_handle Vertex_handle; - typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; - - typedef typename ComplexType::Simplex_handle Simplex_handle; - typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; - - typedef typename ComplexType::Blocker_handle Blocker_handle; - - - typedef typename ComplexType::Simplex_handle::Simplex_vertex_const_iterator Simplex_handle_iterator; - typedef typename ComplexType::Root_simplex_handle::Simplex_vertex_const_iterator Root_simplex_handle_iterator; - - - Skeleton_blocker_link_complex(bool only_superior_vertices=false):only_superior_vertices_(only_superior_vertices){ - } - - /** - * If the parameter only_superior_vertices is true, - * only vertices greater than the one of alpha are added. - */ - Skeleton_blocker_link_complex(const ComplexType & parent_complex,const Simplex_handle& alpha_parent_adress,bool only_superior_vertices = false) - :only_superior_vertices_(only_superior_vertices) { - if(!alpha_parent_adress.empty()) - build_link(parent_complex,alpha_parent_adress); - - } - - /** - * If the parameter only_superior_vertices is true, - * only vertices greater than the one of the vertex are added. - */ - Skeleton_blocker_link_complex(const ComplexType & parent_complex, Vertex_handle a_parent_adress, bool only_superior_vertices = false) - :only_superior_vertices_(only_superior_vertices){ - Simplex_handle alpha_simplex(a_parent_adress); - build_link(parent_complex,alpha_simplex); - } - - - /** - * If the parameter only_superior_vertices is true, - * only vertices greater than the one of the edge are added. - */ - Skeleton_blocker_link_complex(const ComplexType & parent_complex, Edge_handle edge, bool only_superior_vertices = false) - :only_superior_vertices_(only_superior_vertices){ - Simplex_handle alpha_simplex(parent_complex.first_vertex(edge),parent_complex.second_vertex(edge)); - build_link(parent_complex,alpha_simplex); - } - -protected: - - - /** - * @brief compute vertices of the link. - * If the boolean only_superior_vertices is true, then only the vertices - * are greater than vertices of alpha_parent_adress are added. - */ - void compute_link_vertices(const ComplexType & parent_complex,const Simplex_handle& alpha_parent_adress,bool only_superior_vertices,bool is_alpha_blocker = false){ - if(alpha_parent_adress.dimension()==0) - // for a vertex we know exactly the number of vertices of the link (and the size of the corresponding vector) - // thus we call a specific function that will reserve a vector with appropriate size - this->compute_link_vertices(parent_complex,alpha_parent_adress.first_vertex(),only_superior_vertices_); - else{ - // we compute the intersection of neighbors of alpha and store it in link_vertices - Simplex_handle link_vertices_parent; - parent_complex.add_neighbours(alpha_parent_adress,link_vertices_parent,only_superior_vertices); - // For all vertex 'v' in this intersection, we go through all its adjacent blockers. - // If one blocker minus 'v' is included in alpha then the vertex is not in the link complex. - for (auto v_parent : link_vertices_parent){ - bool new_vertex = true; - for (auto beta : parent_complex.const_blocker_range(v_parent)) - { - if(!is_alpha_blocker || *beta!=alpha_parent_adress){ - new_vertex = !(alpha_parent_adress.contains_difference(*beta,v_parent)); - if(!new_vertex) break; - } - } - if (new_vertex) - this->add_vertex(parent_complex.get_id(v_parent)); - } - } - } - - - /** - * @brief compute vertices of the link. - * If the boolean only_superior_vertices is true, then only the vertices - * are greater than vertices of alpha_parent_adress are added. - */ - void compute_link_vertices(const ComplexType & parent_complex, Vertex_handle alpha_parent_adress,bool only_superior_vertices){ - // for a vertex we know exactly the number of vertices of the link (and the size of the corresponding vector - this->skeleton.m_vertices.reserve(parent_complex.degree(alpha_parent_adress)); - - // For all vertex 'v' in this intersection, we go through all its adjacent blockers. - // If one blocker minus 'v' is included in alpha then the vertex is not in the link complex. - for (auto v_parent : parent_complex.vertex_range(alpha_parent_adress)){ - if(!only_superior_vertices || v_parent.vertex> alpha_parent_adress.vertex) - this->add_vertex(parent_complex.get_id(v_parent)); - } - - } - - - void compute_link_edges(const ComplexType & parent_complex,const Simplex_handle& alpha_parent_adress,bool is_alpha_blocker = false){ - Simplex_handle_iterator y_link, x_parent , y_parent; - // ---------------------------- - // Compute edges in the link - // ------------------------- - - if(this->num_vertices()<=1) return; - - for (auto x_link = this->vertex_range().begin() ; - x_link!= this->vertex_range().end(); - ++x_link - ) - { - for (auto y_link = x_link; ++y_link != this->vertex_range().end(); ){ - Vertex_handle x_parent = *parent_complex.get_address(this->get_id(*x_link)); - Vertex_handle y_parent = *parent_complex.get_address(this->get_id(*y_link)); - if (parent_complex.contains_edge(x_parent,y_parent) ){ - // we check that there is no blocker subset of alpha passing trough x and y - bool new_edge=true; - for (auto blocker_parent : parent_complex.const_blocker_range(x_parent)) - { - if(!is_alpha_blocker || *blocker_parent!=alpha_parent_adress){ - if (blocker_parent->contains(y_parent)) - { - new_edge = !(alpha_parent_adress.contains_difference(*blocker_parent,x_parent,y_parent)); - if (!new_edge) break; - } - } - } - if (new_edge) - this->add_edge(*x_link,*y_link); - } - } - } - } - - - /** - * @brief : Given an address in the current complex, it returns the - * corresponding address in 'other_complex'. - * It assumes that other_complex have a vertex 'this.get_id(address)' - */ - boost::optional give_equivalent_vertex(const ComplexType & other_complex, - Vertex_handle address) const{ - Root_vertex_handle id((*this)[address].get_id()); - return other_complex.get_address(id); - } - - /* - * compute the blockers of the link if is_alpha_blocker is false. - * Otherwise, alpha is a blocker, and the link is computed in the complex where - * the blocker is anticollapsed. - */ - void compute_link_blockers(const ComplexType & parent_complex,const Simplex_handle& alpha_parent,bool is_alpha_blocker = false){ - - for (auto x_link : this->vertex_range()){ - - Vertex_handle x_parent = * this->give_equivalent_vertex(parent_complex,x_link); - - for (auto blocker_parent : parent_complex.const_blocker_range(x_parent)){ - if(!is_alpha_blocker || *blocker_parent!=alpha_parent){ - Simplex_handle sigma_parent(*blocker_parent); - - sigma_parent.difference(alpha_parent); - - if (sigma_parent.dimension()>=2 && sigma_parent.first_vertex() == x_parent){ - Root_simplex_handle sigma_id(parent_complex.get_id(sigma_parent)); - auto sigma_link = this->get_simplex_address(sigma_id); - // ie if the vertices of sigma are vertices of the link - if(sigma_link){ - bool is_new_blocker = true; - for(auto a : alpha_parent){ - for(auto eta_parent : parent_complex.const_blocker_range(a)){ - if(!is_alpha_blocker || *eta_parent!=alpha_parent){ - Simplex_handle eta_minus_alpha(*eta_parent); - eta_minus_alpha.difference(alpha_parent); - if(eta_minus_alpha != sigma_parent && - sigma_parent.contains_difference(*eta_parent,alpha_parent)){ - is_new_blocker = false; - break; - } - } - } - if (!is_new_blocker) break; - } - if (is_new_blocker) - this->add_blocker(new Simplex_handle(*sigma_link)); - - } - } - } - } - } - } - - -public: - /** - * @brief compute vertices, edges and blockers of the link. - * @details If the boolean only_superior_vertices is true, then the link is computed only - * with vertices that are greater than vertices of alpha_parent_adress. - */ - void build_link(const ComplexType & parent_complex,const Simplex_handle& alpha_parent_adress,bool is_alpha_blocker = false) - { - assert(is_alpha_blocker || parent_complex.contains(alpha_parent_adress)); - compute_link_vertices(parent_complex,alpha_parent_adress,only_superior_vertices_); - compute_link_edges(parent_complex,alpha_parent_adress,is_alpha_blocker); - compute_link_blockers(parent_complex,alpha_parent_adress,is_alpha_blocker); - } - - - /** - * @brief build the link of a blocker which is the link - * of the blocker's simplex if this simplex had been - * removed from the blockers of the complex. - */ - friend void build_link_of_blocker( - const ComplexType & parent_complex, - Simplex_handle& blocker, - Skeleton_blocker_link_complex & result){ - assert(blocker.dimension()>=2); - assert(parent_complex.contains_blocker(blocker)); - result.clear(); - result.build_link(parent_complex,blocker,true); - } - - +class Skeleton_blocker_link_complex : public Skeleton_blocker_sub_complex< + ComplexType> { + template friend class Skeleton_blocker_link_superior; + typedef typename ComplexType::Edge_handle Edge_handle; + + typedef typename ComplexType::boost_vertex_handle boost_vertex_handle; + + private: + bool only_superior_vertices_; + + public: + typedef typename ComplexType::Vertex_handle Vertex_handle; + typedef typename ComplexType::Root_vertex_handle Root_vertex_handle; + + typedef typename ComplexType::Simplex_handle Simplex_handle; + typedef typename ComplexType::Root_simplex_handle Root_simplex_handle; + + typedef typename ComplexType::Blocker_handle Blocker_handle; + + typedef typename ComplexType::Simplex_handle::Simplex_vertex_const_iterator Simplex_handle_iterator; + typedef typename ComplexType::Root_simplex_handle::Simplex_vertex_const_iterator Root_simplex_handle_iterator; + + explicit Skeleton_blocker_link_complex(bool only_superior_vertices = false) + : only_superior_vertices_(only_superior_vertices) { + } + + /** + * If the parameter only_superior_vertices is true, + * only vertices greater than the one of alpha are added. + */ + Skeleton_blocker_link_complex(const ComplexType & parent_complex, + const Simplex_handle& alpha_parent_adress, + bool only_superior_vertices = false) + : only_superior_vertices_(only_superior_vertices) { + if (!alpha_parent_adress.empty()) + build_link(parent_complex, alpha_parent_adress); + } + + /** + * If the parameter only_superior_vertices is true, + * only vertices greater than the one of the vertex are added. + */ + Skeleton_blocker_link_complex(const ComplexType & parent_complex, + Vertex_handle a_parent_adress, + bool only_superior_vertices = false) + : only_superior_vertices_(only_superior_vertices) { + Simplex_handle alpha_simplex(a_parent_adress); + build_link(parent_complex, alpha_simplex); + } + + /** + * If the parameter only_superior_vertices is true, + * only vertices greater than the one of the edge are added. + */ + Skeleton_blocker_link_complex(const ComplexType & parent_complex, + Edge_handle edge, bool only_superior_vertices = + false) + : only_superior_vertices_(only_superior_vertices) { + Simplex_handle alpha_simplex(parent_complex.first_vertex(edge), + parent_complex.second_vertex(edge)); + build_link(parent_complex, alpha_simplex); + } + + protected: + /** + * @brief compute vertices of the link. + * If the boolean only_superior_vertices is true, then only the vertices + * are greater than vertices of alpha_parent_adress are added. + */ + void compute_link_vertices(const ComplexType & parent_complex, + const Simplex_handle& alpha_parent_adress, + bool only_superior_vertices, + bool is_alpha_blocker = false) { + if (alpha_parent_adress.dimension() == 0) { + // for a vertex we know exactly the number of vertices of the link (and the size of the corresponding vector) + // thus we call a specific function that will reserve a vector with appropriate size + this->compute_link_vertices(parent_complex, + alpha_parent_adress.first_vertex(), + only_superior_vertices_); + } else { + // we compute the intersection of neighbors of alpha and store it in link_vertices + Simplex_handle link_vertices_parent; + parent_complex.add_neighbours(alpha_parent_adress, link_vertices_parent, + only_superior_vertices); + // For all vertex 'v' in this intersection, we go through all its adjacent blockers. + // If one blocker minus 'v' is included in alpha then the vertex is not in the link complex. + for (auto v_parent : link_vertices_parent) { + bool new_vertex = true; + for (auto beta : parent_complex.const_blocker_range(v_parent)) { + if (!is_alpha_blocker || *beta != alpha_parent_adress) { + new_vertex = !(alpha_parent_adress.contains_difference(*beta, + v_parent)); + if (!new_vertex) + break; + } + } + if (new_vertex) + this->add_vertex(parent_complex.get_id(v_parent)); + } + } + } + + /** + * @brief compute vertices of the link. + * If the boolean only_superior_vertices is true, then only the vertices + * are greater than vertices of alpha_parent_adress are added. + */ + void compute_link_vertices(const ComplexType & parent_complex, + Vertex_handle alpha_parent_adress, + bool only_superior_vertices) { + // for a vertex we know exactly the number of vertices of the link (and the size of the corresponding vector + this->skeleton.m_vertices.reserve( + parent_complex.degree(alpha_parent_adress)); + + // For all vertex 'v' in this intersection, we go through all its adjacent blockers. + // If one blocker minus 'v' is included in alpha then the vertex is not in the link complex. + for (auto v_parent : parent_complex.vertex_range(alpha_parent_adress)) { + if (!only_superior_vertices + || v_parent.vertex > alpha_parent_adress.vertex) + this->add_vertex(parent_complex.get_id(v_parent)); + } + } + + void compute_link_edges(const ComplexType & parent_complex, + const Simplex_handle& alpha_parent_adress, + bool is_alpha_blocker = false) { + Simplex_handle_iterator y_link, x_parent, y_parent; + // ---------------------------- + // Compute edges in the link + // ------------------------- + + if (this->num_vertices() <= 1) + return; + + for (auto x_link = this->vertex_range().begin(); + x_link != this->vertex_range().end(); ++x_link) { + for (auto y_link = x_link; ++y_link != this->vertex_range().end();) { + Vertex_handle x_parent = *parent_complex.get_address( + this->get_id(*x_link)); + Vertex_handle y_parent = *parent_complex.get_address( + this->get_id(*y_link)); + if (parent_complex.contains_edge(x_parent, y_parent)) { + // we check that there is no blocker subset of alpha passing trough x and y + bool new_edge = true; + for (auto blocker_parent : parent_complex.const_blocker_range( + x_parent)) { + if (!is_alpha_blocker || *blocker_parent != alpha_parent_adress) { + if (blocker_parent->contains(y_parent)) { + new_edge = !(alpha_parent_adress.contains_difference( + *blocker_parent, x_parent, y_parent)); + if (!new_edge) + break; + } + } + } + if (new_edge) + this->add_edge(*x_link, *y_link); + } + } + } + } + + /** + * @brief : Given an address in the current complex, it returns the + * corresponding address in 'other_complex'. + * It assumes that other_complex have a vertex 'this.get_id(address)' + */ + boost::optional give_equivalent_vertex( + const ComplexType & other_complex, Vertex_handle address) const { + Root_vertex_handle id((*this)[address].get_id()); + return other_complex.get_address(id); + } + + /* + * compute the blockers of the link if is_alpha_blocker is false. + * Otherwise, alpha is a blocker, and the link is computed in the complex where + * the blocker is anticollapsed. + */ + void compute_link_blockers(const ComplexType & parent_complex, + const Simplex_handle& alpha_parent, + bool is_alpha_blocker = false) { + for (auto x_link : this->vertex_range()) { + Vertex_handle x_parent = *this->give_equivalent_vertex(parent_complex, + x_link); + + for (auto blocker_parent : parent_complex.const_blocker_range(x_parent)) { + if (!is_alpha_blocker || *blocker_parent != alpha_parent) { + Simplex_handle sigma_parent(*blocker_parent); + + sigma_parent.difference(alpha_parent); + + if (sigma_parent.dimension() >= 2 + && sigma_parent.first_vertex() == x_parent) { + Root_simplex_handle sigma_id(parent_complex.get_id(sigma_parent)); + auto sigma_link = this->get_simplex_address(sigma_id); + // ie if the vertices of sigma are vertices of the link + if (sigma_link) { + bool is_new_blocker = true; + for (auto a : alpha_parent) { + for (auto eta_parent : parent_complex.const_blocker_range(a)) { + if (!is_alpha_blocker || *eta_parent != alpha_parent) { + Simplex_handle eta_minus_alpha(*eta_parent); + eta_minus_alpha.difference(alpha_parent); + if (eta_minus_alpha != sigma_parent + && sigma_parent.contains_difference(*eta_parent, + alpha_parent)) { + is_new_blocker = false; + break; + } + } + } + if (!is_new_blocker) + break; + } + if (is_new_blocker) + this->add_blocker(new Simplex_handle(*sigma_link)); + } + } + } + } + } + } + + public: + /** + * @brief compute vertices, edges and blockers of the link. + * @details If the boolean only_superior_vertices is true, then the link is computed only + * with vertices that are greater than vertices of alpha_parent_adress. + */ + void build_link(const ComplexType & parent_complex, + const Simplex_handle& alpha_parent_adress, + bool is_alpha_blocker = false) { + assert(is_alpha_blocker || parent_complex.contains(alpha_parent_adress)); + compute_link_vertices(parent_complex, alpha_parent_adress, + only_superior_vertices_); + compute_link_edges(parent_complex, alpha_parent_adress, is_alpha_blocker); + compute_link_blockers(parent_complex, alpha_parent_adress, + is_alpha_blocker); + } + + /** + * @brief build the link of a blocker which is the link + * of the blocker's simplex if this simplex had been + * removed from the blockers of the complex. + */ + friend void build_link_of_blocker(const ComplexType & parent_complex, + Simplex_handle& blocker, + Skeleton_blocker_link_complex & result) { + assert(blocker.dimension() >= 2); + assert(parent_complex.contains_blocker(blocker)); + result.clear(); + result.build_link(parent_complex, blocker, true); + } }; -} +} // namespace skbl -} // namespace GUDHI +} // namespace Gudhi -#endif /* SKELETON_BLOCKERS_LINK_COMPLEX_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_LINK_COMPLEX_H_ 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 f70e585d..305b8829 100644 --- a/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h +++ b/src/Skeleton_blocker/include/gudhi/Skeleton_blocker_simplifiable_complex.h @@ -1,30 +1,34 @@ - /* 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 . - */ -#ifndef GUDHI_SKELETON_BLOCKERS_SIMPLIFIABLE_COMPLEX_H_ -#define GUDHI_SKELETON_BLOCKERS_SIMPLIFIABLE_COMPLEX_H_ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): 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 . + */ +#ifndef SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SIMPLIFIABLE_COMPLEX_H_ +#define SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SIMPLIFIABLE_COMPLEX_H_ #include "gudhi/Skeleton_blocker/Skeleton_blocker_sub_complex.h" -namespace Gudhi{ +#include +#include +#include + +namespace Gudhi { namespace skbl { @@ -35,535 +39,517 @@ namespace skbl { * @extends Skeleton_blocker_complex */ template -class Skeleton_blocker_simplifiable_complex : public Skeleton_blocker_complex -{ - template friend class Skeleton_blocker_sub_complex; - -public: - - typedef Skeleton_blocker_complex SkeletonBlockerComplex; - - typedef typename SkeletonBlockerComplex::Graph_edge Graph_edge; - - typedef typename SkeletonBlockerComplex::boost_adjacency_iterator boost_adjacency_iterator; - typedef typename SkeletonBlockerComplex::Edge_handle Edge_handle; - typedef typename SkeletonBlockerComplex::boost_vertex_handle boost_vertex_handle; - typedef typename SkeletonBlockerComplex::Vertex_handle Vertex_handle; - typedef typename SkeletonBlockerComplex::Root_vertex_handle Root_vertex_handle; - typedef typename SkeletonBlockerComplex::Simplex_handle Simplex_handle; - typedef typename SkeletonBlockerComplex::Root_simplex_handle Root_simplex_handle; - typedef typename SkeletonBlockerComplex::Blocker_handle Blocker_handle; - - - typedef typename SkeletonBlockerComplex::Root_simplex_iterator Root_simplex_iterator; - typedef typename SkeletonBlockerComplex::Simplex_handle_iterator Simplex_handle_iterator; - typedef typename SkeletonBlockerComplex::BlockerMap BlockerMap; - typedef typename SkeletonBlockerComplex::BlockerPair BlockerPair; - typedef typename SkeletonBlockerComplex::BlockerMapIterator BlockerMapIterator; - typedef typename SkeletonBlockerComplex::BlockerMapConstIterator BlockerMapConstIterator; - - typedef typename SkeletonBlockerComplex::Visitor Visitor; - - - /** @name Constructors / Destructors / Initialization - */ - //@{ - Skeleton_blocker_simplifiable_complex(int num_vertices_ = 0,Visitor* visitor_=NULL): - Skeleton_blocker_complex(num_vertices_,visitor_){ } - - - /** - * @brief Constructor with a list of simplices - * @details The list of simplices must be the list - * of simplices of a simplicial complex, sorted with increasing dimension. - * todo take iterator instead - */ - Skeleton_blocker_simplifiable_complex(std::list& simplices,Visitor* visitor_=NULL): - Skeleton_blocker_complex(simplices,visitor_) - {} - - - - virtual ~Skeleton_blocker_simplifiable_complex(){ - } - - //@} - - /** - * Returns true iff the blocker 'sigma' is popable. - * To define popable, let us call 'L' the complex that - * consists in the current complex without the blocker 'sigma'. - * A blocker 'sigma' is then "popable" if the link of 'sigma' - * in L is reducible. - * - */ - virtual bool is_popable_blocker(Blocker_handle sigma) const{ - assert(this->contains_blocker(*sigma)); - Skeleton_blocker_link_complex link_blocker_sigma; - build_link_of_blocker(*this,*sigma,link_blocker_sigma); - - bool res = link_blocker_sigma.is_contractible()==CONTRACTIBLE; - return res; - } - - - - -private: - /** - * @returns the list of blockers of the simplex - * - * @todo a enlever et faire un iterateur sur tous les blockers a la place - */ - std::list get_blockers(){ - std::list res; - for (auto blocker : this->blocker_range()){ - res.push_back(blocker); - } - return res; - } - - - -public: - - /** - * Removes all the popable blockers of the complex and delete them. - * @returns the number of popable blockers deleted - */ - void remove_popable_blockers(){ - std::list vertex_to_check; - for(auto v : this->vertex_range()) - vertex_to_check.push_front(v); - - while(!vertex_to_check.empty()){ - Vertex_handle v = vertex_to_check.front(); - vertex_to_check.pop_front(); - - bool blocker_popable_found=true; - while (blocker_popable_found){ - blocker_popable_found = false; - for(auto block : this->blocker_range(v)){ - if (this->is_popable_blocker(block)) { - for(Vertex_handle nv : *block) - if(nv!=v) vertex_to_check.push_back(nv); - this->delete_blocker(block); - blocker_popable_found = true; - break; - } - } - } - } - } - - - /** - * Removes all the popable blockers of the complex passing through v and delete them. - */ - void remove_popable_blockers(Vertex_handle v){ - bool blocker_popable_found=true; - while (blocker_popable_found){ - blocker_popable_found = false; - for(auto block : this->blocker_range(v)){ - if (is_popable_blocker(block)) { - this->delete_blocker(block); - blocker_popable_found = true; - } - } - } - } - - - /** - * Remove the star of the vertex 'v' - */ - void remove_star(Vertex_handle v){ - // we remove the blockers that are not consistent anymore - - update_blockers_after_remove_star_of_vertex_or_edge(v); - - while (this->degree(v) > 0) - { - Vertex_handle w( * (adjacent_vertices(v.vertex, this->skeleton).first)); - this->remove_edge(v,w); - } - this->remove_vertex(v); - } - -private: - /** - * after removing the star of a simplex, blockers sigma that contains this simplex must be removed. - * Furthermore, all simplices tau of the form sigma \setminus simplex_to_be_removed must be added - * whenever the dimension of tau is at least 2. - */ - void update_blockers_after_remove_star_of_vertex_or_edge(const Simplex_handle& simplex_to_be_removed){ - std::list blockers_to_update; - if(simplex_to_be_removed.empty()) return; - - auto v0 = simplex_to_be_removed.first_vertex(); - for (auto blocker : this->blocker_range(v0)){ - if(blocker->contains(simplex_to_be_removed)) - blockers_to_update.push_back(blocker); - } - - for(auto blocker_to_update : blockers_to_update){ - Simplex_handle sub_blocker_to_be_added; - bool sub_blocker_need_to_be_added = - (blocker_to_update->dimension()-simplex_to_be_removed.dimension()) >= 2; - if(sub_blocker_need_to_be_added){ - sub_blocker_to_be_added = *blocker_to_update; - sub_blocker_to_be_added.difference(simplex_to_be_removed); - } - this->delete_blocker(blocker_to_update); - if(sub_blocker_need_to_be_added) - this->add_blocker(sub_blocker_to_be_added); - } - } - - - -public: - /** - * Remove the star of the edge connecting vertices a and b. - * @returns the number of blocker that have been removed - */ - void remove_star(Vertex_handle a, Vertex_handle b){ - update_blockers_after_remove_star_of_vertex_or_edge(Simplex_handle(a,b)); - // we remove the edge - this->remove_edge(a,b); - } - - - /** - * Remove the star of the edge 'e'. - */ - void remove_star(Edge_handle e){ - return remove_star(this->first_vertex(e),this->second_vertex(e)); - } - - /** - * Remove the star of the simplex 'sigma' which needs to belong to the complex - */ - void remove_star(const Simplex_handle& sigma){ - assert(this->contains(sigma)); - if (sigma.dimension()==0) - remove_star(sigma.first_vertex()); - else - if (sigma.dimension()==1) - remove_star(sigma.first_vertex(),sigma.last_vertex()); - else{ - remove_blocker_containing_simplex(sigma); - this->add_blocker(sigma); - } - } - - /** - * @brief add a maximal simplex plus all its cofaces. - * @details the simplex must have dimension greater than one (otherwise use add_vertex or add_edge). - */ - void add_simplex(const Simplex_handle& sigma){ - assert(!this->contains(sigma)); - assert(sigma.dimension()>1); - remove_blocker_include_in_simplex(sigma); - } - -private: - /** - * remove all blockers that contains sigma - */ - void remove_blocker_containing_simplex(const Simplex_handle& sigma){ - std::vector blockers_to_remove; - for (auto blocker : this->blocker_range(sigma.first_vertex())){ - if(blocker->contains(sigma)) - blockers_to_remove.push_back(blocker); - } - for(auto blocker_to_update : blockers_to_remove) - this->delete_blocker(blocker_to_update); - } - - /** - * remove all blockers that contains sigma - */ - void remove_blocker_include_in_simplex(const Simplex_handle& sigma){ - std::vector blockers_to_remove; - for (auto blocker : this->blocker_range(sigma.first_vertex())){ - if(sigma.contains(*blocker)) - blockers_to_remove.push_back(blocker); - } - for(auto blocker_to_update : blockers_to_remove) - this->delete_blocker(blocker_to_update); - } - - -public: - - enum simplifiable_status{ NOT_HOMOTOPY_EQ,MAYBE_HOMOTOPY_EQ,HOMOTOPY_EQ}; - simplifiable_status is_remove_star_homotopy_preserving(const Simplex_handle& simplex){ - return MAYBE_HOMOTOPY_EQ; - } - - - - enum contractible_status{ NOT_CONTRACTIBLE,MAYBE_CONTRACTIBLE,CONTRACTIBLE}; - /** - * @brief %Test if the complex is reducible using a strategy defined in the class - * (by default it tests if the complex is a cone) - * @details Note that NO could be returned if some invariant ensures that the complex - * is not a point (for instance if the euler characteristic is different from 1). - * This function will surely have to return MAYBE in some case because the - * associated problem is undecidable but it in practice, it can often - * be solved with the help of geometry. - */ - virtual contractible_status is_contractible() const{ - if (this->is_cone()) return CONTRACTIBLE; - else return MAYBE_CONTRACTIBLE; - // return this->is_cone(); - } - - - /** @Edge contraction operations - */ - //@{ - - - - /** - * @return If ignore_popable_blockers is true - * then the result is true iff the link condition at edge ab is satisfied - * or equivalently iff no blocker contains ab. - * If ignore_popable_blockers is false then the - * result is true iff all blocker containing ab are popable. - */ - bool link_condition(Vertex_handle a, Vertex_handle b,bool ignore_popable_blockers = false) const{ - for (auto blocker : this->const_blocker_range(a)) - if ( blocker->contains(b) ){ - // false if ignore_popable_blockers is false - // otherwise the blocker has to be popable - return ignore_popable_blockers && is_popable_blocker(blocker); - } - return true; - } - - /** - * @return If ignore_popable_blockers is true - * then the result is true iff the link condition at edge ab is satisfied - * or equivalently iff no blocker contains ab. - * If ignore_popable_blockers is false then the - * result is true iff all blocker containing ab are popable. - */ - bool link_condition(Edge_handle & e,bool ignore_popable_blockers = false) const{ - const Graph_edge& edge = (*this)[e]; - assert(this->get_address(edge.first())); - assert(this->get_address(edge.second())); - Vertex_handle a(*this->get_address(edge.first())); - Vertex_handle b(*this->get_address(edge.second())); - return link_condition(a,b,ignore_popable_blockers); - } - - -protected: - /** - * Compute simplices beta such that a.beta is an order 0 blocker - * that may be used to construct a new blocker after contracting ab. - * Suppose that the link condition Link(ab) = Link(a) inter Link(b) - * is satisfied. - */ - void tip_blockers(Vertex_handle a, Vertex_handle b, std::vector & buffer) const{ - for (auto const & blocker : this->const_blocker_range(a)) - { - Simplex_handle beta = (*blocker); - beta.remove_vertex(a); - buffer.push_back(beta); - } - - Simplex_handle n; - this->add_neighbours(b,n); - this->remove_neighbours(a,n); - n.remove_vertex(a); - - - for (Vertex_handle y : n) - { - Simplex_handle beta; - beta.add_vertex( y ); - buffer.push_back(beta); - } - } - - -private: - - /** - * @brief "Replace" the edge 'bx' by the edge 'ax'. - * Assume that the edge 'bx' was present whereas 'ax' was not. - * Precisely, it does not replace edges, but remove 'bx' and then add 'ax'. - * The visitor 'on_swaped_edge' is called just after edge 'ax' had been added - * and just before edge 'bx' had been removed. That way, it can - * eventually access to information of 'ax'. - */ - void swap_edge(Vertex_handle a,Vertex_handle b,Vertex_handle x){ - this->add_edge(a,x); - if (this->visitor) this->visitor->on_swaped_edge(a,b,x); - this->remove_edge(b,x); - } - - -private: - /** - * @brief removes all blockers passing through the edge 'ab' - */ - void delete_blockers_around_vertex(Vertex_handle v){ - std::list blockers_to_delete; - for (auto blocker : this->blocker_range(v)){ - blockers_to_delete.push_back(blocker); - } - while (!blockers_to_delete.empty()){ - this->remove_blocker(blockers_to_delete.back()); - blockers_to_delete.pop_back(); - } - - } - /** - * @brief removes all blockers passing through the edge 'ab' - */ - void delete_blockers_around_edge(Vertex_handle a, Vertex_handle b){ - std::list blocker_to_delete; - for (auto blocker : this->blocker_range(a)) - if (blocker->contains(b)) blocker_to_delete.push_back(blocker); - while (!blocker_to_delete.empty()) - { - this->delete_blocker(blocker_to_delete.back()); - blocker_to_delete.pop_back(); - } - } - -public: - - /** - * Contracts the edge. - * @remark If the link condition Link(ab) = Link(a) inter Link(b) is not satisfied, - * it removes first all blockers passing through 'ab' - */ - void contract_edge(Edge_handle edge){ - contract_edge(this->first_vertex(edge),this->second_vertex(edge)); - } - - - /** - * Contracts the edge connecting vertices a and b. - * @remark If the link condition Link(ab) = Link(a) inter Link(b) is not satisfied, - * it removes first all blockers passing through 'ab' - */ - void contract_edge(Vertex_handle a, Vertex_handle b){ - assert(this->contains_vertex(a)); - assert(this->contains_vertex(b)); - assert(this->contains_edge(a,b)); - - // if some blockers passes through 'ab', we remove them. - if (!link_condition(a,b)) - delete_blockers_around_edge(a,b); - - std::set blockers_to_add; - - get_blockers_to_be_added_after_contraction(a,b,blockers_to_add); - - delete_blockers_around_vertices(a,b); - - update_edges_after_contraction(a,b); - - this->remove_vertex(b); - - notify_changed_edges(a); - - for(auto block : blockers_to_add) - this->add_blocker(block); - - assert(this->contains_vertex(a)); - assert(!this->contains_vertex(b)); - } - - -private: - - void get_blockers_to_be_added_after_contraction(Vertex_handle a,Vertex_handle b,std::set& blockers_to_add){ - blockers_to_add.clear(); - - typedef Skeleton_blocker_link_complex > LinkComplexType; - - LinkComplexType link_a(*this,a); - LinkComplexType link_b(*this,b); - - std::vector vector_alpha, vector_beta; - - tip_blockers(a,b,vector_alpha); - tip_blockers(b,a,vector_beta); - - for (auto alpha = vector_alpha.begin(); alpha != vector_alpha.end(); ++alpha){ - for (auto beta = vector_beta.begin(); beta != vector_beta.end(); ++beta) - { - Simplex_handle sigma = *alpha; sigma.union_vertices(*beta); - Root_simplex_handle sigma_id = this->get_id(sigma); - if ( this->contains(sigma) && - proper_faces_in_union(sigma_id,link_a,link_b)) - { - // Blocker_handle blocker = new Simplex_handle(sigma); - sigma.add_vertex(a); - blockers_to_add.insert(sigma); - } - } - } - } - - - /** - * delete all blockers that passes through a or b - */ - void delete_blockers_around_vertices(Vertex_handle a,Vertex_handle b){ - std::vector blocker_to_delete; - for(auto bl : this->blocker_range(a)) - blocker_to_delete.push_back(bl); - for(auto bl : this->blocker_range(b)) - blocker_to_delete.push_back(bl); - while (!blocker_to_delete.empty()) - { - this->delete_blocker(blocker_to_delete.back()); - blocker_to_delete.pop_back(); - } - } - - - void update_edges_after_contraction(Vertex_handle a,Vertex_handle b){ - // We update the set of edges - this->remove_edge(a,b); - - // For all edges {b,x} incident to b, - // we remove {b,x} and add {a,x} if not already there. - while (this->degree(b)> 0) - { - Vertex_handle x(*(adjacent_vertices(b.vertex, this->skeleton).first)); - if(!this->contains_edge(a,x)) - // we 'replace' the edge 'bx' by the edge 'ax' - this->swap_edge(a,b,x); - else - this->remove_edge(b,x); - } - } - - void notify_changed_edges(Vertex_handle a){ - // We notify the visitor that all edges incident to 'a' had changed - boost_adjacency_iterator v, v_end; - - for (tie(v, v_end) = adjacent_vertices(a.vertex, this->skeleton); v != v_end; ++v) - if (this->visitor) this->visitor->on_changed_edge(a,Vertex_handle(*v)); - } - - //@} - - +class Skeleton_blocker_simplifiable_complex : public Skeleton_blocker_complex { + template friend class Skeleton_blocker_sub_complex; + + public: + typedef Skeleton_blocker_complex SkeletonBlockerComplex; + + typedef typename SkeletonBlockerComplex::Graph_edge Graph_edge; + + typedef typename SkeletonBlockerComplex::boost_adjacency_iterator boost_adjacency_iterator; + typedef typename SkeletonBlockerComplex::Edge_handle Edge_handle; + typedef typename SkeletonBlockerComplex::boost_vertex_handle boost_vertex_handle; + typedef typename SkeletonBlockerComplex::Vertex_handle Vertex_handle; + typedef typename SkeletonBlockerComplex::Root_vertex_handle Root_vertex_handle; + typedef typename SkeletonBlockerComplex::Simplex_handle Simplex_handle; + typedef typename SkeletonBlockerComplex::Root_simplex_handle Root_simplex_handle; + typedef typename SkeletonBlockerComplex::Blocker_handle Blocker_handle; + + typedef typename SkeletonBlockerComplex::Root_simplex_iterator Root_simplex_iterator; + typedef typename SkeletonBlockerComplex::Simplex_handle_iterator Simplex_handle_iterator; + typedef typename SkeletonBlockerComplex::BlockerMap BlockerMap; + typedef typename SkeletonBlockerComplex::BlockerPair BlockerPair; + typedef typename SkeletonBlockerComplex::BlockerMapIterator BlockerMapIterator; + typedef typename SkeletonBlockerComplex::BlockerMapConstIterator BlockerMapConstIterator; + + typedef typename SkeletonBlockerComplex::Visitor Visitor; + + /** @name Constructors / Destructors / Initialization + */ + //@{ + Skeleton_blocker_simplifiable_complex(int num_vertices_ = 0, + Visitor* visitor_ = NULL) + : Skeleton_blocker_complex(num_vertices_, visitor_) { + } + + /** + * @brief Constructor with a list of simplices + * @details The list of simplices must be the list + * of simplices of a simplicial complex, sorted with increasing dimension. + * todo take iterator instead + */ + Skeleton_blocker_simplifiable_complex(std::list& simplices, + Visitor* visitor_ = NULL) + : Skeleton_blocker_complex(simplices, visitor_) { + } + + virtual ~Skeleton_blocker_simplifiable_complex() { + } + + //@} + + /** + * Returns true iff the blocker 'sigma' is popable. + * To define popable, let us call 'L' the complex that + * consists in the current complex without the blocker 'sigma'. + * A blocker 'sigma' is then "popable" if the link of 'sigma' + * in L is reducible. + * + */ + virtual bool is_popable_blocker(Blocker_handle sigma) const { + assert(this->contains_blocker(*sigma)); + Skeleton_blocker_link_complex link_blocker_sigma; + build_link_of_blocker(*this, *sigma, link_blocker_sigma); + + bool res = link_blocker_sigma.is_contractible() == CONTRACTIBLE; + return res; + } + + private: + /** + * @returns the list of blockers of the simplex + * + * @todo a enlever et faire un iterateur sur tous les blockers a la place + */ + std::list get_blockers() { + std::list < Blocker_handle > res; + for (auto blocker : this->blocker_range()) { + res.push_back(blocker); + } + return res; + } + + public: + /** + * Removes all the popable blockers of the complex and delete them. + * @returns the number of popable blockers deleted + */ + void remove_popable_blockers() { + std::list < Vertex_handle > vertex_to_check; + for (auto v : this->vertex_range()) + vertex_to_check.push_front(v); + + while (!vertex_to_check.empty()) { + Vertex_handle v = vertex_to_check.front(); + vertex_to_check.pop_front(); + + bool blocker_popable_found = true; + while (blocker_popable_found) { + blocker_popable_found = false; + for (auto block : this->blocker_range(v)) { + if (this->is_popable_blocker(block)) { + for (Vertex_handle nv : *block) + if (nv != v) + vertex_to_check.push_back(nv); + this->delete_blocker(block); + blocker_popable_found = true; + break; + } + } + } + } + } + + /** + * Removes all the popable blockers of the complex passing through v and delete them. + */ + void remove_popable_blockers(Vertex_handle v) { + bool blocker_popable_found = true; + while (blocker_popable_found) { + blocker_popable_found = false; + for (auto block : this->blocker_range(v)) { + if (is_popable_blocker(block)) { + this->delete_blocker(block); + blocker_popable_found = true; + } + } + } + } + + /** + * Remove the star of the vertex 'v' + */ + void remove_star(Vertex_handle v) { + // we remove the blockers that are not consistent anymore + + update_blockers_after_remove_star_of_vertex_or_edge(static_cast(v)); + + while (this->degree(v) > 0) { + Vertex_handle w(*(adjacent_vertices(v.vertex, this->skeleton).first)); + this->remove_edge(v, w); + } + this->remove_vertex(v); + } + + private: + /** + * after removing the star of a simplex, blockers sigma that contains this simplex must be removed. + * Furthermore, all simplices tau of the form sigma \setminus simplex_to_be_removed must be added + * whenever the dimension of tau is at least 2. + */ + void update_blockers_after_remove_star_of_vertex_or_edge( + const Simplex_handle& simplex_to_be_removed) { + std::list < Blocker_handle > blockers_to_update; + if (simplex_to_be_removed.empty()) + return; + + auto v0 = simplex_to_be_removed.first_vertex(); + for (auto blocker : this->blocker_range(v0)) { + if (blocker->contains(simplex_to_be_removed)) + blockers_to_update.push_back(blocker); + } + + for (auto blocker_to_update : blockers_to_update) { + Simplex_handle sub_blocker_to_be_added; + bool sub_blocker_need_to_be_added = (blocker_to_update->dimension() + - simplex_to_be_removed.dimension()) >= 2; + if (sub_blocker_need_to_be_added) { + sub_blocker_to_be_added = *blocker_to_update; + sub_blocker_to_be_added.difference(simplex_to_be_removed); + } + this->delete_blocker(blocker_to_update); + if (sub_blocker_need_to_be_added) + this->add_blocker(sub_blocker_to_be_added); + } + } + + public: + /** + * Remove the star of the edge connecting vertices a and b. + * @returns the number of blocker that have been removed + */ + void remove_star(Vertex_handle a, Vertex_handle b) { + update_blockers_after_remove_star_of_vertex_or_edge(Simplex_handle(a, b)); + // we remove the edge + this->remove_edge(a, b); + } + + /** + * Remove the star of the edge 'e'. + */ + void remove_star(Edge_handle e) { + return remove_star(this->first_vertex(e), this->second_vertex(e)); + } + + /** + * Remove the star of the simplex 'sigma' which needs to belong to the complex + */ + void remove_star(const Simplex_handle& sigma) { + assert(this->contains(sigma)); + if (sigma.dimension() == 0) { + remove_star(sigma.first_vertex()); + } else if (sigma.dimension() == 1) { + remove_star(sigma.first_vertex(), sigma.last_vertex()); + } else { + remove_blocker_containing_simplex(sigma); + this->add_blocker(sigma); + } + } + + /** + * @brief add a maximal simplex plus all its cofaces. + * @details the simplex must have dimension greater than one (otherwise use add_vertex or add_edge). + */ + void add_simplex(const Simplex_handle& sigma) { + assert(!this->contains(sigma)); + assert(sigma.dimension() > 1); + remove_blocker_include_in_simplex(sigma); + } + + private: + /** + * remove all blockers that contains sigma + */ + void remove_blocker_containing_simplex(const Simplex_handle& sigma) { + std::vector < Blocker_handle > blockers_to_remove; + for (auto blocker : this->blocker_range(sigma.first_vertex())) { + if (blocker->contains(sigma)) + blockers_to_remove.push_back(blocker); + } + for (auto blocker_to_update : blockers_to_remove) + this->delete_blocker(blocker_to_update); + } + + /** + * remove all blockers that contains sigma + */ + void remove_blocker_include_in_simplex(const Simplex_handle& sigma) { + std::vector < Blocker_handle > blockers_to_remove; + for (auto blocker : this->blocker_range(sigma.first_vertex())) { + if (sigma.contains(*blocker)) + blockers_to_remove.push_back(blocker); + } + for (auto blocker_to_update : blockers_to_remove) + this->delete_blocker(blocker_to_update); + } + + public: + enum simplifiable_status { + NOT_HOMOTOPY_EQ, + MAYBE_HOMOTOPY_EQ, + HOMOTOPY_EQ + }; + simplifiable_status is_remove_star_homotopy_preserving( + const Simplex_handle& simplex) { + return MAYBE_HOMOTOPY_EQ; + } + + enum contractible_status { + NOT_CONTRACTIBLE, + MAYBE_CONTRACTIBLE, + CONTRACTIBLE + }; + /** + * @brief %Test if the complex is reducible using a strategy defined in the class + * (by default it tests if the complex is a cone) + * @details Note that NO could be returned if some invariant ensures that the complex + * is not a point (for instance if the euler characteristic is different from 1). + * This function will surely have to return MAYBE in some case because the + * associated problem is undecidable but it in practice, it can often + * be solved with the help of geometry. + */ + virtual contractible_status is_contractible() const { + if (this->is_cone()) + return CONTRACTIBLE; + else + return MAYBE_CONTRACTIBLE; + // return this->is_cone(); + } + + /** @Edge contraction operations + */ + //@{ + + /** + * @return If ignore_popable_blockers is true + * then the result is true iff the link condition at edge ab is satisfied + * or equivalently iff no blocker contains ab. + * If ignore_popable_blockers is false then the + * result is true iff all blocker containing ab are popable. + */ + bool link_condition(Vertex_handle a, Vertex_handle b, + bool ignore_popable_blockers = false) const { + for (auto blocker : this->const_blocker_range(a)) + if (blocker->contains(b)) { + // false if ignore_popable_blockers is false + // otherwise the blocker has to be popable + return ignore_popable_blockers && is_popable_blocker(blocker); + } + return true; + } + + /** + * @return If ignore_popable_blockers is true + * then the result is true iff the link condition at edge ab is satisfied + * or equivalently iff no blocker contains ab. + * If ignore_popable_blockers is false then the + * result is true iff all blocker containing ab are popable. + */ + bool link_condition(Edge_handle & e, + bool ignore_popable_blockers = false) const { + const Graph_edge& edge = (*this)[e]; + assert(this->get_address(edge.first())); + assert(this->get_address(edge.second())); + Vertex_handle a(*this->get_address(edge.first())); + Vertex_handle b(*this->get_address(edge.second())); + return link_condition(a, b, ignore_popable_blockers); + } + + protected: + /** + * Compute simplices beta such that a.beta is an order 0 blocker + * that may be used to construct a new blocker after contracting ab. + * Suppose that the link condition Link(ab) = Link(a) inter Link(b) + * is satisfied. + */ + void tip_blockers(Vertex_handle a, Vertex_handle b, + std::vector & buffer) const { + for (auto const & blocker : this->const_blocker_range(a)) { + Simplex_handle beta = (*blocker); + beta.remove_vertex(a); + buffer.push_back(beta); + } + + Simplex_handle n; + this->add_neighbours(b, n); + this->remove_neighbours(a, n); + n.remove_vertex(a); + + for (Vertex_handle y : n) { + Simplex_handle beta; + beta.add_vertex(y); + buffer.push_back(beta); + } + } + + private: + /** + * @brief "Replace" the edge 'bx' by the edge 'ax'. + * Assume that the edge 'bx' was present whereas 'ax' was not. + * Precisely, it does not replace edges, but remove 'bx' and then add 'ax'. + * The visitor 'on_swaped_edge' is called just after edge 'ax' had been added + * and just before edge 'bx' had been removed. That way, it can + * eventually access to information of 'ax'. + */ + void swap_edge(Vertex_handle a, Vertex_handle b, Vertex_handle x) { + this->add_edge(a, x); + if (this->visitor) + this->visitor->on_swaped_edge(a, b, x); + this->remove_edge(b, x); + } + + private: + /** + * @brief removes all blockers passing through the edge 'ab' + */ + void delete_blockers_around_vertex(Vertex_handle v) { + std::list < Blocker_handle > blockers_to_delete; + for (auto blocker : this->blocker_range(v)) { + blockers_to_delete.push_back(blocker); + } + while (!blockers_to_delete.empty()) { + this->remove_blocker(blockers_to_delete.back()); + blockers_to_delete.pop_back(); + } + } + /** + * @brief removes all blockers passing through the edge 'ab' + */ + void delete_blockers_around_edge(Vertex_handle a, Vertex_handle b) { + std::list < Blocker_handle > blocker_to_delete; + for (auto blocker : this->blocker_range(a)) + if (blocker->contains(b)) + blocker_to_delete.push_back(blocker); + while (!blocker_to_delete.empty()) { + this->delete_blocker(blocker_to_delete.back()); + blocker_to_delete.pop_back(); + } + } + + public: + /** + * Contracts the edge. + * @remark If the link condition Link(ab) = Link(a) inter Link(b) is not satisfied, + * it removes first all blockers passing through 'ab' + */ + void contract_edge(Edge_handle edge) { + contract_edge(this->first_vertex(edge), this->second_vertex(edge)); + } + + /** + * Contracts the edge connecting vertices a and b. + * @remark If the link condition Link(ab) = Link(a) inter Link(b) is not satisfied, + * it removes first all blockers passing through 'ab' + */ + void contract_edge(Vertex_handle a, Vertex_handle b) { + assert(this->contains_vertex(a)); + assert(this->contains_vertex(b)); + assert(this->contains_edge(a, b)); + + // if some blockers passes through 'ab', we remove them. + if (!link_condition(a, b)) + delete_blockers_around_edge(a, b); + + std::set < Simplex_handle > blockers_to_add; + + get_blockers_to_be_added_after_contraction(a, b, blockers_to_add); + + delete_blockers_around_vertices(a, b); + + update_edges_after_contraction(a, b); + + this->remove_vertex(b); + + notify_changed_edges(a); + + for (auto block : blockers_to_add) + this->add_blocker(block); + + assert(this->contains_vertex(a)); + assert(!this->contains_vertex(b)); + } + + private: + void get_blockers_to_be_added_after_contraction( + Vertex_handle a, Vertex_handle b, + std::set& blockers_to_add) { + blockers_to_add.clear(); + + typedef Skeleton_blocker_link_complex< + Skeleton_blocker_complex > LinkComplexType; + + LinkComplexType link_a(*this, a); + LinkComplexType link_b(*this, b); + + std::vector < Simplex_handle > vector_alpha, vector_beta; + + tip_blockers(a, b, vector_alpha); + tip_blockers(b, a, vector_beta); + + for (auto alpha = vector_alpha.begin(); alpha != vector_alpha.end(); + ++alpha) { + for (auto beta = vector_beta.begin(); beta != vector_beta.end(); ++beta) { + Simplex_handle sigma = *alpha; + sigma.union_vertices(*beta); + Root_simplex_handle sigma_id = this->get_id(sigma); + if (this->contains(sigma) + && proper_faces_in_union(sigma_id, link_a, + link_b)) { + // Blocker_handle blocker = new Simplex_handle(sigma); + sigma.add_vertex(a); + blockers_to_add.insert(sigma); + } + } + } + } + + /** + * delete all blockers that passes through a or b + */ + void delete_blockers_around_vertices(Vertex_handle a, Vertex_handle b) { + std::vector < Blocker_handle > blocker_to_delete; + for (auto bl : this->blocker_range(a)) + blocker_to_delete.push_back(bl); + for (auto bl : this->blocker_range(b)) + blocker_to_delete.push_back(bl); + while (!blocker_to_delete.empty()) { + this->delete_blocker(blocker_to_delete.back()); + blocker_to_delete.pop_back(); + } + } + + void update_edges_after_contraction(Vertex_handle a, Vertex_handle b) { + // We update the set of edges + this->remove_edge(a, b); + + // For all edges {b,x} incident to b, + // we remove {b,x} and add {a,x} if not already there. + while (this->degree(b) > 0) { + Vertex_handle x(*(adjacent_vertices(b.vertex, this->skeleton).first)); + if (!this->contains_edge(a, x)) + // we 'replace' the edge 'bx' by the edge 'ax' + this->swap_edge(a, b, x); + else + this->remove_edge(b, x); + } + } + + void notify_changed_edges(Vertex_handle a) { + // We notify the visitor that all edges incident to 'a' had changed + boost_adjacency_iterator v, v_end; + + for (tie(v, v_end) = adjacent_vertices(a.vertex, this->skeleton); + v != v_end; ++v) + if (this->visitor) + this->visitor->on_changed_edge(a, Vertex_handle(*v)); + } + + //@} }; -} +} // namespace skbl -} // namespace GUDHI +} // namespace Gudhi -#endif /* GUDHI_SKELETON_BLOCKERS_SIMPLIFIABLE_COMPLEX_H_ */ +#endif // SRC_SKELETON_BLOCKER_INCLUDE_GUDHI_SKELETON_BLOCKER_SIMPLIFIABLE_COMPLEX_H_ diff --git a/src/Skeleton_blocker/test/TestSimplifiable.cpp b/src/Skeleton_blocker/test/TestSimplifiable.cpp index 78412754..8ac9594a 100644 --- a/src/Skeleton_blocker/test/TestSimplifiable.cpp +++ b/src/Skeleton_blocker/test/TestSimplifiable.cpp @@ -73,31 +73,33 @@ bool test_contraction1(){ enum { a, b, x, y, z, n }; Complex complex(n); build_complete(n,complex); - complex.remove_edge(Vertex_handle(b),Vertex_handle(z)); - complex.add_blocker(Simplex_handle(Vertex_handle(a),Vertex_handle(x),Vertex_handle(y))); - complex.add_blocker(Simplex_handle(Vertex_handle(b),Vertex_handle(x),Vertex_handle(y))); + complex.remove_edge(static_cast(b), static_cast(z)); + complex.add_blocker(Simplex_handle(static_cast(a), static_cast(x), + static_cast(y))); + complex.add_blocker(Simplex_handle(static_cast(b), static_cast(x), + static_cast(y))); // Print result cerr << "complex before complex"<< complex.to_string()<(a),static_cast(b)); // Print result cerr << "ContractEdge(0,1)\n"; PRINT(complex.to_string()); // verification for (int i=0;i<5;i++) - if (i!=1) assert_vertex(complex,Vertex_handle(i)); - bool test1 = !complex.contains_edge(Vertex_handle(a),Vertex_handle(b)); + if (i!=1) assert_vertex(complex, static_cast(i)); + bool test1 = !complex.contains_edge(static_cast(a),static_cast(b)); bool test2 = assert_blocker(complex,Root_vertex_handle(a),Root_vertex_handle(x),Root_vertex_handle(y)); bool test3 = complex.num_edges()==6; bool test4 = complex.num_blockers()==1; Simplex_handle sigma; - sigma.add_vertex(Vertex_handle(a)); - sigma.add_vertex(Vertex_handle(x)); - sigma.add_vertex(Vertex_handle(y)); - sigma.add_vertex(Vertex_handle(z)); + sigma.add_vertex(static_cast(a)); + sigma.add_vertex(static_cast(x)); + sigma.add_vertex(static_cast(y)); + sigma.add_vertex(static_cast(z)); bool test5 = !(complex.contains(sigma)); return test1&&test2&&test3&&test4&&test5; } @@ -107,18 +109,18 @@ bool test_contraction2(){ enum { a, b, x, y, z, n }; Complex complex(n); build_complete(n,complex); - complex.remove_edge(Vertex_handle(b),Vertex_handle(x)); + complex.remove_edge(static_cast(b),static_cast(x)); Simplex_handle blocker; - blocker.add_vertex(Vertex_handle(a)); - blocker.add_vertex(Vertex_handle(y)); - blocker.add_vertex(Vertex_handle(z)); + blocker.add_vertex(static_cast(a)); + blocker.add_vertex(static_cast(y)); + blocker.add_vertex(static_cast(z)); complex.add_blocker(blocker); // Print result cerr << "complex complex"<< complex.to_string(); cerr <(a),static_cast(b)); cerr << "complex.ContractEdge(a,b)"<< complex.to_string(); @@ -126,7 +128,8 @@ bool test_contraction2(){ // there should be one blocker (a,c,d,e) in the complex bool test ; - test = complex.contains_blocker(Simplex_handle(a,x,y,z)); + test = complex.contains_blocker(Simplex_handle(static_cast(a), static_cast(x), + static_cast(y),static_cast(z))); test = test && complex.num_blockers()==1; return test; } @@ -136,7 +139,7 @@ bool test_link_condition1(){ Complex complex(0); // Build the complexes build_complete(4,complex); - complex.add_blocker(Simplex_handle(Vertex_handle(0),Vertex_handle(1),Vertex_handle(2))); + complex.add_blocker(Simplex_handle(static_cast(0), static_cast(1), static_cast(2))); // Print result @@ -155,13 +158,13 @@ bool test_collapse0(){ Complex complex(5); build_complete(4,complex); complex.add_vertex(); - complex.add_edge(2,4); - complex.add_edge(3,4); + complex.add_edge(static_cast(2), static_cast(4)); + complex.add_edge(static_cast(3), static_cast(4)); // Print result cerr << "initial complex :\n"<< complex.to_string(); cerr <(1), static_cast(2), static_cast(3)); complex.remove_star(simplex_123); cerr << "complex.remove_star(1,2,3):\n"<< complex.to_string(); cerr <(2)); cerr << "complex after remove star of 2:\n"<< complex.to_string(); bool blocker134_here = complex.contains_blocker(Simplex_handle(Vertex_handle(1),Vertex_handle(3),Vertex_handle(4))); diff --git a/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp b/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp index f7bf289a..6fb42836 100644 --- a/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp +++ b/src/Skeleton_blocker/test/TestSkeletonBlockerComplex.cpp @@ -49,10 +49,10 @@ typedef Complex::Root_simplex_handle Root_simplex_handle; typedef Simplex_handle::Simplex_vertex_const_iterator Simplex_vertex_const_iterator; typedef Complex::Edge_handle Edge_handle; -// true iff v \in complex +// true if v in complex bool assert_vertex(Complex &complex,Vertex_handle v){ //assert(complex.contains(v)); - return complex.contains(v); + return complex.contains(static_cast(v)); } bool assert_simplex(Complex &complex,Root_vertex_handle a,Root_vertex_handle b,Root_vertex_handle c){ @@ -206,7 +206,7 @@ bool test_iterator_triangles(){ int num_triangles_seen=0; //for (auto t : complex.triangle_range(5)){ TEST("triangles around 5 (should be 2 of them):"); - for (auto t : complex.triangle_range(5)){ + for (auto t : complex.triangle_range(Vertex_handle(5))){ PRINT(t); ++num_triangles_seen; } @@ -214,7 +214,7 @@ bool test_iterator_triangles(){ num_triangles_seen=0; TEST("triangles around 0 (should be 6 of them):"); - for (auto t : complex.triangle_range(0)){ + for (auto t : complex.triangle_range(Vertex_handle(0))){ PRINT(t); ++num_triangles_seen; } @@ -732,7 +732,7 @@ list subfaces(Simplex_handle top_face){ bool test_constructor2(){ Simplex_handle simplex; for(int i =0 ; i < 5;++i) - simplex.add_vertex(i); + simplex.add_vertex(static_cast(i)); list simplices(subfaces(simplex)); simplices.remove(simplex); -- cgit v1.2.3