3#ifndef DUNE_GEOGRID_INTERSECTION_HH
4#define DUNE_GEOGRID_INTERSECTION_HH
6#include <dune/grid/geometrygrid/declaration.hh>
7#include <dune/grid/geometrygrid/entitypointer.hh>
8#include <dune/grid/geometrygrid/cornerstorage.hh>
19 template<
class Gr
id,
class HostIntersection >
22 typedef typename HostIntersection::Geometry HostGeometry;
23 typedef typename HostIntersection::LocalGeometry HostLocalGeometry;
25 typedef typename remove_const< Grid >::type::Traits Traits;
28 typedef typename Traits::ctype ctype;
30 static const int dimension = Traits::dimension;
31 static const int dimensionworld = Traits::dimensionworld;
33 typedef typename Traits::template Codim< 0 >::Entity Entity;
34 typedef typename Traits::template Codim< 0 >::EntityPointer EntityPointer;
35 typedef typename Traits::template Codim< 1 >::Geometry Geometry;
36 typedef typename Traits::template Codim< 1 >::LocalGeometry LocalGeometry;
38 typedef typename Traits::template Codim< 0 >::Geometry ElementGeometry;
41 typedef GeoGrid::IntersectionCoordVector< Grid > CoordVector;
43 typedef typename Traits::template Codim< 0 >::EntityPointerImpl EntityPointerImpl;
45 typedef typename Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
46 typedef typename Traits::template Codim< 0 >::GeometryImpl ElementGeometryImpl;
49 explicit Intersection (
const ElementGeometry &insideGeo )
50 : insideGeo_( Grid::getRealImplementation( insideGeo ) ),
51 hostIntersection_( 0 ),
55 Intersection (
const Intersection &other )
56 : insideGeo_( other.insideGeo_ ),
57 hostIntersection_( 0 ),
61 const Intersection &operator= (
const Intersection &other )
63 insideGeo_ = other.insideGeo_;
68 operator bool ()
const {
return bool( hostIntersection_ ); }
70 EntityPointer inside ()
const
72 return EntityPointerImpl( insideGeo_, hostIntersection().inside() );
75 EntityPointer outside ()
const
77 return EntityPointerImpl( grid(), hostIntersection().outside() );
80 bool boundary ()
const {
return hostIntersection().boundary(); }
82 bool conforming ()
const {
return hostIntersection().conforming(); }
84 bool neighbor ()
const {
return hostIntersection().neighbor(); }
86 int boundaryId ()
const {
return hostIntersection().boundaryId(); }
88 size_t boundarySegmentIndex ()
const
90 return hostIntersection().boundarySegmentIndex();
93 LocalGeometry geometryInInside ()
const
95 return hostIntersection().geometryInInside();
98 LocalGeometry geometryInOutside ()
const
100 return hostIntersection().geometryInOutside();
103 Geometry geometry ()
const
107 CoordVector coords( insideGeo_, geometryInInside() );
108 geo_ = GeometryImpl( grid(), type(), coords );
110 return Geometry( geo_ );
113 GeometryType type ()
const {
return hostIntersection().type(); }
115 int indexInInside ()
const
117 return hostIntersection().indexInInside();
120 int indexInOutside ()
const
122 return hostIntersection().indexInOutside();
125 FieldVector< ctype, dimensionworld >
126 integrationOuterNormal (
const FieldVector< ctype, dimension-1 > &local )
const
128 const ReferenceElement< ctype, dimension > &refElement
132 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
135 FieldVector< ctype, dimensionworld > normal;
136 jit.mv( refNormal, normal );
137 normal *= jit.detInv();
142 FieldVector< ctype, dimensionworld >
143 outerNormal (
const FieldVector< ctype, dimension-1 > &local )
const
145 const ReferenceElement< ctype, dimension > &refElement
149 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
152 FieldVector< ctype, dimensionworld > normal;
153 jit.mv( refNormal, normal );
157 FieldVector< ctype, dimensionworld >
158 unitOuterNormal (
const FieldVector< ctype, dimension-1 > &local )
const
160 FieldVector< ctype, dimensionworld > normal = outerNormal( local );
161 normal *= (ctype( 1 ) / normal.two_norm());
165 FieldVector< ctype, dimensionworld > centerUnitOuterNormal ()
const
167 const ReferenceElement< ctype, dimension-1 > &refFace
169 return unitOuterNormal( refFace.position( 0, 0 ) );
172 const HostIntersection &hostIntersection ()
const
175 return *hostIntersection_;
178 const Grid &grid ()
const {
return insideGeo_.grid(); }
182 hostIntersection_ = 0;
183 geo_ = GeometryImpl( grid() );
186 void initialize (
const HostIntersection &hostIntersection )
189 hostIntersection_ = &hostIntersection;
193 ElementGeometryImpl insideGeo_;
194 const HostIntersection *hostIntersection_;
195 mutable GeometryImpl geo_;
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
Dune namespace.
Definition: alignment.hh:14
static const ReferenceElement< ctype, dim > & general(const GeometryType &type)
get general reference elements
Definition: referenceelements.hh:568