dune-grid 3.0-git
geometrygrid/intersection.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_GEOGRID_INTERSECTION_HH
4#define DUNE_GEOGRID_INTERSECTION_HH
5
8
9namespace Dune
10{
11
12 namespace GeoGrid
13 {
14
15 // Intersection
16 // ------------
17
18 template< class Grid, class HostIntersection >
20 {
21 typedef typename HostIntersection::Geometry HostGeometry;
22 typedef typename HostIntersection::LocalGeometry HostLocalGeometry;
23
24 typedef typename std::remove_const< Grid >::type::Traits Traits;
25
26 public:
27 typedef typename Traits::ctype ctype;
28
29 static const int dimension = Traits::dimension;
30 static const int dimensionworld = Traits::dimensionworld;
31
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;
35
36 typedef typename Traits::template Codim< 0 >::Geometry ElementGeometry;
37
38 private:
40
41 typedef typename Traits::template Codim< 0 >::EntityImpl EntityImpl;
42
43 typedef typename Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
44 typedef typename Traits::template Codim< 0 >::GeometryImpl ElementGeometryImpl;
45
46 public:
47
49 {}
50
51 explicit Intersection ( const HostIntersection &hostIntersection, const ElementGeometryImpl &insideGeo )
52 : hostIntersection_( hostIntersection )
53 , insideGeo_ ( insideGeo )
54 , geo_( grid() )
55 {}
56
57 explicit Intersection ( HostIntersection&& hostIntersection, const ElementGeometryImpl &insideGeo )
58 : hostIntersection_( std::move( hostIntersection ) )
59 , insideGeo_ ( insideGeo )
60 , geo_( grid() )
61 {}
62
63 bool equals ( const Intersection &other) const
64 {
65 return hostIntersection_ == other.hostIntersection_;
66 }
67
68 operator bool () const { return bool( hostIntersection_ ); }
69
70 Entity inside () const
71 {
72 return EntityImpl( insideGeo_, hostIntersection().inside() );
73 }
74
75 Entity outside () const
76 {
77 return EntityImpl( grid(), hostIntersection().outside() );
78 }
79
80 bool boundary () const { return hostIntersection().boundary(); }
81
82 bool conforming () const { return hostIntersection().conforming(); }
83
84 bool neighbor () const { return hostIntersection().neighbor(); }
85
86 int boundaryId () const { return hostIntersection().boundaryId(); }
87
88 size_t boundarySegmentIndex () const
89 {
90 return hostIntersection().boundarySegmentIndex();
91 }
92
94 {
95 return hostIntersection().geometryInInside();
96 }
97
99 {
100 return hostIntersection().geometryInOutside();
101 }
102
104 {
105 if( !geo_ )
106 {
107 CoordVector coords( insideGeo_, geometryInInside() );
108 geo_ = GeometryImpl( grid(), type(), coords );
109 }
110 return Geometry( geo_ );
111 }
112
113 GeometryType type () const { return hostIntersection().type(); }
114
115 int indexInInside () const
116 {
117 return hostIntersection().indexInInside();
118 }
119
120 int indexInOutside () const
121 {
122 return hostIntersection().indexInOutside();
123 }
124
125 FieldVector< ctype, dimensionworld >
126 integrationOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
127 {
128 const LocalGeometry geoInInside = geometryInInside();
129 const int idxInInside = indexInInside();
130
131 const ReferenceElement< ctype, dimension > &refElement
132 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
133
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 );
137
138 FieldVector< ctype, dimensionworld > normal;
139 jit.mv( refNormal, normal );
140 if( !conforming() )
141 normal *= geoInInside.volume() / refElement.template geometry< 1 >( idxInInside ).volume();
142 normal *= jit.detInv();
143 //normal *= insideGeo_.integrationElement( x );
144 return normal;
145 }
146
147 FieldVector< ctype, dimensionworld >
148 outerNormal ( const FieldVector< ctype, dimension-1 > &local ) const
149 {
150 const ReferenceElement< ctype, dimension > &refElement
151 = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
152
153 FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
154 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
155 const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( indexInInside() );
156
157 FieldVector< ctype, dimensionworld > normal;
158 jit.mv( refNormal, normal );
159 return normal;
160 }
161
162 FieldVector< ctype, dimensionworld >
163 unitOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
164 {
165 FieldVector< ctype, dimensionworld > normal = outerNormal( local );
166 normal *= (ctype( 1 ) / normal.two_norm());
167 return normal;
168 }
169
170 FieldVector< ctype, dimensionworld > centerUnitOuterNormal () const
171 {
172 const ReferenceElement< ctype, dimension-1 > &refFace
173 = ReferenceElements< ctype, dimension-1 >::general( type() );
174 return unitOuterNormal( refFace.position( 0, 0 ) );
175 }
176
177 const HostIntersection &hostIntersection () const
178 {
179 return hostIntersection_;
180 }
181
182 const Grid &grid () const { return insideGeo_.grid(); }
183
184 private:
185 HostIntersection hostIntersection_;
186 ElementGeometryImpl insideGeo_;
187 mutable GeometryImpl geo_;
188 };
189
190 } // namespace GeoGrid
191
192} // namespace Dune
193
194#endif // #ifndef DUNE_GEOGRID_INTERSECTION_HH
STL namespace.
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