utilities.h
1 /* This file is part of the Gudhi Library. The Gudhi library
2  * (Geometric Understanding in Higher Dimensions) is a generic C++
3  * library for computational topology.
4  *
5  * Author(s): Clement Jamin
6  *
7  * Copyright (C) 2016 Inria
8  *
9  * This program is free software: you can redistribute it and/or modify
10  * it under the terms of the GNU General Public License as published by
11  * the Free Software Foundation, either version 3 of the License, or
12  * (at your option) any later version.
13  *
14  * This program is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  * GNU General Public License for more details.
18  *
19  * You should have received a copy of the GNU General Public License
20  * along with this program. If not, see <http://www.gnu.org/licenses/>.
21  */
22 
23 #ifndef TANGENTIAL_COMPLEX_UTILITIES_H_
24 #define TANGENTIAL_COMPLEX_UTILITIES_H_
25 
26 #include <CGAL/Dimension.h>
27 #include <CGAL/Combination_enumerator.h>
28 #include <CGAL/IO/Triangulation_off_ostream.h>
29 
30 #include <boost/container/flat_set.hpp>
31 
32 #include <Eigen/Core>
33 #include <Eigen/Eigen>
34 
35 #include <set>
36 #include <vector>
37 #include <array>
38 #include <fstream>
39 #include <atomic>
40 #include <cmath> // for std::sqrt
41 
42 namespace Gudhi {
43 namespace tangential_complex {
44 namespace internal {
45 
46 // Provides copy constructors to std::atomic so that
47 // it can be used in a vector
48 template <typename T>
49 struct Atomic_wrapper
50 : public std::atomic<T> {
51  typedef std::atomic<T> Base;
52 
53  Atomic_wrapper() { }
54 
55  Atomic_wrapper(const T &t) : Base(t) { }
56 
57  Atomic_wrapper(const std::atomic<T> &a) : Base(a.load()) { }
58 
59  Atomic_wrapper(const Atomic_wrapper &other) : Base(other.load()) { }
60 
61  Atomic_wrapper &operator=(const T &other) {
62  Base::store(other);
63  return *this;
64  }
65 
66  Atomic_wrapper &operator=(const std::atomic<T> &other) {
67  Base::store(other.load());
68  return *this;
69  }
70 
71  Atomic_wrapper &operator=(const Atomic_wrapper &other) {
72  Base::store(other.load());
73  return *this;
74  }
75 };
76 
77 // Modifies v in-place
78 template <typename K>
79 typename K::Vector_d& normalize_vector(typename K::Vector_d& v,
80  K const& k) {
81  v = k.scaled_vector_d_object()(
82  v, typename K::FT(1) / std::sqrt(k.squared_length_d_object()(v)));
83  return v;
84 }
85 
86 template<typename Kernel>
87 struct Basis {
88  typedef typename Kernel::FT FT;
89  typedef typename Kernel::Point_d Point;
90  typedef typename Kernel::Vector_d Vector;
91  typedef typename std::vector<Vector>::const_iterator const_iterator;
92 
93  std::size_t m_origin;
94  std::vector<Vector> m_vectors;
95 
96  std::size_t origin() const {
97  return m_origin;
98  }
99 
100  void set_origin(std::size_t o) {
101  m_origin = o;
102  }
103 
104  const_iterator begin() const {
105  return m_vectors.begin();
106  }
107 
108  const_iterator end() const {
109  return m_vectors.end();
110  }
111 
112  std::size_t size() const {
113  return m_vectors.size();
114  }
115 
116  Vector& operator[](const std::size_t i) {
117  return m_vectors[i];
118  }
119 
120  const Vector& operator[](const std::size_t i) const {
121  return m_vectors[i];
122  }
123 
124  void push_back(const Vector& v) {
125  m_vectors.push_back(v);
126  }
127 
128  void reserve(const std::size_t s) {
129  m_vectors.reserve(s);
130  }
131 
132  Basis() { }
133 
134  Basis(std::size_t origin) : m_origin(origin) { }
135 
136  Basis(std::size_t origin, const std::vector<Vector>& vectors)
137  : m_origin(origin), m_vectors(vectors) { }
138 
139  int dimension() const {
140  return static_cast<int> (m_vectors.size());
141  }
142 };
143 
144 // 1st line: number of points
145 // Then one point per line
146 template <typename Kernel, typename Point_range>
147 std::ostream &export_point_set(
148  Kernel const& k,
149  Point_range const& points,
150  std::ostream & os,
151  const char *coord_separator = " ") {
152  // Kernel functors
153  typename Kernel::Construct_cartesian_const_iterator_d ccci =
154  k.construct_cartesian_const_iterator_d_object();
155 
156  os << points.size() << "\n";
157 
158  typename Point_range::const_iterator it_p = points.begin();
159  typename Point_range::const_iterator it_p_end = points.end();
160  // For each point p
161  for (; it_p != it_p_end; ++it_p) {
162  for (auto it = ccci(*it_p); it != ccci(*it_p, 0); ++it)
163  os << CGAL::to_double(*it) << coord_separator;
164 
165  os << "\n";
166  }
167 
168  return os;
169 }
170 
171 // Compute all the k-combinations of elements
172 // Output_iterator::value_type must be
173 // boost::container::flat_set<std::size_t>
174 template <typename Elements_container, typename Output_iterator>
175 void combinations(const Elements_container elements, int k,
176  Output_iterator combinations) {
177  std::size_t n = elements.size();
178  std::vector<bool> booleans(n, false);
179  std::fill(booleans.begin() + n - k, booleans.end(), true);
180  do {
181  boost::container::flat_set<std::size_t> combination;
182  typename Elements_container::const_iterator it_elt = elements.begin();
183  for (std::size_t i = 0; i < n; ++i, ++it_elt) {
184  if (booleans[i])
185  combination.insert(*it_elt);
186  }
187  *combinations++ = combination;
188  } while (std::next_permutation(booleans.begin(), booleans.end()));
189 }
190 
191 } // namespace internal
192 } // namespace tangential_complex
193 } // namespace Gudhi
194 
195 #endif // TANGENTIAL_COMPLEX_UTILITIES_H_
Definition: SimplicialComplexForAlpha.h:26
GUDHI  Version 2.3.0  - C++ library for Topological Data Analysis (TDA) and Higher Dimensional Geometry Understanding.  - Copyright : GPL v3 Generated on Fri Oct 5 2018 15:05:03 for GUDHI by Doxygen 1.8.13