summaryrefslogtreecommitdiff
path: root/src/Witness_complex/include
diff options
context:
space:
mode:
authorvrouvrea <vrouvrea@636b058d-ea47-450e-bf9e-a15bfbe3eedb>2016-01-21 16:56:41 +0000
committervrouvrea <vrouvrea@636b058d-ea47-450e-bf9e-a15bfbe3eedb>2016-01-21 16:56:41 +0000
commit6697ad2e4f76e706eae12f5153988691daea1202 (patch)
tree2d1bd3463fad4cc1ff34262d8b6de84bc84c29cd /src/Witness_complex/include
parentc2a938f2de9c5cd995959578a2e6855ec2004ba1 (diff)
parentd13a8867ca368b1f56be3ba151d2042728fb4754 (diff)
cpplint fixes
git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/VR_witness@984 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: c663c383e09719a8c70c1100bb6e21d8fc49e597
Diffstat (limited to 'src/Witness_complex/include')
-rw-r--r--src/Witness_complex/include/gudhi/Landmark_choice_by_furthest_point.h107
-rw-r--r--src/Witness_complex/include/gudhi/Landmark_choice_by_random_point.h96
-rw-r--r--src/Witness_complex/include/gudhi/Witness_complex.h251
-rw-r--r--src/Witness_complex/include/gudhi/Witness_complex_doc.h38
4 files changed, 492 insertions, 0 deletions
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
new file mode 100644
index 00000000..bb7e87f5
--- /dev/null
+++ b/src/Witness_complex/include/gudhi/Landmark_choice_by_furthest_point.h
@@ -0,0 +1,107 @@
+/* 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 <limits> // for numeric_limits<>
+#include <algorithm> // for sort
+#include <vector>
+
+namespace Gudhi {
+
+namespace witness_complex {
+
+/**
+ * \class Landmark_choice_by_furthest_point
+ * \brief The class `Landmark_choice_by_furthest_point` allows to construct the matrix
+ * of closest landmarks per witness by iteratively choosing the furthest witness
+ * from the set of already chosen landmarks as the new landmark.
+ * \ingroup witness_complex
+ */
+
+class Landmark_choice_by_furthest_point {
+ private:
+ typedef std::vector<int> typeVectorVertex;
+
+ public:
+ /**
+ * \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`.
+ */
+
+ template <typename KNearestNeighbours,
+ typename Point_random_access_range>
+ Landmark_choice_by_furthest_point(Point_random_access_range const &points,
+ int nbL,
+ KNearestNeighbours &knn) {
+ int nb_points = points.end() - points.begin();
+ 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
+ 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, *(points.begin() + 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 (unsigned i = 0; i < points.size(); ++i)
+ std::sort(knn[i].begin(),
+ knn[i].end(),
+ [&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
new file mode 100644
index 00000000..215c9e65
--- /dev/null
+++ b/src/Witness_complex/include/gudhi/Landmark_choice_by_random_point.h
@@ -0,0 +1,96 @@
+/* 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 <queue> // for priority_queue<>
+#include <utility> // for pair<>
+#include <vector>
+#include <set>
+
+namespace Gudhi {
+
+namespace witness_complex {
+
+/**
+ * \class Landmark_choice_by_random_point
+ * \brief The class `Landmark_choice_by_random_point` allows to construct the matrix
+ * of closest landmarks per witness by iteratively choosing a random non-chosen witness
+ * as a new landmark.
+ * \ingroup witness_complex
+ */
+
+class Landmark_choice_by_random_point {
+ public:
+ /** \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.
+ */
+
+ template <typename KNearestNeighbours,
+ typename Point_random_access_range>
+ Landmark_choice_by_random_point(Point_random_access_range const &points,
+ int nbL,
+ KNearestNeighbours &knn) {
+ int nbP = points.end() - points.begin();
+ 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 = points.begin()->size();
+ 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/Witness_complex.h b/src/Witness_complex/include/gudhi/Witness_complex.h
new file mode 100644
index 00000000..27f7a9c0
--- /dev/null
+++ b/src/Witness_complex/include/gudhi/Witness_complex.h
@@ -0,0 +1,251 @@
+/* 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 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
+ double density; // Desired density
+ 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 >
+ 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;
+ typeSimplex simplex;
+ typePairSimplexBool returnValue;
+ 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};
+ returnValue = sc.insert_simplex(vv);
+ // TODO(SK) 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);
+ // std::cout << "Successfully added edges" << std::endl;
+ while (!active_w.empty() && k < dim) {
+ // std::cout << "Started the step k=" << k << std::endl;
+ typename ActiveWitnessList::iterator it = active_w.begin();
+ while (it != active_w.end()) {
+ typeVectorVertex simplex_vector;
+ /* THE INSERTION: Checking if all the subfaces are in the simplex tree*/
+ bool ok = all_faces_in(knn, *it, k);
+ if (ok) {
+ for (int i = 0; i != k + 1; ++i)
+ simplex_vector.push_back(knn[*it][i]);
+ returnValue = sc.insert_simplex(simplex_vector, 0.0);
+ it++;
+ } else {
+ active_w.erase(it++); // First increase the iterator and then erase the previous element
+ }
+ }
+ k++;
+ }
+ }
+
+ //@}
+
+ 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>
+ 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/Witness_complex_doc.h b/src/Witness_complex/include/gudhi/Witness_complex_doc.h
new file mode 100644
index 00000000..e9f78170
--- /dev/null
+++ b/src/Witness_complex/include/gudhi/Witness_complex_doc.h
@@ -0,0 +1,38 @@
+#ifndef WITNESS_COMPLEX_DOC_H_
+#define WITNESS_COMPLEX_DOC_H_
+
+/**
+ \defgroup witness_complex Witness complex
+
+ \author Siargey Kachanovich
+
+ \section Definitions
+
+ Witness complex \f$ Wit(W,L) \f$ is a simplicial complex defined on two sets of points in \f$\mathbb{R}^D\f$:
+
+ \li \f$W\f$ set of **witnesses** and
+ \li \f$L \subseteq W\f$ set of **landmarks**.
+
+ The simplices are based on landmarks
+ and a simplex belongs to the witness complex if and only if it is witnessed, that is:
+
+ \f$ \sigma \subset L \f$ is witnessed if there exists a point \f$w \in W\f$ such that
+ w is closer to the vertices of \f$ \sigma \f$ than other points in \f$ L \f$ and all of its faces are witnessed as well.
+
+ \section Implementation
+
+ The principal class of this module is Gudhi::Witness_complex.
+
+ In both cases, the constructor for this class takes a {witness}x{closest_landmarks} table, where each row represents a witness and consists of landmarks sorted by distance to this witness.
+ This table can be constructed by two additional classes Landmark_choice_by_furthest_point and Landmark_choice_by_random_point also included in the module.
+
+ *\image html "bench_Cy8.png" "Running time as function on number of landmarks" width=10cm
+ *\image html "bench_sphere.png" "Running time as function on number of witnesses for |L|=300" width=10cm
+
+
+ \copyright GNU General Public License v3.
+
+
+ */
+
+#endif // WITNESS_COMPLEX_DOC_H_