dune-grid-glue 2.8.0
Loading...
Searching...
No Matches
projection.hh
Go to the documentation of this file.
1#ifndef DUNE_GRIDGLUE_COMMON_PROJECTIONHELPER2_HH
2#define DUNE_GRIDGLUE_COMMON_PROJECTIONHELPER2_HH
3
4#include <array>
5#include <bitset>
6#include <tuple>
7
8namespace Dune {
9namespace GridGlue {
10
17template<typename Coordinate>
19{
20public:
27 {
31 std::array<unsigned, 2> edge;
32
39 std::array<Coordinate, 2> local;
40 };
41
45 constexpr static unsigned dim = Coordinate::dimension;
46
52 constexpr static unsigned maxEdgeIntersections = dim == 3 ? 9 : 0;
53
54 static_assert(dim == 2 || dim == 3, "Projection only implemented for dim=2 or dim=3");
55
59 typedef typename Coordinate::field_type Field;
60
68 typedef std::array<Coordinate, dim> Images;
69
77
78private:
82 const Field m_overlap;
83
92 const Field m_max_normal_product;
93
99 Field m_epsilon = Field(1e-12);
100
102 std::tuple<Images, Preimages> m_images;
103
105 std::tuple<std::bitset<dim>, std::bitset<dim> > m_success;
106
108 unsigned m_number_of_edge_intersections;
109
111 std::array<EdgeIntersection, maxEdgeIntersections> m_edge_intersections;
112
124 bool m_projection_valid;
125
131 template<typename Corners, typename Normals>
132 void doProjection(const std::tuple<Corners&, Corners&>& corners, const std::tuple<Normals&, Normals&>& normals);
133
142 template<typename Corners, typename Normals>
143 void doInverseProjection(const std::tuple<Corners&, Corners&>& corners, const std::tuple<Normals&, Normals&>& normals);
144
153 template<typename Corners, typename Normals>
154 void doEdgeIntersection(const std::tuple<Corners&, Corners&>& corners, const std::tuple<Normals&, Normals&>& normals);
155
181 template<typename Corners, typename Normals>
182 inline bool projectionFeasible(const Coordinate& x, const Coordinate& nx, const Coordinate& px, const Corners& corners, const Normals& normals) const;
183
184public:
189 Projection(const Field overlap = Field(0), const Field max_normal_product = Field(-0.1));
190
196 void epsilon(const Field epsilon);
197
210 template<typename Corners, typename Normals>
211 void project(const std::tuple<Corners&, Corners&>& corners, const std::tuple<Normals&, Normals&>& normals);
212
233 const std::tuple<Images, Preimages>& images() const
234 { return m_images; }
235
250 const std::tuple<std::bitset<dim>, std::bitset<dim> >& success() const
251 { return m_success; }
252
261 { return m_number_of_edge_intersections; }
262
271 const std::array<EdgeIntersection, maxEdgeIntersections>& edgeIntersections() const
272 { return m_edge_intersections; }
273};
274
275} /* namespace GridGlue */
276} /* namespace Dune */
277
278#include "projection_impl.hh"
279
280#endif
Definition: gridglue.hh:35
Projection of a line (triangle) on another line (triangle).
Definition: projection.hh:19
Coordinate::field_type Field
Scalar type.
Definition: projection.hh:59
const std::tuple< std::bitset< dim >, std::bitset< dim > > & success() const
Indicate whether projection (inverse projection) is valid for each corner or not.
Definition: projection.hh:250
static constexpr unsigned maxEdgeIntersections
maximum number of edge-edge intersections
Definition: projection.hh:52
Images Preimages
Definition: projection.hh:76
std::array< Coordinate, dim > Images
List of corner images.
Definition: projection.hh:68
void epsilon(const Field epsilon)
Set epsilon used for floating-point comparisons.
Definition: projection_impl.hh:138
void project(const std::tuple< Corners &, Corners & > &corners, const std::tuple< Normals &, Normals & > &normals)
Do the actual projection.
Definition: projection_impl.hh:467
static constexpr unsigned dim
dimension of coordinates
Definition: projection.hh:45
unsigned numberOfEdgeIntersections() const
Number of edge intersections.
Definition: projection.hh:260
const std::tuple< Images, Preimages > & images() const
Images and preimages of corners.
Definition: projection.hh:233
const std::array< EdgeIntersection, maxEdgeIntersections > & edgeIntersections() const
Edge-edge intersections.
Definition: projection.hh:271
Intersection between two edges of a triangle.
Definition: projection.hh:27
std::array< Coordinate, 2 > local
Local coordinates of intersection and distance along normals.
Definition: projection.hh:39
std::array< unsigned, 2 > edge
Edge numbers in image and preimage triangle.
Definition: projection.hh:31