dune-grid 2.8.0
Loading...
Searching...
No Matches
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 explicit 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 size_t boundarySegmentIndex () const
87 {
88 return hostIntersection().boundarySegmentIndex();
89 }
90
92 {
93 return hostIntersection().geometryInInside();
94 }
95
97 {
98 return hostIntersection().geometryInOutside();
99 }
100
102 {
103 if( !geo_ )
104 {
105 CoordVector coords( insideGeo_, geometryInInside() );
106 geo_ = GeometryImpl( grid(), type(), coords );
107 }
108 return Geometry( geo_ );
109 }
110
111 GeometryType type () const { return hostIntersection().type(); }
112
113 int indexInInside () const
114 {
115 return hostIntersection().indexInInside();
116 }
117
118 int indexInOutside () const
119 {
120 return hostIntersection().indexInOutside();
121 }
122
123 FieldVector< ctype, dimensionworld >
124 integrationOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
125 {
126 const LocalGeometry geoInInside = geometryInInside();
127 const int idxInInside = indexInInside();
128
129 auto refElement = referenceElement< ctype, dimension >( insideGeo_.type() );
130
131 FieldVector< ctype, dimension > x( geoInInside.global( local ) );
132 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
133 FieldVector< ctype, dimension > refNormal = refElement.integrationOuterNormal( idxInInside );
134
135 FieldVector< ctype, dimensionworld > normal;
136 jit.mv( refNormal, normal );
137 if( !conforming() )
138 normal *= geoInInside.volume() / refElement.template geometry< 1 >( idxInInside ).volume();
139 normal *= jit.detInv();
140 //normal *= insideGeo_.integrationElement( x );
141 return normal;
142 }
143
144 FieldVector< ctype, dimensionworld >
145 outerNormal ( const FieldVector< ctype, dimension-1 > &local ) const
146 {
147 auto refElement = referenceElement< ctype, dimension >( insideGeo_.type() );
148
149 FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
150 const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
151 FieldVector< ctype, dimension > refNormal = refElement.integrationOuterNormal( indexInInside() );
152
153 FieldVector< ctype, dimensionworld > normal;
154 jit.mv( refNormal, normal );
155 return normal;
156 }
157
158 FieldVector< ctype, dimensionworld >
159 unitOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
160 {
161 FieldVector< ctype, dimensionworld > normal = outerNormal( local );
162 normal *= (ctype( 1 ) / normal.two_norm());
163 return normal;
164 }
165
166 FieldVector< ctype, dimensionworld > centerUnitOuterNormal () const
167 {
168 auto refFace = referenceElement< ctype, dimension-1 >( type() );
169 return unitOuterNormal( refFace.position( 0, 0 ) );
170 }
171
172 const HostIntersection &hostIntersection () const
173 {
174 return hostIntersection_;
175 }
176
177 const Grid &grid () const { return insideGeo_.grid(); }
178
179 private:
180 HostIntersection hostIntersection_;
181 ElementGeometryImpl insideGeo_;
182 mutable GeometryImpl geo_;
183 };
184
185 } // namespace GeoGrid
186
187} // namespace Dune
188
189#endif // #ifndef DUNE_GEOGRID_INTERSECTION_HH
STL namespace.
Include standard header files.
Definition: agrid.hh:58
auto referenceElement(const Geometry< mydim, cdim, GridImp, GeometryImp > &geo) -> decltype(referenceElement(geo, geo.impl()))
Definition: common/geometry.hh:417
Grid abstract base class.
Definition: common/grid.hh:372
Definition: cornerstorage.hh:121
Definition: geometrygrid/intersection.hh:20
FieldVector< ctype, dimensionworld > integrationOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:124
Traits::ctype ctype
Definition: geometrygrid/intersection.hh:27
GeometryType type() const
Definition: geometrygrid/intersection.hh:111
LocalGeometry geometryInInside() const
Definition: geometrygrid/intersection.hh:91
Geometry geometry() const
Definition: geometrygrid/intersection.hh:101
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:86
FieldVector< ctype, dimensionworld > outerNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:145
FieldVector< ctype, dimensionworld > unitOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:159
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:172
int indexInOutside() const
Definition: geometrygrid/intersection.hh:118
const Grid & grid() const
Definition: geometrygrid/intersection.hh:177
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:113
FieldVector< ctype, dimensionworld > centerUnitOuterNormal() const
Definition: geometrygrid/intersection.hh:166
static const int dimension
Definition: geometrygrid/intersection.hh:29
Entity inside() const
Definition: geometrygrid/intersection.hh:70
LocalGeometry geometryInOutside() const
Definition: geometrygrid/intersection.hh:96