utilities.h
1 /* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
2  * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
3  * Author(s): Clement Jamin
4  *
5  * Copyright (C) 2016 Inria
6  *
7  * Modification(s):
8  * - YYYY/MM Author: Description of the modification
9  */
10 
11 #ifndef TANGENTIAL_COMPLEX_UTILITIES_H_
12 #define TANGENTIAL_COMPLEX_UTILITIES_H_
13 
14 #include <CGAL/Dimension.h>
15 #include <CGAL/Combination_enumerator.h>
16 #include <CGAL/IO/Triangulation_off_ostream.h>
17 
18 #include <boost/container/flat_set.hpp>
19 
20 #include <Eigen/Core>
21 #include <Eigen/Eigen>
22 
23 #include <set>
24 #include <vector>
25 #include <array>
26 #include <fstream>
27 #include <atomic>
28 #include <cmath> // for std::sqrt
29 
30 namespace Gudhi {
31 namespace tangential_complex {
32 namespace internal {
33 
34 // Provides copy constructors to std::atomic so that
35 // it can be used in a vector
36 template <typename T>
37 struct Atomic_wrapper
38 : public std::atomic<T> {
39  typedef std::atomic<T> Base;
40 
41  Atomic_wrapper() { }
42 
43  Atomic_wrapper(const T &t) : Base(t) { }
44 
45  Atomic_wrapper(const std::atomic<T> &a) : Base(a.load()) { }
46 
47  Atomic_wrapper(const Atomic_wrapper &other) : Base(other.load()) { }
48 
49  Atomic_wrapper &operator=(const T &other) {
50  Base::store(other);
51  return *this;
52  }
53 
54  Atomic_wrapper &operator=(const std::atomic<T> &other) {
55  Base::store(other.load());
56  return *this;
57  }
58 
59  Atomic_wrapper &operator=(const Atomic_wrapper &other) {
60  Base::store(other.load());
61  return *this;
62  }
63 };
64 
65 // Modifies v in-place
66 template <typename K>
67 typename K::Vector_d& normalize_vector(typename K::Vector_d& v,
68  K const& k) {
69  v = k.scaled_vector_d_object()(
70  v, typename K::FT(1) / std::sqrt(k.squared_length_d_object()(v)));
71  return v;
72 }
73 
74 template<typename Kernel>
75 struct Basis {
76  typedef typename Kernel::FT FT;
77  typedef typename Kernel::Point_d Point;
78  typedef typename Kernel::Vector_d Vector;
79  typedef typename std::vector<Vector>::const_iterator const_iterator;
80 
81  std::size_t m_origin;
82  std::vector<Vector> m_vectors;
83 
84  std::size_t origin() const {
85  return m_origin;
86  }
87 
88  void set_origin(std::size_t o) {
89  m_origin = o;
90  }
91 
92  const_iterator begin() const {
93  return m_vectors.begin();
94  }
95 
96  const_iterator end() const {
97  return m_vectors.end();
98  }
99 
100  std::size_t size() const {
101  return m_vectors.size();
102  }
103 
104  Vector& operator[](const std::size_t i) {
105  return m_vectors[i];
106  }
107 
108  const Vector& operator[](const std::size_t i) const {
109  return m_vectors[i];
110  }
111 
112  void push_back(const Vector& v) {
113  m_vectors.push_back(v);
114  }
115 
116  void reserve(const std::size_t s) {
117  m_vectors.reserve(s);
118  }
119 
120  Basis() { }
121 
122  Basis(std::size_t origin) : m_origin(origin) { }
123 
124  Basis(std::size_t origin, const std::vector<Vector>& vectors)
125  : m_origin(origin), m_vectors(vectors) { }
126 
127  int dimension() const {
128  return static_cast<int> (m_vectors.size());
129  }
130 };
131 
132 // 1st line: number of points
133 // Then one point per line
134 template <typename Kernel, typename Point_range>
135 std::ostream &export_point_set(
136  Kernel const& k,
137  Point_range const& points,
138  std::ostream & os,
139  const char *coord_separator = " ") {
140  // Kernel functors
141  typename Kernel::Construct_cartesian_const_iterator_d ccci =
142  k.construct_cartesian_const_iterator_d_object();
143 
144  os << points.size() << "\n";
145 
146  typename Point_range::const_iterator it_p = points.begin();
147  typename Point_range::const_iterator it_p_end = points.end();
148  // For each point p
149  for (; it_p != it_p_end; ++it_p) {
150  for (auto it = ccci(*it_p); it != ccci(*it_p, 0); ++it)
151  os << CGAL::to_double(*it) << coord_separator;
152 
153  os << "\n";
154  }
155 
156  return os;
157 }
158 
159 // Compute all the k-combinations of elements
160 // Output_iterator::value_type must be
161 // boost::container::flat_set<std::size_t>
162 template <typename Elements_container, typename Output_iterator>
163 void combinations(const Elements_container elements, int k,
164  Output_iterator combinations) {
165  std::size_t n = elements.size();
166  std::vector<bool> booleans(n, false);
167  std::fill(booleans.begin() + n - k, booleans.end(), true);
168  do {
169  boost::container::flat_set<std::size_t> combination;
170  typename Elements_container::const_iterator it_elt = elements.begin();
171  for (std::size_t i = 0; i < n; ++i, ++it_elt) {
172  if (booleans[i])
173  combination.insert(*it_elt);
174  }
175  *combinations++ = combination;
176  } while (std::next_permutation(booleans.begin(), booleans.end()));
177 }
178 
179 } // namespace internal
180 } // namespace tangential_complex
181 } // namespace Gudhi
182 
183 #endif // TANGENTIAL_COMPLEX_UTILITIES_H_
GUDHI  Version 3.4.0  - C++ library for Topological Data Analysis (TDA) and Higher Dimensional Geometry Understanding.  - Copyright : MIT Generated on Tue Dec 29 2020 10:00:12 for GUDHI by Doxygen 1.8.20