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