dune-grid  2.3.1
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 
9 
10 namespace Dune
11 {
12 
13  namespace GeoGrid
14  {
15 
16  // Intersection
17  // ------------
18 
19  template< class Grid, class HostIntersection >
21  {
22  typedef typename HostIntersection::Geometry HostGeometry;
23  typedef typename HostIntersection::LocalGeometry HostLocalGeometry;
24 
25  typedef typename remove_const< Grid >::type::Traits Traits;
26 
27  public:
28  typedef typename Traits::ctype ctype;
29 
30  static const int dimension = Traits::dimension;
31  static const int dimensionworld = Traits::dimensionworld;
32 
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;
37 
38  typedef typename Traits::template Codim< 0 >::Geometry ElementGeometry;
39 
40  private:
42 
43  typedef typename Traits::template Codim< 0 >::EntityPointerImpl EntityPointerImpl;
44 
45  typedef typename Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
46  typedef typename Traits::template Codim< 0 >::GeometryImpl ElementGeometryImpl;
47 
48  public:
49  explicit Intersection ( const ElementGeometry &insideGeo )
50  : insideGeo_( Grid::getRealImplementation( insideGeo ) ),
51  hostIntersection_( 0 ),
52  geo_( grid() )
53  {}
54 
55  Intersection ( const Intersection &other )
56  : insideGeo_( other.insideGeo_ ),
57  hostIntersection_( 0 ),
58  geo_( grid() )
59  {}
60 
61  const Intersection &operator= ( const Intersection &other )
62  {
63  insideGeo_ = other.insideGeo_;
64  invalidate();
65  return *this;
66  }
67 
68  operator bool () const { return bool( hostIntersection_ ); }
69 
71  {
72  return EntityPointerImpl( insideGeo_, hostIntersection().inside() );
73  }
74 
76  {
77  return EntityPointerImpl( 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 ReferenceElement< ctype, dimension > &refElement
129  = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
130 
131  FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
132  const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
133  const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( indexInInside() );
134 
135  FieldVector< ctype, dimensionworld > normal;
136  jit.mv( refNormal, normal );
137  normal *= jit.detInv();
138  //normal *= insideGeo_.integrationElement( x );
139  return normal;
140  }
141 
142  FieldVector< ctype, dimensionworld >
143  outerNormal ( const FieldVector< ctype, dimension-1 > &local ) const
144  {
145  const ReferenceElement< ctype, dimension > &refElement
146  = ReferenceElements< ctype, dimension>::general( insideGeo_.type() );
147 
148  FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
149  const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
150  const FieldVector< ctype, dimension > &refNormal = refElement.integrationOuterNormal( indexInInside() );
151 
152  FieldVector< ctype, dimensionworld > normal;
153  jit.mv( refNormal, normal );
154  return normal;
155  }
156 
157  FieldVector< ctype, dimensionworld >
158  unitOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
159  {
160  FieldVector< ctype, dimensionworld > normal = outerNormal( local );
161  normal *= (ctype( 1 ) / normal.two_norm());
162  return normal;
163  }
164 
165  FieldVector< ctype, dimensionworld > centerUnitOuterNormal () const
166  {
167  const ReferenceElement< ctype, dimension-1 > &refFace
168  = ReferenceElements< ctype, dimension-1 >::general( type() );
169  return unitOuterNormal( refFace.position( 0, 0 ) );
170  }
171 
172  const HostIntersection &hostIntersection () const
173  {
174  assert( *this );
175  return *hostIntersection_;
176  }
177 
178  const Grid &grid () const { return insideGeo_.grid(); }
179 
180  void invalidate ()
181  {
182  hostIntersection_ = 0;
183  geo_ = GeometryImpl( grid() );
184  }
185 
186  void initialize ( const HostIntersection &hostIntersection )
187  {
188  assert( !(*this) );
189  hostIntersection_ = &hostIntersection;
190  }
191 
192  private:
193  ElementGeometryImpl insideGeo_;
194  const HostIntersection *hostIntersection_;
195  mutable GeometryImpl geo_;
196  };
197 
198  } // namespace GeoGrid
199 
200 } // namespace Dune
201 
202 #endif // #ifndef DUNE_GEOGRID_INTERSECTION_HH
GeometryType type() const
Definition: geometrygrid/intersection.hh:113
Traits::template Codim< 0 >::EntityPointer EntityPointer
Definition: geometrygrid/intersection.hh:34
int indexInOutside() const
Definition: geometrygrid/intersection.hh:120
Traits::ctype ctype
Definition: geometrygrid/intersection.hh:28
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
const Grid & grid() const
Definition: geometrygrid/intersection.hh:178
Geometry geometry() const
Definition: geometrygrid/intersection.hh:103
Intersection(const Intersection &other)
Definition: geometrygrid/intersection.hh:55
bool boundary() const
Definition: geometrygrid/intersection.hh:80
void initialize(const HostIntersection &hostIntersection)
Definition: geometrygrid/intersection.hh:186
static const int dimensionworld
Definition: geometrygrid/intersection.hh:31
LocalGeometry geometryInInside() const
Definition: geometrygrid/intersection.hh:93
LocalGeometry geometryInOutside() const
Definition: geometrygrid/intersection.hh:98
bool neighbor() const
Definition: geometrygrid/intersection.hh:84
FieldVector< ctype, dimensionworld > integrationOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:126
Definition: geometrygrid/intersection.hh:20
Intersection(const ElementGeometry &insideGeo)
Definition: geometrygrid/intersection.hh:49
const HostIntersection & hostIntersection() const
Definition: geometrygrid/intersection.hh:172
Traits::template Codim< 0 >::Geometry ElementGeometry
Definition: geometrygrid/intersection.hh:38
size_t boundarySegmentIndex() const
Definition: geometrygrid/intersection.hh:88
Grid abstract base classThis class is the base class for all grid implementations. Although no virtual functions are used we call it abstract since its methods do not contain an implementation but forward to the methods of the derived class via the Barton-Nackman trick.
Definition: common/grid.hh:386
int indexInInside() const
Definition: geometrygrid/intersection.hh:115
FieldVector< ctype, dimensionworld > outerNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:143
FieldVector< ctype, dimensionworld > centerUnitOuterNormal() const
Definition: geometrygrid/intersection.hh:165
int boundaryId() const
Definition: geometrygrid/intersection.hh:86
Traits::template Codim< 1 >::LocalGeometry LocalGeometry
Definition: geometrygrid/intersection.hh:36
FieldVector< ctype, dimensionworld > unitOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:158
Traits::template Codim< 1 >::Geometry Geometry
Definition: geometrygrid/intersection.hh:35
void invalidate()
Definition: geometrygrid/intersection.hh:180
EntityPointer inside() const
Definition: geometrygrid/intersection.hh:70
bool conforming() const
Definition: geometrygrid/intersection.hh:82
Traits::template Codim< 0 >::Entity Entity
Definition: geometrygrid/intersection.hh:33
static const int dimension
Definition: geometrygrid/intersection.hh:30
const Intersection & operator=(const Intersection &other)
Definition: geometrygrid/intersection.hh:61
EntityPointer outside() const
Definition: geometrygrid/intersection.hh:75
Definition: cornerstorage.hh:121