3#ifndef DUNE_GEOGRID_INTERSECTION_HH
4#define DUNE_GEOGRID_INTERSECTION_HH
18 template<
class Gr
id,
class HostIntersection >
21 typedef typename HostIntersection::Geometry HostGeometry;
22 typedef typename HostIntersection::LocalGeometry HostLocalGeometry;
24 typedef typename std::remove_const< Grid >::type::Traits Traits;
27 typedef typename Traits::ctype
ctype;
32 typedef typename Traits::template Codim< 0 >::Entity
Entity;
33 typedef typename Traits::template Codim< 1 >::Geometry
Geometry;
34 typedef typename Traits::template Codim< 1 >::LocalGeometry
LocalGeometry;
41 typedef typename Traits::template Codim< 0 >::EntityImpl EntityImpl;
43 typedef typename Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
44 typedef typename Traits::template Codim< 0 >::GeometryImpl ElementGeometryImpl;
53 , insideGeo_ ( insideGeo )
59 , insideGeo_ ( insideGeo )
65 return hostIntersection_ == other.hostIntersection_;
68 operator bool ()
const {
return bool( hostIntersection_ ); }
108 geo_ = GeometryImpl(
grid(),
type(), coords );
125 FieldVector< ctype, dimensionworld >
131 const ReferenceElement< ctype, dimension > &refElement
132 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
134 FieldVector< ctype, dimension > x( geoInInside.global( local ) );
135 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
136 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( idxInInside );
138 FieldVector< ctype, dimensionworld > normal;
139 jit.mv( refNormal, normal );
141 normal *= geoInInside.volume() / refElement.template geometry< 1 >( idxInInside ).volume();
142 normal *= jit.detInv();
147 FieldVector< ctype, dimensionworld >
148 outerNormal (
const FieldVector< ctype, dimension-1 > &local )
const
150 const ReferenceElement< ctype, dimension > &refElement
151 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
154 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
155 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal(
indexInInside() );
157 FieldVector< ctype, dimensionworld > normal;
158 jit.mv( refNormal, normal );
162 FieldVector< ctype, dimensionworld >
165 FieldVector< ctype, dimensionworld > normal =
outerNormal( local );
166 normal *= (
ctype( 1 ) / normal.two_norm());
173 = ReferenceElements< ctype, dimension-1 >::general(
type() );
179 return hostIntersection_;
182 const Grid &
grid ()
const {
return insideGeo_.grid(); }
185 HostIntersection hostIntersection_;
186 ElementGeometryImpl insideGeo_;
187 mutable GeometryImpl geo_;
Include standard header files.
Definition agrid.hh:60
Grid abstract base class.
Definition common/grid.hh:373
Definition cornerstorage.hh:122
Definition geometrygrid/intersection.hh:20
FieldVector< ctype, dimensionworld > integrationOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition geometrygrid/intersection.hh:126
Traits::ctype ctype
Definition geometrygrid/intersection.hh:27
GeometryType type() const
Definition geometrygrid/intersection.hh:113
LocalGeometry geometryInInside() const
Definition geometrygrid/intersection.hh:93
Geometry geometry() const
Definition geometrygrid/intersection.hh:103
bool boundary() const
Definition geometrygrid/intersection.hh:80
Traits::template Codim< 1 >::Geometry Geometry
Definition geometrygrid/intersection.hh:33
Intersection()
Definition geometrygrid/intersection.hh:48
bool neighbor() const
Definition geometrygrid/intersection.hh:84
size_t boundarySegmentIndex() const
Definition geometrygrid/intersection.hh:88
FieldVector< ctype, dimensionworld > outerNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition geometrygrid/intersection.hh:148
FieldVector< ctype, dimensionworld > unitOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition geometrygrid/intersection.hh:163
int boundaryId() const
Definition geometrygrid/intersection.hh:86
Intersection(HostIntersection &&hostIntersection, const ElementGeometryImpl &insideGeo)
Definition geometrygrid/intersection.hh:57
Traits::template Codim< 0 >::Geometry ElementGeometry
Definition geometrygrid/intersection.hh:36
Intersection(const HostIntersection &hostIntersection, const ElementGeometryImpl &insideGeo)
Definition geometrygrid/intersection.hh:51
Traits::template Codim< 1 >::LocalGeometry LocalGeometry
Definition geometrygrid/intersection.hh:34
static const int dimensionworld
Definition geometrygrid/intersection.hh:30
const HostIntersection & hostIntersection() const
Definition geometrygrid/intersection.hh:177
int indexInOutside() const
Definition geometrygrid/intersection.hh:120
const Grid & grid() const
Definition geometrygrid/intersection.hh:182
bool equals(const Intersection &other) const
Definition geometrygrid/intersection.hh:63
Traits::template Codim< 0 >::Entity Entity
Definition geometrygrid/intersection.hh:32
bool conforming() const
Definition geometrygrid/intersection.hh:82
Entity outside() const
Definition geometrygrid/intersection.hh:75
int indexInInside() const
Definition geometrygrid/intersection.hh:115
FieldVector< ctype, dimensionworld > centerUnitOuterNormal() const
Definition geometrygrid/intersection.hh:170
static const int dimension
Definition geometrygrid/intersection.hh:29
Entity inside() const
Definition geometrygrid/intersection.hh:70
LocalGeometry geometryInOutside() const
Definition geometrygrid/intersection.hh:98