diff options
author | skachano <skachano@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2016-10-05 15:02:50 +0000 |
---|---|---|
committer | skachano <skachano@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2016-10-05 15:02:50 +0000 |
commit | 93654d5708654a6071c1775580f625da625a08a8 (patch) | |
tree | 50cbe1f4778d6ee88d81a4f66c83782413f9b7a1 /src/Witness_complex/include | |
parent | 1de8ffacd531550f0ce5e871ec0f69924df3ee44 (diff) |
Junk out
git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/relaxed-witness@1646 636b058d-ea47-450e-bf9e-a15bfbe3eedb
Former-commit-id: ce2dfd7033a6a8366f9861ced4ed64ac41dfb82e
Diffstat (limited to 'src/Witness_complex/include')
12 files changed, 1 insertions, 2760 deletions
diff --git a/src/Witness_complex/include/gudhi/A0_complex.h b/src/Witness_complex/include/gudhi/A0_complex.h deleted file mode 100644 index 2018e1c8..00000000 --- a/src/Witness_complex/include/gudhi/A0_complex.h +++ /dev/null @@ -1,337 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef A0_COMPLEX_H_ -#define A0_COMPLEX_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include "gudhi/Simplex_tree.h" -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -// Needed for nearest neighbours -#include <CGAL/Cartesian_d.h> -#include <CGAL/Search_traits.h> -#include <CGAL/Search_traits_adapter.h> -#include <CGAL/property_map.h> -#include <CGAL/Epick_d.h> -#include <CGAL/Orthogonal_k_neighbor_search.h> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplicial_complex > -class A0_complex { - -private: - struct Active_witness { - int witness_id; - int landmark_id; - - Active_witness(int witness_id_, int landmark_id_) - : witness_id(witness_id_), - landmark_id(landmark_id_) { } - }; - -private: - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - typedef typename Simplicial_complex::Filtration_value FT; - - typedef std::vector< double > Point_t; - typedef std::vector< Point_t > Point_Vector; - - // typedef typename Simplicial_complex::Filtration_value Filtration_value; - typedef std::vector< Vertex_handle > typeVectorVertex; - typedef std::pair< typeVectorVertex, Filtration_value> typeSimplex; - typedef std::pair< Simplex_handle, bool > typePairSimplexBool; - - typedef int Witness_id; - typedef int Landmark_id; - typedef std::list< Vertex_handle > ActiveWitnessList; - - private: - int nbL; // Number of landmarks - Simplicial_complex& sc; // Simplicial complex - - public: - /** @name Constructor - */ - - //@{ - - /** - * \brief Iterative construction of the relaxed witness complex. - * \details The witness complex is written in sc_ basing on a matrix knn - * of k nearest neighbours of the form {witnesses}x{landmarks} and - * and a matrix distances of distances to these landmarks from witnesses. - * The parameter alpha defines relaxation and - * limD defines the - * - * The line lengths in one given matrix can differ, - * however both matrices have the same corresponding line lengths. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges. - * - * Constructor takes into account at most (dim+1) - * first landmarks from each landmark range to construct simplices. - * - * Landmarks are supposed to be in [0,nbL_-1] - */ - - template< typename KNearestNeighbours > - A0_complex(std::vector< std::vector<double> > const & distances, - KNearestNeighbours const & knn, - Simplicial_complex & sc_, - int nbL_, - double alpha2, - unsigned limD) : nbL(nbL_), sc(sc_) { - int nbW = knn.size(); - typeVectorVertex vv; - //int counter = 0; - /* The list of still useful witnesses - * it will diminuish in the course of iterations - */ - ActiveWitnessList active_w;// = new ActiveWitnessList(); - for (int i = 0; i != nbL; ++i) { - // initial fill of 0-dimensional simplices - // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore - //counter++; - vv = {i}; - sc.insert_simplex(vv, Filtration_value(0.0)); - /* TODO Error if not inserted : normally no need here though*/ - } - for (int i=0; i != nbW; ++i) { - // int i_end = limD+1; - // if (knn[i].size() < limD+1) - // i_end = knn[i].size(); - // double dist_wL = *(distances[i].begin()); - // while (distances[i][i_end] > dist_wL + alpha2) - // i_end--; - // add_all_witnessed_faces(distances[i].begin(), - // knn[i].begin(), - // knn[i].begin() + i_end + 1); - unsigned j_end = 0; - while (j_end < distances[i].size() && j_end <= limD && distances[i][j_end] <= distances[i][0] + alpha2) { - std::vector<int> simplex; - for (unsigned j = 0; j <= j_end; ++j) - simplex.push_back(knn[i][j]); - assert(distances[i][j_end] - distances[i][0] >= 0); - sc.insert_simplex_and_subfaces(simplex, distances[i][j_end] - distances[i][0]); - j_end++; - } - } - sc.set_dimension(limD); - } - - //@} - -private: - /* \brief Adds recursively all the faces of a certain dimension dim witnessed by the same witness - * Iterator is needed to know until how far we can take landmarks to form simplexes - * simplex is the prefix of the simplexes to insert - * The output value indicates if the witness rests active or not - */ - void add_all_witnessed_faces(std::vector<double>::const_iterator curr_d, - std::vector<int>::const_iterator curr_l, - std::vector<int>::const_iterator end) - { - std::vector<int> simplex; - std::vector<int>::const_iterator l_end = curr_l; - for (; l_end != end; ++l_end) { - std::vector<int>::const_iterator l_it = curr_l; - std::vector<double>::const_iterator d_it = curr_d; - simplex = {}; - for (; l_it != l_end; ++l_it, ++d_it) - simplex.push_back(*l_it); - sc.insert_simplex_and_subfaces(simplex, *(d_it--) - *curr_d); - } - } - - /** \brief Check if the facets of the k-dimensional simplex witnessed - * by witness witness_id are already in the complex. - * inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id - */ - bool all_faces_in(std::vector<int>& simplex, double* filtration_value) - { - std::vector< int > facet; - for (std::vector<int>::iterator not_it = simplex.begin(); not_it != simplex.end(); ++not_it) - { - facet.clear(); - for (std::vector<int>::iterator it = simplex.begin(); it != simplex.end(); ++it) - if (it != not_it) - facet.push_back(*it); - Simplex_handle facet_sh = sc.find(facet); - if (facet_sh == sc.null_simplex()) - return false; - else if (sc.filtration(facet_sh) > *filtration_value) - *filtration_value = sc.filtration(facet_sh); - } - return true; - } - - bool is_face(Simplex_handle face, Simplex_handle coface) - { - // vertex range is sorted in decreasing order - auto fvr = sc.simplex_vertex_range(face); - auto cfvr = sc.simplex_vertex_range(coface); - auto fv_it = fvr.begin(); - auto cfv_it = cfvr.begin(); - while (fv_it != fvr.end() && cfv_it != cfvr.end()) { - if (*fv_it < *cfv_it) - ++cfv_it; - else if (*fv_it == *cfv_it) { - ++cfv_it; - ++fv_it; - } - else - return false; - - } - return (fv_it == fvr.end()); - } - - // void erase_simplex(Simplex_handle sh) - // { - // auto siblings = sc.self_siblings(sh); - // auto oncles = siblings->oncles(); - // int prev_vertex = siblings->parent(); - // siblings->members().erase(sh->first); - // if (siblings->members().empty()) { - // typename typedef Simplicial_complex::Siblings Siblings; - // oncles->members().find(prev_vertex)->second.assign_children(new Siblings(oncles, prev_vertex)); - // assert(!sc.has_children(oncles->members().find(prev_vertex))); - // //delete siblings; - // } - - // } - - void elementary_collapse(Simplex_handle face_sh, Simplex_handle coface_sh) - { - erase_simplex(coface_sh); - erase_simplex(face_sh); - } - -public: - // void collapse(std::vector<Simplex_handle>& simplices) - // { - // // Get a vector of simplex handles ordered by filtration value - // std::cout << sc << std::endl; - // //std::vector<Simplex_handle> simplices; - // for (Simplex_handle sh: sc.filtration_simplex_range()) - // simplices.push_back(sh); - // // std::sort(simplices.begin(), - // // simplices.end(), - // // [&](Simplex_handle sh1, Simplex_handle sh2) - // // { double f1 = sc.filtration(sh1), f2 = sc.filtration(sh2); - // // return (f1 > f2) || (f1 >= f2 && sc.dimension(sh1) > sc.dimension(sh2)); }); - // // Double iteration - // auto face_it = simplices.rbegin(); - // while (face_it != simplices.rend() && sc.filtration(*face_it) != 0) { - // int coface_count = 0; - // auto reduced_coface = simplices.rbegin(); - // for (auto coface_it = simplices.rbegin(); coface_it != simplices.rend() && sc.filtration(*coface_it) != 0; ++coface_it) - // if (face_it != coface_it && is_face(*face_it, *coface_it)) { - // coface_count++; - // if (coface_count == 1) - // reduced_coface = coface_it; - // else - // break; - // } - // if (coface_count == 1) { - // std::cout << "Erase ( "; - // for (auto v: sc.simplex_vertex_range(*(--reduced_coface.base()))) - // std::cout << v << " "; - - // simplices.erase(--(reduced_coface.base())); - // //elementary_collapse(*face_it, *reduced_coface); - // std::cout << ") and then ( "; - // for (auto v: sc.simplex_vertex_range(*(--face_it.base()))) - // std::cout << v << " "; - // std::cout << ")\n"; - // simplices.erase(--((face_it++).base())); - // //face_it = simplices.rbegin(); - // //std::cout << "Size of vector: " << simplices.size() << "\n"; - // } - // else - // face_it++; - // } - // sc.initialize_filtration(); - // //std::cout << sc << std::endl; - // } - - template <class Dim_lists> - void collapse(Dim_lists& dim_lists) - { - dim_lists.collapse(); - } - -private: - - /** Collapse recursively boundary faces of the given simplex - * if its filtration is bigger than alpha_lim. - */ - void rec_boundary_collapse(Simplex_handle sh, FT alpha_lim) - { - for (Simplex_handle face_it : sc.boundary_simplex_range()) { - - } - - } - -}; //class Relaxed_witness_complex - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Construct_closest_landmark_table.h b/src/Witness_complex/include/gudhi/Construct_closest_landmark_table.h deleted file mode 100644 index ef711c34..00000000 --- a/src/Witness_complex/include/gudhi/Construct_closest_landmark_table.h +++ /dev/null @@ -1,90 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef CONSTRUCT_CLOSEST_LANDMARK_TABLE_H_ -#define CONSTRUCT_CLOSEST_LANDMARK_TABLE_H_ - -#include <boost/range/size.hpp> - -#include <queue> // for priority_queue<> -#include <utility> // for pair<> -#include <iterator> -#include <vector> -#include <set> - -namespace Gudhi { - -namespace witness_complex { - - /** - * \ingroup witness_complex - * \brief Construct the closest landmark tables for all witnesses. - * \details Output a table 'knn', each line of which represents a witness and - * consists of landmarks sorted by - * euclidean distance from the corresponding witness. - * - * The type WitnessContainer is a random access range and - * the type LandmarkContainer is a range. - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges and - * Vertex_handle is the label type of a vertex in a simplicial complex. - * Closest_landmark_range needs to have push_back operation. - */ - - template <typename WitnessContainer, - typename LandmarkContainer, - typename KNearestNeighbours> - void construct_closest_landmark_table(WitnessContainer const &points, - LandmarkContainer const &landmarks, - KNearestNeighbours &knn) { - int nbP = boost::size(points); - assert(nbP >= boost::size(landmarks)); - - int dim = boost::size(*std::begin(points)); - typedef std::pair<double, int> dist_i; - typedef bool (*comp)(dist_i, dist_i); - knn = KNearestNeighbours(nbP); - for (int points_i = 0; points_i < nbP; points_i++) { - std::priority_queue<dist_i, std::vector<dist_i>, comp> l_heap([](dist_i j1, dist_i j2) { - return j1.first > j2.first; - }); - typename LandmarkContainer::const_iterator landmarks_it; - int landmarks_i = 0; - for (landmarks_it = landmarks.begin(), landmarks_i = 0; landmarks_it != landmarks.end(); - ++landmarks_it, landmarks_i++) { - dist_i dist = std::make_pair(euclidean_distance(points[points_i], *landmarks_it), landmarks_i); - l_heap.push(dist); - } - for (int i = 0; i < dim + 1; i++) { - dist_i dist = l_heap.top(); - knn[points_i].push_back(dist.second); - l_heap.pop(); - } - } - } - -} // namespace witness_complex - -} // namespace Gudhi - -#endif // CONSTRUCT_CLOSEST_LANDMARK_TABLE_H_ diff --git a/src/Witness_complex/include/gudhi/Dim_list_iterator.h b/src/Witness_complex/include/gudhi/Dim_list_iterator.h deleted file mode 100644 index 3f23e8c9..00000000 --- a/src/Witness_complex/include/gudhi/Dim_list_iterator.h +++ /dev/null @@ -1,155 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef DIM_LISTS_ITERATOR_H_ -#define DIM_LISTS_ITERATOR_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include "gudhi/Simplex_tree.h" -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplex_handle > -class Dim_lists_iterator { - -private: - - typedef typename std::list<Simplex_handle>::iterator List_iterator; - typedef typename std::vector<std::list<Simplex_handle>>::reverse_iterator Vector_iterator; - typedef typename Gudhi::witness_complex::Dim_lists_iterator<Simplex_handle> Iterator; - - - List_iterator sh_; - Vector_iterator curr_list_; - typename std::vector<std::list<Simplex_handle>>& table_; - -public: - - Dim_lists_iterator(List_iterator sh, - Vector_iterator curr_list, - typename std::vector<std::list<Simplex_handle>>& table) - : sh_(sh), curr_list_(curr_list), table_(table) - { - } - - Simplex_handle operator*() - { - return *sh_; - } - - Iterator operator++() - { - increment(); - return (*this); - } - - Iterator operator++(int) - { - Iterator prev_it(sh_, curr_list_, table_); - increment(); - return prev_it; - } - - Iterator dim_begin() - { - return Iterator(curr_list_->begin(), curr_list_, table_); - } - - Iterator dim_end() - { - return Iterator(curr_list_->end(), curr_list_, table_); - } - - Iterator dimp1_begin() - { - return Iterator((curr_list_-1)->begin(), curr_list_-1, table_); - } - - Iterator dimp1_end() - { - return Iterator((curr_list_-1)->end(), curr_list_-1, table_); - } - - bool operator==(const Iterator& it2) const - { - return (sh_ == it2.sh_); - } - - bool operator!=(const Iterator& it2) const - { - return (sh_ != it2.sh_); - } - - void remove_incr() - { - - } - -private: - - void increment() - { - if (++sh_ == curr_list_->end()) - if (++curr_list_ != table_.rend()) - sh_ = curr_list_->begin(); - // The iterator of the end of the table is the end of the last list - } - - -}; //class Dim_lists_iterator - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Dim_lists.h b/src/Witness_complex/include/gudhi/Dim_lists.h deleted file mode 100644 index 1d1b03c3..00000000 --- a/src/Witness_complex/include/gudhi/Dim_lists.h +++ /dev/null @@ -1,195 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef DIM_LISTS_H_ -#define DIM_LISTS_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include <gudhi/Dim_list_iterator.h> -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplicial_complex > -class Dim_lists { - -private: - - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Gudhi::witness_complex::Dim_lists_iterator<Simplex_handle> Iterator; - - std::vector<std::list<Simplex_handle>> table_; - Simplicial_complex& sc_; - -public: - - Dim_lists(Simplicial_complex & sc, int dim, double alpha_max = 100) - : sc_(sc) - { - table_ = std::vector<std::list<Simplex_handle>>(dim+1); - for (auto sh: sc.filtration_simplex_range()) { - if (sc_.filtration(sh) < alpha_max) - table_[sc.dimension(sh)].push_front(sh); - } - auto t_it = table_.rbegin(); - while (t_it->empty()) { - t_it++; - table_.pop_back(); - } - } - - Iterator begin() - { - return Iterator(table_.rbegin()->begin(), table_.rbegin(), table_); - } - - Iterator end() - { - return Iterator(table_[0].end(), table_.rend(), table_); - } - - unsigned size() - { - unsigned curr_size = 0; - for (auto l: table_) - curr_size += l.size(); - return curr_size; - } - - bool is_face(Simplex_handle face, Simplex_handle coface) - { - // vertex range is sorted in decreasing order - auto fvr = sc_.simplex_vertex_range(face); - auto cfvr = sc_.simplex_vertex_range(coface); - auto fv_it = fvr.begin(); - auto cfv_it = cfvr.begin(); - while (fv_it != fvr.end() && cfv_it != cfvr.end()) { - if (*fv_it < *cfv_it) - ++cfv_it; - else if (*fv_it == *cfv_it) { - ++cfv_it; - ++fv_it; - } - else - return false; - - } - return (fv_it == fvr.end()); - } - - void output_simplices() { - std::cout << "Size of vector: " << size() << std::endl; - for (auto line_it = table_.rbegin(); line_it != table_.rend(); ++line_it) - for (auto sh: *line_it) { - std::cout << sc_.dimension(sh) << " "; - for (auto v : sc_.simplex_vertex_range(sh)) - std::cout << v << " "; - std::cout << sc_.filtration(sh) << "\n"; - } - } - - void collapse() - { - auto coface_list_it = table_.rbegin(); - auto face_list_it = table_.rbegin()+1; - for ( ; - face_list_it != table_.rend(); - ++face_list_it) { - auto face_it = face_list_it->begin(); - while (face_it != face_list_it->end() && sc_.filtration(*face_it) != 0) { - int coface_count = 0; - auto reduced_coface = coface_list_it->begin(); - for (auto coface_it = coface_list_it->begin(); coface_it != coface_list_it->end() && sc_.filtration(*coface_it) != 0; ++coface_it) - if (is_face(*face_it, *coface_it)) { - coface_count++; - if (coface_count == 1) - reduced_coface = coface_it; - else - break; - } - if (coface_count == 1) { - /* - std::cout << "Erase ( "; - for (auto v: sc_.simplex_vertex_range(*reduced_coface)) - std::cout << v << " "; - */ - coface_list_it->erase(reduced_coface); - /* - std::cout << ") and then ( "; - for (auto v: sc_.simplex_vertex_range(*face_it)) - std::cout << v << " "; - std::cout << ")\n"; - */ - face_list_it->erase(face_it); - face_it = face_list_it->begin(); - } - else - face_it++; - } - if ((coface_list_it++)->empty()) - table_.pop_back(); - } - } - - bool is_pseudomanifold() - { - - return true; - } - -}; //class Dim_lists - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Good_links.h b/src/Witness_complex/include/gudhi/Good_links.h deleted file mode 100644 index c5e993e5..00000000 --- a/src/Witness_complex/include/gudhi/Good_links.h +++ /dev/null @@ -1,449 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef GOOD_LINKS_H_ -#define GOOD_LINKS_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include <gudhi/Dim_list_iterator.h> -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -template < typename Simplicial_complex > -class Good_links { - - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - typedef std::vector<Vertex_handle> Vertex_vector; - - // graph is an adjacency list - typedef typename boost::adjacency_list<boost::vecS, boost::vecS, boost::undirectedS> Adj_graph; - // map that gives to a certain simplex its node in graph and its dimension - //typedef std::pair<boost::vecS,Color> Reference; - typedef boost::graph_traits<Adj_graph>::vertex_descriptor Vertex_t; - typedef boost::graph_traits<Adj_graph>::edge_descriptor Edge_t; - typedef boost::graph_traits<Adj_graph>::adjacency_iterator Adj_it; - typedef std::pair<Adj_it, Adj_it> Out_edge_it; - - typedef boost::container::flat_map<Simplex_handle, Vertex_t> Graph_map; - typedef boost::container::flat_map<Vertex_t, Simplex_handle> Inv_graph_map; - -public: - Good_links(Simplicial_complex& sc): sc_(sc) - { - int dim = 0; - for (auto sh: sc_.complex_simplex_range()) - if (sc_.dimension(sh) > dim) - dim = sc_.dimension(sh); - dimension = dim; - count_good = Vertex_vector(dim); - count_bad = Vertex_vector(dim); - } - -private: - - Simplicial_complex& sc_; - unsigned dimension; - std::vector<int> count_good; - std::vector<int> count_bad; - - void add_vertices_to_link_graph(Vertex_vector& star_vertices, - typename Vertex_vector::iterator curr_v, - Adj_graph& adj_graph, - Graph_map& d_map, - Graph_map& f_map, - Vertex_vector& curr_simplex, - int curr_d, - int link_dimension) - { - Simplex_handle sh; - Vertex_t vert; - typename Vertex_vector::iterator it; - //std::pair<typename Graph_map::iterator,bool> resPair; - //typename Graph_map::iterator resPair; - //Add vertices - //std::cout << "Entered add vertices\n"; - for (it = curr_v; it != star_vertices.end(); ++it) { - curr_simplex.push_back(*it); //push next vertex in question - curr_simplex.push_back(star_vertices[0]); //push the center of the star - /* - std::cout << "Searching for "; - print_vector(curr_simplex); - std::cout << " curr_dim " << curr_d << " d " << dimension << ""; - */ - Vertex_vector curr_simplex_copy(curr_simplex); - sh = sc_.find(curr_simplex_copy); //a simplex of the star - curr_simplex.pop_back(); //pop the center of the star - curr_simplex_copy = Vertex_vector(curr_simplex); - if (sh != sc_.null_simplex()) { - //std::cout << " added\n"; - if (curr_d == link_dimension) { - sh = sc_.find(curr_simplex_copy); //a simplex of the link - assert(sh != sc_.null_simplex()); //ASSERT! - vert = boost::add_vertex(adj_graph); - d_map.emplace(sh,vert); - } - else { - if (curr_d == link_dimension-1) { - sh = sc_.find(curr_simplex_copy); //a simplex of the link - assert(sh != sc_.null_simplex()); - vert = boost::add_vertex(adj_graph); - f_map.emplace(sh,vert); - } - - //delete (&curr_simplex_copy); //Just so it doesn't stack - add_vertices_to_link_graph(star_vertices, - it+1, - adj_graph, - d_map, - f_map, - curr_simplex, - curr_d+1, link_dimension); - } - } - /* - else - std::cout << "\n"; - */ - curr_simplex.pop_back(); //pop the vertex in question - } - } - - void add_edges_to_link_graph(Adj_graph& adj_graph, - Graph_map& d_map, - Graph_map& f_map) - { - Simplex_handle sh; - // Add edges - //std::cout << "Entered add edges:\n"; - typename Graph_map::iterator map_it; - for (auto d_map_pair : d_map) { - //std::cout << "*"; - sh = d_map_pair.first; - Vertex_t d_vert = d_map_pair.second; - for (auto facet_sh : sc_.boundary_simplex_range(sh)) - //for (auto f_map_it : f_map) - { - //std::cout << "'"; - map_it = f_map.find(facet_sh); - //We must have all the facets in the graph at this point - assert(map_it != f_map.end()); - Vertex_t f_vert = map_it->second; - //std::cout << "Added edge " << sh->first << "-" << map_it->first->first << "\n"; - boost::add_edge(d_vert,f_vert,adj_graph); - } - } - } - - - - /* \brief Verifies if the simplices formed by vertices given by link_vertices - * form a pseudomanifold. - * The idea is to make a bipartite graph, where vertices are the d- and (d-1)-dimensional - * faces and edges represent adjacency between them. - */ - bool link_is_pseudomanifold(Vertex_vector& star_vertices, - int dimension) - { - Adj_graph adj_graph; - Graph_map d_map, f_map; - // d_map = map for d-dimensional simplices - // f_map = map for its facets - Vertex_vector init_vector = {}; - add_vertices_to_link_graph(star_vertices, - star_vertices.begin()+1, - adj_graph, - d_map, - f_map, - init_vector, - 0, dimension); - //std::cout << "DMAP_SIZE: " << d_map.size() << "\n"; - //std::cout << "FMAP_SIZE: " << f_map.size() << "\n"; - add_edges_to_link_graph(adj_graph, - d_map, - f_map); - for (auto f_map_it : f_map) { - //std::cout << "Degree of " << f_map_it.first->first << " is " << boost::out_degree(f_map_it.second, adj_graph) << "\n"; - if (boost::out_degree(f_map_it.second, adj_graph) != 2){ - count_bad[dimension]++; - return false; - } - } - // At this point I know that all (d-1)-simplices are adjacent to exactly 2 d-simplices - // What is left is to check the connexity - //std::vector<int> components(boost::num_vertices(adj_graph)); - return true; //Forget the connexity - //return (boost::connected_components(adj_graph, &components[0]) == 1); - } - - int star_dim(Vertex_vector& star_vertices, - typename Vertex_vector::iterator curr_v, - int curr_d, - Vertex_vector& curr_simplex, - typename std::vector< int >::iterator curr_dc) - { - //std::cout << "Entered star_dim for " << *(curr_v-1) << "\n"; - Simplex_handle sh; - int final_d = curr_d; - typename Vertex_vector::iterator it; - typename Vertex_vector::iterator dc_it; - //std::cout << "Current vertex is " << - for (it = curr_v, dc_it = curr_dc; it != star_vertices.end(); ++it, ++dc_it) - { - curr_simplex.push_back(*it); - Vertex_vector curr_simplex_copy(curr_simplex); - /* - std::cout << "Searching for "; - print_vector(curr_simplex); - std::cout << " curr_dim " << curr_d << " final_dim " << final_d; - */ - sh = sc_.find(curr_simplex_copy); //Need a copy because find sorts the vector and I want star center to be the first - if (sh != sc_.null_simplex()) - { - //std::cout << " -> " << *it << "\n"; - int d = star_dim(star_vertices, - it+1, - curr_d+1, - curr_simplex, - dc_it); - if (d >= final_d) - { - final_d = d; - //std::cout << d << " "; - //print_vector(curr_simplex); - //std::cout << std::endl; - } - if (d >= *dc_it) - *dc_it = d; - } - /* - else - std::cout << "\n"; - */ - curr_simplex.pop_back(); - } - return final_d; - } - -public: - - /** \brief Returns true if the link is good - */ - bool has_good_link(Vertex_handle v) - { - typedef Vertex_vector typeVectorVertex; - typeVectorVertex star_vertices; - // Fill star_vertices - star_vertices.push_back(v); - for (auto u: sc_.complex_vertex_range()) - { - typeVectorVertex edge = {u,v}; - if (u != v && sc_.find(edge) != sc_.null_simplex()) - star_vertices.push_back(u); - } - // Find the dimension - typeVectorVertex init_simplex = {star_vertices[0]}; - bool is_pure = true; - std::vector<int> dim_coface(star_vertices.size(), 1); - int d = star_dim(star_vertices, - star_vertices.begin()+1, - 0, - init_simplex, - dim_coface.begin()+1) - 1; //link_dim = star_dim - 1 - - assert(init_simplex.size() == 1); - // if (!is_pure) - // std::cout << "Found an impure star around " << v << "\n"; - for (int dc: dim_coface) - is_pure = (dc == dim_coface[0]); - /* - if (d == count_good.size()) - { - std::cout << "Found a star of dimension " << (d+1) << " around " << v << "\nThe star is "; - print_vector(star_vertices); std::cout << std::endl; - } - */ - //if (d == -1) count_bad[0]++; - bool b= (is_pure && link_is_pseudomanifold(star_vertices,d)); - if (d != -1) {if (b) count_good[d]++; else count_bad[d]++;} - if (!is_pure) count_bad[0]++; - return (d != -1 && b && is_pure); - } - -private: - void add_max_simplices_to_graph(Vertex_vector& star_vertices, - typename Vertex_vector::iterator curr_v, - Adj_graph& adj_graph, - Graph_map& d_map, - Graph_map& f_map, - Inv_graph_map& inv_d_map, - Vertex_vector& curr_simplex, - int curr_d, - int link_dimension) - { - Simplex_handle sh; - Vertex_t vert; - typename Vertex_vector::iterator it; - //std::pair<typename Graph_map::iterator,bool> resPair; - //typename Graph_map::iterator resPair; - //Add vertices - //std::cout << "Entered add vertices\n"; - for (it = curr_v; it != star_vertices.end(); ++it) { - curr_simplex.push_back(*it); //push next vertex in question - //curr_simplex.push_back(star_vertices[0]); //push the center of the star - /* - std::cout << "Searching for "; - print_vector(curr_simplex); - std::cout << " curr_dim " << curr_d << " d " << dimension << ""; - */ - Vertex_vector curr_simplex_copy(curr_simplex); - sh = sc_.find(curr_simplex_copy); //a simplex of the star - //curr_simplex.pop_back(); //pop the center of the star - curr_simplex_copy = Vertex_vector(curr_simplex); - if (sh != sc_.null_simplex()) { - //std::cout << " added\n"; - if (curr_d == link_dimension) { - sh = sc_.find(curr_simplex_copy); //a simplex of the link - assert(sh != sc_.null_simplex()); //ASSERT! - vert = boost::add_vertex(adj_graph); - d_map.emplace(sh,vert); - inv_d_map.emplace(vert,sh); - } - else { - if (curr_d == link_dimension-1) { - sh = sc_.find(curr_simplex_copy); //a simplex of the link - assert(sh != sc_.null_simplex()); - vert = boost::add_vertex(adj_graph); - f_map.emplace(sh,vert); - } - //delete (&curr_simplex_copy); //Just so it doesn't stack - add_max_simplices_to_graph(star_vertices, - it+1, - adj_graph, - d_map, - f_map, - inv_d_map, - curr_simplex, - curr_d+1, - link_dimension); - } - } - /* - else - std::cout << "\n"; - */ - curr_simplex.pop_back(); //pop the vertex in question - } - } - -public: - bool complex_is_pseudomanifold() - { - Adj_graph adj_graph; - Graph_map d_map, f_map; - // d_map = map for d-dimensional simplices - // f_map = map for its facets - Inv_graph_map inv_d_map; - Vertex_vector init_vector = {}; - std::vector<int> star_vertices; - for (int v: sc_.complex_vertex_range()) - star_vertices.push_back(v); - add_max_simplices_to_graph(star_vertices, - star_vertices.begin(), - adj_graph, - d_map, - f_map, - inv_d_map, - init_vector, - 0, dimension); - //std::cout << "DMAP_SIZE: " << d_map.size() << "\n"; - //std::cout << "FMAP_SIZE: " << f_map.size() << "\n"; - add_edges_to_link_graph(adj_graph, - d_map, - f_map); - for (auto f_map_it : f_map) { - //std::cout << "Degree of " << f_map_it.first->first << " is " << boost::out_degree(f_map_it.second, adj_graph) << "\n"; - if (boost::out_degree(f_map_it.second, adj_graph) != 2) { - // if (boost::out_degree(f_map_it.second, adj_graph) >= 3) { - // // std::cout << "This simplex has 3+ cofaces: "; - // // for(auto v : sc_.simplex_vertex_range(f_map_it.first)) - // // std::cout << v << " "; - // // std::cout << std::endl; - // Adj_it ai, ai_end; - // for (std::tie(ai, ai_end) = boost::adjacent_vertices(f_map_it.second, adj_graph); ai != ai_end; ++ai) { - // auto it = inv_d_map.find(*ai); - // assert (it != inv_d_map.end()); - // typename Simplicial_complex::Simplex_handle sh = it->second; - // for(auto v : sc_.simplex_vertex_range(sh)) - // std::cout << v << " "; - // std::cout << std::endl; - // } - // } - count_bad[dimension]++; - return false; - } - } - // At this point I know that all (d-1)-simplices are adjacent to exactly 2 d-simplices - // What is left is to check the connexity - //std::vector<int> components(boost::num_vertices(adj_graph)); - return true; //Forget the connexity - //return (boost::connected_components(adj_graph, &components[0]) == 1); - } - - int number_good_links(int dim) - { - return count_good[dim]; - } - - int number_bad_links(int dim) - { - return count_bad[dim]; - } - -}; - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Landmark_choice_by_furthest_point.h b/src/Witness_complex/include/gudhi/Landmark_choice_by_furthest_point.h deleted file mode 100644 index df93155b..00000000 --- a/src/Witness_complex/include/gudhi/Landmark_choice_by_furthest_point.h +++ /dev/null @@ -1,105 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef LANDMARK_CHOICE_BY_FURTHEST_POINT_H_ -#define LANDMARK_CHOICE_BY_FURTHEST_POINT_H_ - -#include <boost/range/size.hpp> - -#include <limits> // for numeric_limits<> -#include <iterator> -#include <algorithm> // for sort -#include <vector> - -namespace Gudhi { - -namespace witness_complex { - - typedef std::vector<int> typeVectorVertex; - - /** - * \ingroup witness_complex - * \brief Landmark choice strategy by iteratively adding the furthest witness from the - * current landmark set as the new landmark. - * \details It chooses nbL landmarks from a random access range `points` and - * writes {witness}*{closest landmarks} matrix in `knn`. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges - * - */ - - template <typename KNearestNeighbours, - typename Point_random_access_range> - void landmark_choice_by_furthest_point(Point_random_access_range const &points, - int nbL, - KNearestNeighbours &knn) { - int nb_points = boost::size(points); - assert(nb_points >= nbL); - // distance matrix witness x landmarks - std::vector<std::vector<double>> wit_land_dist(nb_points, std::vector<double>()); - // landmark list - typeVectorVertex chosen_landmarks; - - knn = KNearestNeighbours(nb_points, std::vector<int>()); - int current_number_of_landmarks = 0; // counter for landmarks - double curr_max_dist = 0; // used for defining the furhest point from L - const double infty = std::numeric_limits<double>::infinity(); // infinity (see next entry) - std::vector< double > dist_to_L(nb_points, infty); // vector of current distances to L from points - - // TODO(SK) Consider using rand_r(...) instead of rand(...) for improved thread safety - // or better yet std::uniform_int_distribution - int rand_int = rand() % nb_points; - int curr_max_w = rand_int; // For testing purposes a pseudo-random number is used here - - for (current_number_of_landmarks = 0; current_number_of_landmarks != nbL; current_number_of_landmarks++) { - // curr_max_w at this point is the next landmark - chosen_landmarks.push_back(curr_max_w); - unsigned i = 0; - for (auto& p : points) { - double curr_dist = euclidean_distance(p, *(std::begin(points) + chosen_landmarks[current_number_of_landmarks])); - wit_land_dist[i].push_back(curr_dist); - knn[i].push_back(current_number_of_landmarks); - if (curr_dist < dist_to_L[i]) - dist_to_L[i] = curr_dist; - ++i; - } - curr_max_dist = 0; - for (i = 0; i < dist_to_L.size(); i++) - if (dist_to_L[i] > curr_max_dist) { - curr_max_dist = dist_to_L[i]; - curr_max_w = i; - } - } - for (int i = 0; i < nb_points; ++i) - std::sort(std::begin(knn[i]), - std::end(knn[i]), - [&wit_land_dist, i](int a, int b) { - return wit_land_dist[i][a] < wit_land_dist[i][b]; }); - } - -} // namespace witness_complex - -} // namespace Gudhi - -#endif // LANDMARK_CHOICE_BY_FURTHEST_POINT_H_ diff --git a/src/Witness_complex/include/gudhi/Landmark_choice_by_random_point.h b/src/Witness_complex/include/gudhi/Landmark_choice_by_random_point.h deleted file mode 100644 index ebf6aad1..00000000 --- a/src/Witness_complex/include/gudhi/Landmark_choice_by_random_point.h +++ /dev/null @@ -1,96 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef LANDMARK_CHOICE_BY_RANDOM_POINT_H_ -#define LANDMARK_CHOICE_BY_RANDOM_POINT_H_ - -#include <boost/range/size.hpp> - -#include <queue> // for priority_queue<> -#include <utility> // for pair<> -#include <iterator> -#include <vector> -#include <set> - -namespace Gudhi { - -namespace witness_complex { - - /** - * \ingroup witness_complex - * \brief Landmark choice strategy by taking random vertices for landmarks. - * \details It chooses nbL distinct landmarks from a random access range `points` - * and outputs a matrix {witness}*{closest landmarks} in knn. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges and - * Vertex_handle is the label type of a vertex in a simplicial complex. - * Closest_landmark_range needs to have push_back operation. - */ - - template <typename KNearestNeighbours, - typename Point_random_access_range> - void landmark_choice_by_random_point(Point_random_access_range const &points, - int nbL, - KNearestNeighbours &knn) { - int nbP = boost::size(points); - assert(nbP >= nbL); - std::set<int> landmarks; - int current_number_of_landmarks = 0; // counter for landmarks - - // TODO(SK) Consider using rand_r(...) instead of rand(...) for improved thread safety - int chosen_landmark = rand() % nbP; - for (current_number_of_landmarks = 0; current_number_of_landmarks != nbL; current_number_of_landmarks++) { - while (landmarks.find(chosen_landmark) != landmarks.end()) - chosen_landmark = rand() % nbP; - landmarks.insert(chosen_landmark); - } - - int dim = boost::size(*std::begin(points)); - typedef std::pair<double, int> dist_i; - typedef bool (*comp)(dist_i, dist_i); - knn = KNearestNeighbours(nbP); - for (int points_i = 0; points_i < nbP; points_i++) { - std::priority_queue<dist_i, std::vector<dist_i>, comp> l_heap([](dist_i j1, dist_i j2) { - return j1.first > j2.first; - }); - std::set<int>::iterator landmarks_it; - int landmarks_i = 0; - for (landmarks_it = landmarks.begin(), landmarks_i = 0; landmarks_it != landmarks.end(); - ++landmarks_it, landmarks_i++) { - dist_i dist = std::make_pair(euclidean_distance(points[points_i], points[*landmarks_it]), landmarks_i); - l_heap.push(dist); - } - for (int i = 0; i < dim + 1; i++) { - dist_i dist = l_heap.top(); - knn[points_i].push_back(dist.second); - l_heap.pop(); - } - } - } - -} // namespace witness_complex - -} // namespace Gudhi - -#endif // LANDMARK_CHOICE_BY_RANDOM_POINT_H_ diff --git a/src/Witness_complex/include/gudhi/Relaxed_witness_complex.h b/src/Witness_complex/include/gudhi/Relaxed_witness_complex.h deleted file mode 100644 index 2971149a..00000000 --- a/src/Witness_complex/include/gudhi/Relaxed_witness_complex.h +++ /dev/null @@ -1,379 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef RELAXED_WITNESS_COMPLEX_H_ -#define RELAXED_WITNESS_COMPLEX_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include "gudhi/Simplex_tree.h" -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -// Needed for nearest neighbours -#include <CGAL/Cartesian_d.h> -#include <CGAL/Search_traits.h> -#include <CGAL/Search_traits_adapter.h> -#include <CGAL/property_map.h> -#include <CGAL/Epick_d.h> -#include <CGAL/Orthogonal_k_neighbor_search.h> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplicial_complex > -class Relaxed_witness_complex { - -private: - struct Active_witness { - int witness_id; - int landmark_id; - - Active_witness(int witness_id_, int landmark_id_) - : witness_id(witness_id_), - landmark_id(landmark_id_) { } - }; - -private: - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - typedef typename Simplicial_complex::Filtration_value FT; - - typedef std::vector< double > Point_t; - typedef std::vector< Point_t > Point_Vector; - - // typedef typename Simplicial_complex::Filtration_value Filtration_value; - typedef std::vector< Vertex_handle > typeVectorVertex; - typedef std::pair< typeVectorVertex, Filtration_value> typeSimplex; - typedef std::pair< Simplex_handle, bool > typePairSimplexBool; - - typedef int Witness_id; - typedef int Landmark_id; - typedef std::list< Vertex_handle > ActiveWitnessList; - - private: - int nbL; // Number of landmarks - Simplicial_complex& sc; // Simplicial complex - - public: - /** @name Constructor - */ - - //@{ - - /** - * \brief Iterative construction of the relaxed witness complex. - * \details The witness complex is written in sc_ basing on a matrix knn - * of k nearest neighbours of the form {witnesses}x{landmarks} and - * and a matrix distances of distances to these landmarks from witnesses. - * The parameter alpha defines relaxation and - * limD defines the - * - * The line lengths in one given matrix can differ, - * however both matrices have the same corresponding line lengths. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges. - * - * Constructor takes into account at most (dim+1) - * first landmarks from each landmark range to construct simplices. - * - * Landmarks are supposed to be in [0,nbL_-1] - */ - - template< typename KNearestNeighbours > - Relaxed_witness_complex(std::vector< std::vector<double> > const & distances, - KNearestNeighbours const & knn, - Simplicial_complex & sc_, - int nbL_, - double alpha, - int limD) : nbL(nbL_), sc(sc_) { - int nbW = knn.size(); - typeVectorVertex vv; - //int counter = 0; - /* The list of still useful witnesses - * it will diminuish in the course of iterations - */ - ActiveWitnessList active_w;// = new ActiveWitnessList(); - for (int i = 0; i != nbL; ++i) { - // initial fill of 0-dimensional simplices - // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore - //counter++; - vv = {i}; - sc.insert_simplex(vv, Filtration_value(0.0)); - /* TODO Error if not inserted : normally no need here though*/ - } - int k = 1; /* current dimension in iterative construction */ - for (int i=0; i != nbW; ++i) - active_w.push_back(i); - while (!active_w.empty() && k <= limD && k < nbL ) - { - typename ActiveWitnessList::iterator aw_it = active_w.begin(); - while (aw_it != active_w.end()) - { - std::vector<int> simplex; - bool ok = add_all_faces_of_dimension(k, - alpha, - std::numeric_limits<double>::infinity(), - distances[*aw_it].begin(), - knn[*aw_it].begin(), - simplex, - knn[*aw_it].end()); - if (!ok) - active_w.erase(aw_it++); //First increase the iterator and then erase the previous element - else - aw_it++; - } - k++; - } - sc.set_dimension(limD); - } - - //@} - -private: - /* \brief Adds recursively all the faces of a certain dimension dim witnessed by the same witness - * Iterator is needed to know until how far we can take landmarks to form simplexes - * simplex is the prefix of the simplexes to insert - * The output value indicates if the witness rests active or not - */ - bool add_all_faces_of_dimension(int dim, - double alpha, - double norelax_dist, - std::vector<double>::const_iterator curr_d, - std::vector<int>::const_iterator curr_l, - std::vector<int>& simplex, - std::vector<int>::const_iterator end) - { - if (curr_l == end) - return false; - bool will_be_active = false; - std::vector<int>::const_iterator l_it = curr_l; - std::vector<double>::const_iterator d_it = curr_d; - if (dim > 0) - for (; *d_it - alpha <= norelax_dist && l_it != end; ++l_it, ++d_it) { - simplex.push_back(*l_it); - if (sc.find(simplex) != sc.null_simplex()) - will_be_active = add_all_faces_of_dimension(dim-1, - alpha, - norelax_dist, - d_it+1, - l_it+1, - simplex, - end) || will_be_active; - simplex.pop_back(); - // If norelax_dist is infinity, change to first omitted distance - if (*d_it <= norelax_dist) - norelax_dist = *d_it; - will_be_active = add_all_faces_of_dimension(dim-1, - alpha, - norelax_dist, - d_it+1, - l_it+1, - simplex, - end) || will_be_active; - } - else if (dim == 0) - for (; *d_it - alpha <= norelax_dist && l_it != end; ++l_it, ++d_it) { - simplex.push_back(*l_it); - double filtration_value = 0; - // if norelax_dist is infinite, relaxation is 0. - if (*d_it > norelax_dist) - filtration_value = *d_it - norelax_dist; - if (all_faces_in(simplex, &filtration_value)) { - will_be_active = true; - sc.insert_simplex(simplex, filtration_value); - } - simplex.pop_back(); - // If norelax_dist is infinity, change to first omitted distance - if (*d_it < norelax_dist) - norelax_dist = *d_it; - } - return will_be_active; - } - - /** \brief Check if the facets of the k-dimensional simplex witnessed - * by witness witness_id are already in the complex. - * inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id - */ - bool all_faces_in(std::vector<int>& simplex, double* filtration_value) - { - std::vector< int > facet; - for (std::vector<int>::iterator not_it = simplex.begin(); not_it != simplex.end(); ++not_it) - { - facet.clear(); - for (std::vector<int>::iterator it = simplex.begin(); it != simplex.end(); ++it) - if (it != not_it) - facet.push_back(*it); - Simplex_handle facet_sh = sc.find(facet); - if (facet_sh == sc.null_simplex()) - return false; - else if (sc.filtration(facet_sh) > *filtration_value) - *filtration_value = sc.filtration(facet_sh); - } - return true; - } - - bool is_face(Simplex_handle face, Simplex_handle coface) - { - // vertex range is sorted in decreasing order - auto fvr = sc.simplex_vertex_range(face); - auto cfvr = sc.simplex_vertex_range(coface); - auto fv_it = fvr.begin(); - auto cfv_it = cfvr.begin(); - while (fv_it != fvr.end() && cfv_it != cfvr.end()) { - if (*fv_it < *cfv_it) - ++cfv_it; - else if (*fv_it == *cfv_it) { - ++cfv_it; - ++fv_it; - } - else - return false; - - } - return (fv_it == fvr.end()); - } - - // void erase_simplex(Simplex_handle sh) - // { - // auto siblings = sc.self_siblings(sh); - // auto oncles = siblings->oncles(); - // int prev_vertex = siblings->parent(); - // siblings->members().erase(sh->first); - // if (siblings->members().empty()) { - // typename typedef Simplicial_complex::Siblings Siblings; - // oncles->members().find(prev_vertex)->second.assign_children(new Siblings(oncles, prev_vertex)); - // assert(!sc.has_children(oncles->members().find(prev_vertex))); - // //delete siblings; - // } - - // } - - void elementary_collapse(Simplex_handle face_sh, Simplex_handle coface_sh) - { - erase_simplex(coface_sh); - erase_simplex(face_sh); - } - -public: - // void collapse(std::vector<Simplex_handle>& simplices) - // { - // // Get a vector of simplex handles ordered by filtration value - // std::cout << sc << std::endl; - // //std::vector<Simplex_handle> simplices; - // for (Simplex_handle sh: sc.filtration_simplex_range()) - // simplices.push_back(sh); - // // std::sort(simplices.begin(), - // // simplices.end(), - // // [&](Simplex_handle sh1, Simplex_handle sh2) - // // { double f1 = sc.filtration(sh1), f2 = sc.filtration(sh2); - // // return (f1 > f2) || (f1 >= f2 && sc.dimension(sh1) > sc.dimension(sh2)); }); - // // Double iteration - // auto face_it = simplices.rbegin(); - // while (face_it != simplices.rend() && sc.filtration(*face_it) != 0) { - // int coface_count = 0; - // auto reduced_coface = simplices.rbegin(); - // for (auto coface_it = simplices.rbegin(); coface_it != simplices.rend() && sc.filtration(*coface_it) != 0; ++coface_it) - // if (face_it != coface_it && is_face(*face_it, *coface_it)) { - // coface_count++; - // if (coface_count == 1) - // reduced_coface = coface_it; - // else - // break; - // } - // if (coface_count == 1) { - // std::cout << "Erase ( "; - // for (auto v: sc.simplex_vertex_range(*(--reduced_coface.base()))) - // std::cout << v << " "; - - // simplices.erase(--(reduced_coface.base())); - // //elementary_collapse(*face_it, *reduced_coface); - // std::cout << ") and then ( "; - // for (auto v: sc.simplex_vertex_range(*(--face_it.base()))) - // std::cout << v << " "; - // std::cout << ")\n"; - // simplices.erase(--((face_it++).base())); - // //face_it = simplices.rbegin(); - // //std::cout << "Size of vector: " << simplices.size() << "\n"; - // } - // else - // face_it++; - // } - // sc.initialize_filtration(); - // //std::cout << sc << std::endl; - // } - - template <class Dim_lists> - void collapse(Dim_lists& dim_lists) - { - dim_lists.collapse(); - } - -private: - - /** Collapse recursively boundary faces of the given simplex - * if its filtration is bigger than alpha_lim. - */ - void rec_boundary_collapse(Simplex_handle sh, FT alpha_lim) - { - for (Simplex_handle face_it : sc.boundary_simplex_range()) { - - } - - } - -}; //class Relaxed_witness_complex - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Strange_witness_complex.h b/src/Witness_complex/include/gudhi/Strange_witness_complex.h deleted file mode 100644 index c1c2912b..00000000 --- a/src/Witness_complex/include/gudhi/Strange_witness_complex.h +++ /dev/null @@ -1,232 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef WITNESS_COMPLEX_H_ -#define WITNESS_COMPLEX_H_ - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> - -#include <gudhi/distance_functions.h> - -#include <algorithm> -#include <utility> -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <ctime> -#include <iostream> - -namespace Gudhi { - -namespace witness_complex { - -/** - \class Witness_complex - \brief Constructs the witness complex for the given set of witnesses and landmarks. - \ingroup witness_complex - */ -template< class Simplicial_complex> -class Strange_witness_complex { - private: - struct Active_witness { - int witness_id; - int landmark_id; - - Active_witness(int witness_id_, int landmark_id_) - : witness_id(witness_id_), - landmark_id(landmark_id_) { } - }; - - private: - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - - typedef std::vector< double > Point_t; - typedef std::vector< Point_t > Point_Vector; - - // typedef typename Simplicial_complex::Filtration_value Filtration_value; - typedef std::vector< Vertex_handle > typeVectorVertex; - typedef std::pair< typeVectorVertex, Filtration_value> typeSimplex; - typedef std::pair< Simplex_handle, bool > typePairSimplexBool; - - typedef int Witness_id; - typedef int Landmark_id; - typedef std::list< Vertex_handle > ActiveWitnessList; - - private: - int nbL; // Number of landmarks - Simplicial_complex& sc; // Simplicial complex - - public: - ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////// - /** @name Constructor - */ - - //@{ - - // Witness_range<Closest_landmark_range<Vertex_handle>> - - /** - * \brief Iterative construction of the witness complex. - * \details The witness complex is written in sc_ basing on a matrix knn of - * nearest neighbours of the form {witnesses}x{landmarks}. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges. - * - * Constructor takes into account at most (dim+1) - * first landmarks from each landmark range to construct simplices. - * - * Landmarks are supposed to be in [0,nbL_-1] - */ - template< typename KNearestNeighbors > - Strange_witness_complex(KNearestNeighbors const & knn, - Simplicial_complex & sc_, - int nbL_, - int dim) : nbL(nbL_), sc(sc_) { - // Construction of the active witness list - int nbW = knn.size(); - typeVectorVertex vv; - int counter = 0; - /* The list of still useful witnesses - * it will diminuish in the course of iterations - */ - for (int i = 0; i != nbL; ++i) { - // initial fill of 0-dimensional simplices - // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore - counter++; - vv = {i}; - sc.insert_simplex(vv); - // TODO(SK) Error if not inserted : normally no need here though - } - for (int i = 0; i != nbW; ++i) { - std::vector<int> simplex; - simplex.reserve(dim+1); - for (int j = 0; j <= dim; j++) - simplex.push_back(knn[i][j]); - sc.insert_simplex_and_subfaces(simplex); - } - } - - //@} - - private: - /** \brief Check if the facets of the k-dimensional simplex witnessed - * by witness witness_id are already in the complex. - * inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id - */ - template <typename KNearestNeighbors> - bool all_faces_in(KNearestNeighbors const &knn, int witness_id, int k) { - // std::cout << "All face in with the landmark " << inserted_vertex << std::endl; - std::vector< Vertex_handle > facet; - // Vertex_handle curr_vh = curr_sh->first; - // CHECK ALL THE FACETS - for (int i = 0; i != k + 1; ++i) { - facet = {}; - for (int j = 0; j != k + 1; ++j) { - if (j != i) { - facet.push_back(knn[witness_id][j]); - } - } // endfor - if (sc.find(facet) == sc.null_simplex()) - return false; - // std::cout << "++++ finished loop safely\n"; - } // endfor - return true; - } - - template <typename T> - static void print_vector(const std::vector<T>& v) { - std::cout << "["; - if (!v.empty()) { - std::cout << *(v.begin()); - for (auto it = v.begin() + 1; it != v.end(); ++it) { - std::cout << ","; - std::cout << *it; - } - } - std::cout << "]"; - } - - public: - /** - * \brief Verification if every simplex in the complex is witnessed by witnesses in knn. - * \param print_output =true will print the witnesses for each simplex - * \remark Added for debugging purposes. - */ - template< class KNearestNeighbors > - bool is_witness_complex(KNearestNeighbors const & knn, bool print_output) { - // bool final_result = true; - for (Simplex_handle sh : sc.complex_simplex_range()) { - bool is_witnessed = false; - typeVectorVertex simplex; - int nbV = 0; // number of verticed in the simplex - for (int v : sc.simplex_vertex_range(sh)) - simplex.push_back(v); - nbV = simplex.size(); - for (typeVectorVertex w : knn) { - bool has_vertices = true; - for (int v : simplex) - if (std::find(w.begin(), w.begin() + nbV, v) == w.begin() + nbV) { - has_vertices = false; - // break; - } - if (has_vertices) { - is_witnessed = true; - if (print_output) { - std::cout << "The simplex "; - print_vector(simplex); - std::cout << " is witnessed by the witness "; - print_vector(w); - std::cout << std::endl; - } - break; - } - } - if (!is_witnessed) { - if (print_output) { - std::cout << "The following simplex is not witnessed "; - print_vector(simplex); - std::cout << std::endl; - } - assert(is_witnessed); - return false; - } - } - return true; - } -}; - -} // namespace witness_complex - -} // namespace Gudhi - -#endif // WITNESS_COMPLEX_H_ diff --git a/src/Witness_complex/include/gudhi/Strong_relaxed_witness_complex.h b/src/Witness_complex/include/gudhi/Strong_relaxed_witness_complex.h deleted file mode 100644 index ee863a42..00000000 --- a/src/Witness_complex/include/gudhi/Strong_relaxed_witness_complex.h +++ /dev/null @@ -1,337 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef STRONG_RELAXED_WITNESS_COMPLEX_H_ -#define STRONG_RELAXED_WITNESS_COMPLEX_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include "gudhi/Simplex_tree.h" -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -// Needed for nearest neighbours -#include <CGAL/Cartesian_d.h> -#include <CGAL/Search_traits.h> -#include <CGAL/Search_traits_adapter.h> -#include <CGAL/property_map.h> -#include <CGAL/Epick_d.h> -#include <CGAL/Orthogonal_k_neighbor_search.h> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplicial_complex > -class Strong_relaxed_witness_complex { - -private: - struct Active_witness { - int witness_id; - int landmark_id; - - Active_witness(int witness_id_, int landmark_id_) - : witness_id(witness_id_), - landmark_id(landmark_id_) { } - }; - -private: - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - typedef typename Simplicial_complex::Filtration_value FT; - - typedef std::vector< double > Point_t; - typedef std::vector< Point_t > Point_Vector; - - // typedef typename Simplicial_complex::Filtration_value Filtration_value; - typedef std::vector< Vertex_handle > typeVectorVertex; - typedef std::pair< typeVectorVertex, Filtration_value> typeSimplex; - typedef std::pair< Simplex_handle, bool > typePairSimplexBool; - - typedef int Witness_id; - typedef int Landmark_id; - typedef std::list< Vertex_handle > ActiveWitnessList; - - private: - int nbL; // Number of landmarks - Simplicial_complex& sc; // Simplicial complex - - public: - /** @name Constructor - */ - - //@{ - - /** - * \brief Iterative construction of the relaxed witness complex. - * \details The witness complex is written in sc_ basing on a matrix knn - * of k nearest neighbours of the form {witnesses}x{landmarks} and - * and a matrix distances of distances to these landmarks from witnesses. - * The parameter alpha defines relaxation and - * limD defines the - * - * The line lengths in one given matrix can differ, - * however both matrices have the same corresponding line lengths. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges. - * - * Constructor takes into account at most (dim+1) - * first landmarks from each landmark range to construct simplices. - * - * Landmarks are supposed to be in [0,nbL_-1] - */ - - template< typename KNearestNeighbours > - Strong_relaxed_witness_complex(std::vector< std::vector<double> > const & distances, - KNearestNeighbours const & knn, - Simplicial_complex & sc_, - int nbL_, - double alpha2, - unsigned limD) : nbL(nbL_), sc(sc_) { - int nbW = knn.size(); - typeVectorVertex vv; - //int counter = 0; - /* The list of still useful witnesses - * it will diminuish in the course of iterations - */ - ActiveWitnessList active_w;// = new ActiveWitnessList(); - for (int i = 0; i != nbL; ++i) { - // initial fill of 0-dimensional simplices - // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore - //counter++; - vv = {i}; - sc.insert_simplex(vv, Filtration_value(0.0)); - /* TODO Error if not inserted : normally no need here though*/ - } - for (int i=0; i != nbW; ++i) { - // int i_end = limD+1; - // if (knn[i].size() < limD+1) - // i_end = knn[i].size(); - // double dist_wL = *(distances[i].begin()); - // while (distances[i][i_end] > dist_wL + alpha2) - // i_end--; - // add_all_witnessed_faces(distances[i].begin(), - // knn[i].begin(), - // knn[i].begin() + i_end + 1); - unsigned j_end = 0; - while (j_end < distances[i].size() && j_end <= limD && distances[i][j_end] <= distances[i][0] + alpha2) { - std::vector<int> simplex; - for (unsigned j = 0; j <= j_end; ++j) - simplex.push_back(knn[i][j]); - assert(distances[i][j_end] - distances[i][0] >= 0); - sc.insert_simplex_and_subfaces(simplex, distances[i][j_end] - distances[i][0]); - j_end++; - } - } - sc.set_dimension(limD); - } - - //@} - -private: - /* \brief Adds recursively all the faces of a certain dimension dim witnessed by the same witness - * Iterator is needed to know until how far we can take landmarks to form simplexes - * simplex is the prefix of the simplexes to insert - * The output value indicates if the witness rests active or not - */ - void add_all_witnessed_faces(std::vector<double>::const_iterator curr_d, - std::vector<int>::const_iterator curr_l, - std::vector<int>::const_iterator end) - { - std::vector<int> simplex; - std::vector<int>::const_iterator l_end = curr_l; - for (; l_end != end; ++l_end) { - std::vector<int>::const_iterator l_it = curr_l; - std::vector<double>::const_iterator d_it = curr_d; - simplex = {}; - for (; l_it != l_end; ++l_it, ++d_it) - simplex.push_back(*l_it); - sc.insert_simplex_and_subfaces(simplex, *(d_it--) - *curr_d); - } - } - - /** \brief Check if the facets of the k-dimensional simplex witnessed - * by witness witness_id are already in the complex. - * inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id - */ - bool all_faces_in(std::vector<int>& simplex, double* filtration_value) - { - std::vector< int > facet; - for (std::vector<int>::iterator not_it = simplex.begin(); not_it != simplex.end(); ++not_it) - { - facet.clear(); - for (std::vector<int>::iterator it = simplex.begin(); it != simplex.end(); ++it) - if (it != not_it) - facet.push_back(*it); - Simplex_handle facet_sh = sc.find(facet); - if (facet_sh == sc.null_simplex()) - return false; - else if (sc.filtration(facet_sh) > *filtration_value) - *filtration_value = sc.filtration(facet_sh); - } - return true; - } - - bool is_face(Simplex_handle face, Simplex_handle coface) - { - // vertex range is sorted in decreasing order - auto fvr = sc.simplex_vertex_range(face); - auto cfvr = sc.simplex_vertex_range(coface); - auto fv_it = fvr.begin(); - auto cfv_it = cfvr.begin(); - while (fv_it != fvr.end() && cfv_it != cfvr.end()) { - if (*fv_it < *cfv_it) - ++cfv_it; - else if (*fv_it == *cfv_it) { - ++cfv_it; - ++fv_it; - } - else - return false; - - } - return (fv_it == fvr.end()); - } - - // void erase_simplex(Simplex_handle sh) - // { - // auto siblings = sc.self_siblings(sh); - // auto oncles = siblings->oncles(); - // int prev_vertex = siblings->parent(); - // siblings->members().erase(sh->first); - // if (siblings->members().empty()) { - // typename typedef Simplicial_complex::Siblings Siblings; - // oncles->members().find(prev_vertex)->second.assign_children(new Siblings(oncles, prev_vertex)); - // assert(!sc.has_children(oncles->members().find(prev_vertex))); - // //delete siblings; - // } - - // } - - void elementary_collapse(Simplex_handle face_sh, Simplex_handle coface_sh) - { - erase_simplex(coface_sh); - erase_simplex(face_sh); - } - -public: - // void collapse(std::vector<Simplex_handle>& simplices) - // { - // // Get a vector of simplex handles ordered by filtration value - // std::cout << sc << std::endl; - // //std::vector<Simplex_handle> simplices; - // for (Simplex_handle sh: sc.filtration_simplex_range()) - // simplices.push_back(sh); - // // std::sort(simplices.begin(), - // // simplices.end(), - // // [&](Simplex_handle sh1, Simplex_handle sh2) - // // { double f1 = sc.filtration(sh1), f2 = sc.filtration(sh2); - // // return (f1 > f2) || (f1 >= f2 && sc.dimension(sh1) > sc.dimension(sh2)); }); - // // Double iteration - // auto face_it = simplices.rbegin(); - // while (face_it != simplices.rend() && sc.filtration(*face_it) != 0) { - // int coface_count = 0; - // auto reduced_coface = simplices.rbegin(); - // for (auto coface_it = simplices.rbegin(); coface_it != simplices.rend() && sc.filtration(*coface_it) != 0; ++coface_it) - // if (face_it != coface_it && is_face(*face_it, *coface_it)) { - // coface_count++; - // if (coface_count == 1) - // reduced_coface = coface_it; - // else - // break; - // } - // if (coface_count == 1) { - // std::cout << "Erase ( "; - // for (auto v: sc.simplex_vertex_range(*(--reduced_coface.base()))) - // std::cout << v << " "; - - // simplices.erase(--(reduced_coface.base())); - // //elementary_collapse(*face_it, *reduced_coface); - // std::cout << ") and then ( "; - // for (auto v: sc.simplex_vertex_range(*(--face_it.base()))) - // std::cout << v << " "; - // std::cout << ")\n"; - // simplices.erase(--((face_it++).base())); - // //face_it = simplices.rbegin(); - // //std::cout << "Size of vector: " << simplices.size() << "\n"; - // } - // else - // face_it++; - // } - // sc.initialize_filtration(); - // //std::cout << sc << std::endl; - // } - - template <class Dim_lists> - void collapse(Dim_lists& dim_lists) - { - dim_lists.collapse(); - } - -private: - - /** Collapse recursively boundary faces of the given simplex - * if its filtration is bigger than alpha_lim. - */ - void rec_boundary_collapse(Simplex_handle sh, FT alpha_lim) - { - for (Simplex_handle face_it : sc.boundary_simplex_range()) { - - } - - } - -}; //class Strong_relaxed_witness_complex - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Weak_relaxed_witness_complex.h b/src/Witness_complex/include/gudhi/Weak_relaxed_witness_complex.h deleted file mode 100644 index a1aedd8e..00000000 --- a/src/Witness_complex/include/gudhi/Weak_relaxed_witness_complex.h +++ /dev/null @@ -1,379 +0,0 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ - * library for computational topology. - * - * Author(s): Siargey Kachanovich - * - * Copyright (C) 2015 INRIA Sophia Antipolis-Méditerranée (France) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#ifndef WEAK_RELAXED_WITNESS_COMPLEX_H_ -#define WEAK_RELAXED_WITNESS_COMPLEX_H_ - -#include <boost/container/flat_map.hpp> -#include <boost/iterator/transform_iterator.hpp> -#include <algorithm> -#include <utility> -#include "gudhi/reader_utils.h" -#include "gudhi/distance_functions.h" -#include "gudhi/Simplex_tree.h" -#include <vector> -#include <list> -#include <set> -#include <queue> -#include <limits> -#include <math.h> -#include <ctime> -#include <iostream> - -// Needed for nearest neighbours -#include <CGAL/Cartesian_d.h> -#include <CGAL/Search_traits.h> -#include <CGAL/Search_traits_adapter.h> -#include <CGAL/property_map.h> -#include <CGAL/Epick_d.h> -#include <CGAL/Orthogonal_k_neighbor_search.h> - -#include <boost/tuple/tuple.hpp> -#include <boost/iterator/zip_iterator.hpp> -#include <boost/iterator/counting_iterator.hpp> -#include <boost/range/iterator_range.hpp> - -// Needed for the adjacency graph in bad link search -#include <boost/graph/graph_traits.hpp> -#include <boost/graph/adjacency_list.hpp> -#include <boost/graph/connected_components.hpp> - -namespace Gudhi { - -namespace witness_complex { - - /** \addtogroup simplex_tree - * Witness complex is a simplicial complex defined on two sets of points in \f$\mathbf{R}^D\f$: - * \f$W\f$ set of witnesses and \f$L \subseteq W\f$ set of landmarks. The simplices are based on points in \f$L\f$ - * and a simplex belongs to the witness complex if and only if it is witnessed (there exists a point \f$w \in W\f$ such that - * w is closer to the vertices of this simplex than others) and all of its faces are witnessed as well. - */ -template< class Simplicial_complex > -class Weak_relaxed_witness_complex { - -private: - struct Active_witness { - int witness_id; - int landmark_id; - - Active_witness(int witness_id_, int landmark_id_) - : witness_id(witness_id_), - landmark_id(landmark_id_) { } - }; - -private: - typedef typename Simplicial_complex::Simplex_handle Simplex_handle; - typedef typename Simplicial_complex::Vertex_handle Vertex_handle; - typedef typename Simplicial_complex::Filtration_value FT; - - typedef std::vector< double > Point_t; - typedef std::vector< Point_t > Point_Vector; - - // typedef typename Simplicial_complex::Filtration_value Filtration_value; - typedef std::vector< Vertex_handle > typeVectorVertex; - typedef std::pair< typeVectorVertex, Filtration_value> typeSimplex; - typedef std::pair< Simplex_handle, bool > typePairSimplexBool; - - typedef int Witness_id; - typedef int Landmark_id; - typedef std::list< Vertex_handle > ActiveWitnessList; - - private: - int nbL; // Number of landmarks - Simplicial_complex& sc; // Simplicial complex - - public: - /** @name Constructor - */ - - //@{ - - /** - * \brief Iterative construction of the relaxed witness complex. - * \details The witness complex is written in sc_ basing on a matrix knn - * of k nearest neighbours of the form {witnesses}x{landmarks} and - * and a matrix distances of distances to these landmarks from witnesses. - * The parameter alpha defines relaxation and - * limD defines the - * - * The line lengths in one given matrix can differ, - * however both matrices have the same corresponding line lengths. - * - * The type KNearestNeighbors can be seen as - * Witness_range<Closest_landmark_range<Vertex_handle>>, where - * Witness_range and Closest_landmark_range are random access ranges. - * - * Constructor takes into account at most (dim+1) - * first landmarks from each landmark range to construct simplices. - * - * Landmarks are supposed to be in [0,nbL_-1] - */ - - template< typename KNearestNeighbours > - Weak_relaxed_witness_complex(std::vector< std::vector<double> > const & distances, - KNearestNeighbours const & knn, - Simplicial_complex & sc_, - int nbL_, - double alpha, - int limD) : nbL(nbL_), sc(sc_) { - int nbW = knn.size(); - typeVectorVertex vv; - //int counter = 0; - /* The list of still useful witnesses - * it will diminuish in the course of iterations - */ - ActiveWitnessList active_w;// = new ActiveWitnessList(); - for (int i = 0; i != nbL; ++i) { - // initial fill of 0-dimensional simplices - // by doing it we don't assume that landmarks are necessarily witnesses themselves anymore - //counter++; - vv = {i}; - sc.insert_simplex(vv, Filtration_value(0.0)); - /* TODO Error if not inserted : normally no need here though*/ - } - int k = 1; /* current dimension in iterative construction */ - for (int i=0; i != nbW; ++i) - active_w.push_back(i); - while (!active_w.empty() && k <= limD && k < nbL ) - { - typename ActiveWitnessList::iterator aw_it = active_w.begin(); - while (aw_it != active_w.end()) - { - std::vector<int> simplex; - bool ok = add_all_faces_of_dimension(k, - alpha, - std::numeric_limits<double>::infinity(), - distances[*aw_it].begin(), - knn[*aw_it].begin(), - simplex, - knn[*aw_it].end()); - if (!ok) - active_w.erase(aw_it++); //First increase the iterator and then erase the previous element - else - aw_it++; - } - k++; - } - sc.set_dimension(limD); - } - - //@} - -private: - /* \brief Adds recursively all the faces of a certain dimension dim witnessed by the same witness - * Iterator is needed to know until how far we can take landmarks to form simplexes - * simplex is the prefix of the simplexes to insert - * The output value indicates if the witness rests active or not - */ - bool add_all_faces_of_dimension(int dim, - double alpha, - double norelax_dist, - std::vector<double>::const_iterator curr_d, - std::vector<int>::const_iterator curr_l, - std::vector<int>& simplex, - std::vector<int>::const_iterator end) - { - if (curr_l == end) - return false; - bool will_be_active = false; - std::vector<int>::const_iterator l_it = curr_l; - std::vector<double>::const_iterator d_it = curr_d; - if (dim > 0) - for (; *d_it - alpha <= norelax_dist && l_it != end; ++l_it, ++d_it) { - simplex.push_back(*l_it); - if (sc.find(simplex) != sc.null_simplex()) - will_be_active = add_all_faces_of_dimension(dim-1, - alpha, - norelax_dist, - d_it+1, - l_it+1, - simplex, - end) || will_be_active; - simplex.pop_back(); - // If norelax_dist is infinity, change to first omitted distance - if (*d_it <= norelax_dist) - norelax_dist = *d_it; - will_be_active = add_all_faces_of_dimension(dim-1, - alpha, - norelax_dist, - d_it+1, - l_it+1, - simplex, - end) || will_be_active; - } - else if (dim == 0) - for (; *d_it - alpha <= norelax_dist && l_it != end; ++l_it, ++d_it) { - simplex.push_back(*l_it); - double filtration_value = 0; - // if norelax_dist is infinite, relaxation is 0. - if (*d_it > norelax_dist) - filtration_value = *d_it - norelax_dist; - if (all_faces_in(simplex, &filtration_value)) { - will_be_active = true; - sc.insert_simplex(simplex, filtration_value); - } - simplex.pop_back(); - // If norelax_dist is infinity, change to first omitted distance - if (*d_it < norelax_dist) - norelax_dist = *d_it; - } - return will_be_active; - } - - /** \brief Check if the facets of the k-dimensional simplex witnessed - * by witness witness_id are already in the complex. - * inserted_vertex is the handle of the (k+1)-th vertex witnessed by witness_id - */ - bool all_faces_in(std::vector<int>& simplex, double* filtration_value) - { - std::vector< int > facet; - for (std::vector<int>::iterator not_it = simplex.begin(); not_it != simplex.end(); ++not_it) - { - facet.clear(); - for (std::vector<int>::iterator it = simplex.begin(); it != simplex.end(); ++it) - if (it != not_it) - facet.push_back(*it); - Simplex_handle facet_sh = sc.find(facet); - if (facet_sh == sc.null_simplex()) - return false; - else if (sc.filtration(facet_sh) > *filtration_value) - *filtration_value = sc.filtration(facet_sh); - } - return true; - } - - bool is_face(Simplex_handle face, Simplex_handle coface) - { - // vertex range is sorted in decreasing order - auto fvr = sc.simplex_vertex_range(face); - auto cfvr = sc.simplex_vertex_range(coface); - auto fv_it = fvr.begin(); - auto cfv_it = cfvr.begin(); - while (fv_it != fvr.end() && cfv_it != cfvr.end()) { - if (*fv_it < *cfv_it) - ++cfv_it; - else if (*fv_it == *cfv_it) { - ++cfv_it; - ++fv_it; - } - else - return false; - - } - return (fv_it == fvr.end()); - } - - // void erase_simplex(Simplex_handle sh) - // { - // auto siblings = sc.self_siblings(sh); - // auto oncles = siblings->oncles(); - // int prev_vertex = siblings->parent(); - // siblings->members().erase(sh->first); - // if (siblings->members().empty()) { - // typename typedef Simplicial_complex::Siblings Siblings; - // oncles->members().find(prev_vertex)->second.assign_children(new Siblings(oncles, prev_vertex)); - // assert(!sc.has_children(oncles->members().find(prev_vertex))); - // //delete siblings; - // } - - // } - - void elementary_collapse(Simplex_handle face_sh, Simplex_handle coface_sh) - { - erase_simplex(coface_sh); - erase_simplex(face_sh); - } - -public: - // void collapse(std::vector<Simplex_handle>& simplices) - // { - // // Get a vector of simplex handles ordered by filtration value - // std::cout << sc << std::endl; - // //std::vector<Simplex_handle> simplices; - // for (Simplex_handle sh: sc.filtration_simplex_range()) - // simplices.push_back(sh); - // // std::sort(simplices.begin(), - // // simplices.end(), - // // [&](Simplex_handle sh1, Simplex_handle sh2) - // // { double f1 = sc.filtration(sh1), f2 = sc.filtration(sh2); - // // return (f1 > f2) || (f1 >= f2 && sc.dimension(sh1) > sc.dimension(sh2)); }); - // // Double iteration - // auto face_it = simplices.rbegin(); - // while (face_it != simplices.rend() && sc.filtration(*face_it) != 0) { - // int coface_count = 0; - // auto reduced_coface = simplices.rbegin(); - // for (auto coface_it = simplices.rbegin(); coface_it != simplices.rend() && sc.filtration(*coface_it) != 0; ++coface_it) - // if (face_it != coface_it && is_face(*face_it, *coface_it)) { - // coface_count++; - // if (coface_count == 1) - // reduced_coface = coface_it; - // else - // break; - // } - // if (coface_count == 1) { - // std::cout << "Erase ( "; - // for (auto v: sc.simplex_vertex_range(*(--reduced_coface.base()))) - // std::cout << v << " "; - - // simplices.erase(--(reduced_coface.base())); - // //elementary_collapse(*face_it, *reduced_coface); - // std::cout << ") and then ( "; - // for (auto v: sc.simplex_vertex_range(*(--face_it.base()))) - // std::cout << v << " "; - // std::cout << ")\n"; - // simplices.erase(--((face_it++).base())); - // //face_it = simplices.rbegin(); - // //std::cout << "Size of vector: " << simplices.size() << "\n"; - // } - // else - // face_it++; - // } - // sc.initialize_filtration(); - // //std::cout << sc << std::endl; - // } - - template <class Dim_lists> - void collapse(Dim_lists& dim_lists) - { - dim_lists.collapse(); - } - -private: - - /** Collapse recursively boundary faces of the given simplex - * if its filtration is bigger than alpha_lim. - */ - void rec_boundary_collapse(Simplex_handle sh, FT alpha_lim) - { - for (Simplex_handle face_it : sc.boundary_simplex_range()) { - - } - - } - -}; //class Weak_relaxed_witness_complex - -} // namespace witness_complex - -} // namespace Gudhi - -#endif diff --git a/src/Witness_complex/include/gudhi/Witness_complex.h b/src/Witness_complex/include/gudhi/Witness_complex.h index e732cb18..a16f9270 100644 --- a/src/Witness_complex/include/gudhi/Witness_complex.h +++ b/src/Witness_complex/include/gudhi/Witness_complex.h @@ -173,8 +173,6 @@ private: else aw_it++; } - std::cout << "Active witnesses after dim=" << k << " is finished: " << active_witnesses.size() << "\n"; - std::cout << complex << "\n"; k++; } return true; @@ -233,11 +231,8 @@ private: simplex.push_back(l_it->first); double filtration_value = 0; // if norelax_dist is infinite, relaxation is 0. - //std::cout << "landmark_id=" << l_it->first << " distance=" << l_it->second << "\n"; - // std::size_t landmark_id = l_it->first; - // double distance = l_it->second; if (l_it->second > norelax_dist2) - filtration_value = l_it->second - norelax_dist2; + filtration_value = l_it->second - norelax_dist2; if (all_faces_in(simplex, &filtration_value, sc)) { will_be_active = true; sc.insert_simplex(simplex, filtration_value); |