1#ifndef DUNE_FEM_GRIDPART_GEOMETRYGRIDPART_INTERSECTION_HH
2#define DUNE_FEM_GRIDPART_GEOMETRYGRIDPART_INTERSECTION_HH
6#include <dune/common/version.hh>
9#include <dune/geometry/affinegeometry.hh>
21 template<
class Gr
idFamily >
25 typedef typename std::remove_const< GridFamily >::type::Traits Traits;
28 typedef typename std::remove_const< GridFamily >::type::ctype
ctype;
30 static const int dimension = std::remove_const< GridFamily >::type::dimension;
31 static const int dimensionworld = std::remove_const< GridFamily >::type::dimensionworld;
33 typedef typename Traits::template Codim< 0 >::Entity
Entity;
35 typedef typename Traits::template Codim< 1 >::Geometry
Geometry;
36 typedef typename Traits::template Codim< 1 >::LocalGeometry
LocalGeometry;
44 typedef typename Geometry::Implementation GeometryImplType;
46 typedef typename Traits::HostGridPartType HostGridPartType;
47 typedef typename HostGridPartType::IntersectionType HostIntersectionType;
55 operator bool ()
const {
return bool( hostIntersection_ ); }
81 typedef typename Geometry::Implementation Impl;
94 const auto &refElement = ReferenceElements< ctype, dimension >::general( insideGeo_.type() );
95 const auto &refNormal = refElement.integrationOuterNormal(
indexInInside() );
97 const auto jit = insideGeo_.jacobianInverseTransposed(
geometryInInside().global( local ) );
100 jit.mv( refNormal, normal );
103 return normal *=
geometry().integrationElement( local ) / normal.two_norm();
108FieldVector< ctype, dimensionworld > normal, tau, nu;
109FieldMatrix< ctype, 1, dimensionworld > tauT =
geometry().jacobianTransposed(local);
114tau /= tau.two_norm();
116FieldMatrix< ctype, dimensionworld, dimensionworld > localFnGrad, localFnGradT;
121 localFnGradT[i][j] = localFnGrad[j][i];
123crossProduct(localFnGradT[0], localFnGradT[1], nu );
131crossProduct(nu, tau, normal);
139normal *=
geometry().integrationElement(local)/normal.two_norm();
146 void crossProduct(
const FieldVector< ctype, dimensionworld > &vec1,
const FieldVector< ctype, dimensionworld > &vec2, FieldVector< ctype, dimensionworld > &ret)
const
150 ret[0] = vec1[1]*vec2[2] - vec1[2]*vec2[1];
151 ret[1] = vec1[2]*vec2[0] - vec1[0]*vec2[2];
152 ret[2] = vec1[0]*vec2[1] - vec1[1]*vec2[0];
158 const auto &refElement = Dune::ReferenceElements< ctype, dimension >::general( insideGeo_.type() );
159 const auto &refNormal = refElement.integrationOuterNormal(
indexInInside() );
162 insideGeo_.jacobianInverseTransposed(
geometryInInside().global( local ) ).mv( refNormal, normal );
169 return normal *= (
ctype( 1 ) / normal.two_norm());
174 const auto &refElement = Dune::ReferenceElements< ctype, dimension >::general( insideGeo_.type() );
175 const auto &refNormal = refElement.integrationOuterNormal(
indexInInside() );
178 insideGeo_.jacobianInverseTransposed(
geometryInInside().center() ).mv( refNormal, normal );
179 return normal *= (
ctype( 1 ) / normal.two_norm());
186 assert( gridFunction_ );
187 return *gridFunction_;
191 HostIntersectionType hostIntersection_;
193 typename ElementGeometry::Implementation insideGeo_;
Definition: bindguard.hh:11
Definition: geometrygridpart/intersection.hh:23
Traits::template Codim< 1 >::LocalGeometry LocalGeometry
Definition: geometrygridpart/intersection.hh:36
Geometry geometry() const
Definition: geometrygridpart/intersection.hh:79
GeometryGridPartIntersection(const GridFunctionType &gridFunction, const typename ElementGeometry::Implementation &insideGeo, HostIntersectionType hostIntersection)
Definition: geometrygridpart/intersection.hh:51
int indexInOutside() const
Definition: geometrygridpart/intersection.hh:90
bool boundary() const
Definition: geometrygridpart/intersection.hh:64
FieldVector< ctype, dimensionworld > GlobalCoordinate
Definition: geometrygridpart/intersection.hh:40
FieldVector< ctype, dimension-1 > LocalCoordinate
Definition: geometrygridpart/intersection.hh:41
int twistInNeighbor() const
Definition: geometrygridpart/intersection.hh:70
static const int dimension
Definition: geometrygridpart/intersection.hh:30
Traits::template Codim< 0 >::Geometry ElementGeometry
Definition: geometrygridpart/intersection.hh:34
LocalGeometry geometryInOutside() const
Definition: geometrygridpart/intersection.hh:77
Traits::template Codim< 0 >::Entity Entity
Definition: geometrygridpart/intersection.hh:33
int indexInInside() const
Definition: geometrygridpart/intersection.hh:89
GlobalCoordinate integrationOuterNormal(const LocalCoordinate &local) const
Definition: geometrygridpart/intersection.hh:92
static const int dimensionworld
Definition: geometrygridpart/intersection.hh:31
GlobalCoordinate centerUnitOuterNormal() const
Definition: geometrygridpart/intersection.hh:172
bool equals(const This &other) const
Definition: geometrygridpart/intersection.hh:85
LocalGeometry geometryInInside() const
Definition: geometrygridpart/intersection.hh:76
bool conforming() const
Definition: geometrygridpart/intersection.hh:66
const GridFunctionType & gridFunction() const
Definition: geometrygridpart/intersection.hh:184
std::size_t boundarySegmentIndex() const
Definition: geometrygridpart/intersection.hh:74
std::remove_const< GridFamily >::type::ctype ctype
Definition: geometrygridpart/intersection.hh:28
Entity inside() const
Definition: geometrygridpart/intersection.hh:57
GlobalCoordinate unitOuterNormal(const LocalCoordinate &local) const
Definition: geometrygridpart/intersection.hh:166
GeometryType type() const
Definition: geometrygridpart/intersection.hh:87
Traits::template Codim< 1 >::Geometry Geometry
Definition: geometrygridpart/intersection.hh:35
GlobalCoordinate outerNormal(const LocalCoordinate &local) const
Definition: geometrygridpart/intersection.hh:156
Entity outside() const
Definition: geometrygridpart/intersection.hh:58
Traits::GridFunctionType GridFunctionType
Definition: geometrygridpart/intersection.hh:38
const HostIntersectionType & hostIntersection() const
Definition: geometrygridpart/intersection.hh:182
GeometryGridPartIntersection()=default
bool neighbor() const
Definition: geometrygridpart/intersection.hh:72
int twistInSelf() const
Definition: geometrygridpart/intersection.hh:68
int boundaryId() const
Definition: geometrygridpart/intersection.hh:60
Definition: boundaryidprovider.hh:36