preCICE v3.1.2
Loading...
Searching...
No Matches
RTreeAdapter.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <Eigen/Core>
4#include <boost/geometry.hpp>
5#include "mesh/Edge.hpp"
6#include "mesh/Mesh.hpp"
8#include "mesh/Vertex.hpp"
9#include "utils/assertion.hpp"
10
11namespace precice {
12namespace mesh {
13class Triangle;
14} // namespace mesh
15} // namespace precice
16
17namespace pm = precice::mesh;
18namespace bg = boost::geometry;
19
20namespace boost {
21namespace geometry {
22namespace traits {
23
25/*
26 * This adapts every VectorXd to a 3d point. For non-existing dimensions, zero is returned.
27 */
28template <>
29struct tag<Eigen::VectorXd> {
30 using type = point_tag;
31};
32template <>
33struct coordinate_type<Eigen::VectorXd> {
34 using type = double;
35};
36template <>
37struct coordinate_system<Eigen::VectorXd> {
38 using type = cs::cartesian;
39};
40template <>
41struct dimension<Eigen::VectorXd> : boost::mpl::int_<3> {
42};
43
44template <size_t Dimension>
45struct access<Eigen::VectorXd, Dimension> {
46 static double get(Eigen::VectorXd const &p)
47 {
48 if (Dimension >= static_cast<size_t>(p.rows())) {
49 return 0;
50 } else {
51 return p(Dimension);
52 }
53 }
54
55 static void set(Eigen::VectorXd &p, double const &value)
56 {
57 if (Dimension >= static_cast<size_t>(p.rows())) {
58 Eigen::VectorXd tmp(Dimension);
59 std::copy_n(p.data(), Dimension - 1, tmp.data());
60 p = std::move(tmp);
61 }
62 p(Dimension) = value;
63 }
64};
65
67template <>
68struct tag<pm::Vertex::RawCoords> {
69 using type = point_tag;
70};
71template <>
72struct coordinate_type<pm::Vertex::RawCoords> {
73 using type = double;
74};
75template <>
76struct coordinate_system<pm::Vertex::RawCoords> {
77 using type = cs::cartesian;
78};
79template <>
80struct dimension<pm::Vertex::RawCoords> : boost::mpl::int_<3> {
81};
82
83template <size_t Dimension>
84struct access<pm::Vertex::RawCoords, Dimension> {
85 static double get(pm::Vertex::RawCoords const &p)
86 {
87 return p[Dimension];
88 }
89
90 static void set(pm::Vertex::RawCoords &p, double const &value)
91 {
92 p[Dimension] = value;
93 }
94};
95
96BOOST_CONCEPT_ASSERT((bg::concepts::Point<pm::Vertex::RawCoords>) );
97
99/*
100 * This adapts every Vertex to a 3d point. For non-existing dimensions, zero is returned.
101 */
102template <>
103struct tag<pm::Vertex> {
104 using type = point_tag;
105};
106template <>
107struct coordinate_type<pm::Vertex> {
108 using type = double;
109};
110template <>
111struct coordinate_system<pm::Vertex> {
112 using type = cs::cartesian;
113};
114template <>
115struct dimension<pm::Vertex> : boost::mpl::int_<3> {
116};
117
118template <size_t Dimension>
119struct access<pm::Vertex, Dimension> {
120 static double get(pm::Vertex const &p)
121 {
122 return p.rawCoords()[Dimension];
123 }
124
125 static void set(pm::Vertex &p, double const &value)
126 {
127 PRECICE_UNREACHABLE("Boost.Geometry is not allowed to write to mesh::Vertex.");
128 }
129};
130
131BOOST_CONCEPT_ASSERT((concepts::Point<pm::Vertex>) );
132
138template <>
139struct tag<pm::Edge> {
140 using type = segment_tag;
141};
142template <>
143struct point_type<pm::Edge> {
144 using type = pm::Vertex::RawCoords;
145};
146
147template <size_t Index, size_t Dimension>
148struct indexed_access<pm::Edge, Index, Dimension> {
149 static_assert((Index <= 1), "Valid Indices are {0, 1}");
150 static_assert((Dimension <= 2), "Valid Dimensions are {0, 1, 2}");
151
152 static double get(pm::Edge const &e)
153 {
154 return access<point_type<pm::Edge>::type, Dimension>::get(e.vertex(Index).rawCoords());
155 }
156
157 static void set(pm::Edge &e, double const &value)
158 {
159 PRECICE_UNREACHABLE("Boost.Geometry is not allowed to write to mesh::Edge.");
160 }
161};
162
168template <>
169struct tag<pm::Triangle> {
170 using type = ring_tag;
171};
172template <>
173struct point_order<pm::Triangle> {
174 static const order_selector value = clockwise;
175};
176template <>
177struct closure<pm::Triangle> {
178 static const closure_selector value = open;
179};
180
181// BOOST_CONCEPT_ASSERT( (bg::concepts::Ring<pm::Triangle>));
182
183} // namespace traits
184} // namespace geometry
185} // namespace boost
186
187namespace precice {
188namespace query {
189
191using RTreeBox = boost::geometry::model::box<pm::Vertex::RawCoords>;
192
193inline RTreeBox makeBox(const pm::Vertex::RawCoords &min, const pm::Vertex::RawCoords &max)
194{
195 return {min, max};
196}
197
198inline pm::Vertex::RawCoords eigenToRaw(const Eigen::VectorXd &v)
199{
200 const auto size = v.size();
201 PRECICE_ASSERT(size == 2 || size == 3, size);
202 pm::Vertex::RawCoords r{0.0, 0.0, 0.0};
203 std::copy_n(v.data(), size, r.data());
204 return r;
205}
206
207inline Eigen::VectorXd rawToEigen(const pm::Vertex::RawCoords &v)
208{
209 Eigen::VectorXd r(3);
210 std::copy(v.begin(), v.end(), r.data());
211 return r;
212}
213
214inline RTreeBox makeBox(const Eigen::VectorXd &min, const Eigen::VectorXd &max)
215{
216 return {eigenToRaw(min), eigenToRaw(max)};
217}
218
219// Overload for a tetrahedron
221{
222
224 for (int i = 0; i < 4; ++i) {
225 box.expandBy(tetra.vertex(i));
226 }
227
228 // Convert to Boost type
229 return makeBox(box.minCorner(), box.maxCorner());
230}
231
232namespace impl {
233
235using RTreeParameters = boost::geometry::index::rstar<16>;
236
238template <class T>
240
241template <>
245
246template <>
250
251template <>
255
256template <>
257struct PrimitiveTraits<mesh::Tetrahedron> {
259};
260
262template <typename Container>
264 using size_type = typename Container::container::size_type;
265 using cref = const typename Container::value_type &;
266 Container const &container;
267
268public:
270
271 explicit PtrVectorIndexable(Container const &c)
272 : container(c)
273 {
274 }
275
277 {
278 return container[i];
279 }
280};
281
283template <typename Container>
285 using size_type = typename Container::size_type;
286 using cref = const typename Container::value_type &;
287 Container const &container;
288
289public:
291
292 explicit VectorIndexable(Container const &c)
293 : container(c)
294 {
295 }
296
298 {
299 return container[i];
300 }
301};
302
303template <typename Primitive>
305private:
306 template <typename T, typename = typename std::enable_if<
308 typename boost::geometry::traits::tag<T>::type,
309 boost::geometry::point_tag>::value,
311 static std::true_type test(char *);
312 template <typename T, typename = typename std::enable_if<
314 typename boost::geometry::traits::tag<T>::type,
315 boost::geometry::segment_tag>::value,
317 static std::true_type test(int *);
318 template <typename T, typename = typename std::enable_if<
320 typename boost::geometry::traits::tag<T>::type,
321 boost::geometry::box_tag>::value,
323 static std::true_type test(void *);
324
325 template <typename T>
327
328public:
329 using type = decltype(test<Primitive>(nullptr));
330};
331
332template <class Primitive>
334};
335
337template <class Primitive>
340 using MeshContainerIndex = typename MeshContainer::size_type;
341
342 using IndexType = typename std::conditional<
346
350 boost::geometry::index::indexable<IndexType>>::type;
351
352 using RTree = boost::geometry::index::rtree<IndexType, RTreeParameters, IndexGetter>;
354};
355
356} // namespace impl
357} // namespace query
358} // namespace precice
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
#define PRECICE_UNREACHABLE(...)
Definition assertion.hpp:95
An axis-aligned bounding box around a (partition of a) mesh.
Eigen::VectorXd maxCorner() const
the max corner of the bounding box
Eigen::VectorXd minCorner() const
the min corner of the bounding box
void expandBy(const BoundingBox &otherBB)
Expand bounding box using another bounding box.
Linear edge of a mesh, defined by two Vertex objects.
Definition Edge.hpp:16
Vertex & vertex(int i)
Returns the edge's vertex with index 0 or 1.
Definition Edge.hpp:77
std::deque< Triangle > TriangleContainer
Definition Mesh.hpp:43
std::deque< Tetrahedron > TetraContainer
Definition Mesh.hpp:44
std::deque< Edge > EdgeContainer
Definition Mesh.hpp:42
std::deque< Vertex > VertexContainer
Definition Mesh.hpp:41
Tetrahedron of a mesh, defined by 4 vertices.
Vertex & vertex(int i)
Returns tetrahedron vertex with index 0, 1, 2 or 3.
int getDimensions() const
Returns dimensionalty of space the Tetrahedron is embedded in.
Vertex of a mesh.
Definition Vertex.hpp:16
const RawCoords & rawCoords() const
Direct access to the coordinates.
Definition Vertex.hpp:123
static std::true_type test(void *)
static std::true_type test(char *)
decltype(test< Primitive >(nullptr)) type
Makes a utils::PtrVector indexable and thus be usable in boost::geometry::rtree.
typename Container::container::size_type size_type
const typename Container::value_type & cref
result_type operator()(size_type i) const
Makes a std::vector indexable and thus be usable in boost::geometry::rtree.
const typename Container::value_type & cref
result_type operator()(size_type i) const
typename Container::size_type size_type
T copy(T... args)
T copy_n(T... args)
BOOST_CONCEPT_ASSERT((bg::concepts::Point< pm::Vertex::RawCoords >))
provides Mesh, Data and primitives.
boost::geometry::index::rstar< 16 > RTreeParameters
The general rtree parameter type used in precice.
RTreeBox makeBox(const pm::Vertex::RawCoords &min, const pm::Vertex::RawCoords &max)
pm::Vertex::RawCoords eigenToRaw(const Eigen::VectorXd &v)
Eigen::VectorXd rawToEigen(const pm::Vertex::RawCoords &v)
boost::geometry::model::box< pm::Vertex::RawCoords > RTreeBox
The RTree box type.
Main namespace of the precice library.
static std::unique_ptr< precice::Participant > impl
Definition preciceC.cpp:21
static void set(Eigen::VectorXd &p, double const &value)
static void set(pm::Vertex &p, double const &value)
static void set(pm::Vertex::RawCoords &p, double const &value)
Type trait to extract information based on the type of a Primitive.
The type traits of a rtree based on a Primitive.
typename PrimitiveTraits< Primitive >::MeshContainer MeshContainer
typename std::conditional< IsDirectIndexable< Primitive >::value, MeshContainerIndex, std::pair< RTreeBox, MeshContainerIndex > >::type IndexType
typename MeshContainer::size_type MeshContainerIndex
typename std::conditional< IsDirectIndexable< Primitive >::value, impl::VectorIndexable< MeshContainer >, boost::geometry::index::indexable< IndexType > >::type IndexGetter
boost::geometry::index::rtree< IndexType, RTreeParameters, IndexGetter > RTree