1 #ifndef DUNE_FEM_GRIDPART_GEOGRIDPART_INTERSECTION_HH 2 #define DUNE_FEM_GRIDPART_GEOGRIDPART_INTERSECTION_HH 19 template<
class Gr
idFamily >
22 typedef typename std::remove_const< GridFamily >::type::Traits Traits;
25 typedef typename std::remove_const< GridFamily >::type::ctype
ctype;
27 static const int dimension = std::remove_const< GridFamily >::type::dimension;
28 static const int dimensionworld = std::remove_const< GridFamily >::type::dimensionworld;
30 typedef typename Traits::template Codim< 0 >::Entity
Entity;
32 typedef typename Traits::template Codim< 1 >::Geometry
Geometry;
33 typedef typename Traits::template Codim< 1 >::LocalGeometry
LocalGeometry;
38 typedef typename Entity::Implementation EntityImplType;
39 typedef typename ElementGeometry::Implementation ElementGeometryImplType;
40 typedef typename Geometry::Implementation GeometryImplType;
42 typedef typename Traits::HostGridPartType HostGridPartType;
43 typedef typename HostGridPartType::IntersectionType HostIntersectionType;
49 : coordFunction_( &coordFunction ),
50 insideGeo_( insideGeo.impl() ),
51 hostIntersection_(
std::
move( hostIntersection ) )
102 CoordVectorType coords( insideGeo_, localGeo );
121 FieldVector< ctype, dimensionworld >
124 const ReferenceElement< ctype, dimension > &refElement
125 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
128 typedef typename ElementGeometryImplType::JacobianInverseTransposed JacobianInverseTransposed;
129 const JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
130 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal(
indexInInside() );
132 FieldVector< ctype, dimensionworld > normal;
133 jit.mv( refNormal, normal );
134 normal *=
ctype( 1 ) / jit.det();
138 FieldVector< ctype, dimensionworld >
139 outerNormal (
const FieldVector< ctype, dimension-1 > &local )
const 141 const ReferenceElement< ctype, dimension > &refElement
142 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
145 typedef typename ElementGeometryImplType::JacobianInverseTransposed JacobianInverseTransposed;
146 const JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
147 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal(
indexInInside() );
149 FieldVector< ctype, dimensionworld > normal;
150 jit.mv( refNormal, normal );
154 FieldVector< ctype, dimensionworld >
157 FieldVector< ctype, dimensionworld > normal =
outerNormal( local );
158 normal *= (
ctype( 1 ) / normal.two_norm());
164 const ReferenceElement<
ctype, dimension-1 > &refFace
165 = ReferenceElements< ctype, dimension-1 >::general(
type() );
171 assert( coordFunction_ );
172 return *coordFunction_;
177 return hostIntersection_;
181 const CoordFunctionType *coordFunction_ =
nullptr;
182 ElementGeometryImplType insideGeo_;
183 HostIntersectionType hostIntersection_;
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
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