dune-fem 2.8.0
Loading...
Searching...
No Matches
geogridpart/intersection.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_GRIDPART_GEOGRIDPART_INTERSECTION_HH
2#define DUNE_FEM_GRIDPART_GEOGRIDPART_INTERSECTION_HH
3
4#include <type_traits>
5#include <utility>
6
8
9namespace Dune
10{
11
12 namespace Fem
13 {
14
15 // GeoIntersection
16 // --------------
17
18 template< class GridFamily >
20 {
21 typedef typename std::remove_const< GridFamily >::type::Traits Traits;
22
23 public:
24 typedef typename std::remove_const< GridFamily >::type::ctype ctype;
25
26 static const int dimension = std::remove_const< GridFamily >::type::dimension;
27 static const int dimensionworld = std::remove_const< GridFamily >::type::dimensionworld;
28
29 typedef typename Traits::template Codim< 0 >::Entity Entity;
30 typedef typename Traits::template Codim< 0 >::Geometry ElementGeometry;
31 typedef typename Traits::template Codim< 1 >::Geometry Geometry;
32 typedef typename Traits::template Codim< 1 >::LocalGeometry LocalGeometry;
33
34 typedef typename Traits::CoordFunctionType CoordFunctionType;
35
36 private:
37 typedef typename Entity::Implementation EntityImplType;
38 typedef typename ElementGeometry::Implementation ElementGeometryImplType;
39 typedef typename Geometry::Implementation GeometryImplType;
40
41 typedef typename Traits::HostGridPartType HostGridPartType;
42 typedef typename HostGridPartType::IntersectionType HostIntersectionType;
43
45
46 public:
47 GeoIntersection ( const CoordFunctionType &coordFunction, const ElementGeometry &insideGeo, HostIntersectionType hostIntersection )
48 : coordFunction_( &coordFunction ),
49 insideGeo_( insideGeo.impl() ),
50 hostIntersection_( std::move( hostIntersection ) )
51 {}
52
53 Entity inside () const
54 {
55 return EntityImplType( coordFunction(), hostIntersection().inside() );
56 }
57
58 Entity outside () const
59 {
60 return EntityImplType( coordFunction(), hostIntersection().outside() );
61 }
62
63 bool boundary () const
64 {
65 return hostIntersection().boundary();
66 }
67
68 bool conforming () const
69 {
70 return hostIntersection().conforming();
71 }
72
73 bool neighbor () const
74 {
75 return hostIntersection().neighbor();
76 }
77
78 int boundaryId () const
79 {
80 return hostIntersection().boundaryId();
81 }
82
83 size_t boundarySegmentIndex () const
84 {
85 return hostIntersection().boundarySegmentIndex();
86 }
87
89 {
90 return hostIntersection().geometryInInside();
91 }
92
94 {
95 return hostIntersection().geometryInOutside();
96 }
97
99 {
100 const LocalGeometry &localGeo = geometryInInside();
101 CoordVectorType coords( insideGeo_, localGeo );
102 return Geometry( GeometryImplType( type(), coords ) );
103 }
104
105 GeometryType type () const
106 {
107 return hostIntersection().type();
108 }
109
110 int indexInInside () const
111 {
112 return hostIntersection().indexInInside();
113 }
114
115 int indexInOutside () const
116 {
117 return hostIntersection().indexInOutside();
118 }
119
120 FieldVector< ctype, dimensionworld >
121 integrationOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
122 {
123 auto refElement = referenceElement< ctype, dimension>( insideGeo_.type() );
124
125 FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
126 typedef typename ElementGeometryImplType::JacobianInverseTransposed JacobianInverseTransposed;
127 const JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
128 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( indexInInside() );
129
130 FieldVector< ctype, dimensionworld > normal;
131 jit.mv( refNormal, normal );
132 normal *= ctype( 1 ) / jit.det();
133 return normal;
134 }
135
136 FieldVector< ctype, dimensionworld >
137 outerNormal ( const FieldVector< ctype, dimension-1 > &local ) const
138 {
139 auto refElement = referenceElement< ctype, dimension>( insideGeo_.type() );
140
141 FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
142 typedef typename ElementGeometryImplType::JacobianInverseTransposed JacobianInverseTransposed;
143 const JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
144 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( indexInInside() );
145
146 FieldVector< ctype, dimensionworld > normal;
147 jit.mv( refNormal, normal );
148 return normal;
149 }
150
151 FieldVector< ctype, dimensionworld >
152 unitOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
153 {
154 FieldVector< ctype, dimensionworld > normal = outerNormal( local );
155 normal *= (ctype( 1 ) / normal.two_norm());
156 return normal;
157 }
158
159 FieldVector< ctype, dimensionworld > centerUnitOuterNormal () const
160 {
161 auto refFace = referenceElement< ctype, dimension-1 >( type() );
162 return unitOuterNormal( refFace.position( 0, 0 ) );
163 }
164
166 {
167 assert( coordFunction_ );
168 return *coordFunction_;
169 }
170
171 const HostIntersectionType &hostIntersection () const
172 {
173 return hostIntersection_;
174 }
175
176 private:
177 const CoordFunctionType *coordFunction_ = nullptr;
178 ElementGeometryImplType insideGeo_;
179 HostIntersectionType hostIntersection_;
180 };
181
182 } // namespace Fem
183
184} // namespace Dune
185
186#endif // #ifndef DUNE_FEM_GRIDPART_GEOGRIDPART_INTERSECTION_HH
STL namespace.
Definition: bindguard.hh:11
Definition: cornerstorage.hh:234
Definition: geogridpart/intersection.hh:20
FieldVector< ctype, dimensionworld > unitOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geogridpart/intersection.hh:152
bool neighbor() const
Definition: geogridpart/intersection.hh:73
static const int dimensionworld
Definition: geogridpart/intersection.hh:27
std::remove_const< GridFamily >::type::ctype ctype
Definition: geogridpart/intersection.hh:24
bool conforming() const
Definition: geogridpart/intersection.hh:68
LocalGeometry geometryInInside() const
Definition: geogridpart/intersection.hh:88
Entity outside() const
Definition: geogridpart/intersection.hh:58
int indexInInside() const
Definition: geogridpart/intersection.hh:110
int boundaryId() const
Definition: geogridpart/intersection.hh:78
GeoIntersection(const CoordFunctionType &coordFunction, const ElementGeometry &insideGeo, HostIntersectionType hostIntersection)
Definition: geogridpart/intersection.hh:47
FieldVector< ctype, dimensionworld > centerUnitOuterNormal() const
Definition: geogridpart/intersection.hh:159
FieldVector< ctype, dimensionworld > outerNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geogridpart/intersection.hh:137
const CoordFunctionType & coordFunction() const
Definition: geogridpart/intersection.hh:165
LocalGeometry geometryInOutside() const
Definition: geogridpart/intersection.hh:93
Traits::template Codim< 0 >::Entity Entity
Definition: geogridpart/intersection.hh:29
bool boundary() const
Definition: geogridpart/intersection.hh:63
Traits::template Codim< 1 >::LocalGeometry LocalGeometry
Definition: geogridpart/intersection.hh:32
Geometry geometry() const
Definition: geogridpart/intersection.hh:98
const HostIntersectionType & hostIntersection() const
Definition: geogridpart/intersection.hh:171
Traits::CoordFunctionType CoordFunctionType
Definition: geogridpart/intersection.hh:34
GeometryType type() const
Definition: geogridpart/intersection.hh:105
Entity inside() const
Definition: geogridpart/intersection.hh:53
int indexInOutside() const
Definition: geogridpart/intersection.hh:115
FieldVector< ctype, dimensionworld > integrationOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geogridpart/intersection.hh:121
Traits::template Codim< 0 >::Geometry ElementGeometry
Definition: geogridpart/intersection.hh:30
static const int dimension
Definition: geogridpart/intersection.hh:26
Traits::template Codim< 1 >::Geometry Geometry
Definition: geogridpart/intersection.hh:31
size_t boundarySegmentIndex() const
Definition: geogridpart/intersection.hh:83