diff options
Diffstat (limited to 'src/Coxeter_triangulation/include/gudhi')
38 files changed, 1418 insertions, 1916 deletions
diff --git a/src/Coxeter_triangulation/include/gudhi/Cell_complex.h b/src/Coxeter_triangulation/include/gudhi/Cell_complex.h index 3bc60a50..9cac58d9 100644 --- a/src/Coxeter_triangulation/include/gudhi/Cell_complex.h +++ b/src/Coxeter_triangulation/include/gudhi/Cell_complex.h @@ -11,14 +11,14 @@ #ifndef CELL_COMPLEX_H_ #define CELL_COMPLEX_H_ -#include <Eigen/Dense> +#include <Eigen/Dense> #include <vector> #include <map> #include <utility> // for std::make_pair #include <gudhi/Permutahedral_representation/Simplex_comparator.h> -#include <gudhi/Cell_complex/Hasse_diagram_cell.h> // for Hasse_cell +#include <gudhi/Cell_complex/Hasse_diagram_cell.h> // for Hasse_cell namespace Gudhi { @@ -28,8 +28,8 @@ namespace coxeter_triangulation { * \ingroup coxeter_triangulation */ -/** \class Cell_complex - * \brief A class that constructs the cell complex from the output provided by the class +/** \class Cell_complex + * \brief A class that constructs the cell complex from the output provided by the class * Gudhi::Manifold_tracing<Triangulation_>. * * \tparam Out_simplex_map_ The type of a map from a simplex type that is a @@ -39,13 +39,12 @@ namespace coxeter_triangulation { */ template <class Out_simplex_map_> class Cell_complex { -public: - + public: /** \brief Type of a simplex in the ambient triangulation. * Is a model of the concept SimplexInCoxeterTriangulation. */ using Simplex_handle = typename Out_simplex_map_::key_type; - /** \brief Type of a cell in the cell complex. + /** \brief Type of a cell in the cell complex. * Always is Gudhi::Hasse_cell from the Hasse diagram module. * The additional information is the boolean that is true if and only if the cell lies * on the boundary. @@ -60,7 +59,7 @@ public: * ambient triangulation to the corresponding cells in the cell complex of various dimensions. */ using Simplex_cell_maps = std::vector<Simplex_cell_map>; - + /** \brief Type of a map from cells in the cell complex to the permutahedral representations * of the corresponding simplices in the ambient triangulation. */ @@ -71,21 +70,19 @@ public: */ using Cell_point_map = std::map<Hasse_cell*, Eigen::VectorXd>; -private: - + private: Hasse_cell* insert_cell(const Simplex_handle& simplex, std::size_t cell_d, bool is_boundary) { - Simplex_cell_maps& simplex_cell_maps - = (is_boundary? boundary_simplex_cell_maps_ : interior_simplex_cell_maps_); + Simplex_cell_maps& simplex_cell_maps = (is_boundary ? boundary_simplex_cell_maps_ : interior_simplex_cell_maps_); #ifdef GUDHI_COX_OUTPUT_TO_HTML CC_detail_list& cc_detail_list = - (is_boundary? cc_boundary_detail_lists[cell_d] : cc_interior_detail_lists[cell_d]); + (is_boundary ? cc_boundary_detail_lists[cell_d] : cc_interior_detail_lists[cell_d]); cc_detail_list.emplace_back(CC_detail_info(simplex)); #endif Simplex_cell_map& simplex_cell_map = simplex_cell_maps[cell_d]; auto map_it = simplex_cell_map.find(simplex); if (map_it == simplex_cell_map.end()) { hasse_cells_.push_back(new Hasse_cell(is_boundary, cell_d)); - Hasse_cell* new_cell = hasse_cells_.back(); + Hasse_cell* new_cell = hasse_cells_.back(); simplex_cell_map.emplace(std::make_pair(simplex, new_cell)); cell_simplex_map_.emplace(std::make_pair(new_cell, simplex)); #ifdef GUDHI_COX_OUTPUT_TO_HTML @@ -103,21 +100,21 @@ private: void expand_level(std::size_t cell_d) { bool is_manifold_with_boundary = boundary_simplex_cell_maps_.size() > 0; - for (auto& sc_pair: interior_simplex_cell_maps_[cell_d - 1]) { + for (auto& sc_pair : interior_simplex_cell_maps_[cell_d - 1]) { const Simplex_handle& simplex = sc_pair.first; Hasse_cell* cell = sc_pair.second; - for (Simplex_handle coface: simplex.coface_range(cod_d_ + cell_d)) { + for (Simplex_handle coface : simplex.coface_range(cod_d_ + cell_d)) { Hasse_cell* new_cell = insert_cell(coface, cell_d, false); new_cell->get_boundary().emplace_back(std::make_pair(cell, 1)); } } - + if (is_manifold_with_boundary) { - for (auto& sc_pair: boundary_simplex_cell_maps_[cell_d - 1]) { + for (auto& sc_pair : boundary_simplex_cell_maps_[cell_d - 1]) { const Simplex_handle& simplex = sc_pair.first; Hasse_cell* cell = sc_pair.second; if (cell_d != intr_d_) - for (Simplex_handle coface: simplex.coface_range(cod_d_ + cell_d + 1)) { + for (Simplex_handle coface : simplex.coface_range(cod_d_ + cell_d + 1)) { Hasse_cell* new_cell = insert_cell(coface, cell_d, true); new_cell->get_boundary().emplace_back(std::make_pair(cell, 1)); } @@ -131,29 +128,26 @@ private: } } } - + void construct_complex_(const Out_simplex_map_& out_simplex_map) { #ifdef GUDHI_COX_OUTPUT_TO_HTML cc_interior_summary_lists.resize(interior_simplex_cell_maps_.size()); cc_interior_prejoin_lists.resize(interior_simplex_cell_maps_.size()); cc_interior_detail_lists.resize(interior_simplex_cell_maps_.size()); -#endif - for (auto& os_pair: out_simplex_map) { +#endif + for (auto& os_pair : out_simplex_map) { const Simplex_handle& simplex = os_pair.first; const Eigen::VectorXd& point = os_pair.second; Hasse_cell* new_cell = insert_cell(simplex, 0, false); cell_point_map_.emplace(std::make_pair(new_cell, point)); } for (std::size_t cell_d = 1; - cell_d < interior_simplex_cell_maps_.size() && - !interior_simplex_cell_maps_[cell_d - 1].empty(); - ++cell_d) { + cell_d < interior_simplex_cell_maps_.size() && !interior_simplex_cell_maps_[cell_d - 1].empty(); ++cell_d) { expand_level(cell_d); } } - - void construct_complex_(const Out_simplex_map_& interior_simplex_map, - const Out_simplex_map_& boundary_simplex_map) { + + void construct_complex_(const Out_simplex_map_& interior_simplex_map, const Out_simplex_map_& boundary_simplex_map) { #ifdef GUDHI_COX_OUTPUT_TO_HTML cc_interior_summary_lists.resize(interior_simplex_cell_maps_.size()); cc_interior_prejoin_lists.resize(interior_simplex_cell_maps_.size()); @@ -161,76 +155,70 @@ private: cc_boundary_summary_lists.resize(boundary_simplex_cell_maps_.size()); cc_boundary_prejoin_lists.resize(boundary_simplex_cell_maps_.size()); cc_boundary_detail_lists.resize(boundary_simplex_cell_maps_.size()); -#endif - for (auto& os_pair: boundary_simplex_map) { +#endif + for (auto& os_pair : boundary_simplex_map) { const Simplex_handle& simplex = os_pair.first; const Eigen::VectorXd& point = os_pair.second; Hasse_cell* new_cell = insert_cell(simplex, 0, true); cell_point_map_.emplace(std::make_pair(new_cell, point)); } - for (auto& os_pair: interior_simplex_map) { + for (auto& os_pair : interior_simplex_map) { const Simplex_handle& simplex = os_pair.first; const Eigen::VectorXd& point = os_pair.second; Hasse_cell* new_cell = insert_cell(simplex, 0, false); cell_point_map_.emplace(std::make_pair(new_cell, point)); } #ifdef GUDHI_COX_OUTPUT_TO_HTML - for (const auto& sc_pair: interior_simplex_cell_maps_[0]) + for (const auto& sc_pair : interior_simplex_cell_maps_[0]) cc_interior_summary_lists[0].push_back(CC_summary_info(sc_pair)); - for (const auto& sc_pair: boundary_simplex_cell_maps_[0]) + for (const auto& sc_pair : boundary_simplex_cell_maps_[0]) cc_boundary_summary_lists[0].push_back(CC_summary_info(sc_pair)); -#endif +#endif for (std::size_t cell_d = 1; - cell_d < interior_simplex_cell_maps_.size() && - !interior_simplex_cell_maps_[cell_d - 1].empty(); - ++cell_d) { + cell_d < interior_simplex_cell_maps_.size() && !interior_simplex_cell_maps_[cell_d - 1].empty(); ++cell_d) { expand_level(cell_d); - + #ifdef GUDHI_COX_OUTPUT_TO_HTML - for (const auto& sc_pair: interior_simplex_cell_maps_[cell_d]) + for (const auto& sc_pair : interior_simplex_cell_maps_[cell_d]) cc_interior_summary_lists[cell_d].push_back(CC_summary_info(sc_pair)); if (cell_d < boundary_simplex_cell_maps_.size()) - for (const auto& sc_pair: boundary_simplex_cell_maps_[cell_d]) + for (const auto& sc_pair : boundary_simplex_cell_maps_[cell_d]) cc_boundary_summary_lists[cell_d].push_back(CC_summary_info(sc_pair)); #endif } } - -public: - + + public: /** * \brief Constructs the the cell complex that approximates an \f$m\f$-dimensional manifold * without boundary embedded in the \f$ d \f$-dimensional Euclidean space * from the output of the class Gudhi::Manifold_tracing. * - * \param[in] out_simplex_map A map from simplices of dimension \f$(d-m)\f$ - * in the ambient triangulation that intersect the relative interior of the manifold + * \param[in] out_simplex_map A map from simplices of dimension \f$(d-m)\f$ + * in the ambient triangulation that intersect the relative interior of the manifold * to the intersection points. */ void construct_complex(const Out_simplex_map_& out_simplex_map) { interior_simplex_cell_maps_.resize(intr_d_ + 1); - if (!out_simplex_map.empty()) - cod_d_ = out_simplex_map.begin()->first.dimension(); + if (!out_simplex_map.empty()) cod_d_ = out_simplex_map.begin()->first.dimension(); construct_complex_(out_simplex_map); } - + /** * \brief Constructs the skeleton of the cell complex that approximates - * an \f$m\f$-dimensional manifold without boundary embedded + * an \f$m\f$-dimensional manifold without boundary embedded * in the \f$d\f$-dimensional Euclidean space * up to a limit dimension from the output of the class Gudhi::Manifold_tracing. * - * \param[in] out_simplex_map A map from simplices of dimension \f$(d-m)\f$ - * in the ambient triangulation that intersect the relative interior of the manifold + * \param[in] out_simplex_map A map from simplices of dimension \f$(d-m)\f$ + * in the ambient triangulation that intersect the relative interior of the manifold * to the intersection points. * \param[in] limit_dimension The dimension of the constructed skeleton. */ - void construct_complex(const Out_simplex_map_& out_simplex_map, - std::size_t limit_dimension) { + void construct_complex(const Out_simplex_map_& out_simplex_map, std::size_t limit_dimension) { interior_simplex_cell_maps_.resize(limit_dimension + 1); - if (!out_simplex_map.empty()) - cod_d_ = out_simplex_map.begin()->first.dimension(); + if (!out_simplex_map.empty()) cod_d_ = out_simplex_map.begin()->first.dimension(); construct_complex_(out_simplex_map); } @@ -239,74 +227,64 @@ public: * with boundary embedded in the \f$ d \f$-dimensional Euclidean space * from the output of the class Gudhi::Manifold_tracing. * - * \param[in] interior_simplex_map A map from simplices of dimension \f$(d-m)\f$ - * in the ambient triangulation that intersect the relative interior of the manifold + * \param[in] interior_simplex_map A map from simplices of dimension \f$(d-m)\f$ + * in the ambient triangulation that intersect the relative interior of the manifold * to the intersection points. - * \param[in] boundary_simplex_map A map from simplices of dimension \f$(d-m+1)\f$ - * in the ambient triangulation that intersect the boundary of the manifold + * \param[in] boundary_simplex_map A map from simplices of dimension \f$(d-m+1)\f$ + * in the ambient triangulation that intersect the boundary of the manifold * to the intersection points. */ - void construct_complex(const Out_simplex_map_& interior_simplex_map, - const Out_simplex_map_& boundary_simplex_map) { + void construct_complex(const Out_simplex_map_& interior_simplex_map, const Out_simplex_map_& boundary_simplex_map) { interior_simplex_cell_maps_.resize(intr_d_ + 1); boundary_simplex_cell_maps_.resize(intr_d_); - if (!interior_simplex_map.empty()) - cod_d_ = interior_simplex_map.begin()->first.dimension(); + if (!interior_simplex_map.empty()) cod_d_ = interior_simplex_map.begin()->first.dimension(); construct_complex_(interior_simplex_map, boundary_simplex_map); } /** * \brief Constructs the skeleton of the cell complex that approximates - * an \f$m\f$-dimensional manifold with boundary embedded + * an \f$m\f$-dimensional manifold with boundary embedded * in the \f$d\f$-dimensional Euclidean space * up to a limit dimension from the output of the class Gudhi::Manifold_tracing. * - * \param[in] interior_simplex_map A map from simplices of dimension \f$(d-m)\f$ - * in the ambient triangulation that intersect the relative interior of the manifold + * \param[in] interior_simplex_map A map from simplices of dimension \f$(d-m)\f$ + * in the ambient triangulation that intersect the relative interior of the manifold * to the intersection points. - * \param[in] boundary_simplex_map A map from simplices of dimension \f$(d-m+1)\f$ - * in the ambient triangulation that intersect the boundary of the manifold + * \param[in] boundary_simplex_map A map from simplices of dimension \f$(d-m+1)\f$ + * in the ambient triangulation that intersect the boundary of the manifold * to the intersection points. * \param[in] limit_dimension The dimension of the constructed skeleton. */ - void construct_complex(const Out_simplex_map_& interior_simplex_map, - const Out_simplex_map_& boundary_simplex_map, + void construct_complex(const Out_simplex_map_& interior_simplex_map, const Out_simplex_map_& boundary_simplex_map, std::size_t limit_dimension) { interior_simplex_cell_maps_.resize(limit_dimension + 1); boundary_simplex_cell_maps_.resize(limit_dimension); - if (!interior_simplex_map.empty()) - cod_d_ = interior_simplex_map.begin()->first.dimension(); + if (!interior_simplex_map.empty()) cod_d_ = interior_simplex_map.begin()->first.dimension(); construct_complex_(interior_simplex_map, boundary_simplex_map); } /** * \brief Returns the dimension of the cell complex. */ - std::size_t intrinsic_dimension() const { - return intr_d_; - } + std::size_t intrinsic_dimension() const { return intr_d_; } /** * \brief Returns a vector of maps from the cells of various dimensions in the interior - * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations + * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations * of the corresponding simplices in the ambient triangulation. */ - const Simplex_cell_maps& interior_simplex_cell_maps() const { - return interior_simplex_cell_maps_; - } + const Simplex_cell_maps& interior_simplex_cell_maps() const { return interior_simplex_cell_maps_; } /** * \brief Returns a vector of maps from the cells of various dimensions on the boundary - * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations + * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations * of the corresponding simplices in the ambient triangulation. */ - const Simplex_cell_maps& boundary_simplex_cell_maps() const { - return boundary_simplex_cell_maps_; - } - + const Simplex_cell_maps& boundary_simplex_cell_maps() const { return boundary_simplex_cell_maps_; } + /** - * \brief Returns a map from the cells of a given dimension in the interior - * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations + * \brief Returns a map from the cells of a given dimension in the interior + * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations * of the corresponding simplices in the ambient triangulation. * * \param[in] cell_d The dimension of the cells. @@ -316,8 +294,8 @@ public: } /** - * \brief Returns a map from the cells of a given dimension on the boundary - * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations + * \brief Returns a map from the cells of a given dimension on the boundary + * of the cell complex of type Gudhi::Hasse_cell to the permutahedral representations * of the corresponding simplices in the ambient triangulation. * * \param[in] cell_d The dimension of the cells. @@ -328,35 +306,29 @@ public: /** * \brief Returns a map from the cells in the cell complex of type Gudhi::Hasse_cell - * to the permutahedral representations of the corresponding simplices in the + * to the permutahedral representations of the corresponding simplices in the * ambient triangulation. */ - const Cell_simplex_map& cell_simplex_map() const { - return cell_simplex_map_; - } + const Cell_simplex_map& cell_simplex_map() const { return cell_simplex_map_; } /** * \brief Returns a map from the vertex cells in the cell complex of type Gudhi::Hasse_cell * to their Cartesian coordinates. */ - const Cell_point_map& cell_point_map() const { - return cell_point_map_; - } + const Cell_point_map& cell_point_map() const { return cell_point_map_; } /** * \brief Conxtructor for the class Cell_complex. * * \param[in] intrinsic_dimension The dimension of the cell complex. */ - Cell_complex(std::size_t intrinsic_dimension) - : intr_d_(intrinsic_dimension) {} + Cell_complex(std::size_t intrinsic_dimension) : intr_d_(intrinsic_dimension) {} ~Cell_complex() { - for (Hasse_cell* hs_ptr: hasse_cells_) - delete hs_ptr; + for (Hasse_cell* hs_ptr : hasse_cells_) delete hs_ptr; } - -private: + + private: std::size_t intr_d_, cod_d_; Simplex_cell_maps interior_simplex_cell_maps_, boundary_simplex_cell_maps_; Cell_simplex_map cell_simplex_map_; @@ -364,8 +336,8 @@ private: std::vector<Hasse_cell*> hasse_cells_; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Cell_complex/Hasse_diagram_cell.h b/src/Coxeter_triangulation/include/gudhi/Cell_complex/Hasse_diagram_cell.h index a59026fa..9847a6d3 100644 --- a/src/Coxeter_triangulation/include/gudhi/Cell_complex/Hasse_diagram_cell.h +++ b/src/Coxeter_triangulation/include/gudhi/Cell_complex/Hasse_diagram_cell.h @@ -20,10 +20,8 @@ namespace Gudhi { namespace Hasse_diagram { - -template < typename Cell_type > class Hasse_diagram; - - +template <typename Cell_type> +class Hasse_diagram; /** * \class Hasse_diagram_cell @@ -32,270 +30,255 @@ template < typename Cell_type > class Hasse_diagram; * \ingroup Hasse_diagram * * \details - * This is a data structure to store a cell in a general Hasse diagram data structure. It store the following information about the cell: - * References to boundary and coBoundary elements, dimension of a cell and its filtration. It also allow to store any additional - * information of a type Additional_information which is a template parameter of the class (set by default to void). + * This is a data structure to store a cell in a general Hasse diagram data structure. It store the following + * information about the cell: References to boundary and coBoundary elements, dimension of a cell and its filtration. + * It also allow to store any additional information of a type Additional_information which is a template parameter of + * the class (set by default to void). * * Please refer to \ref Hasse_diagram for examples. * * The complex is a template class requiring the following parameters: - * Incidence_type_ - determine the type of incidence coefficients. Use integers in most general case. + * Incidence_type_ - determine the type of incidence coefficients. Use integers in most general case. * Filtration_type_ - type of filtration of cells. - * Additional_information_ (set by default to void) - allows to store any + * Additional_information_ (set by default to void) - allows to store any * additional information in the cells of Hasse diagrams. * */ -template <typename Incidence_type_, typename Filtration_type_ , typename Additional_information_ = void> -class Hasse_diagram_cell -{ -public: - typedef Incidence_type_ Incidence_type; - typedef Filtration_type_ Filtration_type; - typedef Additional_information_ Additional_information; - using Cell_range = std::vector< std::pair<Hasse_diagram_cell*,Incidence_type> >; - - /** - * Default constructor. - **/ - Hasse_diagram_cell():dimension(0),position(0),deleted_(false){} - - /** - * Constructor of a cell of dimension dim. - **/ - Hasse_diagram_cell( int dim ):dimension(dim),position(0),deleted_(false){} - - /** - * Constructor of a cell of dimension dim. - **/ - Hasse_diagram_cell( int dim , Filtration_type filt_ ):dimension(dim),position(0),deleted_(false),filtration(filt_){} - - /** - * Constructor of a cell of dimension dim with a given boundary. - **/ - Hasse_diagram_cell( const Cell_range& boundary_ , int dim ): - dimension(dim),boundary(boundary_),position(0),deleted_(false){} - - /** - * Constructor of a cell of dimension dim with a given boundary and coboundary. - **/ - Hasse_diagram_cell( const Cell_range& boundary_ , const Cell_range& coboundary_, - int dim ):dimension(dim),boundary(boundary_),coBoundary(coboundary_), - position(0),deleted_(false){} - - /** - * Constructor of a cell of dimension dim with a given boundary, coboundary and - * additional information. - **/ - Hasse_diagram_cell( const Cell_range& boundary_ , const Cell_range& coboundary_, - const Additional_information& ai, int dim ): - dimension(dim),boundary(boundary_),coBoundary(coboundary_),additional_info(ai), - position(0),deleted_(false){} - - /** - * Construcor of a cell of dimension dim having given additional information. - **/ - Hasse_diagram_cell(Additional_information ai, int dim ): - dimension(dim),additional_info(ai),position(0),deleted_(false){} - - /** - * Procedure to get the boundary of a fiven cell. The output format - * is a vector of pairs of pointers to boundary elements and incidence - * coefficients. - **/ - inline Cell_range& get_boundary(){return this->boundary;} - - /** - * Procedure to get the coboundary of a fiven cell. The output format - * is a vector of pairs of pointers to coboundary elements and incidence - * coefficients. - **/ - inline Cell_range& get_coBoundary(){return this->coBoundary;} - - /** - * Procedure to get the dimension of a cell. - **/ - inline int& get_dimension(){return this->dimension;} - - /** - * Procedure to get additional information about the cell.s - **/ - inline Additional_information& get_additional_information(){return this->additional_info;} - - /** - * Procedure to retrive position of the cell in the structure. It is used in - * the implementation of Hasse diagram and set by it. Note that removal of - * cell and subsequent call of clean_up_the_structure will change those - * positions. - **/ - inline unsigned& get_position(){return this->position;} - - /** - * Accessing the filtration of the cell. - **/ - inline Filtration_type& get_filtration() - { - //std::cout << "Accessing the filtration of a cell : " << *this << std::endl; - return this->filtration; - } - - /** - * A procedure used to check if the cell is deleted. It is used by the - * subsequent implementation of Hasse diagram that is absed on lazy - * delete. - **/ - inline bool deleted(){return this->deleted_;} - - template < typename Cell_type > - friend class Hasse_diagram; - - template < typename Cell_type > - friend class is_before_in_filtration; - - - template <typename Complex_type , typename Cell_type > - friend std::vector<Cell_type*> convert_to_vector_of_Cell_type( Complex_type& cmplx ); - - /** - * Procedure to remove deleted boundary and coboundary elements from the - * vectors of boundary and coboundary elements of this cell. - **/ - void remove_deleted_elements_from_boundary_and_coboundary() - { - Cell_range new_boundary; - new_boundary.reserve( this->boundary.size() ); - for ( size_t bd = 0 ; bd != this->boundary.size() ; ++bd ) - { - if ( !this->boundary[bd].first->deleted() ) - { - new_boundary.push_back( this->boundary[bd] ); - } - } - this->boundary.swap( new_boundary ); - - Cell_range new_coBoundary; - new_coBoundary.reserve( this->coBoundary.size() ); - for ( size_t cbd = 0 ; cbd != this->coBoundary.size() ; ++cbd ) - { - if ( !this->coBoundary[cbd].first->deleted() ) - { - new_coBoundary.push_back( this->coBoundary[cbd] ); - } - } - this->coBoundary.swap( new_coBoundary ); - } - - /** - * Writing to a stream operator. - **/ - friend std::ostream& operator<<( std::ostream& out, const Hasse_diagram_cell<Incidence_type,Filtration_type,Additional_information>& c ) - { - //cout << "position : " << c.position << ", dimension : " << c.dimension << ", filtration: " << c.filtration << ", size of boudary : " << c.boundary.size() << "\n"; - out << c.position << " " << c.dimension << " " << c.filtration << std::endl; - for ( size_t bd = 0 ; bd != c.boundary.size() ; ++bd ) - { - //do not write out the cells that has been deleted - if ( c.boundary[bd].first->deleted() )continue; - out << c.boundary[bd].first->position << " " << c.boundary[bd].second << " "; - } - out << std::endl; - return out; - } - - - /** - * Procedure that return vector of pointers to boundary elements of a given cell. - **/ - inline std::vector< Hasse_diagram_cell* > get_list_of_boundary_elements() - { - std::vector< Hasse_diagram_cell* > result; - size_t size_of_boundary = this->boundary.size(); - result.reserve( size_of_boundary ); - for ( size_t bd = 0 ; bd != size_of_boundary ; ++bd ) - { - result.push_back( this->boundary[bd].first ); - } - return result; - } - - /** - * Procedure that return vector of positios of boundary elements of a given cell. - **/ - inline std::vector< unsigned > get_list_of_positions_of_boundary_elements() - { - std::vector< unsigned > result; - size_t size_of_boundary = this->boundary.size(); - result.reserve( size_of_boundary ); - for ( size_t bd = 0 ; bd != size_of_boundary ; ++bd ) - { - result.push_back( this->boundary[bd].first->position ); - } - return result; - } - - /** - * Function that display a string being a signature of a structure. - * Used mainly for debugging purposes. - **/ - std::string full_signature_of_the_structure() - { - std::string result; - result += "dimension: "; - result += std::to_string(this->dimension); - result += " filtration: "; - result += std::to_string(this->filtration); - result += " position: "; - result += std::to_string(this->position); - result += " deleted_: "; - result += std::to_string(this->deleted_); - - //if the Additional_information is not void, add them to - //the signature as well. - if ( std::is_same<Additional_information, void>::value ) - { - result += " Additional_information: "; - result += std::to_string(this->additional_info); - } - result += " boundary "; - for ( size_t bd = 0 ; bd != this->boundary.size() ; ++bd ) - { - result += "( " + std::to_string(this->boundary[bd].first->position); - result += " " + std::to_string(this->boundary[bd].second); - result += ") "; - } - - result += " coBoundary "; - for ( size_t cbd = 0 ; cbd != this->coBoundary.size() ; ++cbd ) - { - result += "( " + std::to_string(this->coBoundary[cbd].first->position); - result += " " + std::to_string(this->coBoundary[cbd].second); - result += ") "; - } - - return result; - } - - -protected: - Cell_range boundary; - Cell_range coBoundary; - int dimension; - Additional_information additional_info; - unsigned position; - bool deleted_; - Filtration_type filtration; - - /** - * A procedure to delete a cell. It is a private function of the Hasse_diagram_cell - * class, since in the Hasse_diagram class I want to have a control - * of removal of cells. Therefore, to remove cell please use - * remove_cell in the Hasse_diagram structure. - **/ - void delete_cell(){ this->deleted_ = true; } -};//Hasse_diagram_cell - - - -}//namespace Hasse_diagram -}//namespace Gudhi - -#endif //CELL_H +template <typename Incidence_type_, typename Filtration_type_, typename Additional_information_ = void> +class Hasse_diagram_cell { + public: + typedef Incidence_type_ Incidence_type; + typedef Filtration_type_ Filtration_type; + typedef Additional_information_ Additional_information; + using Cell_range = std::vector<std::pair<Hasse_diagram_cell*, Incidence_type> >; + + /** + * Default constructor. + **/ + Hasse_diagram_cell() : dimension(0), position(0), deleted_(false) {} + + /** + * Constructor of a cell of dimension dim. + **/ + Hasse_diagram_cell(int dim) : dimension(dim), position(0), deleted_(false) {} + + /** + * Constructor of a cell of dimension dim. + **/ + Hasse_diagram_cell(int dim, Filtration_type filt_) + : dimension(dim), position(0), deleted_(false), filtration(filt_) {} + + /** + * Constructor of a cell of dimension dim with a given boundary. + **/ + Hasse_diagram_cell(const Cell_range& boundary_, int dim) + : dimension(dim), boundary(boundary_), position(0), deleted_(false) {} + + /** + * Constructor of a cell of dimension dim with a given boundary and coboundary. + **/ + Hasse_diagram_cell(const Cell_range& boundary_, const Cell_range& coboundary_, int dim) + : dimension(dim), boundary(boundary_), coBoundary(coboundary_), position(0), deleted_(false) {} + + /** + * Constructor of a cell of dimension dim with a given boundary, coboundary and + * additional information. + **/ + Hasse_diagram_cell(const Cell_range& boundary_, const Cell_range& coboundary_, const Additional_information& ai, + int dim) + : dimension(dim), + boundary(boundary_), + coBoundary(coboundary_), + additional_info(ai), + position(0), + deleted_(false) {} + + /** + * Construcor of a cell of dimension dim having given additional information. + **/ + Hasse_diagram_cell(Additional_information ai, int dim) + : dimension(dim), additional_info(ai), position(0), deleted_(false) {} + + /** + * Procedure to get the boundary of a fiven cell. The output format + * is a vector of pairs of pointers to boundary elements and incidence + * coefficients. + **/ + inline Cell_range& get_boundary() { return this->boundary; } + + /** + * Procedure to get the coboundary of a fiven cell. The output format + * is a vector of pairs of pointers to coboundary elements and incidence + * coefficients. + **/ + inline Cell_range& get_coBoundary() { return this->coBoundary; } + + /** + * Procedure to get the dimension of a cell. + **/ + inline int& get_dimension() { return this->dimension; } + + /** + * Procedure to get additional information about the cell.s + **/ + inline Additional_information& get_additional_information() { return this->additional_info; } + + /** + * Procedure to retrive position of the cell in the structure. It is used in + * the implementation of Hasse diagram and set by it. Note that removal of + * cell and subsequent call of clean_up_the_structure will change those + * positions. + **/ + inline unsigned& get_position() { return this->position; } + + /** + * Accessing the filtration of the cell. + **/ + inline Filtration_type& get_filtration() { + // std::cout << "Accessing the filtration of a cell : " << *this << std::endl; + return this->filtration; + } + + /** + * A procedure used to check if the cell is deleted. It is used by the + * subsequent implementation of Hasse diagram that is absed on lazy + * delete. + **/ + inline bool deleted() { return this->deleted_; } + + template <typename Cell_type> + friend class Hasse_diagram; + + template <typename Cell_type> + friend class is_before_in_filtration; + + template <typename Complex_type, typename Cell_type> + friend std::vector<Cell_type*> convert_to_vector_of_Cell_type(Complex_type& cmplx); + + /** + * Procedure to remove deleted boundary and coboundary elements from the + * vectors of boundary and coboundary elements of this cell. + **/ + void remove_deleted_elements_from_boundary_and_coboundary() { + Cell_range new_boundary; + new_boundary.reserve(this->boundary.size()); + for (size_t bd = 0; bd != this->boundary.size(); ++bd) { + if (!this->boundary[bd].first->deleted()) { + new_boundary.push_back(this->boundary[bd]); + } + } + this->boundary.swap(new_boundary); + + Cell_range new_coBoundary; + new_coBoundary.reserve(this->coBoundary.size()); + for (size_t cbd = 0; cbd != this->coBoundary.size(); ++cbd) { + if (!this->coBoundary[cbd].first->deleted()) { + new_coBoundary.push_back(this->coBoundary[cbd]); + } + } + this->coBoundary.swap(new_coBoundary); + } + + /** + * Writing to a stream operator. + **/ + friend std::ostream& operator<<( + std::ostream& out, const Hasse_diagram_cell<Incidence_type, Filtration_type, Additional_information>& c) { + // cout << "position : " << c.position << ", dimension : " << c.dimension << ", filtration: " << c.filtration << ", + // size of boudary : " << c.boundary.size() << "\n"; + out << c.position << " " << c.dimension << " " << c.filtration << std::endl; + for (size_t bd = 0; bd != c.boundary.size(); ++bd) { + // do not write out the cells that has been deleted + if (c.boundary[bd].first->deleted()) continue; + out << c.boundary[bd].first->position << " " << c.boundary[bd].second << " "; + } + out << std::endl; + return out; + } + + /** + * Procedure that return vector of pointers to boundary elements of a given cell. + **/ + inline std::vector<Hasse_diagram_cell*> get_list_of_boundary_elements() { + std::vector<Hasse_diagram_cell*> result; + size_t size_of_boundary = this->boundary.size(); + result.reserve(size_of_boundary); + for (size_t bd = 0; bd != size_of_boundary; ++bd) { + result.push_back(this->boundary[bd].first); + } + return result; + } + + /** + * Procedure that return vector of positios of boundary elements of a given cell. + **/ + inline std::vector<unsigned> get_list_of_positions_of_boundary_elements() { + std::vector<unsigned> result; + size_t size_of_boundary = this->boundary.size(); + result.reserve(size_of_boundary); + for (size_t bd = 0; bd != size_of_boundary; ++bd) { + result.push_back(this->boundary[bd].first->position); + } + return result; + } + + /** + * Function that display a string being a signature of a structure. + * Used mainly for debugging purposes. + **/ + std::string full_signature_of_the_structure() { + std::string result; + result += "dimension: "; + result += std::to_string(this->dimension); + result += " filtration: "; + result += std::to_string(this->filtration); + result += " position: "; + result += std::to_string(this->position); + result += " deleted_: "; + result += std::to_string(this->deleted_); + + // if the Additional_information is not void, add them to + // the signature as well. + if (std::is_same<Additional_information, void>::value) { + result += " Additional_information: "; + result += std::to_string(this->additional_info); + } + result += " boundary "; + for (size_t bd = 0; bd != this->boundary.size(); ++bd) { + result += "( " + std::to_string(this->boundary[bd].first->position); + result += " " + std::to_string(this->boundary[bd].second); + result += ") "; + } + + result += " coBoundary "; + for (size_t cbd = 0; cbd != this->coBoundary.size(); ++cbd) { + result += "( " + std::to_string(this->coBoundary[cbd].first->position); + result += " " + std::to_string(this->coBoundary[cbd].second); + result += ") "; + } + + return result; + } + + protected: + Cell_range boundary; + Cell_range coBoundary; + int dimension; + Additional_information additional_info; + unsigned position; + bool deleted_; + Filtration_type filtration; + + /** + * A procedure to delete a cell. It is a private function of the Hasse_diagram_cell + * class, since in the Hasse_diagram class I want to have a control + * of removal of cells. Therefore, to remove cell please use + * remove_cell in the Hasse_diagram structure. + **/ + void delete_cell() { this->deleted_ = true; } +}; // Hasse_diagram_cell + +} // namespace Hasse_diagram +} // namespace Gudhi + +#endif // CELL_H diff --git a/src/Coxeter_triangulation/include/gudhi/Coxeter_triangulation.h b/src/Coxeter_triangulation/include/gudhi/Coxeter_triangulation.h index 19ceb007..2a57666d 100644 --- a/src/Coxeter_triangulation/include/gudhi/Coxeter_triangulation.h +++ b/src/Coxeter_triangulation/include/gudhi/Coxeter_triangulation.h @@ -26,10 +26,10 @@ #include <gudhi/Permutahedral_representation.h> namespace Gudhi { - + namespace coxeter_triangulation { -/** +/** * \class Coxeter_triangulation * \brief A class that stores Coxeter triangulation of type \f$\tilde{A}_d\f$. * This triangulation has the greatest simplex quality out of all linear transformations @@ -40,50 +40,48 @@ namespace coxeter_triangulation { * \tparam Permutahedral_representation_ Type of a simplex given by a permutahedral representation. * Needs to be a model of SimplexInCoxeterTriangulation. */ -template <class Permutahedral_representation_ - = Permutahedral_representation<std::vector<int>, std::vector<std::vector<std::size_t> > > > +template <class Permutahedral_representation_ = + Permutahedral_representation<std::vector<int>, std::vector<std::vector<std::size_t> > > > class Coxeter_triangulation : public Freudenthal_triangulation<Permutahedral_representation_> { using Matrix = Eigen::MatrixXd; Matrix root_matrix(unsigned d) { - Matrix cartan(d,d); + Matrix cartan(d, d); for (unsigned i = 0; i < d; i++) { - cartan(i,i) = 1.0; + cartan(i, i) = 1.0; } for (unsigned i = 1; i < d; i++) { - cartan(i-1,i) = -0.5; - cartan(i,i-1) = -0.5; + cartan(i - 1, i) = -0.5; + cartan(i, i - 1) = -0.5; } for (unsigned i = 0; i < d; i++) for (unsigned j = 0; j < d; j++) - if (j+1 < i || j > i+1) - cartan(i,j) = 0; + if (j + 1 < i || j > i + 1) cartan(i, j) = 0; Eigen::SelfAdjointEigenSolver<Matrix> saes(cartan); Eigen::VectorXd sqrt_diag(d); - for (unsigned i = 0; i < d; ++i) - sqrt_diag(i) = std::sqrt(saes.eigenvalues()[i]); + for (unsigned i = 0; i < d; ++i) sqrt_diag(i) = std::sqrt(saes.eigenvalues()[i]); - Matrix lower(d,d); + Matrix lower(d, d); for (unsigned i = 0; i < d; i++) for (unsigned j = 0; j < d; j++) if (i < j) - lower(i,j) = 0; + lower(i, j) = 0; else - lower(i,j) = 1; - Matrix result = (lower * saes.eigenvectors()*sqrt_diag.asDiagonal()).inverse(); + lower(i, j) = 1; + Matrix result = (lower * saes.eigenvectors() * sqrt_diag.asDiagonal()).inverse(); return result; } public: - /** \brief Constructor of Coxeter triangulation of a given dimension. + /** \brief Constructor of Coxeter triangulation of a given dimension. * @param[in] dimension The dimension of the triangulation. */ Coxeter_triangulation(std::size_t dimension) - : Freudenthal_triangulation<Permutahedral_representation_>(dimension, root_matrix(dimension)) {} + : Freudenthal_triangulation<Permutahedral_representation_>(dimension, root_matrix(dimension)) {} }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Freudenthal_triangulation.h b/src/Coxeter_triangulation/include/gudhi/Freudenthal_triangulation.h index 6b33a00c..43b33d1f 100644 --- a/src/Coxeter_triangulation/include/gudhi/Freudenthal_triangulation.h +++ b/src/Coxeter_triangulation/include/gudhi/Freudenthal_triangulation.h @@ -13,9 +13,9 @@ #include <vector> #include <algorithm> // for std::sort -#include <cmath> // for std::floor -#include <numeric> // for std::iota -#include <cstdlib> // for std::size_t +#include <cmath> // for std::floor +#include <numeric> // for std::iota +#include <cstdlib> // for std::size_t #include <Eigen/Eigenvalues> #include <Eigen/Sparse> @@ -27,7 +27,7 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** * \class Freudenthal_triangulation * \brief A class that stores any affine transformation of the Freudenthal-Kuhn * triangulation. @@ -41,32 +41,30 @@ namespace coxeter_triangulation { * \tparam Permutahedral_representation_ Type of a simplex given by a permutahedral representation. * Needs to be a model of SimplexInCoxeterTriangulation. */ -template <class Permutahedral_representation_ - = Permutahedral_representation<std::vector<int>, std::vector<std::vector<std::size_t> > > > +template <class Permutahedral_representation_ = + Permutahedral_representation<std::vector<int>, std::vector<std::vector<std::size_t> > > > class Freudenthal_triangulation { using Matrix = Eigen::MatrixXd; using Vector = Eigen::VectorXd; using SparseMatrix = Eigen::SparseMatrix<double>; using Triplet = Eigen::Triplet<double>; - + public: /** \brief Type of the simplices in the triangulation. */ using Simplex_handle = Permutahedral_representation_; /** \brief Type of the vertices in the triangulation. */ using Vertex_handle = typename Permutahedral_representation_::Vertex; - - - /** \brief Constructor of the Freudenthal-Kuhn triangulation of a given dimension. + + /** \brief Constructor of the Freudenthal-Kuhn triangulation of a given dimension. * @param[in] dimension The dimension of the triangulation. */ Freudenthal_triangulation(std::size_t dimension) - : dimension_(dimension), - matrix_(Matrix::Identity(dimension, dimension)), - offset_(Vector::Zero(dimension)), - colpivhouseholderqr_(matrix_.colPivHouseholderQr()), - is_freudenthal(true) { - } + : dimension_(dimension), + matrix_(Matrix::Identity(dimension, dimension)), + offset_(Vector::Zero(dimension)), + colpivhouseholderqr_(matrix_.colPivHouseholderQr()), + is_freudenthal(true) {} /** \brief Constructor of the Freudenthal-Kuhn triangulation of a given dimension under * a linear transformation by a given matrix. @@ -75,12 +73,11 @@ class Freudenthal_triangulation { * Needs to be invertible. */ Freudenthal_triangulation(std::size_t dimension, const Matrix& matrix) - : dimension_(dimension), - matrix_(matrix), - offset_(Vector::Zero(dimension)), - colpivhouseholderqr_(matrix_.colPivHouseholderQr()), - is_freudenthal(false) { - } + : dimension_(dimension), + matrix_(matrix), + offset_(Vector::Zero(dimension)), + colpivhouseholderqr_(matrix_.colPivHouseholderQr()), + is_freudenthal(false) {} /** \brief Constructor of the Freudenthal-Kuhn triangulation of a given dimension under * an affine transformation by a given matrix and a translation vector. @@ -90,64 +87,56 @@ class Freudenthal_triangulation { * @param[in] offset The offset vector. */ Freudenthal_triangulation(unsigned dimension, const Matrix& matrix, const Vector& offset) - : dimension_(dimension), - matrix_(matrix), - offset_(offset), - colpivhouseholderqr_(matrix_.colPivHouseholderQr()), - is_freudenthal(false) { - } - + : dimension_(dimension), + matrix_(matrix), + offset_(offset), + colpivhouseholderqr_(matrix_.colPivHouseholderQr()), + is_freudenthal(false) {} /** \brief Dimension of the triangulation. */ - unsigned dimension() const { - return dimension_; - } + unsigned dimension() const { return dimension_; } /** \brief Matrix that defines the linear transformation of the triangulation. */ - const Matrix& matrix() const { - return matrix_; - } + const Matrix& matrix() const { return matrix_; } /** \brief Vector that defines the offset of the triangulation. */ - const Vector& offset() const { - return offset_; - } + const Vector& offset() const { return offset_; } - /** \brief Change the linear transformation matrix to a given value. + /** \brief Change the linear transformation matrix to a given value. * @param[in] matrix New value of the linear transformation matrix. */ void change_matrix(const Eigen::MatrixXd& matrix) { matrix_ = matrix; colpivhouseholderqr_ = matrix.colPivHouseholderQr(); is_freudenthal = false; - } + } - /** \brief Change the offset vector to a given value. + /** \brief Change the offset vector to a given value. * @param[in] offset New value of the offset vector. */ void change_offset(const Eigen::VectorXd& offset) { offset_ = offset; is_freudenthal = false; } - + /** \brief Returns the permutahedral representation of the simplex in the * triangulation that contains a given query point. - * \details Using the additional parameter scale, the search can be done in a + * \details Using the additional parameter scale, the search can be done in a * triangulation that shares the origin, but is scaled by a given factor. * This parameter can be useful to simulate the point location in a subdivided * triangulation. * The returned simplex is always minimal by inclusion. - * + * * \tparam Point_d A class that represents a point in d-dimensional Euclidean space. * The coordinates should be random-accessible. Needs to provide the method size(). * * @param[in] point The query point. - * @param[in] scale The scale of the triangulation. + * @param[in] scale The scale of the triangulation. */ template <class Point_d> - Simplex_handle locate_point(const Point_d& point, double scale=1) const { - using Ordered_set_partition = typename Simplex_handle::OrderedSetPartition; - using Part = typename Ordered_set_partition::value_type; + Simplex_handle locate_point(const Point_d& point, double scale = 1) const { + using Ordered_set_partition = typename Simplex_handle::OrderedSetPartition; + using Part = typename Ordered_set_partition::value_type; unsigned d = point.size(); assert(d == dimension_); double error = 1e-9; @@ -160,11 +149,9 @@ class Freudenthal_triangulation { output.vertex().push_back(y_i); z.push_back(x_i - y_i); } - } - else { + } else { Eigen::VectorXd p_vect(d); - for (std::size_t i = 0; i < d; i++) - p_vect(i) = point[i]; + for (std::size_t i = 0; i < d; i++) p_vect(i) = point[i]; Eigen::VectorXd x_vect = colpivhouseholderqr_.solve(p_vect - offset_); for (std::size_t i = 0; i < d; i++) { double x_i = scale * x_vect(i); @@ -174,15 +161,13 @@ class Freudenthal_triangulation { } } z.push_back(0); - Part indices(d+1); + Part indices(d + 1); std::iota(indices.begin(), indices.end(), 0); - std::sort(indices.begin(), - indices.end(), - [&z](std::size_t i1, std::size_t i2) {return z[i1] > z[i2];}); + std::sort(indices.begin(), indices.end(), [&z](std::size_t i1, std::size_t i2) { return z[i1] > z[i2]; }); output.partition().push_back(Part(1, indices[0])); for (std::size_t i = 1; i <= d; ++i) - if (z[indices[i-1]] > z[indices[i]] + error) + if (z[indices[i - 1]] > z[indices[i]] + error) output.partition().push_back(Part(1, indices[i])); else output.partition().back().push_back(indices[i]); @@ -190,39 +175,37 @@ class Freudenthal_triangulation { } /** \brief Returns the Cartesian coordinates of the given vertex. - * \details Using the additional parameter scale, the search can be done in a + * \details Using the additional parameter scale, the search can be done in a * triangulation that shares the origin, but is scaled by a given factor. - * This parameter can be useful to simulate the computation of Cartesian coordinates + * This parameter can be useful to simulate the computation of Cartesian coordinates * of a vertex in a subdivided triangulation. * @param[in] vertex The query vertex. - * @param[in] scale The scale of the triangulation. + * @param[in] scale The scale of the triangulation. */ Eigen::VectorXd cartesian_coordinates(const Vertex_handle& vertex, double scale = 1) const { Eigen::VectorXd v_vect(dimension_); - for (std::size_t j = 0; j < dimension_; j++) - v_vect(j) = vertex[j] / scale; + for (std::size_t j = 0; j < dimension_; j++) v_vect(j) = vertex[j] / scale; return matrix_ * v_vect + offset_; } /** \brief Returns the Cartesian coordinates of the barycenter of a given simplex. - * \details Using the additional parameter scale, the search can be done in a + * \details Using the additional parameter scale, the search can be done in a * triangulation that shares the origin, but is scaled by a given factor. - * This parameter can be useful to simulate the computation of Cartesian coordinates + * This parameter can be useful to simulate the computation of Cartesian coordinates * of the barycenter of a simplex in a subdivided triangulation. * @param[in] simplex The query simplex. - * @param[in] scale The scale of the triangulation. + * @param[in] scale The scale of the triangulation. */ Eigen::VectorXd barycenter(const Simplex_handle& simplex, double scale = 1) const { Eigen::VectorXd res_vector(dimension_); - for (size_t i = 0; i < dimension_; ++i) - res_vector(i) = 0; - for (auto v: simplex.vertex_range()) { + for (size_t i = 0; i < dimension_; ++i) res_vector(i) = 0; + for (auto v : simplex.vertex_range()) { res_vector += cartesian_coordinates(v, scale); } - return (1./(simplex.dimension()+1)) * res_vector; + return (1. / (simplex.dimension() + 1)) * res_vector; } - -protected: + + protected: unsigned dimension_; Matrix matrix_; Vector offset_; @@ -230,8 +213,8 @@ protected: bool is_freudenthal; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Cartesian_product.h b/src/Coxeter_triangulation/include/gudhi/Functions/Cartesian_product.h index 0a6f264d..c0d6aec4 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Cartesian_product.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Cartesian_product.h @@ -14,7 +14,7 @@ #include <cstdlib> #include <tuple> #include <type_traits> // for std::enable_if -#include <cstdlib> // for std::size_t +#include <cstdlib> // for std::size_t #include <gudhi/Functions/Function.h> @@ -27,86 +27,79 @@ namespace coxeter_triangulation { /* Get the domain dimension of the tuple of functions. */ template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I == sizeof... (T), std::size_t>::type -get_amb_d (const std::tuple<T...>& tuple) { +inline typename std::enable_if<I == sizeof...(T), std::size_t>::type get_amb_d(const std::tuple<T...>& tuple) { return 0; } template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I != sizeof... (T), std::size_t>::type -get_amb_d (const std::tuple<T...>& tuple) { - return std::get<I>(tuple).amb_d() + get_amb_d<I+1, T...>(tuple); +inline typename std::enable_if<I != sizeof...(T), std::size_t>::type get_amb_d(const std::tuple<T...>& tuple) { + return std::get<I>(tuple).amb_d() + get_amb_d<I + 1, T...>(tuple); } /* Get the codomain dimension of the tuple of functions. */ template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I == sizeof... (T), std::size_t>::type -get_cod_d (const std::tuple<T...>& tuple) { +inline typename std::enable_if<I == sizeof...(T), std::size_t>::type get_cod_d(const std::tuple<T...>& tuple) { return 0; } template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I != sizeof... (T), std::size_t>::type -get_cod_d (const std::tuple<T...>& tuple) { - return std::get<I>(tuple).cod_d() + get_cod_d<I+1, T...>(tuple); +inline typename std::enable_if<I != sizeof...(T), std::size_t>::type get_cod_d(const std::tuple<T...>& tuple) { + return std::get<I>(tuple).cod_d() + get_cod_d<I + 1, T...>(tuple); } /* Get the seed of the tuple of functions. */ template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I == sizeof... (T), void>::type -get_seed (const std::tuple<T...>& tuple, Eigen::VectorXd& point, std::size_t i = 0) { -} +inline typename std::enable_if<I == sizeof...(T), void>::type get_seed(const std::tuple<T...>& tuple, + Eigen::VectorXd& point, std::size_t i = 0) {} template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I != sizeof... (T), void>::type -get_seed (const std::tuple<T...>& tuple, Eigen::VectorXd& point, std::size_t i = 0) { +inline typename std::enable_if<I != sizeof...(T), void>::type get_seed(const std::tuple<T...>& tuple, + Eigen::VectorXd& point, std::size_t i = 0) { const auto& f = std::get<I>(tuple); std::size_t n = f.amb_d(); Eigen::VectorXd seed = f.seed(); - for (std::size_t j = 0; j < n; ++j) - point(i+j) = seed(j); - get_seed<I+1, T...>(tuple, point, i+n); + for (std::size_t j = 0; j < n; ++j) point(i + j) = seed(j); + get_seed<I + 1, T...>(tuple, point, i + n); } /* Get the seed of the tuple of functions. */ template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I == sizeof... (T), void>::type -get_value (const std::tuple<T...>& tuple, const Eigen::VectorXd& x, Eigen::VectorXd& point, std::size_t i = 0, std::size_t j = 0) { -} +inline typename std::enable_if<I == sizeof...(T), void>::type get_value(const std::tuple<T...>& tuple, + const Eigen::VectorXd& x, + Eigen::VectorXd& point, std::size_t i = 0, + std::size_t j = 0) {} template <std::size_t I = 0, typename... T> -inline typename std::enable_if<I != sizeof... (T), void>::type -get_value (const std::tuple<T...>& tuple, const Eigen::VectorXd& x, Eigen::VectorXd& point, std::size_t i = 0, std::size_t j = 0) { +inline typename std::enable_if<I != sizeof...(T), void>::type get_value(const std::tuple<T...>& tuple, + const Eigen::VectorXd& x, + Eigen::VectorXd& point, std::size_t i = 0, + std::size_t j = 0) { const auto& f = std::get<I>(tuple); std::size_t n = f.amb_d(); std::size_t k = f.cod_d(); Eigen::VectorXd x_i(n); - for (std::size_t l = 0; l < n; ++l) - x_i(l) = x(i+l); + for (std::size_t l = 0; l < n; ++l) x_i(l) = x(i + l); Eigen::VectorXd res = f(x_i); - for (std::size_t l = 0; l < k; ++l) - point(j+l) = res(l); - get_value<I+1, T...>(tuple, x, point, i+n, j+k); + for (std::size_t l = 0; l < k; ++l) point(j + l) = res(l); + get_value<I + 1, T...>(tuple, x, point, i + n, j + k); } - -/** +/** * \class Cartesian_product * \brief Constructs the function the zero-set of which is the Cartesian product * of the zero-sets of some given functions. * - * \tparam Functions A pack template parameter for functions. All functions should be models of + * \tparam Functions A pack template parameter for functions. All functions should be models of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class... Functions> struct Cartesian_product : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -117,10 +110,10 @@ struct Cartesian_product : public Function { } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return amb_d_;} + virtual std::size_t amb_d() const override { return amb_d_; } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return cod_d_;} + virtual std::size_t cod_d() const override { return cod_d_; } /** \brief Returns a point on the zero-set. */ virtual Eigen::VectorXd seed() const override { @@ -128,15 +121,14 @@ struct Cartesian_product : public Function { get_seed(function_tuple_, result, 0); return result; } - - /** + + /** * \brief Constructor of the Cartesian product function. * * @param[in] functions The functions the zero-sets of which are factors in the * Cartesian product of the resulting function. */ - Cartesian_product(const Functions&... functions) - : function_tuple_(std::make_tuple(functions...)) { + Cartesian_product(const Functions&... functions) : function_tuple_(std::make_tuple(functions...)) { amb_d_ = get_amb_d(function_tuple_); cod_d_ = get_cod_d(function_tuple_); } @@ -146,14 +138,13 @@ struct Cartesian_product : public Function { std::size_t amb_d_, cod_d_; }; - -/** +/** * \brief Static constructor of a Cartesian product function. * * @param[in] functions The functions the zero-sets of which are factors in the * Cartesian product of the resulting function. * - * \tparam Functions A pack template parameter for functions. All functions should be models of + * \tparam Functions A pack template parameter for functions. All functions should be models of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation @@ -163,9 +154,8 @@ Cartesian_product<Functions...> make_product_function(const Functions&... functi return Cartesian_product<Functions...>(functions...); } +} // namespace coxeter_triangulation -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Constant_function.h b/src/Coxeter_triangulation/include/gudhi/Functions/Constant_function.h index c03a2a24..424bcbff 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Constant_function.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Constant_function.h @@ -21,15 +21,14 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Constant_function +/** + * \class Constant_function * \brief A class that encodes a constant function from R^d to R^k. * This class does not have any implicit manifold in correspondence. * * \ingroup coxeter_triangulation */ struct Constant_function : public Function { - /** \brief Value of the function at a specified point. The value is constant. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -37,37 +36,34 @@ struct Constant_function : public Function { Eigen::VectorXd result = value_; return result; } - + /** \brief Returns the domain dimension. Same as the ambient dimension of the sphere. */ - virtual std::size_t amb_d() const override {return d_;}; + virtual std::size_t amb_d() const override { return d_; }; /** \brief Returns the codomain dimension. Same as the codimension of the sphere. */ - virtual std::size_t cod_d() const override {return k_;}; + virtual std::size_t cod_d() const override { return k_; }; /** \brief No seed point is available. Throws an exception on evocation. */ - virtual Eigen::VectorXd seed() const override { - throw "Seed invoked on a constant function.\n"; - } + virtual Eigen::VectorXd seed() const override { throw "Seed invoked on a constant function.\n"; } Constant_function() {} - - /** + + /** * \brief Constructor of a constant function from R^d to R^m. * * @param[in] d The domain dimension. * @param[in] k The codomain dimension. * @param[in] value The constant value of the function. */ - Constant_function(std::size_t d, std::size_t k, const Eigen::VectorXd& value) - : d_(d), k_(k), value_(value) {} + Constant_function(std::size_t d, std::size_t k, const Eigen::VectorXd& value) : d_(d), k_(k), value_(value) {} private: std::size_t d_, k_; Eigen::VectorXd value_; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Embed_in_Rd.h b/src/Coxeter_triangulation/include/gudhi/Functions/Embed_in_Rd.h index 5ebeac75..69bd269a 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Embed_in_Rd.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Embed_in_Rd.h @@ -17,7 +17,6 @@ #include <Eigen/Dense> - namespace Gudhi { namespace coxeter_triangulation { @@ -26,70 +25,62 @@ namespace coxeter_triangulation { * \class Embed_in_Rd * \brief Embedding of an implicit manifold in a higher dimension. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class Function_> struct Embed_in_Rd : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { Eigen::VectorXd x = p; Eigen::VectorXd x_k(fun_.amb_d()), x_rest(d_ - fun_.amb_d()); - for (std::size_t i = 0; i < fun_.amb_d(); ++i) - x_k(i) = x(i); - for (std::size_t i = fun_.amb_d(); i < d_; ++i) - x_rest(i - fun_.amb_d()) = x(i); + for (std::size_t i = 0; i < fun_.amb_d(); ++i) x_k(i) = x(i); + for (std::size_t i = fun_.amb_d(); i < d_; ++i) x_rest(i - fun_.amb_d()) = x(i); Eigen::VectorXd result = fun_(x_k); result.conservativeResize(this->cod_d()); - for (std::size_t i = fun_.cod_d(); i < this->cod_d(); ++i) - result(i) = x_rest(i - fun_.cod_d()); + for (std::size_t i = fun_.cod_d(); i < this->cod_d(); ++i) result(i) = x_rest(i - fun_.cod_d()); return result; } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return d_;} + virtual std::size_t amb_d() const override { return d_; } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return d_-(fun_.amb_d() - fun_.cod_d());} + virtual std::size_t cod_d() const override { return d_ - (fun_.amb_d() - fun_.cod_d()); } /** \brief Returns a point on the zero-set of the embedded function. */ virtual Eigen::VectorXd seed() const override { Eigen::VectorXd result = fun_.seed(); result.conservativeResize(d_); - for (std::size_t l = fun_.amb_d(); l < d_; ++l) - result(l) = 0; + for (std::size_t l = fun_.amb_d(); l < d_; ++l) result(l) = 0; return result; } - /** + /** * \brief Constructor of the embedding function. * * @param[in] function The function to be embedded in higher dimension. * @param[in] d Embedding dimension. */ - Embed_in_Rd(const Function_& function, std::size_t d) : - fun_(function), d_(d) { - } + Embed_in_Rd(const Function_& function, std::size_t d) : fun_(function), d_(d) {} private: Function_ fun_; std::size_t d_; }; - -/** +/** * \brief Static constructor of an embedding function. * * @param[in] function The function to be embedded in higher dimension. * @param[in] d Embedding dimension. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation @@ -99,8 +90,8 @@ Embed_in_Rd<Function_> make_embedding(const Function_& function, std::size_t d) return Embed_in_Rd<Function_>(function, d); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function.h index 29fadcbf..eddfedf5 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function.h @@ -17,7 +17,7 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** * \class Function * \brief The parent class for all functions implemented in the module. * Contains virtual methods needed to be a model of the concept FunctionForImplicitManifold. @@ -25,7 +25,6 @@ namespace coxeter_triangulation { * \ingroup coxeter_triangulation */ struct Function { - /** \brief Virtual method for the value of the function at a specified point. * @param[in] p The input point. */ @@ -44,9 +43,8 @@ struct Function { virtual ~Function() {} }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_Sm_in_Rd.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_Sm_in_Rd.h index 7ebb7fc8..25b9b5fc 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_Sm_in_Rd.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_Sm_in_Rd.h @@ -21,47 +21,42 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_Sm_in_Rd +/** + * \class Function_Sm_in_Rd * \brief A class for the function that defines an m-dimensional implicit sphere embedded * in the d-dimensional Euclidean space. * * \ingroup coxeter_triangulation */ -struct Function_Sm_in_Rd: public Function { - -/** \brief Value of the function at a specified point. - * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. - */ +struct Function_Sm_in_Rd : public Function { + /** \brief Value of the function at a specified point. + * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. + */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { Eigen::VectorXd x = p; - for (std::size_t i = 0; i < d_; ++i) - x(i) -= center_[i]; - Eigen::VectorXd result = Eigen::VectorXd::Zero(k_); - for (std::size_t i = 0; i < m_+1; ++i) - result(0) += x(i)*x(i); - result(0) -= r_*r_; - for (std::size_t j = 1; j < k_; ++j) - result(j) = x(m_+j); + for (std::size_t i = 0; i < d_; ++i) x(i) -= center_[i]; + Eigen::VectorXd result = Eigen::VectorXd::Zero(k_); + for (std::size_t i = 0; i < m_ + 1; ++i) result(0) += x(i) * x(i); + result(0) -= r_ * r_; + for (std::size_t j = 1; j < k_; ++j) result(j) = x(m_ + j); return result; } - + /** \brief Returns the domain dimension. Same as the ambient dimension of the sphere. */ - virtual std::size_t amb_d() const override {return d_;}; + virtual std::size_t amb_d() const override { return d_; }; /** \brief Returns the codomain dimension. Same as the codimension of the sphere. */ - virtual std::size_t cod_d() const override {return k_;}; + virtual std::size_t cod_d() const override { return k_; }; /** \brief Returns a point on the sphere. */ virtual Eigen::VectorXd seed() const override { Eigen::VectorXd result = Eigen::VectorXd::Zero(d_); result(0) += r_; - for (std::size_t i = 0; i < d_; ++i) - result(i) += center_[i]; + for (std::size_t i = 0; i < d_; ++i) result(i) += center_[i]; return result; } - - /** + + /** * \brief Constructor of the function that defines an m-dimensional implicit sphere embedded * in the d-dimensional Euclidean space. * @@ -70,13 +65,10 @@ struct Function_Sm_in_Rd: public Function { * @param[in] d The ambient dimension of the sphere. * @param[in] center The center of the sphere. */ - Function_Sm_in_Rd(double r, - std::size_t m, - std::size_t d, - Eigen::VectorXd center) - : m_(m), k_(d-m), d_(d), r_(r), center_(center) {} + Function_Sm_in_Rd(double r, std::size_t m, std::size_t d, Eigen::VectorXd center) + : m_(m), k_(d - m), d_(d), r_(r), center_(center) {} - /** + /** * \brief Constructor of the function that defines an m-dimensional implicit sphere embedded * in the d-dimensional Euclidean space centered at the origin. * @@ -84,13 +76,10 @@ struct Function_Sm_in_Rd: public Function { * @param[in] m The dimension of the sphere. * @param[in] d The ambient dimension of the sphere. */ - Function_Sm_in_Rd(double r, - std::size_t m, - std::size_t d) - : m_(m), k_(d-m), d_(d), r_(r), center_(Eigen::VectorXd::Zero(d_)) {} + Function_Sm_in_Rd(double r, std::size_t m, std::size_t d) + : m_(m), k_(d - m), d_(d), r_(r), center_(Eigen::VectorXd::Zero(d_)) {} - - /** + /** * \brief Constructor of the function that defines an m-dimensional implicit sphere embedded * in the (m+1)-dimensional Euclidean space. * @@ -98,35 +87,28 @@ struct Function_Sm_in_Rd: public Function { * @param[in] m The dimension of the sphere. * @param[in] center The center of the sphere. */ - Function_Sm_in_Rd(double r, - std::size_t m, - Eigen::VectorXd center) - : m_(m), k_(1), d_(m_+1), r_(r), center_(center) {} + Function_Sm_in_Rd(double r, std::size_t m, Eigen::VectorXd center) + : m_(m), k_(1), d_(m_ + 1), r_(r), center_(center) {} - /** + /** * \brief Constructor of the function that defines an m-dimensional implicit sphere embedded * in the (m+1)-dimensional Euclidean space centered at the origin. * * @param[in] r The radius of the sphere. * @param[in] m The dimension of the sphere. */ - Function_Sm_in_Rd(double r, - std::size_t m) - : m_(m), k_(1), d_(m_+1), r_(r), center_(Eigen::VectorXd::Zero(d_)) {} + Function_Sm_in_Rd(double r, std::size_t m) : m_(m), k_(1), d_(m_ + 1), r_(r), center_(Eigen::VectorXd::Zero(d_)) {} - Function_Sm_in_Rd(const Function_Sm_in_Rd& rhs) - : Function_Sm_in_Rd(rhs.r_, rhs.m_, rhs.d_, rhs.center_) {} + Function_Sm_in_Rd(const Function_Sm_in_Rd& rhs) : Function_Sm_in_Rd(rhs.r_, rhs.m_, rhs.d_, rhs.center_) {} - private: std::size_t m_, k_, d_; double r_; Eigen::VectorXd center_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_affine_plane_in_Rd.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_affine_plane_in_Rd.h index 1e950d4e..d0089149 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_affine_plane_in_Rd.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_affine_plane_in_Rd.h @@ -21,16 +21,15 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_affine_plane_in_Rd - * \brief A class for the function that defines an m-dimensional implicit affine plane +/** + * \class Function_affine_plane_in_Rd + * \brief A class for the function that defines an m-dimensional implicit affine plane * embedded in d-dimensional Euclidean space. * * \ingroup coxeter_triangulation */ struct Function_affine_plane_in_Rd : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -38,65 +37,59 @@ struct Function_affine_plane_in_Rd : public Function { Eigen::VectorXd result = normal_matrix_.transpose() * (p - off_); return result; } - + /** \brief Returns the domain dimension. Same as the ambient dimension of the sphere. */ - virtual std::size_t amb_d() const override {return d_;}; + virtual std::size_t amb_d() const override { return d_; }; /** \brief Returns the codomain dimension. Same as the codimension of the sphere. */ - virtual std::size_t cod_d() const override {return k_;}; + virtual std::size_t cod_d() const override { return k_; }; /** \brief Returns a point on the affine plane. */ virtual Eigen::VectorXd seed() const override { Eigen::VectorXd result = off_; return result; } - - /** + + /** * \brief Constructor of the function that defines an m-dimensional implicit affine * plane in the d-dimensional Euclidean space. * * @param[in] normal_matrix A normal matrix of the affine plane. The number of rows should - * correspond to the ambient dimension, the number of columns should corespond to - * the size of the normal basis (codimension). + * correspond to the ambient dimension, the number of columns should corespond to + * the size of the normal basis (codimension). * @param[in] offset The offset vector of the affine plane. * The dimension of the vector should be the ambient dimension of the manifold. */ - Function_affine_plane_in_Rd(const Eigen::MatrixXd& normal_matrix, - const Eigen::VectorXd& offset) - : normal_matrix_(normal_matrix), - d_(normal_matrix.rows()), - k_(normal_matrix.cols()), - m_(d_ - k_), - off_(offset) { + Function_affine_plane_in_Rd(const Eigen::MatrixXd& normal_matrix, const Eigen::VectorXd& offset) + : normal_matrix_(normal_matrix), d_(normal_matrix.rows()), k_(normal_matrix.cols()), m_(d_ - k_), off_(offset) { normal_matrix_.colwise().normalize(); } - /** + /** * \brief Constructor of the function that defines an m-dimensional implicit affine * plane in the d-dimensional Euclidean space that passes through origin. * * @param[in] normal_matrix A normal matrix of the affine plane. The number of rows should - * correspond to the ambient dimension, the number of columns should corespond to - * the size of the normal basis (codimension). + * correspond to the ambient dimension, the number of columns should corespond to + * the size of the normal basis (codimension). */ Function_affine_plane_in_Rd(const Eigen::MatrixXd& normal_matrix) - : normal_matrix_(normal_matrix), - d_(normal_matrix.rows()), - k_(normal_matrix.cols()), - m_(d_ - k_), - off_(Eigen::VectorXd::Zero(d_)) { + : normal_matrix_(normal_matrix), + d_(normal_matrix.rows()), + k_(normal_matrix.cols()), + m_(d_ - k_), + off_(Eigen::VectorXd::Zero(d_)) { normal_matrix_.colwise().normalize(); } -private: + private: Eigen::MatrixXd normal_matrix_; std::size_t d_, k_, m_; - Eigen::VectorXd off_; + Eigen::VectorXd off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_chair_in_R3.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_chair_in_R3.h index fe16d37b..e2f65144 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_chair_in_R3.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_chair_in_R3.h @@ -12,7 +12,7 @@ #define FUNCTIONS_FUNCTION_CHAIR_IN_R3_H_ #include <cstdlib> // for std::size_t -#include <cmath> // for std::pow +#include <cmath> // for std::pow #include <gudhi/Functions/Function.h> @@ -22,42 +22,42 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_chair_in_R3 +/** + * \class Function_chair_in_R3 * \brief A class that encodes the function, the zero-set of which is a so-called * "chair" surface embedded in R^3. * * \ingroup coxeter_triangulation */ struct Function_chair_in_R3 : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { - double x = p(0)-off_[0], y = p(1)-off_[1], z = p(2)-off_[2]; + double x = p(0) - off_[0], y = p(1) - off_[1], z = p(2) - off_[2]; Eigen::VectorXd result(cod_d()); - result(0) = std::pow(x*x + y*y + z*z - a_*k_*k_, 2) - b_*((z-k_)*(z-k_) - 2*x*x)*((z+k_)*(z+k_) - 2*y*y); + result(0) = std::pow(x * x + y * y + z * z - a_ * k_ * k_, 2) - + b_ * ((z - k_) * (z - k_) - 2 * x * x) * ((z + k_) * (z + k_) - 2 * y * y); return result; } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return 3;} + virtual std::size_t amb_d() const override { return 3; } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return 1;} + virtual std::size_t cod_d() const override { return 1; } /** \brief Returns a point on the surface. */ virtual Eigen::VectorXd seed() const override { - double t1 = a_-b_; - double discr = t1*t1 - (1.0 - b_)*(a_*a_ - b_); - double z0 = k_*std::sqrt((t1+std::sqrt(discr))/(1-b_)); - Eigen::Vector3d result(off_[0], off_[1], z0+off_[2]); + double t1 = a_ - b_; + double discr = t1 * t1 - (1.0 - b_) * (a_ * a_ - b_); + double z0 = k_ * std::sqrt((t1 + std::sqrt(discr)) / (1 - b_)); + Eigen::Vector3d result(off_[0], off_[1], z0 + off_[2]); return result; } - /** + /** * \brief Constructor of the function that defines the 'chair' surface * embedded in R^3. * @@ -66,20 +66,17 @@ struct Function_chair_in_R3 : public Function { * @param[in] k A numerical parameter. * @param[in] off Offset vector. */ - Function_chair_in_R3(double a = 0.8, - double b = 0.4, - double k = 1.0, - Eigen::Vector3d off = Eigen::Vector3d::Zero()) : - a_(a), b_(b), k_(k), off_(off) {} - -protected: + Function_chair_in_R3(double a = 0.8, double b = 0.4, double k = 1.0, Eigen::Vector3d off = Eigen::Vector3d::Zero()) + : a_(a), b_(b), k_(k), off_(off) {} + + protected: double a_, b_, k_; Eigen::Vector3d off_; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_iron_in_R3.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_iron_in_R3.h index fa5e4e66..5842c58f 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_iron_in_R3.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_iron_in_R3.h @@ -12,7 +12,7 @@ #define FUNCTIONS_FUNCTION_IRON_IN_R3_H_ #include <cstdlib> // for std::size_t -#include <cmath> // for std::pow +#include <cmath> // for std::pow #include <gudhi/Functions/Function.h> @@ -22,54 +22,52 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_iron_in_R3 +/** + * \class Function_iron_in_R3 * \brief A class that encodes the function, the zero-set of which is a surface * embedded in R^3 that ressembles an iron. * * \ingroup coxeter_triangulation */ struct Function_iron_in_R3 : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { double x = p(0), y = p(1), z = p(2); Eigen::VectorXd result(cod_d()); - result(0) = - std::pow(x,6)/300. - std::pow(y,6)/300. - std::pow(z,6)/300. + x*y*y*z/2.1 + y*y + std::pow(z-2, 4) - 1; + result(0) = -std::pow(x, 6) / 300. - std::pow(y, 6) / 300. - std::pow(z, 6) / 300. + x * y * y * z / 2.1 + y * y + + std::pow(z - 2, 4) - 1; return result; } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return 3;}; + virtual std::size_t amb_d() const override { return 3; }; /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return 1;}; + virtual std::size_t cod_d() const override { return 1; }; /** \brief Returns a point on the surface. */ virtual Eigen::VectorXd seed() const override { - Eigen::Vector3d result(std::pow(4500, 1./6), 0, 0); + Eigen::Vector3d result(std::pow(4500, 1. / 6), 0, 0); return result; } - - /** + + /** * \brief Constructor of the function that defines a surface embedded in R^3 * that ressembles an iron. * * @param[in] off Offset vector. */ - Function_iron_in_R3(Eigen::Vector3d off = Eigen::Vector3d::Zero()) : - off_(off) {} - + Function_iron_in_R3(Eigen::Vector3d off = Eigen::Vector3d::Zero()) : off_(off) {} + private: Eigen::Vector3d off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_lemniscate_revolution_in_R3.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_lemniscate_revolution_in_R3.h index b715932b..b5d54b24 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_lemniscate_revolution_in_R3.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_lemniscate_revolution_in_R3.h @@ -12,7 +12,7 @@ #define FUNCTIONS_FUNCTION_LEMNISCATE_REVOLUTION_IN_R3_H_ #include <cstdlib> // for std::size_t -#include <cmath> // for std::sqrt +#include <cmath> // for std::sqrt #include <gudhi/Functions/Function.h> @@ -22,71 +22,68 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_lemniscate_revolution_in_R3 +/** + * \class Function_lemniscate_revolution_in_R3 * \brief A class that encodes the function, the zero-set of which is a surface of revolution * around the x axis based on the lemniscate of Bernoulli embedded in R^3. * * \ingroup coxeter_triangulation */ struct Function_lemniscate_revolution_in_R3 : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { - double x = p(0)-off_[0], y = p(1)-off_[1], z = p(2)-off_[2]; + double x = p(0) - off_[0], y = p(1) - off_[1], z = p(2) - off_[2]; Eigen::VectorXd result(cod_d()); - double x2 = x*x, y2 = y*y, z2 = z*z, a2 = a_*a_; + double x2 = x * x, y2 = y * y, z2 = z * z, a2 = a_ * a_; double t1 = x2 + y2 + z2; - result(0) = t1*t1 - 2*a2*(x2 - y2 - z2); + result(0) = t1 * t1 - 2 * a2 * (x2 - y2 - z2); return result; } /** \brief Returns the (ambient) domain dimension.*/ - virtual std::size_t amb_d() const override {return 3;}; + virtual std::size_t amb_d() const override { return 3; }; /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return 1;}; + virtual std::size_t cod_d() const override { return 1; }; - /** \brief Returns a point on the surface. This seed point is only one of + /** \brief Returns a point on the surface. This seed point is only one of * two necessary seed points for the manifold tracing algorithm. - * See the method seed2() for the other point. + * See the method seed2() for the other point. */ virtual Eigen::VectorXd seed() const override { - Eigen::Vector3d result(std::sqrt(2*a_)+off_[0], off_[1], off_[2]); + Eigen::Vector3d result(std::sqrt(2 * a_) + off_[0], off_[1], off_[2]); return result; } - /** \brief Returns a point on the surface. This seed point is only one of - * two necessary seed points for the manifold tracing algorithm. - * See the method seed() for the other point. + /** \brief Returns a point on the surface. This seed point is only one of + * two necessary seed points for the manifold tracing algorithm. + * See the method seed() for the other point. */ Eigen::VectorXd seed2() const { - Eigen::Vector3d result(-std::sqrt(2*a_)+off_[0], off_[1], off_[2]); + Eigen::Vector3d result(-std::sqrt(2 * a_) + off_[0], off_[1], off_[2]); return result; } - - /** + + /** * \brief Constructor of the function that defines a surface of revolution * around the x axis based on the lemniscate of Bernoulli embedded in R^3. * * @param[in] a A numerical parameter. * @param[in] off Offset vector. */ - Function_lemniscate_revolution_in_R3(double a = 1, - Eigen::Vector3d off = Eigen::Vector3d::Zero()) - : a_(a), off_(off) {} + Function_lemniscate_revolution_in_R3(double a = 1, Eigen::Vector3d off = Eigen::Vector3d::Zero()) + : a_(a), off_(off) {} private: double a_; Eigen::Vector3d off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_moment_curve_in_Rd.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_moment_curve_in_Rd.h index eecbc914..3ea5cf12 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_moment_curve_in_Rd.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_moment_curve_in_Rd.h @@ -21,49 +21,45 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_moment_curve_in_Rd +/** + * \class Function_moment_curve_in_Rd * \brief A class for the function that defines an implicit moment curve * in the d-dimensional Euclidean space. * * \ingroup coxeter_triangulation */ struct Function_moment_curve_in_Rd : public Function { - -/** \brief Value of the function at a specified point. - * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. - */ + /** \brief Value of the function at a specified point. + * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. + */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { Eigen::VectorXd result(k_); - for (std::size_t i = 1; i < d_; ++i) - result(i-1) = p(i) - p(0) * p(i-1); + for (std::size_t i = 1; i < d_; ++i) result(i - 1) = p(i) - p(0) * p(i - 1); return result; } - + /** \brief Returns the domain (ambient) dimension.. */ - virtual std::size_t amb_d() const override {return d_;}; + virtual std::size_t amb_d() const override { return d_; }; /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return k_;}; + virtual std::size_t cod_d() const override { return k_; }; /** \brief Returns a point on the moment curve. */ virtual Eigen::VectorXd seed() const override { Eigen::VectorXd result = Eigen::VectorXd::Zero(d_); return result; } - - /** + + /** * \brief Constructor of the function that defines an implicit moment curve * in the d-dimensional Euclidean space. * * @param[in] r Numerical parameter. * @param[in] d The ambient dimension. */ - Function_moment_curve_in_Rd(double r, - std::size_t d) - : m_(1), k_(d-1), d_(d), r_(r) {} + Function_moment_curve_in_Rd(double r, std::size_t d) : m_(1), k_(d - 1), d_(d), r_(r) {} - /** + /** * \brief Constructor of the function that defines an implicit moment curve * in the d-dimensional Euclidean space. * @@ -71,20 +67,17 @@ struct Function_moment_curve_in_Rd : public Function { * @param[in] d The ambient dimension. * @param[in] offset The offset of the moment curve. */ - Function_moment_curve_in_Rd(double r, - std::size_t d, - Eigen::VectorXd& offset) - : m_(1), k_(d-1), d_(d), r_(r), off_(offset) {} - + Function_moment_curve_in_Rd(double r, std::size_t d, Eigen::VectorXd& offset) + : m_(1), k_(d - 1), d_(d), r_(r), off_(offset) {} + private: std::size_t m_, k_, d_; double r_; Eigen::VectorXd off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_torus_in_R3.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_torus_in_R3.h index c756c1a6..cdc39199 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_torus_in_R3.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_torus_in_R3.h @@ -12,7 +12,7 @@ #define FUNCTIONS_FUNCTION_TORUS_IN_R3_H_ #include <cstdlib> // for std::size_t -#include <cmath> // for std::sqrt +#include <cmath> // for std::sqrt #include <gudhi/Functions/Function.h> @@ -22,58 +22,54 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_torus_in_R3 +/** + * \class Function_torus_in_R3 * \brief A class that encodes the function, the zero-set of which is a torus * surface embedded in R^3. * * \ingroup coxeter_triangulation */ struct Function_torus_in_R3 : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { - double x = p(0)-off_[0], y = p(1)-off_[1], z = p(2)-off_[2]; + double x = p(0) - off_[0], y = p(1) - off_[1], z = p(2) - off_[2]; Eigen::VectorXd result(cod_d()); - result(0) = (z*z + (std::sqrt(x*x + y*y) - r_)*(std::sqrt(x*x + y*y) - r_) - R_*R_); + result(0) = (z * z + (std::sqrt(x * x + y * y) - r_) * (std::sqrt(x * x + y * y) - r_) - R_ * R_); return result; } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return 3;}; + virtual std::size_t amb_d() const override { return 3; }; /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return 1;}; + virtual std::size_t cod_d() const override { return 1; }; /** \brief Returns a point on the surface. */ virtual Eigen::VectorXd seed() const override { - Eigen::Vector3d result(R_ + r_ +off_[0], off_[1], off_[2]); + Eigen::Vector3d result(R_ + r_ + off_[0], off_[1], off_[2]); return result; } - - /** + + /** * \brief Constructor of the function that defines a torus embedded in R^3. * * @param[in] R The outer radius of the torus. * @param[in] r The inner radius of the torus. * @param[in] off Offset vector. */ - Function_torus_in_R3(double R = 1, - double r = 0.5, - Eigen::Vector3d off = Eigen::Vector3d::Zero()) : - R_(R), r_(r), off_(off) {} - + Function_torus_in_R3(double R = 1, double r = 0.5, Eigen::Vector3d off = Eigen::Vector3d::Zero()) + : R_(R), r_(r), off_(off) {} + private: double R_, r_; Eigen::Vector3d off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Function_whitney_umbrella_in_R3.h b/src/Coxeter_triangulation/include/gudhi/Functions/Function_whitney_umbrella_in_R3.h index 09306561..4415ab87 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Function_whitney_umbrella_in_R3.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Function_whitney_umbrella_in_R3.h @@ -21,65 +21,62 @@ namespace Gudhi { namespace coxeter_triangulation { -/** - * \class Function_whitney_umbrella_in_R3 +/** + * \class Function_whitney_umbrella_in_R3 * \brief A class that encodes the function, the zero-set of which is the Whitney umbrella * surface embedded in R^3. * * \ingroup coxeter_triangulation */ struct Function_whitney_umbrella_in_R3 : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ virtual Eigen::VectorXd operator()(const Eigen::VectorXd& p) const override { - double x = p(0)-off_[0], y = p(1)-off_[1], z = p(2)-off_[2]; + double x = p(0) - off_[0], y = p(1) - off_[1], z = p(2) - off_[2]; Eigen::VectorXd result(cod_d()); - result(0) = x*x - y*y*z; + result(0) = x * x - y * y * z; return result; } /** \brief Returns the (ambient) domain dimension.*/ - virtual std::size_t amb_d() const override {return 3;}; + virtual std::size_t amb_d() const override { return 3; }; /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return 1;}; + virtual std::size_t cod_d() const override { return 1; }; - /** \brief Returns a point on the surface. This seed point is only one of + /** \brief Returns a point on the surface. This seed point is only one of * two necessary seed points for the manifold tracing algorithm. - * See the method seed2() for the other point. + * See the method seed2() for the other point. */ virtual Eigen::VectorXd seed() const override { - Eigen::Vector3d result(1+off_[0], 1+off_[1], 1+off_[2]); + Eigen::Vector3d result(1 + off_[0], 1 + off_[1], 1 + off_[2]); return result; } - /** \brief Returns a point on the surface. This seed point is only one of - * two necessary seed points for the manifold tracing algorithm. - * See the method seed() for the other point. + /** \brief Returns a point on the surface. This seed point is only one of + * two necessary seed points for the manifold tracing algorithm. + * See the method seed() for the other point. */ Eigen::VectorXd seed2() const { - Eigen::Vector3d result(-1+off_[0], -1+off_[1], 1+off_[2]); + Eigen::Vector3d result(-1 + off_[0], -1 + off_[1], 1 + off_[2]); return result; } - - /** + + /** * \brief Constructor of the function that defines the Whitney umbrella in R^3. * * @param[in] off Offset vector. */ - Function_whitney_umbrella_in_R3(Eigen::Vector3d off = Eigen::Vector3d::Zero()) - : off_(off) {} + Function_whitney_umbrella_in_R3(Eigen::Vector3d off = Eigen::Vector3d::Zero()) : off_(off) {} private: Eigen::Vector3d off_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Linear_transformation.h b/src/Coxeter_triangulation/include/gudhi/Functions/Linear_transformation.h index a5d17fe0..1a82984b 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Linear_transformation.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Linear_transformation.h @@ -25,15 +25,14 @@ namespace coxeter_triangulation { * \brief Transforms the zero-set of the function by a given linear transformation. * The underlying function corresponds to f(M*x), where M is the transformation matrix. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class Function_> struct Linear_transformation : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -43,10 +42,10 @@ struct Linear_transformation : public Function { } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return fun_.amb_d();} + virtual std::size_t amb_d() const override { return fun_.amb_d(); } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return fun_.cod_d();} + virtual std::size_t cod_d() const override { return fun_.cod_d(); } /** \brief Returns a point on the zero-set. */ virtual Eigen::VectorXd seed() const override { @@ -55,43 +54,39 @@ struct Linear_transformation : public Function { return result; } - /** + /** * \brief Constructor of a linearly transformed function. * * @param[in] function The function to be linearly transformed. * @param[in] matrix The transformation matrix. Its dimension should be d*d, * where d is the domain (ambient) dimension of 'function'. */ - Linear_transformation(const Function_& function, const Eigen::MatrixXd& matrix) : - fun_(function), matrix_(matrix) { - } + Linear_transformation(const Function_& function, const Eigen::MatrixXd& matrix) : fun_(function), matrix_(matrix) {} private: Function_ fun_; Eigen::MatrixXd matrix_; }; - -/** +/** * \brief Static constructor of a linearly transformed function. * * @param[in] function The function to be linearly transformed. * @param[in] matrix The transformation matrix. Its dimension should be d*d, * where d is the domain (ambient) dimension of 'function'. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class Function_> -Linear_transformation<Function_> make_linear_transformation(const Function_& function, - const Eigen::MatrixXd& matrix) { - return Linear_transformation<Function_>(function, matrix); +Linear_transformation<Function_> make_linear_transformation(const Function_& function, const Eigen::MatrixXd& matrix) { + return Linear_transformation<Function_>(function, matrix); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Negation.h b/src/Coxeter_triangulation/include/gudhi/Functions/Negation.h index 3439dbad..58d8ebe0 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Negation.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Negation.h @@ -21,20 +21,19 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** *\class Negation * \brief Constructs the "minus" function. The zero-set is the same, but * the values at other points are the negative of their original value. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class Function_> struct Negation : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -44,10 +43,10 @@ struct Negation : public Function { } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return fun_.amb_d();} + virtual std::size_t amb_d() const override { return fun_.amb_d(); } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return fun_.cod_d();} + virtual std::size_t cod_d() const override { return fun_.cod_d(); } /** \brief Returns a point on the zero-set. */ virtual Eigen::VectorXd seed() const override { @@ -55,27 +54,24 @@ struct Negation : public Function { return result; } - /** + /** * \brief Constructor of the negative function. * * @param[in] function The function to be negated. */ - Negation(const Function_& function) : - fun_(function) { - } + Negation(const Function_& function) : fun_(function) {} private: Function_ fun_; }; - -/** +/** * \brief Static constructor of the negative function. * * @param[in] function The function to be translated. * domain (ambient) dimension of 'function'. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation @@ -85,8 +81,8 @@ Negation<Function_> negation(const Function_& function) { return Negation<Function_>(function); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/PL_approximation.h b/src/Coxeter_triangulation/include/gudhi/Functions/PL_approximation.h index ab2b9294..3d1f3564 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/PL_approximation.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/PL_approximation.h @@ -21,23 +21,21 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** * \class PL_approximation * \brief Constructs a piecewise-linear approximation of a function induced by * an ambient triangulation. * - * \tparam Function The function template parameter. Should be a model of + * \tparam Function The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * \tparam Triangulation The triangulation template parameter. Should be a model of * the concept TriangulationForManifoldTracing. * * \ingroup coxeter_triangulation */ -template <class Function_, - class Triangulation_> +template <class Function_, class Triangulation_> struct PL_approximation : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -47,33 +45,29 @@ struct PL_approximation : public Function { auto s = tr_.locate_point(p); Eigen::MatrixXd matrix(cod_d, s.dimension() + 1); Eigen::MatrixXd vertex_matrix(amb_d + 1, s.dimension() + 1); - for (std::size_t i = 0; i < s.dimension() + 1; ++i) - vertex_matrix(0, i) = 1; + for (std::size_t i = 0; i < s.dimension() + 1; ++i) vertex_matrix(0, i) = 1; std::size_t j = 0; - for (auto v: s.vertex_range()) { + for (auto v : s.vertex_range()) { Eigen::VectorXd pt_v = tr_.cartesian_coordinates(v); Eigen::VectorXd fun_v = fun_(pt_v); - for (std::size_t i = 1; i < amb_d + 1; ++i) - vertex_matrix(i, j) = pt_v(i-1); - for (std::size_t i = 0; i < cod_d; ++i) - matrix(i, j) = fun_v(i); + for (std::size_t i = 1; i < amb_d + 1; ++i) vertex_matrix(i, j) = pt_v(i - 1); + for (std::size_t i = 0; i < cod_d; ++i) matrix(i, j) = fun_v(i); j++; } - assert(j == s.dimension()+1); + assert(j == s.dimension() + 1); Eigen::VectorXd z(amb_d + 1); z(0) = 1; - for (std::size_t i = 1; i < amb_d + 1; ++i) - z(i) = p(i-1); - Eigen::VectorXd lambda = vertex_matrix.colPivHouseholderQr().solve(z); + for (std::size_t i = 1; i < amb_d + 1; ++i) z(i) = p(i - 1); + Eigen::VectorXd lambda = vertex_matrix.colPivHouseholderQr().solve(z); Eigen::VectorXd result = matrix * lambda; return result; } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return fun_.amb_d();} + virtual std::size_t amb_d() const override { return fun_.amb_d(); } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return fun_.cod_d();} + virtual std::size_t cod_d() const override { return fun_.cod_d(); } /** \brief Returns a point on the zero-set. */ virtual Eigen::VectorXd seed() const override { @@ -81,44 +75,41 @@ struct PL_approximation : public Function { return Eigen::VectorXd(amb_d()); } - /** - * \brief Constructor of the piecewise-linear approximation of a function + /** + * \brief Constructor of the piecewise-linear approximation of a function * induced by an ambient triangulation. * * @param[in] function The function. * @param[in] triangulation The ambient triangulation. */ PL_approximation(const Function_& function, const Triangulation_& triangulation) - : fun_(function), tr_(triangulation) {} + : fun_(function), tr_(triangulation) {} - private: + private: Function_ fun_; Triangulation_ tr_; }; - -/** - * \brief Static constructor of the piecewise-linear approximation of a function +/** + * \brief Static constructor of the piecewise-linear approximation of a function * induced by an ambient triangulation. * * @param[in] function The function. * @param[in] triangulation The ambient triangulation. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ -template <class Function_, - class Triangulation_> -PL_approximation<Function_, Triangulation_> -make_pl_approximation(const Function_& function, - const Triangulation_& triangulation) { +template <class Function_, class Triangulation_> +PL_approximation<Function_, Triangulation_> make_pl_approximation(const Function_& function, + const Triangulation_& triangulation) { return PL_approximation<Function_, Triangulation_>(function, triangulation); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/Translate.h b/src/Coxeter_triangulation/include/gudhi/Functions/Translate.h index b3f53fe0..bd91c707 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/Translate.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/Translate.h @@ -21,20 +21,19 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** * \class Translate * \brief Translates the zero-set of the function by a vector. * The underlying function corresponds to f(x-off), where off is the offset vector. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ template <class Function_> struct Translate : public Function { - - /** + /** * \brief Value of the function at a specified point. * @param[in] p The input point. The dimension needs to coincide with the ambient dimension. */ @@ -44,10 +43,10 @@ struct Translate : public Function { } /** \brief Returns the domain (ambient) dimension. */ - virtual std::size_t amb_d() const override {return fun_.amb_d();} + virtual std::size_t amb_d() const override { return fun_.amb_d(); } /** \brief Returns the codomain dimension. */ - virtual std::size_t cod_d() const override {return fun_.cod_d();} + virtual std::size_t cod_d() const override { return fun_.cod_d(); } /** \brief Returns a point on the zero-set. */ virtual Eigen::VectorXd seed() const override { @@ -56,31 +55,28 @@ struct Translate : public Function { return result; } - /** + /** * \brief Constructor of the translated function. * * @param[in] function The function to be translated. - * @param[in] off The offset vector. The dimension should correspond to the + * @param[in] off The offset vector. The dimension should correspond to the * domain (ambient) dimension of 'function'. */ - Translate(const Function_& function, const Eigen::VectorXd& off) : - fun_(function), off_(off) { - } + Translate(const Function_& function, const Eigen::VectorXd& off) : fun_(function), off_(off) {} private: Function_ fun_; Eigen::VectorXd off_; }; - -/** +/** * \brief Static constructor of a translated function. * * @param[in] function The function to be translated. - * @param[in] off The offset vector. The dimension should correspond to the + * @param[in] off The offset vector. The dimension should correspond to the * domain (ambient) dimension of 'function'. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation @@ -90,8 +86,8 @@ Translate<Function_> translate(const Function_& function, Eigen::VectorXd off) { return Translate<Function_>(function, off); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Functions/random_orthogonal_matrix.h b/src/Coxeter_triangulation/include/gudhi/Functions/random_orthogonal_matrix.h index 34fc1a67..fbc1c895 100644 --- a/src/Coxeter_triangulation/include/gudhi/Functions/random_orthogonal_matrix.h +++ b/src/Coxeter_triangulation/include/gudhi/Functions/random_orthogonal_matrix.h @@ -12,8 +12,8 @@ #define FUNCTIONS_RANDOM_ORTHOGONAL_MATRIX_H_ #include <cstdlib> // for std::size_t -#include <cmath> // for std::cos, std::sin -#include <random> // for std::uniform_real_distribution, std::random_device +#include <cmath> // for std::cos, std::sin +#include <random> // for std::uniform_real_distribution, std::random_device #include <gudhi/math.h> @@ -28,10 +28,10 @@ namespace Gudhi { namespace coxeter_triangulation { -/** \brief Generates a uniform random orthogonal matrix using the "subgroup algorithm" by - * Diaconis & Shashahani. +/** \brief Generates a uniform random orthogonal matrix using the "subgroup algorithm" by + * Diaconis & Shashahani. * \details Taken from https://en.wikipedia.org/wiki/Rotation_matrix#Uniform_random_rotation_matrices. - * The idea: take a random rotation matrix of dimension d-1, embed it + * The idea: take a random rotation matrix of dimension d-1, embed it * as a d*d matrix M with the last column (0,...,0,1). * Pick a random vector v on a sphere S^d. rotate the matrix M so that its last column is v. * The determinant of the matrix can be either 1 or -1 @@ -41,8 +41,7 @@ namespace coxeter_triangulation { Eigen::MatrixXd random_orthogonal_matrix(std::size_t d) { typedef CGAL::Epick_d<CGAL::Dynamic_dimension_tag> Kernel; typedef typename Kernel::Point_d Point_d; - if (d == 1) - return Eigen::VectorXd::Constant(1, 1.0); + if (d == 1) return Eigen::VectorXd::Constant(1, 1.0); if (d == 2) { // 0. < alpha < 2 Pi std::uniform_real_distribution<double> unif(0., 2 * Gudhi::PI); @@ -54,22 +53,20 @@ Eigen::MatrixXd random_orthogonal_matrix(std::size_t d) { rot << std::cos(alpha), -std::sin(alpha), std::sin(alpha), cos(alpha); return rot; } - Eigen::MatrixXd low_dim_rot = random_orthogonal_matrix(d-1); - Eigen::MatrixXd rot(d,d); + Eigen::MatrixXd low_dim_rot = random_orthogonal_matrix(d - 1); + Eigen::MatrixXd rot(d, d); Point_d v = *CGAL::Random_points_on_sphere_d<Point_d>(d, 1); - for (std::size_t i = 0; i < d; ++i) - rot(i,0) = v[i]; - for (std::size_t i = 0; i < d-1; ++i) - for (std::size_t j = 1; j < d-1; ++j) - rot(i,j) = low_dim_rot(i,j-1); - for (std::size_t j = 1; j < d; ++j) - rot(d-1,j) = 0; - rot = rot.householderQr().householderQ(); // a way to do Gram-Schmidt, see https://forum.kde.org/viewtopic.php?f=74&t=118568#p297246 + for (std::size_t i = 0; i < d; ++i) rot(i, 0) = v[i]; + for (std::size_t i = 0; i < d - 1; ++i) + for (std::size_t j = 1; j < d - 1; ++j) rot(i, j) = low_dim_rot(i, j - 1); + for (std::size_t j = 1; j < d; ++j) rot(d - 1, j) = 0; + rot = rot.householderQr() + .householderQ(); // a way to do Gram-Schmidt, see https://forum.kde.org/viewtopic.php?f=74&t=118568#p297246 return rot; } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/IO/Mesh_medit.h b/src/Coxeter_triangulation/include/gudhi/IO/Mesh_medit.h index 28276bf7..203ed6c5 100644 --- a/src/Coxeter_triangulation/include/gudhi/IO/Mesh_medit.h +++ b/src/Coxeter_triangulation/include/gudhi/IO/Mesh_medit.h @@ -11,15 +11,15 @@ #ifndef IO_MESH_MEDIT_H_ #define IO_MESH_MEDIT_H_ -namespace Gudhi { - -namespace coxeter_triangulation { - #include <Eigen/Dense> #include <vector> #include <utility> // for std::pair +namespace Gudhi { + +namespace coxeter_triangulation { + /* \class Mesh_medit * \brief Structure to store a mesh that can be output in Medit .mesh file format * using the output_meshes_to_medit method. @@ -53,8 +53,8 @@ struct Mesh_medit { Scalar_field_range tetrahedra_scalar_range; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/IO/build_mesh_from_cell_complex.h b/src/Coxeter_triangulation/include/gudhi/IO/build_mesh_from_cell_complex.h index abe6cdbf..0dc98f95 100644 --- a/src/Coxeter_triangulation/include/gudhi/IO/build_mesh_from_cell_complex.h +++ b/src/Coxeter_triangulation/include/gudhi/IO/build_mesh_from_cell_complex.h @@ -19,7 +19,7 @@ #include <map> #include <set> #include <string> -#include <utility> // for std::make_pair +#include <utility> // for std::make_pair #include <algorithm> // for std::min namespace Gudhi { @@ -27,104 +27,90 @@ namespace Gudhi { namespace coxeter_triangulation { struct Configuration { - Configuration(bool t_edges, bool t_triangles, bool t_tetrahedra, - std::size_t r_edges, std::size_t r_triangles, std::size_t r_tetrahedra) - : toggle_edges(t_edges), toggle_triangles(t_triangles), toggle_tetrahedra(t_tetrahedra), - ref_edges(r_edges), ref_triangles(r_triangles), ref_tetrahedra(r_tetrahedra) {} + Configuration(bool t_edges, bool t_triangles, bool t_tetrahedra, std::size_t r_edges, std::size_t r_triangles, + std::size_t r_tetrahedra) + : toggle_edges(t_edges), + toggle_triangles(t_triangles), + toggle_tetrahedra(t_tetrahedra), + ref_edges(r_edges), + ref_triangles(r_triangles), + ref_tetrahedra(r_tetrahedra) {} Configuration() {} - bool toggle_edges = true, - toggle_triangles = true, - toggle_tetrahedra = true; - std::size_t ref_edges = 1, - ref_triangles = 1, - ref_tetrahedra = 1; + bool toggle_edges = true, toggle_triangles = true, toggle_tetrahedra = true; + std::size_t ref_edges = 1, ref_triangles = 1, ref_tetrahedra = 1; }; - template <class Hasse_cell, class Simplex_cell_map> -void populate_mesh(Mesh_medit& output, - Simplex_cell_map& sc_map, - Configuration configuration, - std::size_t amb_d, +void populate_mesh(Mesh_medit& output, Simplex_cell_map& sc_map, Configuration configuration, std::size_t amb_d, std::map<Hasse_cell*, std::size_t> vi_map) { using Mesh_element_vertices = Mesh_medit::Mesh_elements::value_type::first_type; std::map<Hasse_cell*, std::size_t> ci_map; - std::size_t index = vi_map.size() + 1; // current size of output.vertex_points - if (sc_map.size() >= 3) - for (const auto& sc_pair: sc_map[2]) { + std::size_t index = vi_map.size() + 1; // current size of output.vertex_points + if (sc_map.size() >= 3) + for (const auto& sc_pair : sc_map[2]) { Eigen::VectorXd barycenter = Eigen::VectorXd::Zero(amb_d); std::set<std::size_t> vertex_indices; Hasse_cell* cell = sc_pair.second; - for (const auto& ei_pair: cell->get_boundary()) - for (const auto& vi_pair: ei_pair.first->get_boundary()) - vertex_indices.emplace(vi_map[vi_pair.first]); - for (const std::size_t& v: vertex_indices) - barycenter += output.vertex_points[v-1]; + for (const auto& ei_pair : cell->get_boundary()) + for (const auto& vi_pair : ei_pair.first->get_boundary()) vertex_indices.emplace(vi_map[vi_pair.first]); + for (const std::size_t& v : vertex_indices) barycenter += output.vertex_points[v - 1]; ci_map.emplace(std::make_pair(cell, index++)); - output.vertex_points.emplace_back((1./vertex_indices.size()) * barycenter); + output.vertex_points.emplace_back((1. / vertex_indices.size()) * barycenter); #ifdef GUDHI_COX_OUTPUT_TO_HTML - std::string vlist = " (" + std::to_string(index-1) + ")"; - for (const std::size_t& v: vertex_indices) - vlist += " " + std::to_string(v); + std::string vlist = " (" + std::to_string(index - 1) + ")"; + for (const std::size_t& v : vertex_indices) vlist += " " + std::to_string(v); cell_vlist_map.emplace(std::make_pair(to_string(cell), vlist)); #endif } if (configuration.toggle_edges && sc_map.size() >= 2) - for (const auto& sc_map: sc_map[1]) { + for (const auto& sc_map : sc_map[1]) { Hasse_cell* edge_cell = sc_map.second; Mesh_element_vertices edge; - for (const auto& vi_pair: edge_cell->get_boundary()) - edge.push_back(vi_map[vi_pair.first]); + for (const auto& vi_pair : edge_cell->get_boundary()) edge.push_back(vi_map[vi_pair.first]); output.edges.emplace_back(std::make_pair(edge, configuration.ref_edges)); #ifdef GUDHI_COX_OUTPUT_TO_HTML std::string vlist; - for (const std::size_t& v: edge) - vlist += " " + std::to_string(v); + for (const std::size_t& v : edge) vlist += " " + std::to_string(v); cell_vlist_map.emplace(std::make_pair(to_string(edge_cell), vlist)); #endif } - + if (configuration.toggle_triangles && sc_map.size() >= 3) - for (const auto& sc_pair: sc_map[2]) { - for (const auto& ei_pair: sc_pair.second->get_boundary()) { + for (const auto& sc_pair : sc_map[2]) { + for (const auto& ei_pair : sc_pair.second->get_boundary()) { Mesh_element_vertices triangle(1, ci_map[sc_pair.second]); - for (const auto& vi_pair: ei_pair.first->get_boundary()) - triangle.push_back(vi_map[vi_pair.first]); + for (const auto& vi_pair : ei_pair.first->get_boundary()) triangle.push_back(vi_map[vi_pair.first]); output.triangles.emplace_back(std::make_pair(triangle, configuration.ref_triangles)); } } - + if (configuration.toggle_tetrahedra && sc_map.size() >= 4) - for (const auto& sc_pair: sc_map[3]) { + for (const auto& sc_pair : sc_map[3]) { Eigen::VectorXd barycenter = Eigen::VectorXd::Zero(amb_d); std::set<std::size_t> vertex_indices; Hasse_cell* cell = sc_pair.second; - for (const auto& ci_pair: cell->get_boundary()) - for (const auto& ei_pair: ci_pair.first->get_boundary()) - for (const auto& vi_pair: ei_pair.first->get_boundary()) - vertex_indices.emplace(vi_map[vi_pair.first]); - for (const std::size_t& v: vertex_indices) - barycenter += output.vertex_points[v-1]; - output.vertex_points.emplace_back((1./vertex_indices.size()) * barycenter); + for (const auto& ci_pair : cell->get_boundary()) + for (const auto& ei_pair : ci_pair.first->get_boundary()) + for (const auto& vi_pair : ei_pair.first->get_boundary()) vertex_indices.emplace(vi_map[vi_pair.first]); + for (const std::size_t& v : vertex_indices) barycenter += output.vertex_points[v - 1]; + output.vertex_points.emplace_back((1. / vertex_indices.size()) * barycenter); #ifdef GUDHI_COX_OUTPUT_TO_HTML std::string vlist = " (" + std::to_string(index) + ")"; - for (const std::size_t& v: vertex_indices) - vlist += " " + std::to_string(v); + for (const std::size_t& v : vertex_indices) vlist += " " + std::to_string(v); cell_vlist_map.emplace(std::make_pair(to_string(cell), vlist)); #endif - for (const auto& ci_pair: cell->get_boundary()) - for (const auto& ei_pair: ci_pair.first->get_boundary()) { + for (const auto& ci_pair : cell->get_boundary()) + for (const auto& ei_pair : ci_pair.first->get_boundary()) { Mesh_element_vertices tetrahedron = {index, ci_map[sc_pair.second]}; - for (const auto& vi_pair: ei_pair.first->get_boundary()) - tetrahedron.push_back(vi_map[vi_pair.first]); + for (const auto& vi_pair : ei_pair.first->get_boundary()) tetrahedron.push_back(vi_map[vi_pair.first]); output.tetrahedra.emplace_back(std::make_pair(tetrahedron, configuration.ref_tetrahedra)); } index++; - } + } } template <class Cell_complex> @@ -133,14 +119,13 @@ Mesh_medit build_mesh_from_cell_complex(const Cell_complex& cell_complex, Configuration b_configuration = Configuration()) { using Hasse_cell = typename Cell_complex::Hasse_cell; Mesh_medit output; - std::map<Hasse_cell*, std::size_t> vi_map; // one for vertices, other for 2d-cells - std::size_t index = 1; // current size of output.vertex_points - - if (cell_complex.cell_point_map().empty()) - return output; - std::size_t amb_d = std::min((int) cell_complex.cell_point_map().begin()->second.size(), 3); - - for (const auto& cp_pair: cell_complex.cell_point_map()) { + std::map<Hasse_cell*, std::size_t> vi_map; // one for vertices, other for 2d-cells + std::size_t index = 1; // current size of output.vertex_points + + if (cell_complex.cell_point_map().empty()) return output; + std::size_t amb_d = std::min((int)cell_complex.cell_point_map().begin()->second.size(), 3); + + for (const auto& cp_pair : cell_complex.cell_point_map()) { #ifdef GUDHI_COX_OUTPUT_TO_HTML std::string vlist; vlist += " " + std::to_string(index); @@ -153,29 +138,29 @@ Mesh_medit build_mesh_from_cell_complex(const Cell_complex& cell_complex, populate_mesh(output, cell_complex.interior_simplex_cell_maps(), i_configuration, amb_d, vi_map); #ifdef GUDHI_COX_OUTPUT_TO_HTML - for (const auto& sc_map: cell_complex.interior_simplex_cell_maps()) - for (const auto& sc_pair: sc_map) { + for (const auto& sc_map : cell_complex.interior_simplex_cell_maps()) + for (const auto& sc_pair : sc_map) { std::string simplex = "I" + to_string(sc_pair.first); std::string cell = to_string(sc_pair.second); std::string vlist = cell_vlist_map.at(cell).substr(1); simplex_vlist_map.emplace(std::make_pair(simplex, vlist)); } -#endif - populate_mesh(output, cell_complex.boundary_simplex_cell_maps(), b_configuration, amb_d, vi_map); +#endif + populate_mesh(output, cell_complex.boundary_simplex_cell_maps(), b_configuration, amb_d, vi_map); #ifdef GUDHI_COX_OUTPUT_TO_HTML - for (const auto& sc_map: cell_complex.boundary_simplex_cell_maps()) - for (const auto& sc_pair: sc_map) { + for (const auto& sc_map : cell_complex.boundary_simplex_cell_maps()) + for (const auto& sc_pair : sc_map) { std::string simplex = "B" + to_string(sc_pair.first); std::string cell = to_string(sc_pair.second); std::string vlist = cell_vlist_map.at(cell).substr(1); simplex_vlist_map.emplace(std::make_pair(simplex, vlist)); } -#endif +#endif return output; } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/IO/output_debug_traces_to_html.h b/src/Coxeter_triangulation/include/gudhi/IO/output_debug_traces_to_html.h index 43cca967..147ab908 100644 --- a/src/Coxeter_triangulation/include/gudhi/IO/output_debug_traces_to_html.h +++ b/src/Coxeter_triangulation/include/gudhi/IO/output_debug_traces_to_html.h @@ -16,22 +16,21 @@ namespace coxeter_triangulation { template <class T> std::ostream& operator<<(std::ostream& os, const std::vector<T>& vector) { - os << "("; + os << "("; if (vector.empty()) { os << ")"; return os; } auto v_it = vector.begin(); os << *v_it++; - for (; v_it != vector.end(); ++v_it) - os << ", " << *v_it; + for (; v_it != vector.end(); ++v_it) os << ", " << *v_it; os << ")"; return os; } /* A class to make the vector horizontal instead of vertical */ struct Straighten { - Straighten (const Eigen::VectorXd& vector) : vector_(vector) {} + Straighten(const Eigen::VectorXd& vector) : vector_(vector) {} const Eigen::VectorXd& vector_; }; @@ -42,8 +41,7 @@ std::ostream& operator<<(std::ostream& os, const Straighten& str) { os << ")"; return os; } - for (std::size_t i = 1; i < size; ++i) - os << ", " << str.vector_(i); + for (std::size_t i = 1; i < size; ++i) os << ", " << str.vector_(i); os << ")"; return os; } @@ -68,61 +66,57 @@ std::string to_string(const T& t) { struct MT_inserted_info { std::string qr_face_, init_face_, qr_intersection_; bool qr_success_, is_boundary_; - template <class Query_result, - class Simplex_handle> + template <class Query_result, class Simplex_handle> MT_inserted_info(const Query_result& qr, const Simplex_handle& face, bool is_boundary) - : qr_face_(to_string(face)), init_face_(to_string(face)), - qr_intersection_(to_string(qr.intersection)), - qr_success_(qr.success), is_boundary_(is_boundary) {} + : qr_face_(to_string(face)), + init_face_(to_string(face)), + qr_intersection_(to_string(qr.intersection)), + qr_success_(qr.success), + is_boundary_(is_boundary) {} }; std::list<MT_inserted_info> mt_seed_inserted_list, mt_inserted_list; struct CC_summary_info { std::string face_, cell_; template <class SC_pair> - CC_summary_info(const SC_pair& sc_pair) - : face_(to_string(sc_pair.first)), cell_(to_string(sc_pair.second)) {} + CC_summary_info(const SC_pair& sc_pair) : face_(to_string(sc_pair.first)), cell_(to_string(sc_pair.second)) {} }; using CC_summary_list = std::list<CC_summary_info>; std::vector<CC_summary_list> cc_interior_summary_lists, cc_boundary_summary_lists; struct CC_detail_info { - enum class Result_type {self, face, coface, inserted, join_single, join_is_face}; + enum class Result_type { self, face, coface, inserted, join_single, join_is_face }; std::string simplex_, trigger_, init_simplex_; Result_type status_; bool join_trigger_ = false; std::list<std::string> faces_, post_faces_, cofaces_; template <class Simplex_handle> - CC_detail_info(const Simplex_handle& simplex) - : simplex_(to_string(simplex)) {} + CC_detail_info(const Simplex_handle& simplex) : simplex_(to_string(simplex)) {} }; using CC_detail_list = std::list<CC_detail_info>; std::vector<CC_detail_list> cc_interior_detail_lists, cc_boundary_detail_lists; std::vector<CC_detail_list> cc_interior_insert_detail_lists, cc_boundary_insert_detail_lists; struct CC_prejoin_info { - enum class Result_type {join_single, join_is_face, join_different, join_same}; + enum class Result_type { join_single, join_is_face, join_different, join_same }; std::string simplex_, join_; std::vector<std::string> faces_; std::size_t dimension_; Result_type status_; template <class Simplex_handle> - CC_prejoin_info(const Simplex_handle& simplex) - : simplex_(to_string(simplex)), dimension_(simplex.dimension()) {} + CC_prejoin_info(const Simplex_handle& simplex) : simplex_(to_string(simplex)), dimension_(simplex.dimension()) {} }; using CC_prejoin_list = std::list<CC_prejoin_info>; std::vector<CC_prejoin_list> cc_interior_prejoin_lists, cc_boundary_prejoin_lists; - struct CC_join_info { - enum class Result_type {self, face, coface, inserted, join_single, join_is_face}; + enum class Result_type { self, face, coface, inserted, join_single, join_is_face }; std::string simplex_, join_, trigger_; Result_type status_; std::list<std::string> boundary_faces_; std::list<std::string> faces_, post_faces_, cofaces_; template <class Simplex_handle> - CC_join_info(const Simplex_handle& simplex) - : simplex_(to_string(simplex)) {} + CC_join_info(const Simplex_handle& simplex) : simplex_(to_string(simplex)) {} }; bool join_switch = false; std::vector<CC_detail_list> cc_interior_join_detail_lists, cc_boundary_join_detail_lists; @@ -134,16 +128,15 @@ std::ostringstream mt_ostream, vis_ostream; std::vector<std::ostringstream> cc_summary_ostream, cc_traces_ostream; std::string simplex_format(const std::string& simplex, bool is_boundary) { - std::string b_simplex = (is_boundary? "B": "I") + simplex; + std::string b_simplex = (is_boundary ? "B" : "I") + simplex; std::string tooltiptext; auto it = simplex_vlist_map.find(b_simplex); if (it == simplex_vlist_map.end()) tooltiptext = "deleted"; else tooltiptext = simplex_vlist_map.at(b_simplex); - return (std::string)"<a class=\"" + (is_boundary? "boundary": "interior") - + "\" href=\"#" + id_from_simplex(b_simplex) + "\">" + b_simplex - + "<span class=\"tooltiptext\">" + tooltiptext + "</span></a>"; + return (std::string) "<a class=\"" + (is_boundary ? "boundary" : "interior") + "\" href=\"#" + + id_from_simplex(b_simplex) + "\">" + b_simplex + "<span class=\"tooltiptext\">" + tooltiptext + "</span></a>"; } std::string simplex_format(const std::string& b_simplex) { @@ -154,12 +147,10 @@ std::string simplex_format(const std::string& b_simplex) { tooltiptext = "deleted"; else tooltiptext = simplex_vlist_map.at(b_simplex); - return (std::string)"<a class=\"" + (is_boundary? "boundary": "interior") - + "\" href=\"#" + id_from_simplex(b_simplex) + "\">" + b_simplex - + "<span class=\"tooltiptext\">" + tooltiptext + "</span></a>"; + return (std::string) "<a class=\"" + (is_boundary ? "boundary" : "interior") + "\" href=\"#" + + id_from_simplex(b_simplex) + "\">" + b_simplex + "<span class=\"tooltiptext\">" + tooltiptext + "</span></a>"; } - void write_head(std::ofstream& ofs) { ofs << " <head>\n" << " <title>Cell complex debug trace</title>\n" @@ -215,11 +206,11 @@ void write_head(std::ofstream& ofs) { << " .boundary:hover .tooltiptext {\n" << " visibility: visible;\n" << " opacity: 1;\n" - << " }\n" + << " }\n" << " .interior:hover .tooltiptext {\n" << " visibility: visible;\n" << " opacity: 1;\n" - << " }\n" + << " }\n" << " ul.nav {\n" << " list-style-type: none;\n" << " margin: 0;\n" @@ -284,30 +275,26 @@ void write_mt(std::ofstream& ofs) { ofs << " <h2> Manifold debug trace </h2>\n"; ofs << " <h3> Simplices inserted during the seed phase </h3>\n"; ofs << " <ul>\n"; - for (const MT_inserted_info& mt_info: mt_seed_inserted_list) { + for (const MT_inserted_info& mt_info : mt_seed_inserted_list) { if (mt_info.qr_success_) { ofs << " <li>Inserted " << simplex_format(mt_info.qr_face_, mt_info.is_boundary_); if (mt_info.qr_face_ != mt_info.init_face_) - ofs << " (initially " << simplex_format(mt_info.init_face_, mt_info.is_boundary_) << ")"; + ofs << " (initially " << simplex_format(mt_info.init_face_, mt_info.is_boundary_) << ")"; ofs << " intersection point is " << mt_info.qr_intersection_ << "</li>\n"; - } - else - ofs << " <li>Failed to insert " << mt_info.init_face_ - << "</li>\n"; + } else + ofs << " <li>Failed to insert " << mt_info.init_face_ << "</li>\n"; } ofs << " </ul>\n"; ofs << " <h3> Simplices inserted during the while loop phase </h3>\n"; ofs << " <ul>\n"; - for (const MT_inserted_info& mt_info: mt_inserted_list) { + for (const MT_inserted_info& mt_info : mt_inserted_list) { if (mt_info.qr_success_) { ofs << " <li>Inserted " << simplex_format(mt_info.qr_face_, mt_info.is_boundary_); if (mt_info.qr_face_ != mt_info.init_face_) - ofs << " (initially " << simplex_format(mt_info.init_face_, mt_info.is_boundary_) << ")"; + ofs << " (initially " << simplex_format(mt_info.init_face_, mt_info.is_boundary_) << ")"; ofs << " intersection point is " << mt_info.qr_intersection_ << "</li>\n"; - } - else - ofs << " <li>Failed to insert " << mt_info.init_face_ - << ")</li>\n"; + } else + ofs << " <li>Failed to insert " << mt_info.init_face_ << ")</li>\n"; } ofs << " </ul>\n"; ofs << " </div>\n"; @@ -325,13 +312,11 @@ void write_cc(std::ofstream& ofs) { for (std::size_t i = 0; i < cc_interior_summary_lists.size(); ++i) { ofs << " <h3 id=\"dim" << i << "\"> Dimension " << i << "</h3>\n"; ofs << " <h4 id=\"dim" << i << "i\"> Summary for interior simplices</h4>\n"; - if (i < cc_boundary_summary_lists.size()) - ofs << " <p><a href=\"#dim" << i << "b\">Go to boundary</a></p>\n"; + if (i < cc_boundary_summary_lists.size()) ofs << " <p><a href=\"#dim" << i << "b\">Go to boundary</a></p>\n"; ofs << " <ul>\n"; - for (const CC_summary_info& cc_info: cc_interior_summary_lists[i]) + for (const CC_summary_info& cc_info : cc_interior_summary_lists[i]) ofs << " <li id = \"" << id_from_simplex("I" + cc_info.face_) << "\">" - << simplex_format(cc_info.face_, false) - << " cell =" << cc_info.cell_ << "</li>\n"; + << simplex_format(cc_info.face_, false) << " cell =" << cc_info.cell_ << "</li>\n"; ofs << " </ul>\n"; ofs << " <h4> Prejoin state of the interior cells of dimension " << i << "</h4>\n"; auto prejoin_it = cc_interior_prejoin_lists[i].begin(); @@ -340,57 +325,48 @@ void write_cc(std::ofstream& ofs) { ofs << " <h5>" << j << "-dimensional ambient simplices</h5>\n"; ofs << " <ul>\n"; for (; prejoin_it->dimension_ == j; ++prejoin_it) { - ofs << " <li>" << simplex_format(prejoin_it->simplex_, false) - << " join = " << simplex_format(prejoin_it->join_, false) - << " boundary:\n" - << " <ul>\n"; - for (const auto& face: prejoin_it->faces_) - ofs << " <li>" << simplex_format(face) << "</li>"; - ofs << " </ul>\n"; - switch (prejoin_it->status_) { - case (CC_prejoin_info::Result_type::join_single): - ofs << " <p style=\"color: red\">Deleted " - << simplex_format(prejoin_it->simplex_, false) - << " as it has a single face.</p>"; - break; - case (CC_prejoin_info::Result_type::join_is_face): - ofs << " <p style=\"color: red\">Deleted " - << simplex_format(prejoin_it->simplex_, false) - << " as its join " << simplex_format(prejoin_it->join_, false) - << " is one of the faces.</p>"; - break; - case (CC_prejoin_info::Result_type::join_different): - ofs << " <p style=\"color: magenta\">Deleted " << simplex_format(prejoin_it->simplex_, false) - << " and replaced by its join " << simplex_format(prejoin_it->join_, false) - << ".</p>"; - break; - case (CC_prejoin_info::Result_type::join_same): - ofs << " <p style=\"color: green\">Kept " << simplex_format(prejoin_it->simplex_, false) - << ".</p>"; - } - ofs << " </li>"; + ofs << " <li>" << simplex_format(prejoin_it->simplex_, false) + << " join = " << simplex_format(prejoin_it->join_, false) << " boundary:\n" + << " <ul>\n"; + for (const auto& face : prejoin_it->faces_) ofs << " <li>" << simplex_format(face) << "</li>"; + ofs << " </ul>\n"; + switch (prejoin_it->status_) { + case (CC_prejoin_info::Result_type::join_single): + ofs << " <p style=\"color: red\">Deleted " << simplex_format(prejoin_it->simplex_, false) + << " as it has a single face.</p>"; + break; + case (CC_prejoin_info::Result_type::join_is_face): + ofs << " <p style=\"color: red\">Deleted " << simplex_format(prejoin_it->simplex_, false) + << " as its join " << simplex_format(prejoin_it->join_, false) << " is one of the faces.</p>"; + break; + case (CC_prejoin_info::Result_type::join_different): + ofs << " <p style=\"color: magenta\">Deleted " << simplex_format(prejoin_it->simplex_, false) + << " and replaced by its join " << simplex_format(prejoin_it->join_, false) << ".</p>"; + break; + case (CC_prejoin_info::Result_type::join_same): + ofs << " <p style=\"color: green\">Kept " << simplex_format(prejoin_it->simplex_, false) + << ".</p>"; + } + ofs << " </li>"; } ofs << " </ul>\n"; } ofs << " <h4> Details for interior simplices</h4>\n"; ofs << " <ul>\n"; - for (const CC_detail_info& cc_info: cc_interior_detail_lists[i]) { + for (const CC_detail_info& cc_info : cc_interior_detail_lists[i]) { if (cc_info.status_ == CC_detail_info::Result_type::join_single) { - ofs << " <li style=\"color:magenta\" id = \"" - << id_from_simplex("I" + cc_info.simplex_) << "\"> Simplex " - << simplex_format(cc_info.simplex_, false) << " has only one face (" - << simplex_format(cc_info.trigger_, false) << ") and is deleted."; - continue; + ofs << " <li style=\"color:magenta\" id = \"" << id_from_simplex("I" + cc_info.simplex_) + << "\"> Simplex " << simplex_format(cc_info.simplex_, false) << " has only one face (" + << simplex_format(cc_info.trigger_, false) << ") and is deleted."; + continue; } if (cc_info.status_ == CC_detail_info::Result_type::join_single) { - ofs << " <li style=\"color:darkmagenta\" id = \"" - << id_from_simplex("I" + cc_info.simplex_) << "\"> The join of the simplex " - << simplex_format(cc_info.simplex_, false) << " is one of its faces (" - << simplex_format(cc_info.trigger_, false) << "), hence it is is deleted."; - continue; + ofs << " <li style=\"color:darkmagenta\" id = \"" << id_from_simplex("I" + cc_info.simplex_) + << "\"> The join of the simplex " << simplex_format(cc_info.simplex_, false) << " is one of its faces (" + << simplex_format(cc_info.trigger_, false) << "), hence it is is deleted."; + continue; } - ofs << " <li> Insert_cell called for " << simplex_format(cc_info.simplex_, false) - << "\n"; + ofs << " <li> Insert_cell called for " << simplex_format(cc_info.simplex_, false) << "\n"; ofs << " <ul>\n"; // for (const std::string& cof: cc_info.faces_) // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, false) @@ -398,34 +374,31 @@ void write_cc(std::ofstream& ofs) { ofs << " </ul>\n"; ofs << " <ul>\n"; if (cc_info.status_ == CC_detail_info::Result_type::self) { - ofs << " <p><span style=\"color:blue\">The simplex " - << simplex_format(cc_info.simplex_, false) - << " already exists in the cell complex!</span></p>\n"; + ofs << " <p><span style=\"color:blue\">The simplex " << simplex_format(cc_info.simplex_, false) + << " already exists in the cell complex!</span></p>\n"; } if (cc_info.status_ == CC_detail_info::Result_type::face) { - ofs << " <p><span style=\"color:red\">The simplex " - << simplex_format(cc_info.simplex_, false) << " is a face of the simplex " - << simplex_format(cc_info.trigger_, false) << "!</span><br>\n"; - ofs << " <ul>\n"; - for (const std::string post_face: cc_info.post_faces_) - ofs << " <li id = \"" << id_from_simplex("I" + post_face) << "\">" - << "Post deleting " << simplex_format(post_face, false) << "</li>\n"; - ofs << " </ul>\n"; - ofs << " </p>\n"; - ofs << " <p id = \"" << id_from_simplex("I" + cc_info.trigger_) << "\">" - << "Deleting " << simplex_format(cc_info.trigger_, false) << "</p>\n"; + ofs << " <p><span style=\"color:red\">The simplex " << simplex_format(cc_info.simplex_, false) + << " is a face of the simplex " << simplex_format(cc_info.trigger_, false) << "!</span><br>\n"; + ofs << " <ul>\n"; + for (const std::string post_face : cc_info.post_faces_) + ofs << " <li id = \"" << id_from_simplex("I" + post_face) << "\">" + << "Post deleting " << simplex_format(post_face, false) << "</li>\n"; + ofs << " </ul>\n"; + ofs << " </p>\n"; + ofs << " <p id = \"" << id_from_simplex("I" + cc_info.trigger_) << "\">" + << "Deleting " << simplex_format(cc_info.trigger_, false) << "</p>\n"; } // for (const std::string& fac: cc_info.cofaces_) // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, false) // << " is a coface of " << simplex_format(fac, false) << "\n"; if (cc_info.status_ == CC_detail_info::Result_type::coface) { - ofs << " <p><span style=\"color:darkorange\">The simplex " - << simplex_format(cc_info.simplex_, false) << " is a coface of the simplex " - << simplex_format(cc_info.trigger_, false) << "!</span><p>\n"; + ofs << " <p><span style=\"color:darkorange\">The simplex " << simplex_format(cc_info.simplex_, false) + << " is a coface of the simplex " << simplex_format(cc_info.trigger_, false) << "!</span><p>\n"; } if (cc_info.status_ == CC_detail_info::Result_type::inserted) { - ofs << " <p><span style=\"color:green\">Successfully inserted " - << simplex_format(cc_info.simplex_, false) << "!</span><p>\n"; + ofs << " <p><span style=\"color:green\">Successfully inserted " + << simplex_format(cc_info.simplex_, false) << "!</span><p>\n"; } ofs << " </ul>\n"; ofs << " </li>\n"; @@ -436,113 +409,101 @@ void write_cc(std::ofstream& ofs) { ofs << " <h4 id=\"dim" << i << "b\"> Summary for boundary simplices</h4>\n"; ofs << " <p><a href=\"#dim" << i << "i\">Go to interior</a></p>\n"; ofs << " <ul>\n"; - for (const CC_summary_info& cc_info: cc_boundary_summary_lists[i]) - ofs << " <li id = \"" << id_from_simplex("B" + cc_info.face_) << "\">" - << simplex_format(cc_info.face_, true) - << " cell =" << cc_info.cell_ << "</li>\n"; + for (const CC_summary_info& cc_info : cc_boundary_summary_lists[i]) + ofs << " <li id = \"" << id_from_simplex("B" + cc_info.face_) << "\">" + << simplex_format(cc_info.face_, true) << " cell =" << cc_info.cell_ << "</li>\n"; ofs << " </ul>\n"; ofs << " <h4> Prejoin state of the boundary cells of dimension " << i << "</h4>\n"; auto prejoin_it = cc_boundary_prejoin_lists[i].begin(); while (prejoin_it != cc_boundary_prejoin_lists[i].end()) { - std::size_t j = prejoin_it->dimension_; - ofs << " <h5>" << j << "-dimensional ambient simplices</h5>\n"; - ofs << " <ul>\n"; - for (; prejoin_it->dimension_ == j; ++prejoin_it) { - ofs << " <li>" << simplex_format(prejoin_it->simplex_, true) - << " join = " << simplex_format(prejoin_it->join_, true) - << " boundary:\n" - << " <ul>\n"; - for (const auto& face: prejoin_it->faces_) - ofs << " <li>" << simplex_format(face) << "</li>"; - ofs << " </ul>\n"; - switch (prejoin_it->status_) { - case (CC_prejoin_info::Result_type::join_single): - ofs << " <p style=\"color: red\">Deleted " - << simplex_format(prejoin_it->simplex_, true) - << " as it has a single face.</p>"; - break; - case (CC_prejoin_info::Result_type::join_is_face): - ofs << " <p style=\"color: red\">Deleted " - << simplex_format(prejoin_it->simplex_, true) - << " as its join " << simplex_format(prejoin_it->join_, true) - << " is one of the faces.</p>"; - break; - case (CC_prejoin_info::Result_type::join_different): - ofs << " <p style=\"color: magenta\">Deleted " << simplex_format(prejoin_it->simplex_, true) - << " and replaced by its join " << simplex_format(prejoin_it->join_, true) - << ".</p>"; - break; - case (CC_prejoin_info::Result_type::join_same): - ofs << " <p style=\"color: green\">Kept " << simplex_format(prejoin_it->simplex_, true) - << ".</p>"; - } - ofs << " </li>"; - } - ofs << " </ul>\n"; + std::size_t j = prejoin_it->dimension_; + ofs << " <h5>" << j << "-dimensional ambient simplices</h5>\n"; + ofs << " <ul>\n"; + for (; prejoin_it->dimension_ == j; ++prejoin_it) { + ofs << " <li>" << simplex_format(prejoin_it->simplex_, true) + << " join = " << simplex_format(prejoin_it->join_, true) << " boundary:\n" + << " <ul>\n"; + for (const auto& face : prejoin_it->faces_) ofs << " <li>" << simplex_format(face) << "</li>"; + ofs << " </ul>\n"; + switch (prejoin_it->status_) { + case (CC_prejoin_info::Result_type::join_single): + ofs << " <p style=\"color: red\">Deleted " << simplex_format(prejoin_it->simplex_, true) + << " as it has a single face.</p>"; + break; + case (CC_prejoin_info::Result_type::join_is_face): + ofs << " <p style=\"color: red\">Deleted " << simplex_format(prejoin_it->simplex_, true) + << " as its join " << simplex_format(prejoin_it->join_, true) << " is one of the faces.</p>"; + break; + case (CC_prejoin_info::Result_type::join_different): + ofs << " <p style=\"color: magenta\">Deleted " << simplex_format(prejoin_it->simplex_, true) + << " and replaced by its join " << simplex_format(prejoin_it->join_, true) << ".</p>"; + break; + case (CC_prejoin_info::Result_type::join_same): + ofs << " <p style=\"color: green\">Kept " << simplex_format(prejoin_it->simplex_, true) + << ".</p>"; + } + ofs << " </li>"; + } + ofs << " </ul>\n"; } } if (i < cc_boundary_detail_lists.size()) { ofs << " <h4> Details for boundary simplices</h4>\n" - << " <ul>\n"; - for (const CC_detail_info& cc_info: cc_boundary_detail_lists[i]) { - if (cc_info.status_ == CC_detail_info::Result_type::join_single) { - ofs << " <li style=\"color:magenta\" id = \"" - << id_from_simplex("B" + cc_info.simplex_) << "\"> Simplex " - << simplex_format(cc_info.simplex_, true) << " has only one face (" - << simplex_format(cc_info.trigger_, true) << ") and is deleted."; - continue; - } - if (cc_info.status_ == CC_detail_info::Result_type::join_single) { - ofs << " <li style=\"color:darkmagenta\" id = \"" - << id_from_simplex("B" + cc_info.simplex_) << "\"> The join of the simplex " - << simplex_format(cc_info.simplex_, true) << " is one of its faces (" - << simplex_format(cc_info.trigger_, true) << "), hence it is is deleted."; - continue; - } - ofs << " <li> Insert_simplex called on " << simplex_format(cc_info.simplex_, true); - ofs << " <ul>\n"; - // for (const std::string& cof: cc_info.faces_) - // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, true) - // << " is a face of " << simplex_format(cof, true) << "\n"; - ofs << " </ul>\n"; - ofs << " <ul>\n"; - if (cc_info.status_ == CC_detail_info::Result_type::self) { - ofs << " <p><span style=\"color:blue\">The simplex " - << simplex_format(cc_info.simplex_, true) - << " already exists in the cell complex!</span></p>\n"; - } - if (cc_info.status_ == CC_detail_info::Result_type::face) { - ofs << " <p><span style=\"color:red\">The simplex " - << simplex_format(cc_info.simplex_, true) << " is a face of the simplex " - << simplex_format(cc_info.trigger_, true) << "!</span><br>\n"; - ofs << " <ul>\n"; - for (const std::string post_face: cc_info.post_faces_) - ofs << " <li id=\"" << id_from_simplex("B" + post_face) - << "\">Post deleting " << simplex_format(post_face, true) - << "</li>\n"; - ofs << " </ul>\n"; - ofs << " </p>\n"; - ofs << " <p id=\"" << id_from_simplex(cc_info.trigger_) - << "\">Deleting " << simplex_format(cc_info.trigger_, true) << "</p>\n"; - } - // for (const std::string& fac: cc_info.cofaces_) - // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, true) - // << " is a coface of " << simplex_format(fac, true) << "\n"; - ofs << " </ul>\n"; - ofs << " </li>\n"; - if (cc_info.status_ == CC_detail_info::Result_type::coface) { - ofs << " <p><span style=\"color:darkorange\">The simplex " - << simplex_format(cc_info.simplex_, true) << " is a coface of the simplex " - << simplex_format(cc_info.trigger_, true) << "!</span><p>\n"; - } - if (cc_info.status_ == CC_detail_info::Result_type::inserted) { - ofs << " <p><span style=\"color:green\">Successfully inserted " - << simplex_format(cc_info.simplex_, true) << "!</span><p>\n"; - } + << " <ul>\n"; + for (const CC_detail_info& cc_info : cc_boundary_detail_lists[i]) { + if (cc_info.status_ == CC_detail_info::Result_type::join_single) { + ofs << " <li style=\"color:magenta\" id = \"" << id_from_simplex("B" + cc_info.simplex_) + << "\"> Simplex " << simplex_format(cc_info.simplex_, true) << " has only one face (" + << simplex_format(cc_info.trigger_, true) << ") and is deleted."; + continue; + } + if (cc_info.status_ == CC_detail_info::Result_type::join_single) { + ofs << " <li style=\"color:darkmagenta\" id = \"" << id_from_simplex("B" + cc_info.simplex_) + << "\"> The join of the simplex " << simplex_format(cc_info.simplex_, true) << " is one of its faces (" + << simplex_format(cc_info.trigger_, true) << "), hence it is is deleted."; + continue; + } + ofs << " <li> Insert_simplex called on " << simplex_format(cc_info.simplex_, true); + ofs << " <ul>\n"; + // for (const std::string& cof: cc_info.faces_) + // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, true) + // << " is a face of " << simplex_format(cof, true) << "\n"; + ofs << " </ul>\n"; + ofs << " <ul>\n"; + if (cc_info.status_ == CC_detail_info::Result_type::self) { + ofs << " <p><span style=\"color:blue\">The simplex " << simplex_format(cc_info.simplex_, true) + << " already exists in the cell complex!</span></p>\n"; + } + if (cc_info.status_ == CC_detail_info::Result_type::face) { + ofs << " <p><span style=\"color:red\">The simplex " << simplex_format(cc_info.simplex_, true) + << " is a face of the simplex " << simplex_format(cc_info.trigger_, true) << "!</span><br>\n"; + ofs << " <ul>\n"; + for (const std::string post_face : cc_info.post_faces_) + ofs << " <li id=\"" << id_from_simplex("B" + post_face) << "\">Post deleting " + << simplex_format(post_face, true) << "</li>\n"; + ofs << " </ul>\n"; + ofs << " </p>\n"; + ofs << " <p id=\"" << id_from_simplex(cc_info.trigger_) << "\">Deleting " + << simplex_format(cc_info.trigger_, true) << "</p>\n"; + } + // for (const std::string& fac: cc_info.cofaces_) + // ofs << " <li>Checking if " << simplex_format(cc_info.simplex_, true) + // << " is a coface of " << simplex_format(fac, true) << "\n"; + ofs << " </ul>\n"; + ofs << " </li>\n"; + if (cc_info.status_ == CC_detail_info::Result_type::coface) { + ofs << " <p><span style=\"color:darkorange\">The simplex " + << simplex_format(cc_info.simplex_, true) << " is a coface of the simplex " + << simplex_format(cc_info.trigger_, true) << "!</span><p>\n"; + } + if (cc_info.status_ == CC_detail_info::Result_type::inserted) { + ofs << " <p><span style=\"color:green\">Successfully inserted " + << simplex_format(cc_info.simplex_, true) << "!</span><p>\n"; + } } ofs << " </ul>\n"; } - } + } ofs << " </div>\n"; } @@ -551,14 +512,12 @@ void write_visu(std::ofstream& ofs) { << " <h2> Visualization details debug trace </h2>\n"; // std::vector<std::map<std::string, std::string> > vs_maps(cc_interior_summary_lists.size()); std::map<std::string, std::string> vs_map; - for (const auto& sv_pair: simplex_vlist_map) - vs_map.emplace(std::make_pair(sv_pair.second, sv_pair.first)); + for (const auto& sv_pair : simplex_vlist_map) vs_map.emplace(std::make_pair(sv_pair.second, sv_pair.first)); ofs << " <ul>\n"; - for (const auto& vs_pair: vs_map) { + for (const auto& vs_pair : vs_map) { std::string w_simplex = vs_pair.second.substr(1); bool is_boundary = vs_pair.second[0] == 'B'; - ofs << " <li><b>" << vs_pair.first << "</b>: " - << simplex_format(w_simplex, is_boundary) << "</li>\n"; + ofs << " <li><b>" << vs_pair.first << "</b>: " << simplex_format(w_simplex, is_boundary) << "</li>\n"; } ofs << " </ul>\n"; ofs << " </div>\n"; @@ -574,16 +533,15 @@ void write_to_html(std::string file_name) { ofs << " <h1> Debug traces for " << file_name << " </h1>\n"; write_mt(ofs); write_cc(ofs); - write_visu(ofs); + write_visu(ofs); ofs << " </body>\n"; ofs << "</html>\n"; ofs.close(); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif - diff --git a/src/Coxeter_triangulation/include/gudhi/IO/output_meshes_to_medit.h b/src/Coxeter_triangulation/include/gudhi/IO/output_meshes_to_medit.h index 850736e9..4b454373 100644 --- a/src/Coxeter_triangulation/include/gudhi/IO/output_meshes_to_medit.h +++ b/src/Coxeter_triangulation/include/gudhi/IO/output_meshes_to_medit.h @@ -19,8 +19,8 @@ #include <fstream> // for std::ofstream #include <vector> #include <type_traits> // for std::enable_if -#include <tuple> // for std::get -#include <utility> // for std::make_pair +#include <tuple> // for std::get +#include <utility> // for std::make_pair namespace Gudhi { @@ -30,67 +30,46 @@ using Vertex_points = Mesh_medit::Vertex_points; using Mesh_elements = Mesh_medit::Mesh_elements; using Scalar_field_range = Mesh_medit::Scalar_field_range; -template <std::size_t I = 0, - typename... Meshes> -typename std::enable_if<I == sizeof... (Meshes), void>::type -fill_meshes (Vertex_points& vertex_points, - Mesh_elements& edges, - Mesh_elements& triangles, - Mesh_elements& tetrahedra, - Scalar_field_range& triangles_scalar_range, - Scalar_field_range& tetrahedra_scalar_range, - std::size_t index, - const Meshes&... meshes) { -} - -template <std::size_t I = 0, - typename... Meshes> -typename std::enable_if<I != sizeof... (Meshes), void>::type -fill_meshes (Vertex_points& vertex_points, - Mesh_elements& edges, - Mesh_elements& triangles, - Mesh_elements& tetrahedra, - Scalar_field_range& triangles_scalar_range, - Scalar_field_range& tetrahedra_scalar_range, - std::size_t index, - const Meshes&... meshes) { +template <std::size_t I = 0, typename... Meshes> +typename std::enable_if<I == sizeof...(Meshes), void>::type fill_meshes(Vertex_points& vertex_points, + Mesh_elements& edges, Mesh_elements& triangles, + Mesh_elements& tetrahedra, + Scalar_field_range& triangles_scalar_range, + Scalar_field_range& tetrahedra_scalar_range, + std::size_t index, const Meshes&... meshes) {} + +template <std::size_t I = 0, typename... Meshes> +typename std::enable_if<I != sizeof...(Meshes), void>::type fill_meshes(Vertex_points& vertex_points, + Mesh_elements& edges, Mesh_elements& triangles, + Mesh_elements& tetrahedra, + Scalar_field_range& triangles_scalar_range, + Scalar_field_range& tetrahedra_scalar_range, + std::size_t index, const Meshes&... meshes) { auto mesh = std::get<I>(std::forward_as_tuple(meshes...)); - for (const auto& v: mesh.vertex_points) - vertex_points.push_back(v); - for (const auto& e: mesh.edges) { + for (const auto& v : mesh.vertex_points) vertex_points.push_back(v); + for (const auto& e : mesh.edges) { std::vector<std::size_t> edge; - for (const auto& v_i: e.first) - edge.push_back(v_i + index); + for (const auto& v_i : e.first) edge.push_back(v_i + index); edges.emplace_back(std::make_pair(edge, e.second)); } - for (const auto& t: mesh.triangles) { + for (const auto& t : mesh.triangles) { std::vector<std::size_t> triangle; - for (const auto& v_i: t.first) - triangle.push_back(v_i + index); + for (const auto& v_i : t.first) triangle.push_back(v_i + index); triangles.emplace_back(std::make_pair(triangle, t.second)); } - for (const auto& t: mesh.tetrahedra) { + for (const auto& t : mesh.tetrahedra) { std::vector<std::size_t> tetrahedron; - for (const auto& v_i: t.first) - tetrahedron.push_back(v_i + index); + for (const auto& v_i : t.first) tetrahedron.push_back(v_i + index); tetrahedra.emplace_back(std::make_pair(tetrahedron, t.second)); } - for (const auto& b: mesh.triangles_scalar_range) - triangles_scalar_range.push_back(b); - for (const auto& b: mesh.tetrahedra_scalar_range) - tetrahedra_scalar_range.push_back(b); - fill_meshes<I+1, Meshes...>(vertex_points, - edges, - triangles, - tetrahedra, - triangles_scalar_range, - tetrahedra_scalar_range, - index + mesh.vertex_points.size(), - meshes...); + for (const auto& b : mesh.triangles_scalar_range) triangles_scalar_range.push_back(b); + for (const auto& b : mesh.tetrahedra_scalar_range) tetrahedra_scalar_range.push_back(b); + fill_meshes<I + 1, Meshes...>(vertex_points, edges, triangles, tetrahedra, triangles_scalar_range, + tetrahedra_scalar_range, index + mesh.vertex_points.size(), meshes...); } /** \brief Outputs a text file with specified meshes that can be visualized in Medit. - * + * * @param[in] amb_d Ambient dimension. Can be 2 or 3. * @param[in] file_name The name of the output file. * @param[in] meshes A pack of meshes to be specified separated by commas. @@ -100,87 +79,73 @@ void output_meshes_to_medit(std::size_t amb_d, std::string file_name, const Mesh Vertex_points vertex_points; Mesh_elements edges, triangles, tetrahedra; Scalar_field_range triangles_scalar_range, tetrahedra_scalar_range; - fill_meshes(vertex_points, - edges, - triangles, - tetrahedra, - triangles_scalar_range, - tetrahedra_scalar_range, - 0, - meshes...); - - std::ofstream ofs (file_name + ".mesh", std::ofstream::out); - std::ofstream ofs_bb (file_name + ".bb", std::ofstream::out); - + fill_meshes(vertex_points, edges, triangles, tetrahedra, triangles_scalar_range, tetrahedra_scalar_range, 0, + meshes...); + + std::ofstream ofs(file_name + ".mesh", std::ofstream::out); + std::ofstream ofs_bb(file_name + ".bb", std::ofstream::out); + if (amb_d == 2) { ofs << "MeshVersionFormatted 1\nDimension 2\n"; ofs_bb << "2 1 "; ofs << "Vertices\n" << vertex_points.size() << "\n"; - for (auto p: vertex_points) { + for (auto p : vertex_points) { ofs << p[0] << " " << p[1] << " 2\n"; } ofs << "Edges " << edges.size() << "\n"; - for (auto e: edges) { - for (auto v: e.first) - ofs << v << " "; + for (auto e : edges) { + for (auto v : e.first) ofs << v << " "; ofs << e.second << std::endl; } ofs << "Triangles " << triangles.size() << "\n"; - for (auto s: triangles) { - for (auto v: s.first) { - ofs << v << " "; + for (auto s : triangles) { + for (auto v : s.first) { + ofs << v << " "; } ofs << s.second << std::endl; } ofs_bb << triangles_scalar_range.size() << " 1\n"; - for (auto& b: triangles_scalar_range) - ofs_bb << b << "\n"; + for (auto& b : triangles_scalar_range) ofs_bb << b << "\n"; - } - else { + } else { ofs << "MeshVersionFormatted 1\nDimension 3\n"; ofs_bb << "3 1 "; ofs << "Vertices\n" << vertex_points.size() << "\n"; - for (auto p: vertex_points) { + for (auto p : vertex_points) { ofs << p[0] << " " << p[1] << " " << p[2] << " 2\n"; } ofs << "Edges " << edges.size() << "\n"; - for (auto e: edges) { - for (auto v: e.first) - ofs << v << " "; + for (auto e : edges) { + for (auto v : e.first) ofs << v << " "; ofs << e.second << std::endl; } ofs << "Triangles " << triangles.size() << "\n"; - for (auto s: triangles) { - for (auto v: s.first) { - ofs << v << " "; + for (auto s : triangles) { + for (auto v : s.first) { + ofs << v << " "; } ofs << s.second << std::endl; } ofs << "Tetrahedra " << tetrahedra.size() << "\n"; - for (auto s: tetrahedra) { - for (auto v: s.first) { - ofs << v << " "; + for (auto s : tetrahedra) { + for (auto v : s.first) { + ofs << v << " "; } ofs << s.second << std::endl; } ofs_bb << triangles_scalar_range.size() + tetrahedra_scalar_range.size() << " 1\n"; - for (auto& b: triangles_scalar_range) - ofs_bb << b << "\n"; - for (auto& b: tetrahedra_scalar_range) - ofs_bb << b << "\n"; + for (auto& b : triangles_scalar_range) ofs_bb << b << "\n"; + for (auto& b : tetrahedra_scalar_range) ofs_bb << b << "\n"; } - ofs.close(); ofs_bb.close(); - } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Implicit_manifold_intersection_oracle.h b/src/Coxeter_triangulation/include/gudhi/Implicit_manifold_intersection_oracle.h index 51d84274..4f44df7c 100644 --- a/src/Coxeter_triangulation/include/gudhi/Implicit_manifold_intersection_oracle.h +++ b/src/Coxeter_triangulation/include/gudhi/Implicit_manifold_intersection_oracle.h @@ -24,169 +24,142 @@ namespace Gudhi { namespace coxeter_triangulation { - - /** \class Implicit_manifold_intersection_oracle * \brief An oracle that supports the intersection query on an implicit manifold. * - * \tparam Function_ The function template parameter. Should be a model of + * \tparam Function_ The function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * \tparam Domain_function_ The domain function template parameter. Should be a model of * the concept FunctionForImplicitManifold. * * \ingroup coxeter_triangulation */ -template<class Function_, - class Domain_function_ = Constant_function> +template <class Function_, class Domain_function_ = Constant_function> class Implicit_manifold_intersection_oracle { - /* Computes the affine coordinates of the intersection point of the implicit manifold * and the affine hull of the simplex. */ - template <class Simplex_handle, - class Triangulation> - Eigen::VectorXd compute_lambda(const Simplex_handle& simplex, - const Triangulation& triangulation) const { + template <class Simplex_handle, class Triangulation> + Eigen::VectorXd compute_lambda(const Simplex_handle& simplex, const Triangulation& triangulation) const { std::size_t cod_d = this->cod_d(); Eigen::MatrixXd matrix(cod_d + 1, cod_d + 1); - for (std::size_t i = 0; i < cod_d + 1; ++i) - matrix(0, i) = 1; + for (std::size_t i = 0; i < cod_d + 1; ++i) matrix(0, i) = 1; std::size_t j = 0; - for (auto v: simplex.vertex_range()) { + for (auto v : simplex.vertex_range()) { Eigen::VectorXd v_coords = fun_(triangulation.cartesian_coordinates(v)); - for (std::size_t i = 1; i < cod_d + 1; ++i) - matrix(i, j) = v_coords(i-1); + for (std::size_t i = 1; i < cod_d + 1; ++i) matrix(i, j) = v_coords(i - 1); j++; } Eigen::VectorXd z(cod_d + 1); z(0) = 1; - for (std::size_t i = 1; i < cod_d + 1; ++i) - z(i) = 0; + for (std::size_t i = 1; i < cod_d + 1; ++i) z(i) = 0; Eigen::VectorXd lambda = matrix.colPivHouseholderQr().solve(z); return lambda; } /* Computes the affine coordinates of the intersection point of the boundary * of the implicit manifold and the affine hull of the simplex. */ - template <class Simplex_handle, - class Triangulation> - Eigen::VectorXd compute_boundary_lambda(const Simplex_handle& simplex, - const Triangulation& triangulation) const { + template <class Simplex_handle, class Triangulation> + Eigen::VectorXd compute_boundary_lambda(const Simplex_handle& simplex, const Triangulation& triangulation) const { std::size_t cod_d = this->cod_d(); Eigen::MatrixXd matrix(cod_d + 2, cod_d + 2); - for (std::size_t i = 0; i < cod_d + 2; ++i) - matrix(0, i) = 1; + for (std::size_t i = 0; i < cod_d + 2; ++i) matrix(0, i) = 1; std::size_t j = 0; - for (auto v: simplex.vertex_range()) { + for (auto v : simplex.vertex_range()) { Eigen::VectorXd v_coords = fun_(triangulation.cartesian_coordinates(v)); - for (std::size_t i = 1; i < cod_d + 1; ++i) - matrix(i, j) = v_coords(i-1); + for (std::size_t i = 1; i < cod_d + 1; ++i) matrix(i, j) = v_coords(i - 1); Eigen::VectorXd bv_coords = domain_fun_(triangulation.cartesian_coordinates(v)); matrix(cod_d + 1, j) = bv_coords(0); j++; } Eigen::VectorXd z(cod_d + 2); z(0) = 1; - for (std::size_t i = 1; i < cod_d + 2; ++i) - z(i) = 0; + for (std::size_t i = 1; i < cod_d + 2; ++i) z(i) = 0; Eigen::VectorXd lambda = matrix.colPivHouseholderQr().solve(z); return lambda; } /* Computes the intersection result for a given simplex in a triangulation. */ - template <class Simplex_handle, - class Triangulation> - Query_result<Simplex_handle> intersection_result(const Eigen::VectorXd& lambda, - const Simplex_handle& simplex, + template <class Simplex_handle, class Triangulation> + Query_result<Simplex_handle> intersection_result(const Eigen::VectorXd& lambda, const Simplex_handle& simplex, const Triangulation& triangulation) const { using QR = Query_result<Simplex_handle>; std::size_t amb_d = triangulation.dimension(); std::size_t cod_d = simplex.dimension(); for (std::size_t i = 0; i < (std::size_t)lambda.size(); ++i) - if (lambda(i) < 0 || lambda(i) > 1) - return QR({Eigen::VectorXd(), false}); + if (lambda(i) < 0 || lambda(i) > 1) return QR({Eigen::VectorXd(), false}); Eigen::MatrixXd vertex_matrix(cod_d + 1, amb_d); auto v_range = simplex.vertex_range(); auto v_it = v_range.begin(); for (std::size_t i = 0; i < cod_d + 1 && v_it != v_range.end(); ++v_it, ++i) { Eigen::VectorXd v_coords = triangulation.cartesian_coordinates(*v_it); - for (std::size_t j = 0; j < amb_d; ++j) - vertex_matrix(i, j) = v_coords(j); + for (std::size_t j = 0; j < amb_d; ++j) vertex_matrix(i, j) = v_coords(j); } - Eigen::VectorXd intersection = lambda.transpose()*vertex_matrix; + Eigen::VectorXd intersection = lambda.transpose() * vertex_matrix; return QR({intersection, true}); } - -public: + public: /** \brief Ambient dimension of the implicit manifold. */ - std::size_t amb_d() const { - return fun_.amb_d(); - } - + std::size_t amb_d() const { return fun_.amb_d(); } + /** \brief Codimension of the implicit manifold. */ - std::size_t cod_d() const { - return fun_.cod_d(); - } + std::size_t cod_d() const { return fun_.cod_d(); } /** \brief Intersection query with the relative interior of the manifold. - * + * * \details The returned structure Query_result contains the boolean value * that is true only if the intersection point of the query simplex and * the relative interior of the manifold exists, the intersection point - * and the face of the query simplex that contains + * and the face of the query simplex that contains * the intersection point. - * + * * \tparam Simplex_handle The class of the query simplex. * Needs to be a model of the concept SimplexInCoxeterTriangulation. * \tparam Triangulation The class of the triangulation. * Needs to be a model of the concept TriangulationForManifoldTracing. * * @param[in] simplex The query simplex. The dimension of the simplex - * should be the same as the codimension of the manifold + * should be the same as the codimension of the manifold * (the codomain dimension of the function). - * @param[in] triangulation The ambient triangulation. The dimension of - * the triangulation should be the same as the ambient dimension of the manifold + * @param[in] triangulation The ambient triangulation. The dimension of + * the triangulation should be the same as the ambient dimension of the manifold * (the domain dimension of the function). */ - template <class Simplex_handle, - class Triangulation> - Query_result<Simplex_handle> intersects(const Simplex_handle& simplex, - const Triangulation& triangulation) const { + template <class Simplex_handle, class Triangulation> + Query_result<Simplex_handle> intersects(const Simplex_handle& simplex, const Triangulation& triangulation) const { Eigen::VectorXd lambda = compute_lambda(simplex, triangulation); return intersection_result(lambda, simplex, triangulation); } /** \brief Intersection query with the boundary of the manifold. - * + * * \details The returned structure Query_result contains the boolean value * that is true only if the intersection point of the query simplex and * the boundary of the manifold exists, the intersection point - * and the face of the query simplex that contains + * and the face of the query simplex that contains * the intersection point. - * + * * \tparam Simplex_handle The class of the query simplex. * Needs to be a model of the concept SimplexInCoxeterTriangulation. * \tparam Triangulation The class of the triangulation. * Needs to be a model of the concept TriangulationForManifoldTracing. * * @param[in] simplex The query simplex. The dimension of the simplex - * should be the same as the codimension of the boundary of the manifold + * should be the same as the codimension of the boundary of the manifold * (the codomain dimension of the function + 1). - * @param[in] triangulation The ambient triangulation. The dimension of - * the triangulation should be the same as the ambient dimension of the manifold + * @param[in] triangulation The ambient triangulation. The dimension of + * the triangulation should be the same as the ambient dimension of the manifold * (the domain dimension of the function). */ - template <class Simplex_handle, - class Triangulation> + template <class Simplex_handle, class Triangulation> Query_result<Simplex_handle> intersects_boundary(const Simplex_handle& simplex, const Triangulation& triangulation) const { Eigen::VectorXd lambda = compute_boundary_lambda(simplex, triangulation); return intersection_result(lambda, simplex, triangulation); } - /** \brief Returns true if the input point lies inside the piecewise-linear * domain induced by the given ambient triangulation that defines the relative * interior of the piecewise-linear approximation of the manifold. @@ -194,22 +167,19 @@ public: * @param p The input point. Needs to have the same dimension as the ambient * dimension of the manifold (the domain dimension of the function). * @param triangulation The ambient triangulation. Needs to have the same - * dimension as the ambient dimension of the manifold + * dimension as the ambient dimension of the manifold * (the domain dimension of the function). */ template <class Triangulation> - bool lies_in_domain(const Eigen::VectorXd& p, - const Triangulation& triangulation) const { + bool lies_in_domain(const Eigen::VectorXd& p, const Triangulation& triangulation) const { Eigen::VectorXd pl_p = make_pl_approximation(domain_fun_, triangulation)(p); return pl_p(0) < 0; } /** \brief Returns the function that defines the interior of the manifold */ - const Function_& function() const { - return fun_; - } - - /** \brief Constructs an intersection oracle for an implicit manifold potentially + const Function_& function() const { return fun_; } + + /** \brief Constructs an intersection oracle for an implicit manifold potentially * with boundary from given function and domain. * * @param function The input function that represents the implicit manifold @@ -217,24 +187,22 @@ public: * @param domain_function The input domain function that can be used to define an implicit * manifold with boundary. */ - Implicit_manifold_intersection_oracle(const Function_& function, - const Domain_function_& domain_function) - : fun_(function), domain_fun_(domain_function) {} + Implicit_manifold_intersection_oracle(const Function_& function, const Domain_function_& domain_function) + : fun_(function), domain_fun_(domain_function) {} - /** \brief Constructs an intersection oracle for an implicit manifold + /** \brief Constructs an intersection oracle for an implicit manifold * without boundary from a given function. * - * \details To use this constructor, the template Domain_function_ needs to be left + * \details To use this constructor, the template Domain_function_ needs to be left * at its default value (Gudhi::coxeter_triangulation::Constant_function). * * @param function The input function that represents the implicit manifold * without boundary. */ Implicit_manifold_intersection_oracle(const Function_& function) - : fun_(function), - domain_fun_(function.amb_d(), 1, Eigen::VectorXd::Constant(1,-1)) {} - -private: + : fun_(function), domain_fun_(function.amb_d(), 1, Eigen::VectorXd::Constant(1, -1)) {} + + private: Function_ fun_; Domain_function_ domain_fun_; }; @@ -248,16 +216,12 @@ private: * * \ingroup coxeter_triangulation */ -template<class Function_, - class Domain_function_> -Implicit_manifold_intersection_oracle<Function_, Domain_function_> -make_oracle(const Function_& function, - const Domain_function_& domain_function){ - return Implicit_manifold_intersection_oracle<Function_, Domain_function_>(function, - domain_function); +template <class Function_, class Domain_function_> +Implicit_manifold_intersection_oracle<Function_, Domain_function_> make_oracle( + const Function_& function, const Domain_function_& domain_function) { + return Implicit_manifold_intersection_oracle<Function_, Domain_function_>(function, domain_function); } - /** \brief Static constructor of an intersection oracle from a function without a domain. * * @param function The input function that represents the implicit manifold @@ -265,13 +229,13 @@ make_oracle(const Function_& function, * * \ingroup coxeter_triangulation */ -template<class Function_> -Implicit_manifold_intersection_oracle<Function_> make_oracle(const Function_& function){ +template <class Function_> +Implicit_manifold_intersection_oracle<Function_> make_oracle(const Function_& function) { return Implicit_manifold_intersection_oracle<Function_>(function); } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Manifold_tracing.h b/src/Coxeter_triangulation/include/gudhi/Manifold_tracing.h index 25b664eb..b678566a 100644 --- a/src/Coxeter_triangulation/include/gudhi/Manifold_tracing.h +++ b/src/Coxeter_triangulation/include/gudhi/Manifold_tracing.h @@ -28,7 +28,7 @@ namespace coxeter_triangulation { * \ingroup coxeter_triangulation */ -/** \class Manifold_tracing +/** \class Manifold_tracing * \brief A class that assembles methods for manifold tracing algorithm. * * \tparam Triangulation_ The type of the ambient triangulation. @@ -36,9 +36,8 @@ namespace coxeter_triangulation { */ template <class Triangulation_> class Manifold_tracing { - typedef typename Triangulation_::Simplex_handle Simplex_handle; - + struct Simplex_hash { typedef Simplex_handle argument_type; typedef std::size_t result_type; @@ -46,8 +45,8 @@ class Manifold_tracing { return boost::hash<typename Simplex_handle::Vertex>()(s.vertex()); } }; - -public: + + public: /** \brief Type of the output simplex map with keys of type Triangulation_::Simplex_handle * and values of type Eigen::VectorXd. * This type should be used for the output in the method manifold_tracing_algorithm. @@ -58,38 +57,34 @@ public: * \brief Computes the set of k-simplices that intersect * a boundaryless implicit manifold given by an intersection oracle, where k * is the codimension of the manifold. - * The computation is based on the seed propagation --- it starts at the + * The computation is based on the seed propagation --- it starts at the * given seed points and then propagates along the manifold. * * \tparam Point_range Range of points of type Eigen::VectorXd. * \tparam Intersection_oracle Intersection oracle that represents the manifold. * Needs to be a model of the concept IntersectionOracle. * - * \param[in] seed_points The range of points on the manifold from which + * \param[in] seed_points The range of points on the manifold from which * the computation begins. * \param[in] triangulation The ambient triangulation. * \param[in] oracle The intersection oracle for the manifold. * The ambient dimension needs to match the dimension of the * triangulation. * \param[out] out_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the input manifold and the mapped values + * the input triangulation that intersect the input manifold and the mapped values * are the intersection points. */ - template <class Point_range, - class Intersection_oracle> - void manifold_tracing_algorithm(const Point_range& seed_points, - const Triangulation_& triangulation, - const Intersection_oracle& oracle, - Out_simplex_map& out_simplex_map) { + template <class Point_range, class Intersection_oracle> + void manifold_tracing_algorithm(const Point_range& seed_points, const Triangulation_& triangulation, + const Intersection_oracle& oracle, Out_simplex_map& out_simplex_map) { std::size_t cod_d = oracle.cod_d(); std::queue<Simplex_handle> queue; - for (const auto& p: seed_points) { + for (const auto& p : seed_points) { Simplex_handle full_simplex = triangulation.locate_point(p); - for (Simplex_handle face: full_simplex.face_range(cod_d)) { + for (Simplex_handle face : full_simplex.face_range(cod_d)) { Query_result<Simplex_handle> qr = oracle.intersects(face, triangulation); - if (qr.success && - out_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) { + if (qr.success && out_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) { #ifdef GUDHI_COX_OUTPUT_TO_HTML mt_seed_inserted_list.push_back(MT_inserted_info(qr, face, false)); #endif @@ -98,17 +93,14 @@ public: } } } - while (!queue.empty()) { Simplex_handle s = queue.front(); queue.pop(); - for (auto cof: s.coface_range(cod_d+1)) { - for (auto face: cof.face_range(cod_d)) { + for (auto cof : s.coface_range(cod_d + 1)) { + for (auto face : cof.face_range(cod_d)) { Query_result<Simplex_handle> qr = oracle.intersects(face, triangulation); - if (qr.success && - out_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) - queue.emplace(face); + if (qr.success && out_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) queue.emplace(face); } } } @@ -118,85 +110,76 @@ public: * \brief Computes the set of k-simplices that intersect * the dimensional manifold given by an intersection oracle, where k * is the codimension of the manifold. - * The computation is based on the seed propagation --- it starts at the + * The computation is based on the seed propagation --- it starts at the * given seed points and then propagates along the manifold. * * \tparam Point_range Range of points of type Eigen::VectorXd. * \tparam Intersection_oracle Intersection oracle that represents the manifold. * Needs to be a model of the concept IntersectionOracle. * - * \param[in] seed_points The range of points on the manifold from which + * \param[in] seed_points The range of points on the manifold from which * the computation begins. * \param[in] triangulation The ambient triangulation. * \param[in] oracle The intersection oracle for the manifold. * The ambient dimension needs to match the dimension of the * triangulation. * \param[out] interior_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the relative interior of the input manifold + * the input triangulation that intersect the relative interior of the input manifold * and the mapped values are the intersection points. * \param[out] boundary_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the boundary of the input manifold + * the input triangulation that intersect the boundary of the input manifold * and the mapped values are the intersection points. */ - template <class Point_range, - class Intersection_oracle> - void manifold_tracing_algorithm(const Point_range& seed_points, - const Triangulation_& triangulation, - const Intersection_oracle& oracle, - Out_simplex_map& interior_simplex_map, + template <class Point_range, class Intersection_oracle> + void manifold_tracing_algorithm(const Point_range& seed_points, const Triangulation_& triangulation, + const Intersection_oracle& oracle, Out_simplex_map& interior_simplex_map, Out_simplex_map& boundary_simplex_map) { std::size_t cod_d = oracle.cod_d(); std::queue<Simplex_handle> queue; - for (const auto& p: seed_points) { + for (const auto& p : seed_points) { Simplex_handle full_simplex = triangulation.locate_point(p); - for (Simplex_handle face: full_simplex.face_range(cod_d)) { + for (Simplex_handle face : full_simplex.face_range(cod_d)) { auto qr = oracle.intersects(face, triangulation); #ifdef GUDHI_COX_OUTPUT_TO_HTML mt_seed_inserted_list.push_back(MT_inserted_info(qr, face, false)); -#endif +#endif if (qr.success) { if (oracle.lies_in_domain(qr.intersection, triangulation)) { - if (interior_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) - queue.emplace(face); - } - else { - for (Simplex_handle cof: face.coface_range(cod_d+1)) { + if (interior_simplex_map.emplace(std::make_pair(face, qr.intersection)).second) queue.emplace(face); + } else { + for (Simplex_handle cof : face.coface_range(cod_d + 1)) { auto qrb = oracle.intersects_boundary(cof, triangulation); #ifdef GUDHI_COX_OUTPUT_TO_HTML mt_seed_inserted_list.push_back(MT_inserted_info(qrb, cof, true)); -#endif - if (qrb.success) - boundary_simplex_map.emplace(cof, qrb.intersection); +#endif + if (qrb.success) boundary_simplex_map.emplace(cof, qrb.intersection); } } // break; } } } - + while (!queue.empty()) { Simplex_handle s = queue.front(); queue.pop(); - for (auto cof: s.coface_range(cod_d+1)) { - for (auto face: cof.face_range(cod_d)) { + for (auto cof : s.coface_range(cod_d + 1)) { + for (auto face : cof.face_range(cod_d)) { auto qr = oracle.intersects(face, triangulation); #ifdef GUDHI_COX_OUTPUT_TO_HTML mt_inserted_list.push_back(MT_inserted_info(qr, face, false)); -#endif +#endif if (qr.success) { if (oracle.lies_in_domain(qr.intersection, triangulation)) { - if (interior_simplex_map.emplace(face, qr.intersection).second) - queue.emplace(face); - } - else { + if (interior_simplex_map.emplace(face, qr.intersection).second) queue.emplace(face); + } else { auto qrb = oracle.intersects_boundary(cof, triangulation); #ifdef GUDHI_COX_OUTPUT_TO_HTML mt_inserted_list.push_back(MT_inserted_info(qrb, cof, true)); -#endif +#endif // assert (qrb.success); // always a success - if (qrb.success) - boundary_simplex_map.emplace(cof, qrb.intersection); + if (qrb.success) boundary_simplex_map.emplace(cof, qrb.intersection); } } } @@ -206,7 +189,6 @@ public: /** \brief Empty constructor */ Manifold_tracing() {} - }; /** @@ -214,7 +196,7 @@ public: * that computes the set of k-simplices that intersect * a boundaryless implicit manifold given by an intersection oracle, where k * is the codimension of the manifold. - * The computation is based on the seed propagation --- it starts at the + * The computation is based on the seed propagation --- it starts at the * given seed points and then propagates along the manifold. * * \tparam Point_range Range of points of type Eigen::VectorXd. @@ -224,38 +206,30 @@ public: * Needs to be a model of the concept IntersectionOracle. * \tparam Out_simplex_map Needs to be Manifold_tracing<Triangulation_>::Out_simplex_map. * - * \param[in] seed_points The range of points on the manifold from which + * \param[in] seed_points The range of points on the manifold from which * the computation begins. * \param[in] triangulation The ambient triangulation. * \param[in] oracle The intersection oracle for the manifold. * The ambient dimension needs to match the dimension of the * triangulation. * \param[out] out_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the input manifold and the mapped values + * the input triangulation that intersect the input manifold and the mapped values * are the intersection points. * * \ingroup coxeter_triangulation */ -template <class Point_range, - class Triangulation, - class Intersection_oracle, - class Out_simplex_map> -void manifold_tracing_algorithm(const Point_range& seed_points, - const Triangulation& triangulation, - const Intersection_oracle& oracle, - Out_simplex_map& out_simplex_map) { +template <class Point_range, class Triangulation, class Intersection_oracle, class Out_simplex_map> +void manifold_tracing_algorithm(const Point_range& seed_points, const Triangulation& triangulation, + const Intersection_oracle& oracle, Out_simplex_map& out_simplex_map) { Manifold_tracing<Triangulation> mt; - mt.manifold_tracing_algorithm(seed_points, - triangulation, - oracle, - out_simplex_map); + mt.manifold_tracing_algorithm(seed_points, triangulation, oracle, out_simplex_map); } /** * \brief Static method for Manifold_tracing<Triangulation_>::manifold_tracing_algorithm * the dimensional manifold given by an intersection oracle, where k * is the codimension of the manifold. - * The computation is based on the seed propagation --- it starts at the + * The computation is based on the seed propagation --- it starts at the * given seed points and then propagates along the manifold. * * \tparam Point_range Range of points of type Eigen::VectorXd. @@ -265,41 +239,31 @@ void manifold_tracing_algorithm(const Point_range& seed_points, * Needs to be a model of the concept IntersectionOracle. * \tparam Out_simplex_map Needs to be Manifold_tracing<Triangulation_>::Out_simplex_map. * - * \param[in] seed_points The range of points on the manifold from which + * \param[in] seed_points The range of points on the manifold from which * the computation begins. * \param[in] triangulation The ambient triangulation. * \param[in] oracle The intersection oracle for the manifold. * The ambient dimension needs to match the dimension of the * triangulation. * \param[out] interior_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the relative interior of the input manifold + * the input triangulation that intersect the relative interior of the input manifold * and the mapped values are the intersection points. * \param[out] boundary_simplex_map The output map, where the keys are k-simplices in - * the input triangulation that intersect the boundary of the input manifold + * the input triangulation that intersect the boundary of the input manifold * and the mapped values are the intersection points. * * \ingroup coxeter_triangulation */ -template <class Point_range, - class Triangulation, - class Intersection_oracle, - class Out_simplex_map> -void manifold_tracing_algorithm(const Point_range& seed_points, - const Triangulation& triangulation, - const Intersection_oracle& oracle, - Out_simplex_map& interior_simplex_map, +template <class Point_range, class Triangulation, class Intersection_oracle, class Out_simplex_map> +void manifold_tracing_algorithm(const Point_range& seed_points, const Triangulation& triangulation, + const Intersection_oracle& oracle, Out_simplex_map& interior_simplex_map, Out_simplex_map& boundary_simplex_map) { Manifold_tracing<Triangulation> mt; - mt.manifold_tracing_algorithm(seed_points, - triangulation, - oracle, - interior_simplex_map, - boundary_simplex_map); + mt.manifold_tracing_algorithm(seed_points, triangulation, oracle, interior_simplex_map, boundary_simplex_map); } +} // namespace coxeter_triangulation -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation.h index b1f57b08..76438c91 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation.h @@ -19,7 +19,7 @@ namespace Gudhi { namespace coxeter_triangulation { -/** +/** * \class Permutahedral_representation * \brief A class that stores the permutahedral representation of a simplex * in a Coxeter triangulation or a Freudenthal-Kuhn triangulation. @@ -32,93 +32,72 @@ namespace coxeter_triangulation { * * \tparam Vertex_ needs to be a random-access range. * \tparam Ordered_set_partition_ needs to be a a random-access range that consists of - * random-access ranges. + * random-access ranges. */ -template <class Vertex_, - class Ordered_set_partition_> +template <class Vertex_, class Ordered_set_partition_> class Permutahedral_representation { - typedef Permutahedral_representation<Vertex_, Ordered_set_partition_> Self; -public: + public: /** \brief Type of the vertex. */ typedef Vertex_ Vertex; - + /** \brief Type of the ordered partition. */ typedef Ordered_set_partition_ OrderedSetPartition; - /** \brief Permutahedral_representation constructor from a vertex and an ordered set partition. * * @param[in] vertex Vertex. * @param[in] partition Ordered set partition. - * + * * \details If the size of vertex is d, the ranges in partition must consist * of the integers 0,...,d without repetition or collision between the ranges. */ Permutahedral_representation(const Vertex& vertex, const OrderedSetPartition& partition) - : vertex_(vertex), partition_(partition) {} + : vertex_(vertex), partition_(partition) {} /** \brief Constructor for an empty permutahedral representation that does not correspond * to any simplex. */ Permutahedral_representation() {} - + /** \brief Dimension of the simplex. */ - std::size_t dimension() const { - return partition_.size() - 1; - } + std::size_t dimension() const { return partition_.size() - 1; } /** \brief Lexicographically-minimal vertex. */ - Vertex& vertex() { - return vertex_; - } + Vertex& vertex() { return vertex_; } /** \brief Lexicographically-minimal vertex. */ - const Vertex& vertex() const { - return vertex_; - } + const Vertex& vertex() const { return vertex_; } /** \brief Ordered set partition. */ - OrderedSetPartition& partition() { - return partition_; - } + OrderedSetPartition& partition() { return partition_; } /** \brief Identifying vertex. */ - const OrderedSetPartition& partition() const { - return partition_; - } + const OrderedSetPartition& partition() const { return partition_; } /** \brief Equality operator. * Returns true if an only if both vertex and the ordered set partition coincide. */ bool operator==(const Permutahedral_representation& other) const { - if (dimension() != other.dimension()) - return false; - if (vertex_ != other.vertex_) - return false; + if (dimension() != other.dimension()) return false; + if (vertex_ != other.vertex_) return false; for (std::size_t k = 0; k < partition_.size(); ++k) - if (partition_[k] != other.partition_[k]) - return false; + if (partition_[k] != other.partition_[k]) return false; return true; } /** \brief Inequality operator. * Returns true if an only if either vertex or the ordered set partition are different. */ - bool operator!=(const Permutahedral_representation& other) const { - return !(*this == other); - } + bool operator!=(const Permutahedral_representation& other) const { return !(*this == other); } typedef Gudhi::coxeter_triangulation::Vertex_iterator<Self> Vertex_iterator; typedef boost::iterator_range<Vertex_iterator> Vertex_range; /** \brief Returns a range of vertices of the simplex. * The type of vertices is Vertex. */ - Vertex_range vertex_range() const { - return Vertex_range(Vertex_iterator(*this), - Vertex_iterator()); - } + Vertex_range vertex_range() const { return Vertex_range(Vertex_iterator(*this), Vertex_iterator()); } typedef Gudhi::coxeter_triangulation::Face_iterator<Self> Face_iterator; typedef boost::iterator_range<Face_iterator> Face_range; @@ -126,34 +105,29 @@ public: * @param[in] value_dim The dimension of the faces. Must be between 0 and the dimension of the simplex. */ Face_range face_range(std::size_t value_dim) const { - return Face_range(Face_iterator(*this, value_dim), - Face_iterator()); + return Face_range(Face_iterator(*this, value_dim), Face_iterator()); } /** \brief Returns a range of permutahedral representations of facets of the simplex. * The dimension of the simplex must be strictly positive. */ - Face_range facet_range() const { - return Face_range(Face_iterator(*this, dimension()-1), - Face_iterator()); - } - + Face_range facet_range() const { return Face_range(Face_iterator(*this, dimension() - 1), Face_iterator()); } + typedef Gudhi::coxeter_triangulation::Coface_iterator<Self> Coface_iterator; typedef boost::iterator_range<Coface_iterator> Coface_range; /** \brief Returns a range of permutahedral representations of cofaces of the simplex. - * @param[in] value_dim The dimension of the cofaces. Must be between the dimension of the simplex and the ambient dimension (the size of the vertex). + * @param[in] value_dim The dimension of the cofaces. Must be between the dimension of the simplex and the ambient + * dimension (the size of the vertex). */ Coface_range coface_range(std::size_t value_dim) const { - return Coface_range(Coface_iterator(*this, value_dim), - Coface_iterator()); + return Coface_range(Coface_iterator(*this, value_dim), Coface_iterator()); } /** \brief Returns a range of permutahedral representations of cofacets of the simplex. * The dimension of the simplex must be strictly different from the ambient dimension (the size of the vertex). */ Coface_range cofacet_range() const { - return Coface_range(Coface_iterator(*this, dimension()+1), - Coface_iterator()); + return Coface_range(Coface_iterator(*this, dimension() + 1), Coface_iterator()); } /** \brief Returns true, if the simplex is a face of other simplex. @@ -162,38 +136,31 @@ public: */ bool is_face_of(const Permutahedral_representation& other) const { using Part = typename OrderedSetPartition::value_type; - - if (other.dimension() < dimension()) - return false; + + if (other.dimension() < dimension()) return false; if (other.vertex_.size() != vertex_.size()) std::cerr << "Error: Permutahedral_representation::is_face_of: incompatible ambient dimensions.\n"; - + Vertex v_self = vertex_, v_other = other.vertex_; auto self_partition_it = partition_.begin(); auto other_partition_it = other.partition_.begin(); while (self_partition_it != partition_.end()) { while (other_partition_it != other.partition_.end() && v_self != v_other) { const Part& other_part = *other_partition_it++; - if (other_partition_it == other.partition_.end()) - return false; - for (const auto& k: other_part) - v_other[k]++; + if (other_partition_it == other.partition_.end()) return false; + for (const auto& k : other_part) v_other[k]++; } - if (other_partition_it == other.partition_.end()) - return false; + if (other_partition_it == other.partition_.end()) return false; const Part& self_part = *self_partition_it++; - if (self_partition_it == partition_.end()) - return true; - for (const auto& k: self_part) - v_self[k]++; + if (self_partition_it == partition_.end()) return true; + for (const auto& k : self_part) v_self[k]++; } return true; } -private: + private: Vertex vertex_; OrderedSetPartition partition_; - }; /** \brief Print a permutahedral representation to a stream. @@ -202,10 +169,8 @@ private: * @param[in] os The output stream. * @param[in] simplex A simplex represented by its permutahedral representation. */ -template <class Vertex, - class OrderedSetPartition> -std::ostream& operator<<(std::ostream& os, - const Permutahedral_representation<Vertex, OrderedSetPartition>& simplex) { +template <class Vertex, class OrderedSetPartition> +std::ostream& operator<<(std::ostream& os, const Permutahedral_representation<Vertex, OrderedSetPartition>& simplex) { // vertex part os << "("; if (simplex.vertex().empty()) { @@ -214,24 +179,21 @@ std::ostream& operator<<(std::ostream& os, } auto v_it = simplex.vertex().begin(); os << *v_it++; - for (; v_it != simplex.vertex().end(); ++v_it) - os << ", " << *v_it; + for (; v_it != simplex.vertex().end(); ++v_it) os << ", " << *v_it; os << ")"; - + // ordered partition part using Part = typename OrderedSetPartition::value_type; - auto print_part = - [&os](const Part& p) { - os << "{"; - if (p.empty()) { - os << "}"; - } - auto p_it = p.begin(); - os << *p_it++; - for (; p_it != p.end(); ++p_it) - os << ", " << *p_it; + auto print_part = [&os](const Part& p) { + os << "{"; + if (p.empty()) { os << "}"; - }; + } + auto p_it = p.begin(); + os << *p_it++; + for (; p_it != p.end(); ++p_it) os << ", " << *p_it; + os << "}"; + }; os << " ["; if (simplex.partition().empty()) { os << "]"; @@ -247,8 +209,8 @@ std::ostream& operator<<(std::ostream& os, return os; } -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi -#endif // PERMUTAHEDRAL_REPRESENTATION_H_ +#endif // PERMUTAHEDRAL_REPRESENTATION_H_ diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Combination_iterator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Combination_iterator.h index ce6a34ec..5f382e31 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Combination_iterator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Combination_iterator.h @@ -20,25 +20,20 @@ namespace coxeter_triangulation { typedef unsigned uint; -/** \brief Class that allows the user to generate combinations of +/** \brief Class that allows the user to generate combinations of * k elements in a set of n elements. * Based on the algorithm by Mifsud. -*/ -class Combination_iterator : public boost::iterator_facade< Combination_iterator, - std::vector<uint> const, - boost::forward_traversal_tag> { + */ +class Combination_iterator + : public boost::iterator_facade<Combination_iterator, std::vector<uint> const, boost::forward_traversal_tag> { typedef std::vector<uint> value_t; - -protected: + + protected: friend class boost::iterator_core_access; - - bool equal(Combination_iterator const& other) const { - return (is_end_ && other.is_end_); - } - value_t const& dereference() const { - return value_; - } + bool equal(Combination_iterator const& other) const { return (is_end_ && other.is_end_); } + + value_t const& dereference() const { return value_; } void increment() { if (value_[0] == n_ - k_) { @@ -51,25 +46,16 @@ protected: return; } for (; j > 0; --j) - if (value_[j-1] < n_ - k_ + j-1) { - value_[j-1]++; - for (uint s = j; s < k_; s++) - value_[s] = value_[j-1] + s - (j-1); - return; + if (value_[j - 1] < n_ - k_ + j - 1) { + value_[j - 1]++; + for (uint s = j; s < k_; s++) value_[s] = value_[j - 1] + s - (j - 1); + return; } } -public: - - Combination_iterator(const uint& n, const uint& k) - : - value_(k), - is_end_(n == 0), - n_(n), - k_(k) - { - for (uint i = 0; i < k; ++i) - value_[i] = i; + public: + Combination_iterator(const uint& n, const uint& k) : value_(k), is_end_(n == 0), n_(n), k_(k) { + for (uint i = 0; i < k; ++i) value_[i] = i; } // Used for the creating an end iterator @@ -78,22 +64,20 @@ public: void reinitialize() { if (n_ > 0) { is_end_ = false; - for (uint i = 0; i < n_; ++i) - value_[i] = i; + for (uint i = 0; i < n_; ++i) value_[i] = i; } } - + private: - value_t value_; // the dereference value - bool is_end_; // is true when the current permutation is the final one + value_t value_; // the dereference value + bool is_end_; // is true when the current permutation is the final one uint n_; uint k_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Integer_combination_iterator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Integer_combination_iterator.h index c4e86a36..3ee73754 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Integer_combination_iterator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Integer_combination_iterator.h @@ -20,37 +20,32 @@ namespace coxeter_triangulation { typedef unsigned uint; -/** \brief Class that allows the user to generate combinations of +/** \brief Class that allows the user to generate combinations of * k elements in a set of n elements. * Based on the algorithm by Mifsud. -*/ -class Integer_combination_iterator : public boost::iterator_facade< Integer_combination_iterator, - std::vector<uint> const, - boost::forward_traversal_tag> { + */ +class Integer_combination_iterator + : public boost::iterator_facade<Integer_combination_iterator, std::vector<uint> const, + boost::forward_traversal_tag> { using value_t = std::vector<uint>; - + private: friend class boost::iterator_core_access; - - bool equal(Integer_combination_iterator const& other) const { - return (is_end_ && other.is_end_); - } - value_t const& dereference() const { - return value_; - } + bool equal(Integer_combination_iterator const& other) const { return (is_end_ && other.is_end_); } + + value_t const& dereference() const { return value_; } void increment() { uint j1 = 0; uint s = 0; - while (value_[j1] == 0 && j1 < k_) - j1++; - uint j2 = j1+1; + while (value_[j1] == 0 && j1 < k_) j1++; + uint j2 = j1 + 1; while (value_[j2] == bounds_[j2]) { if (bounds_[j2] != 0) { - s += value_[j1]; - value_[j1] = 0; - j1 = j2; + s += value_[j1]; + value_[j1] = 0; + j1 = j2; } j2++; } @@ -70,20 +65,15 @@ class Integer_combination_iterator : public boost::iterator_facade< Integer_comb value_[i++] = s; } -public: + public: template <class Bound_range> Integer_combination_iterator(const uint& n, const uint& k, const Bound_range& bounds) - : - value_(k+2), - is_end_(n == 0 || k == 0), - n_(n), - k_(k) - { - bounds_.reserve(k+2); + : value_(k + 2), is_end_(n == 0 || k == 0), n_(n), k_(k) { + bounds_.reserve(k + 2); uint sum_radices = 0; - for (auto b: bounds) { + for (auto b : bounds) { bounds_.push_back(b); - sum_radices += b; + sum_radices += b; } bounds_.push_back(2); bounds_.push_back(1); @@ -100,26 +90,25 @@ public: } value_[i++] = s; - while (i < k_) - value_[i++] = 0; + while (i < k_) value_[i++] = 0; value_[k] = 1; - value_[k+1] = 0; + value_[k + 1] = 0; } // Used for the creating an end iterator Integer_combination_iterator() : is_end_(true), n_(0), k_(0) {} private: - value_t value_; // the dereference value - bool is_end_; // is true when the current integer combination is the final one + value_t value_; // the dereference value + bool is_end_; // is true when the current integer combination is the final one uint n_; uint k_; std::vector<uint> bounds_; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Ordered_set_partition_iterator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Ordered_set_partition_iterator.h index d6f9f121..866079fa 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Ordered_set_partition_iterator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Ordered_set_partition_iterator.h @@ -35,52 +35,40 @@ struct Ordered_set_partition { // Ordered_set_partition(const Set_partition_iterator& s_it, const Permutation_iterator& p_it) // : s_it_(s_it), p_it_(p_it) {} - - const std::vector<uint> operator[](const uint& i) const { - return (*s_it_)[(*p_it_)[i]]; - } - - std::size_t size() const { - return s_it_->size(); - } - + + const std::vector<uint> operator[](const uint& i) const { return (*s_it_)[(*p_it_)[i]]; } + + std::size_t size() const { return s_it_->size(); } }; /** \brief Class that allows the user to generate set partitions of a set {0,...,n-1} in k parts. - * -*/ -class Ordered_set_partition_iterator : public boost::iterator_facade< Ordered_set_partition_iterator, - Ordered_set_partition const, - boost::forward_traversal_tag> { + * + */ +class Ordered_set_partition_iterator + : public boost::iterator_facade<Ordered_set_partition_iterator, Ordered_set_partition const, + boost::forward_traversal_tag> { using value_t = Ordered_set_partition; - + private: friend class boost::iterator_core_access; - bool equal(Ordered_set_partition_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Ordered_set_partition_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void increment() { if (++value_.p_it_ == p_end_) { if (++value_.s_it_ == s_end_) { is_end_ = true; return; - } - else + } else value_.p_it_.reinitialize(); } } - + public: Ordered_set_partition_iterator(const uint& n, const uint& k) - : - value_({Set_partition_iterator(n,k), Permutation_iterator(k)}), - is_end_(n == 0) {} + : value_({Set_partition_iterator(n, k), Permutation_iterator(k)}), is_end_(n == 0) {} // Used for the creating an end iterator Ordered_set_partition_iterator() : is_end_(true) {} @@ -90,16 +78,16 @@ class Ordered_set_partition_iterator : public boost::iterator_facade< Ordered_se value_.p_it_.reinitialize(); value_.s_it_.reinitialize(); } - + private: - Set_partition_iterator s_end_; // Set partition iterator and the corresponding end iterator - Permutation_iterator p_end_; // Permutation iterator and the corresponding end iterator - value_t value_; // the dereference value - bool is_end_; // is true when the current permutation is the final one + Set_partition_iterator s_end_; // Set partition iterator and the corresponding end iterator + Permutation_iterator p_end_; // Permutation iterator and the corresponding end iterator + value_t value_; // the dereference value + bool is_end_; // is true when the current permutation is the final one }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutahedral_representation_iterators.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutahedral_representation_iterators.h index d42e892a..9263c67e 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutahedral_representation_iterators.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutahedral_representation_iterators.h @@ -36,9 +36,9 @@ namespace coxeter_triangulation { * * Forward iterator, 'value_type' is Permutahedral_representation::Vertex.*/ template <class Permutahedral_representation> -class Vertex_iterator : public boost::iterator_facade< Vertex_iterator<Permutahedral_representation>, - typename Permutahedral_representation::Vertex const, - boost::forward_traversal_tag> { +class Vertex_iterator + : public boost::iterator_facade<Vertex_iterator<Permutahedral_representation>, + typename Permutahedral_representation::Vertex const, boost::forward_traversal_tag> { private: friend class boost::iterator_core_access; @@ -47,50 +47,40 @@ class Vertex_iterator : public boost::iterator_facade< Vertex_iterator<Permutahe using value_t = Vertex; - - bool equal(Vertex_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Vertex_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void update_value() { std::size_t d = value_.size(); - for (auto i: *o_it_) + for (auto i : *o_it_) if (i != d) value_[i]++; else - for (std::size_t j = 0; j < d; j++) - value_[j]--; + for (std::size_t j = 0; j < d; j++) value_[j]--; } - + void increment() { - if (is_end_) - return; + if (is_end_) return; update_value(); - if (++o_it_ == o_end_) - is_end_ = true; + if (++o_it_ == o_end_) is_end_ = true; } -public: + public: Vertex_iterator(const Permutahedral_representation& simplex) - : - o_it_(simplex.partition().begin()), - o_end_(simplex.partition().end()), - value_(simplex.vertex()), - is_end_(o_it_ == o_end_) - {} + : o_it_(simplex.partition().begin()), + o_end_(simplex.partition().end()), + value_(simplex.vertex()), + is_end_(o_it_ == o_end_) {} Vertex_iterator() : is_end_(true) {} private: - typename Ordered_partition::const_iterator o_it_, o_end_; + typename Ordered_partition::const_iterator o_it_, o_end_; value_t value_; bool is_end_; -}; // Vertex_iterator +}; // Vertex_iterator /*---------------------------------------------------------------------------*/ /* \brief Iterator over the k-faces of a simplex @@ -98,24 +88,19 @@ public: * * Forward iterator, value_type is Permutahedral_representation. */ template <class Permutahedral_representation> -class Face_iterator : public boost::iterator_facade< Face_iterator<Permutahedral_representation>, - Permutahedral_representation const, - boost::forward_traversal_tag> { +class Face_iterator : public boost::iterator_facade<Face_iterator<Permutahedral_representation>, + Permutahedral_representation const, boost::forward_traversal_tag> { using value_t = Permutahedral_representation; - + private: friend class boost::iterator_core_access; using Vertex = typename Permutahedral_representation::Vertex; using Ordered_partition = typename Permutahedral_representation::OrderedSetPartition; - bool equal(Face_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Face_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void increment() { if (++c_it_ == c_end_) { @@ -130,31 +115,30 @@ class Face_iterator : public boost::iterator_facade< Face_iterator<Permutahedral value_ = face_from_indices<Permutahedral_representation>(simplex_, *c_it_); } -public: + public: Face_iterator(const Permutahedral_representation& simplex, const uint& k) - : simplex_(simplex), - k_(k), - l_(simplex.dimension()), - c_it_(l_+1, k_+1), - is_end_(k_ > l_), - value_({Vertex(simplex.vertex().size()), Ordered_partition(k+1)}) - { + : simplex_(simplex), + k_(k), + l_(simplex.dimension()), + c_it_(l_ + 1, k_ + 1), + is_end_(k_ > l_), + value_({Vertex(simplex.vertex().size()), Ordered_partition(k + 1)}) { update_value(); } // Used for the creating an end iterator Face_iterator() : is_end_(true) {} - + private: - Permutahedral_representation simplex_; // Input simplex + Permutahedral_representation simplex_; // Input simplex uint k_; - uint l_; // Dimension of the input simplex - Combination_iterator c_it_, c_end_; // indicates the vertices in the current face + uint l_; // Dimension of the input simplex + Combination_iterator c_it_, c_end_; // indicates the vertices in the current face - bool is_end_; // is true when the current permutation is the final one - value_t value_; // the dereference value + bool is_end_; // is true when the current permutation is the final one + value_t value_; // the dereference value -}; // Face_iterator +}; // Face_iterator /*---------------------------------------------------------------------------*/ /* \brief Iterator over the k-cofaces of a simplex @@ -162,9 +146,9 @@ public: * * Forward iterator, value_type is Permutahedral_representation. */ template <class Permutahedral_representation> -class Coface_iterator : public boost::iterator_facade< Coface_iterator<Permutahedral_representation>, - Permutahedral_representation const, - boost::forward_traversal_tag> { +class Coface_iterator + : public boost::iterator_facade<Coface_iterator<Permutahedral_representation>, Permutahedral_representation const, + boost::forward_traversal_tag> { using value_t = Permutahedral_representation; private: @@ -173,75 +157,62 @@ class Coface_iterator : public boost::iterator_facade< Coface_iterator<Permutahe using Vertex = typename Permutahedral_representation::Vertex; using Ordered_partition = typename Permutahedral_representation::OrderedSetPartition; - bool equal(Coface_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Coface_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void increment() { uint i = 0; - for (; i < k_+1; i++) { - if (++(o_its_[i]) != o_end_) - break; + for (; i < k_ + 1; i++) { + if (++(o_its_[i]) != o_end_) break; } - if (i == k_+1) { + if (i == k_ + 1) { if (++i_it_ == i_end_) { is_end_ = true; return; } o_its_.clear(); for (uint j = 0; j < k_ + 1; j++) - o_its_.emplace_back(Ordered_set_partition_iterator(simplex_.partition()[j].size(), (*i_it_)[j]+1)); - } - else - for (uint j = 0; j < i; j++) - o_its_[j].reinitialize(); + o_its_.emplace_back(Ordered_set_partition_iterator(simplex_.partition()[j].size(), (*i_it_)[j] + 1)); + } else + for (uint j = 0; j < i; j++) o_its_[j].reinitialize(); update_value(); } void update_value() { value_.vertex() = simplex_.vertex(); - for (auto& p: value_.partition()) - p.clear(); - uint u_ = 0; // the part in o_its_[k_] that contains t_ + for (auto& p : value_.partition()) p.clear(); + uint u_ = 0; // the part in o_its_[k_] that contains t_ for (; u_ <= (*i_it_)[k_]; u_++) { auto range = (*o_its_[k_])[u_]; - if (std::find(range.begin(), range.end(), t_) != range.end()) - break; + if (std::find(range.begin(), range.end(), t_) != range.end()) break; } uint i = 0; - for (uint j = u_+1; j <= (*i_it_)[k_]; j++, i++) - for (uint b: (*o_its_[k_])[j]) { + for (uint j = u_ + 1; j <= (*i_it_)[k_]; j++, i++) + for (uint b : (*o_its_[k_])[j]) { uint c = simplex_.partition()[k_][b]; value_.partition()[i].push_back(c); value_.vertex()[c]--; } for (uint h = 0; h < k_; h++) for (uint j = 0; j <= (*i_it_)[h]; j++, i++) { - for (uint b: (*o_its_[h])[j]) - value_.partition()[i].push_back(simplex_.partition()[h][b]); + for (uint b : (*o_its_[h])[j]) value_.partition()[i].push_back(simplex_.partition()[h][b]); } for (uint j = 0; j <= u_; j++, i++) - for (uint b: (*o_its_[k_])[j]) - value_.partition()[i].push_back(simplex_.partition()[k_][b]); + for (uint b : (*o_its_[k_])[j]) value_.partition()[i].push_back(simplex_.partition()[k_][b]); // sort the values in each part (probably not needed) - for (auto& part: value_.partition()) - std::sort(part.begin(), part.end()); + for (auto& part : value_.partition()) std::sort(part.begin(), part.end()); } public: Coface_iterator(const Permutahedral_representation& simplex, const uint& l) - : simplex_(simplex), - d_(simplex.vertex().size()), - l_(l), - k_(simplex.dimension()), - i_it_(l_-k_ , k_+1, Size_range<Ordered_partition>(simplex.partition())), - is_end_(k_ > l_), - value_({Vertex(d_), Ordered_partition(l_+1)}) - { + : simplex_(simplex), + d_(simplex.vertex().size()), + l_(l), + k_(simplex.dimension()), + i_it_(l_ - k_, k_ + 1, Size_range<Ordered_partition>(simplex.partition())), + is_end_(k_ > l_), + value_({Vertex(d_), Ordered_partition(l_ + 1)}) { uint j = 0; for (; j < simplex_.partition()[k_].size(); j++) if (simplex_.partition()[k_][j] == d_) { @@ -253,8 +224,8 @@ class Coface_iterator : public boost::iterator_facade< Coface_iterator<Permutahe is_end_ = true; return; } - for (uint i = 0; i < k_+1; i++) - o_its_.emplace_back(Ordered_set_partition_iterator(simplex_.partition()[i].size(), (*i_it_)[i]+1)); + for (uint i = 0; i < k_ + 1; i++) + o_its_.emplace_back(Ordered_set_partition_iterator(simplex_.partition()[i].size(), (*i_it_)[i] + 1)); update_value(); } @@ -262,22 +233,22 @@ class Coface_iterator : public boost::iterator_facade< Coface_iterator<Permutahe Coface_iterator() : is_end_(true) {} private: - Permutahedral_representation simplex_; // Input simplex - uint d_; // Ambient dimension - uint l_; // Dimension of the coface - uint k_; // Dimension of the input simplex - uint t_; // The position of d in simplex_.partition()[k_] - Integer_combination_iterator i_it_, i_end_; // indicates in how many parts each simplex_[i] is subdivided - std::vector<Ordered_set_partition_iterator> o_its_; // indicates subdivision for each simplex_[i] - Ordered_set_partition_iterator o_end_; // one end for all o_its_ + Permutahedral_representation simplex_; // Input simplex + uint d_; // Ambient dimension + uint l_; // Dimension of the coface + uint k_; // Dimension of the input simplex + uint t_; // The position of d in simplex_.partition()[k_] + Integer_combination_iterator i_it_, i_end_; // indicates in how many parts each simplex_[i] is subdivided + std::vector<Ordered_set_partition_iterator> o_its_; // indicates subdivision for each simplex_[i] + Ordered_set_partition_iterator o_end_; // one end for all o_its_ + + bool is_end_; // is true when the current permutation is the final one + value_t value_; // the dereference value - bool is_end_; // is true when the current permutation is the final one - value_t value_; // the dereference value +}; // Coface_iterator -}; // Coface_iterator +} // namespace coxeter_triangulation -} // namespace coxeter_triangulation +} // namespace Gudhi -} // namespace Gudhi - #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutation_iterator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutation_iterator.h index e0142bf4..0f91d41c 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutation_iterator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Permutation_iterator.h @@ -24,32 +24,27 @@ typedef unsigned uint; /** \brief Class that allows the user to generate permutations. * Based on the optimization of the Heap's algorithm by Sedgewick. -*/ -class Permutation_iterator : public boost::iterator_facade< Permutation_iterator, - std::vector<uint> const, - boost::forward_traversal_tag> { + */ +class Permutation_iterator + : public boost::iterator_facade<Permutation_iterator, std::vector<uint> const, boost::forward_traversal_tag> { using value_t = std::vector<uint>; - + private: friend class boost::iterator_core_access; - bool equal(Permutation_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Permutation_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void swap_two_indices(std::size_t i, std::size_t j) { uint t = value_[i]; value_[i] = value_[j]; value_[j] = t; } - + void elementary_increment() { uint j = 0; - while (d_[j] == j+1) { + while (d_[j] == j + 1) { d_[j] = 0; ++j; } @@ -57,8 +52,8 @@ class Permutation_iterator : public boost::iterator_facade< Permutation_iterator is_end_ = true; return; } - uint k = j+1; - uint x = (k%2 ? d_[j] : 0); + uint k = j + 1; + uint x = (k % 2 ? d_[j] : 0); swap_two_indices(k, x); ++d_[j]; } @@ -66,12 +61,11 @@ class Permutation_iterator : public boost::iterator_facade< Permutation_iterator void elementary_increment_optim_3() { if (ct_ != 0) { --ct_; - swap_two_indices(1 + (ct_%2), 0); - } - else { + swap_two_indices(1 + (ct_ % 2), 0); + } else { ct_ = 5; uint j = 2; - while (d_[j] == j+1) { + while (d_[j] == j + 1) { d_[j] = 0; ++j; } @@ -79,59 +73,48 @@ class Permutation_iterator : public boost::iterator_facade< Permutation_iterator is_end_ = true; return; } - uint k = j+1; - uint x = (k%2 ? d_[j] : 0); + uint k = j + 1; + uint x = (k % 2 ? d_[j] : 0); swap_two_indices(k, x); ++d_[j]; } } - + void increment() { if (optim_3_) elementary_increment_optim_3(); else - elementary_increment(); + elementary_increment(); } -public: - - Permutation_iterator(const uint& n) - : - value_(n), - is_end_(n == 0), - optim_3_(n >= 3), - n_(n), - d_(n), - ct_(5) - { + public: + Permutation_iterator(const uint& n) : value_(n), is_end_(n == 0), optim_3_(n >= 3), n_(n), d_(n), ct_(5) { for (uint i = 0; i < n; ++i) { value_[i] = i; d_[i] = 0; } - if (n > 0) - d_[n-1] = -1; + if (n > 0) d_[n - 1] = -1; } // Used for the creating an end iterator Permutation_iterator() : is_end_(true), n_(0) {} void reinitialize() { - if (n_ > 0) - is_end_ = false; + if (n_ > 0) is_end_ = false; } private: - value_t value_; // the dereference value - bool is_end_; // is true when the current permutation is the final one - bool optim_3_; // true if n>=3. for n >= 3, the algorithm is optimized + value_t value_; // the dereference value + bool is_end_; // is true when the current permutation is the final one + bool optim_3_; // true if n>=3. for n >= 3, the algorithm is optimized uint n_; - std::vector<uint> d_; // mix radix digits with radix [2 3 4 ... n-1 (sentinel=-1)] - uint ct_; // counter with values in {0,...,5} used in the n>=3 optimization. + std::vector<uint> d_; // mix radix digits with radix [2 3 4 ... n-1 (sentinel=-1)] + uint ct_; // counter with values in {0,...,5} used in the n>=3 optimization. }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Set_partition_iterator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Set_partition_iterator.h index bd1770bc..94ac10c2 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Set_partition_iterator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Set_partition_iterator.h @@ -22,40 +22,32 @@ namespace coxeter_triangulation { typedef unsigned uint; /** \brief Class that allows the user to generate set partitions of a set {0,...,n-1} in k parts. - * -*/ -class Set_partition_iterator : public boost::iterator_facade< Set_partition_iterator, - std::vector<std::vector<uint> > const, - boost::forward_traversal_tag> { + * + */ +class Set_partition_iterator + : public boost::iterator_facade<Set_partition_iterator, std::vector<std::vector<uint>> const, + boost::forward_traversal_tag> { using value_t = std::vector<std::vector<uint>>; private: friend class boost::iterator_core_access; - bool equal(Set_partition_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Set_partition_iterator const& other) const { return (is_end_ && other.is_end_); } - value_t const& dereference() const { - return value_; - } + value_t const& dereference() const { return value_; } void update_value() { - for (uint i = 0; i < k_; i++) - value_[i].clear(); - for (uint i = 0; i < n_; i++) - value_[rgs_[i]].push_back(i); + for (uint i = 0; i < k_; i++) value_[i].clear(); + for (uint i = 0; i < n_; i++) value_[rgs_[i]].push_back(i); } - + void increment() { if (k_ <= 1) { is_end_ = true; return; } uint i = n_ - 1; - while (rgs_[i] + 1 > max_[i] || - rgs_[i] + 1 >= k_) - i--; + while (rgs_[i] + 1 > max_[i] || rgs_[i] + 1 >= k_) i--; if (i == 0) { is_end_ = true; return; @@ -63,14 +55,13 @@ class Set_partition_iterator : public boost::iterator_facade< Set_partition_iter rgs_[i]++; uint mm = max_[i]; mm += (rgs_[i] >= mm); - max_[i+1] = mm; + max_[i + 1] = mm; while (++i < n_) { rgs_[i] = 0; - max_[i+1] = mm; + max_[i + 1] = mm; } uint p = k_; - if (mm < p) - do { + if (mm < p) do { max_[i] = p; --i; --p; @@ -78,26 +69,17 @@ class Set_partition_iterator : public boost::iterator_facade< Set_partition_iter } while (max_[i] < p); update_value(); } - + public: Set_partition_iterator(const uint& n, const uint& k) - : - value_(k), - rgs_(n, 0), - max_(n+1), - is_end_(n == 0), - n_(n), - k_(k) - { + : value_(k), rgs_(n, 0), max_(n + 1), is_end_(n == 0), n_(n), k_(k) { max_[0] = std::numeric_limits<uint>::max(); - for (uint i = 0; i <= n-k; ++i) - value_[0].push_back(i); - for (uint i = n-k+1, j = 1; i < n; ++i, ++j) { + for (uint i = 0; i <= n - k; ++i) value_[0].push_back(i); + for (uint i = n - k + 1, j = 1; i < n; ++i, ++j) { rgs_[i] = j; value_[j].push_back(i); } - for (uint i = 1; i <= n; i++) - max_[i] = rgs_[i-1] + 1; + for (uint i = 1; i <= n; i++) max_[i] = rgs_[i - 1] + 1; update_value(); } @@ -105,30 +87,25 @@ class Set_partition_iterator : public boost::iterator_facade< Set_partition_iter Set_partition_iterator() : is_end_(true), n_(0), k_(0) {} void reinitialize() { - if (n_ > 0) - is_end_ = false; - for (uint i = 0; i <= n_-k_; ++i) - rgs_[i] = 0; - for (uint i = n_-k_+1, j = 1; i < n_; ++i, ++j) - rgs_[i] = j; - for (uint i = 1; i <= n_; i++) - max_[i] = rgs_[i-1] + 1; + if (n_ > 0) is_end_ = false; + for (uint i = 0; i <= n_ - k_; ++i) rgs_[i] = 0; + for (uint i = n_ - k_ + 1, j = 1; i < n_; ++i, ++j) rgs_[i] = j; + for (uint i = 1; i <= n_; i++) max_[i] = rgs_[i - 1] + 1; update_value(); } private: - value_t value_; // the dereference value - std::vector<uint> rgs_; // restricted growth string - std::vector<uint> max_; // max_[i] = max(rgs_[0],...,rgs[i-1]) + 1 - bool is_end_; // is true when the current permutation is the final one + value_t value_; // the dereference value + std::vector<uint> rgs_; // restricted growth string + std::vector<uint> max_; // max_[i] = max(rgs_[0],...,rgs[i-1]) + 1 + bool is_end_; // is true when the current permutation is the final one uint n_; uint k_; }; -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Simplex_comparator.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Simplex_comparator.h index b8925d96..905d68d5 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Simplex_comparator.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Simplex_comparator.h @@ -15,48 +15,40 @@ namespace Gudhi { namespace coxeter_triangulation { -/** \class Simplex_comparator +/** \class Simplex_comparator * \brief A comparator class for Permutahedral_representation. * The comparison is in lexicographic order first on * vertices and then on ordered partitions with sorted parts. * The lexicographic order forces that any face is larger than * a coface. * - * \tparam Permutahdral_representation_ Needs to be + * \tparam Permutahdral_representation_ Needs to be * Permutahedral_representation<Vertex_, Ordered_set_partition_> * * \ingroup coxeter_triangulation */ template <class Permutahedral_representation_> struct Simplex_comparator { - /** \brief Comparison between two permutahedral representations. * Both permutahedral representations need to be valid and * the vertices of both permutahedral representations need to be of the same size. */ - bool operator() (const Permutahedral_representation_& lhs, - const Permutahedral_representation_& rhs) const { - if (lhs.vertex() < rhs.vertex()) - return true; - if (lhs.vertex() > rhs.vertex()) - return false; - - if (lhs.partition().size() > rhs.partition().size()) - return true; - if (lhs.partition().size() < rhs.partition().size()) - return false; - - if (lhs.partition() < rhs.partition()) - return true; - if (lhs.partition() > rhs.partition()) - return false; - + bool operator()(const Permutahedral_representation_& lhs, const Permutahedral_representation_& rhs) const { + if (lhs.vertex() < rhs.vertex()) return true; + if (lhs.vertex() > rhs.vertex()) return false; + + if (lhs.partition().size() > rhs.partition().size()) return true; + if (lhs.partition().size() < rhs.partition().size()) return false; + + if (lhs.partition() < rhs.partition()) return true; + if (lhs.partition() > rhs.partition()) return false; + return false; } }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Size_range.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Size_range.h index f41335e9..c43effc8 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Size_range.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/Size_range.h @@ -22,33 +22,26 @@ namespace coxeter_triangulation { /** \brief Auxillary iterator class for sizes of parts in an ordered set partition. */ template <class T_it> -class Size_iterator : public boost::iterator_facade< Size_iterator<T_it>, - std::size_t const, - boost::forward_traversal_tag> { +class Size_iterator + : public boost::iterator_facade<Size_iterator<T_it>, std::size_t const, boost::forward_traversal_tag> { friend class boost::iterator_core_access; private: - bool equal(Size_iterator const& other) const { - return (is_end_ && other.is_end_); - } + bool equal(Size_iterator const& other) const { return (is_end_ && other.is_end_); } - std::size_t const& dereference() const { - return value_; - } + std::size_t const& dereference() const { return value_; } void increment() { if (++t_it_ == t_end_) { is_end_ = true; return; } - value_ = t_it_->size()-1; + value_ = t_it_->size() - 1; } public: - Size_iterator(const T_it& t_begin, const T_it& t_end) - : t_it_(t_begin), t_end_(t_end), is_end_(t_begin == t_end) { - if (!is_end_) - value_ = t_it_->size()-1; + Size_iterator(const T_it& t_begin, const T_it& t_end) : t_it_(t_begin), t_end_(t_end), is_end_(t_begin == t_end) { + if (!is_end_) value_ = t_it_->size() - 1; } private: @@ -63,25 +56,18 @@ class Size_range { public: typedef Size_iterator<typename T::const_iterator> iterator; - + Size_range(const T& t) : t_(t) {} - - std::size_t operator[](std::size_t i) const { - return t_[i].size()-1; - } - iterator begin() const { - return iterator(t_.begin(), t_.end()); - } + std::size_t operator[](std::size_t i) const { return t_[i].size() - 1; } - iterator end() const { - return iterator(t_.end(), t_.end()); - } + iterator begin() const { return iterator(t_.begin(), t_.end()); } + iterator end() const { return iterator(t_.end(), t_.end()); } }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/face_from_indices.h b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/face_from_indices.h index 21ee3c8b..47120689 100644 --- a/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/face_from_indices.h +++ b/src/Coxeter_triangulation/include/gudhi/Permutahedral_representation/face_from_indices.h @@ -11,23 +11,23 @@ #ifndef PERMUTAHEDRAL_REPRESENTATION_FACE_FROM_INDICES_H_ #define PERMUTAHEDRAL_REPRESENTATION_FACE_FROM_INDICES_H_ +#include <cstdlib> // for std::size_t +#include <algorithm> + namespace Gudhi { namespace coxeter_triangulation { -#include <cstdlib> // for std::size_t -#include <algorithm> - - /** \brief Computes the permutahedral representation of a face of a given simplex - * and a range of the vertex indices that compose the face. - * - * \tparam Permutahedral_representation has to be Permutahedral_representation - * \tparam Index_range is a range of unsigned integers taking values in 0,...,k, - * where k is the dimension of the simplex simplex. - * - * @param[in] simplex Input simplex. - * @param[in] indices Input range of indices. - */ +/** \brief Computes the permutahedral representation of a face of a given simplex + * and a range of the vertex indices that compose the face. + * + * \tparam Permutahedral_representation has to be Permutahedral_representation + * \tparam Index_range is a range of unsigned integers taking values in 0,...,k, + * where k is the dimension of the simplex simplex. + * + * @param[in] simplex Input simplex. + * @param[in] indices Input range of indices. + */ template <class Permutahedral_representation, class Index_range> Permutahedral_representation face_from_indices(const Permutahedral_representation& simplex, const Index_range& indices) { @@ -38,34 +38,29 @@ Permutahedral_representation face_from_indices(const Permutahedral_representatio Permutahedral_representation value; std::size_t d = simplex.vertex().size(); value.vertex() = simplex.vertex(); - std::size_t k = indices.size()-1; - value.partition().resize(k+1); - std::size_t l = simplex.partition().size()-1; - for (std::size_t h = 1; h < k+1; h++) - for (range_index i = indices[h-1]; i < indices[h]; i++) - for (part_index j: simplex.partition()[i]) - value.partition()[h-1].push_back(j); - for (range_index i = indices[k]; i < l+1; i++) - for (part_index j: simplex.partition()[i]) - value.partition()[k].push_back(j); + std::size_t k = indices.size() - 1; + value.partition().resize(k + 1); + std::size_t l = simplex.partition().size() - 1; + for (std::size_t h = 1; h < k + 1; h++) + for (range_index i = indices[h - 1]; i < indices[h]; i++) + for (part_index j : simplex.partition()[i]) value.partition()[h - 1].push_back(j); + for (range_index i = indices[k]; i < l + 1; i++) + for (part_index j : simplex.partition()[i]) value.partition()[k].push_back(j); for (range_index i = 0; i < indices[0]; i++) - for (part_index j: simplex.partition()[i]) { + for (part_index j : simplex.partition()[i]) { if (j != d) - value.vertex()[j]++; + value.vertex()[j]++; else - for (std::size_t l = 0; l < d; l++) - value.vertex()[l]--; + for (std::size_t l = 0; l < d; l++) value.vertex()[l]--; value.partition()[k].push_back(j); } // sort the values in each part (probably not needed) - for (auto& part: value.partition()) - std::sort(part.begin(), part.end()); + for (auto& part : value.partition()) std::sort(part.begin(), part.end()); return value; } -} // namespace coxeter_triangulation - -} // namespace Gudhi +} // namespace coxeter_triangulation +} // namespace Gudhi #endif diff --git a/src/Coxeter_triangulation/include/gudhi/Query_result.h b/src/Coxeter_triangulation/include/gudhi/Query_result.h index ecc40c25..5543c2fb 100644 --- a/src/Coxeter_triangulation/include/gudhi/Query_result.h +++ b/src/Coxeter_triangulation/include/gudhi/Query_result.h @@ -16,7 +16,7 @@ namespace Gudhi { namespace coxeter_triangulation { /** \class Query_result - * \brief The result of a query by an oracle such as Implicit_manifold_intersection_oracle. + * \brief The result of a query by an oracle such as Implicit_manifold_intersection_oracle. * * \tparam Simplex_handle The class of the query simplex. * @@ -33,8 +33,8 @@ struct Query_result { bool success; }; -} // namespace coxeter_triangulation +} // namespace coxeter_triangulation -} // namespace Gudhi +} // namespace Gudhi #endif |