dune-fem  2.8-git
elementpointlistbase.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_ELEMENTPOINTLISTBASE_HH
2 #define DUNE_FEM_ELEMENTPOINTLISTBASE_HH
3 
4 #include <dune/geometry/referenceelements.hh>
5 
7 
9 
10 namespace Dune
11 {
12 
13  namespace Fem
14  {
15 
17  template< class GridPartImp, int codim, class IntegrationTraits >
18  class ElementPointListBase;
19 
20 
21  template< class GridPartImp, class IntegrationTraits >
22  class ElementPointListBase< GridPartImp, 0, IntegrationTraits >
23  {
25 
26  public:
28  typedef GridPartImp GridPartType;
29 
31  enum Side { INSIDE, OUTSIDE };
32 
34  static const int codimension = 0;
35 
37  typedef typename GridPartType::ctype RealType;
38 
40  static const int dimension = GridPartType::dimension;
41 
43  typedef typename IntegrationTraits::IntegrationPointListType IntegrationPointListType;
44 
45  typedef typename IntegrationTraits::CoordinateType CoordinateType;
46  typedef typename IntegrationPointListType::CoordinateType LocalCoordinateType;
47 
48  typedef typename IntegrationPointListType :: QuadratureKeyType QuadratureKeyType;
49 
55  ElementPointListBase ( const GeometryType &geometry, const QuadratureKeyType& quadKey )
56  : quad_( geometry, quadKey )
57  {}
58 
65  : quad_( ipList )
66  {}
67 
69  size_t nop () const
70  {
71  return quadImp().nop();
72  }
73 
84  const LocalCoordinateType &localPoint( size_t i ) const
85  {
86  return quadImp().point( i );
87  }
88 
91  size_t id () const
92  {
93  return quadImp().id();
94  }
95 
98  int order () const
99  {
100  return quadImp().order();
101  }
102 
105  GeometryType geometry () const
106  {
107  return quadImp().geometryType();
108  }
109 
112  GeometryType type () const
113  {
114  return quadImp().geometryType();
115  }
116 
119  GeometryType geometryType () const
120  {
121  return quadImp().geometryType();
122  }
123 
137  GeometryType elementGeometry () const
138  {
139  return quadImp().geometry();
140  }
141 
143  size_t cachingPoint( const size_t quadraturePoint ) const
144  {
145  return quadraturePoint;
146  }
147 
149  size_t localCachingPoint( const size_t quadraturePoint ) const
150  {
151  return quadraturePoint;
152  }
153 
155  inline bool twisted () const { return false; }
156 
158  inline int twistId () const { return 0; }
159 
160  int localFaceIndex () const
161  {
162  return 0;
163  }
164 
165  inline int nCachingPoints () const { return nop(); }
166  inline int cachingPointStart () const { return 0; }
167 
168  protected:
176  {
177  return quad_;
178  }
179 
180  protected:
182  };
183 
184 
185 
187  template< class GridPartImp, int codim, class IntegrationTraits >
189  {
191 
192  public:
194  typedef GridPartImp GridPartType;
195 
197  enum Side { INSIDE, OUTSIDE };
198 
200  static const int codimension = codim;
201 
203  typedef typename GridPartType::ctype RealType;
204 
206  static const int dimension = GridPartType::dimension;
207 
209  typedef typename IntegrationTraits::IntegrationPointListType IntegrationPointListType;
210 
211  typedef typename IntegrationTraits::CoordinateType CoordinateType;
212  typedef typename IntegrationPointListType::CoordinateType LocalCoordinateType;
213 
214  typedef typename IntegrationPointListType :: QuadratureKeyType QuadratureKeyType;
215 
223  ElementPointListBase ( const GeometryType &elementGeo,
224  const GeometryType &faceGeo,
225  const int localFaceIndex,
226  const QuadratureKeyType& quadKey )
227  : quad_( faceGeo, quadKey ),
228  elementGeometry_( elementGeo ),
229  localFaceIndex_( localFaceIndex )
230  {}
231 
238  ElementPointListBase ( const GeometryType &elementGeo,
239  const int localFaceIndex,
240  const QuadratureKeyType& quadKey )
241  : quad_( getFaceGeometry( elementGeo, localFaceIndex ), quadKey ),
242  elementGeometry_( elementGeo ),
243  localFaceIndex_( localFaceIndex )
244  {}
245 
248  size_t nop () const
249  {
250  return quadImp().nop();
251  }
252 
263  const LocalCoordinateType &localPoint ( size_t i ) const
264  {
265  return quad_.point( i );
266  }
267 
270  size_t id () const
271  {
272  return quadImp().id();
273  }
274 
277  int order () const
278  {
279  return quadImp().order();
280  }
281 
284  GeometryType geometry () const
285  {
286  return quadImp().geo();
287  }
288 
302  GeometryType elementGeometry () const
303  {
304  return elementGeometry_;
305  }
306 
307  size_t cachingPoint( const size_t quadraturePoint ) const
308  {
309  return quadraturePoint;
310  }
311 
312  size_t localCachingPoint( const size_t quadraturePoint ) const
313  {
314  return quadraturePoint;
315  }
316 
318  inline bool twisted () const { return false; }
319 
321  inline int twistId () const { return 0; }
322 
323  inline int nCachingPoints () const { return nop(); }
324  inline int cachingPointStart () const { return 0; }
325 
326  int localFaceIndex () const
327  {
328  return localFaceIndex_;
329  }
330 
331  protected:
339  {
340  return quad_;
341  }
342 
343  static GeometryType
344  getFaceGeometry ( const GeometryType &elementGeo, const int face )
345  {
346  // for cube and simplex geom types the dim-1 geom type
347  // is also cube or simplex
348 
349  static const bool isCube =
352 
353  static const bool isSimplex =
356 
357  if( isCube || isSimplex )
358  {
359  assert( elementGeo.dim() == dimension );
360  if( isCube )
361  {
362  return Dune::GeometryTypes::cube( dimension-1 );
363  }
364  else
365  {
366  assert( isSimplex );
367  return Dune::GeometryTypes::simplex( dimension-1 );
368  }
369  }
370  else if( elementGeo.isNone() )
371  {
372  // if cell geometry is none and dim is 2 then the
373  // face is a normal edge which is of type cube
374  if( elementGeo.dim() == 2 )
375  {
376  return Dune::GeometryTypes::cube( 1 );
377  }
378  else
379  {
380  return Dune::GeometryTypes::none( elementGeo.dim()-1 );
381  }
382  }
383  else // use reference element to determine type
384  {
385  assert( ! elementGeo.isNone() );
386  typedef Dune::ReferenceElements< RealType, dimension > RefElements;
387  return RefElements::general( elementGeo ).type( face, codimension );
388  }
389  }
390 
391  private:
393  GeometryType elementGeometry_;
394  int localFaceIndex_;
395  };
396 
397  } // namespace Fem
398 
399 } // namespace Dune
400 
401 #endif // #ifndef DUNE_FEM_ELEMENTPOINTLISTBASE_HH
Definition: bindguard.hh:11
specialize with 'true' for if the codimension 0 entity of the grid part has only one possible geometr...
Definition: gridpart/common/capabilities.hh:29
ElementPointListBase.
Definition: elementpointlistbase.hh:189
IntegrationPointListType::CoordinateType LocalCoordinateType
Definition: elementpointlistbase.hh:212
ElementPointListBase(const GeometryType &elementGeo, const GeometryType &faceGeo, const int localFaceIndex, const QuadratureKeyType &quadKey)
constructor
Definition: elementpointlistbase.hh:223
Side
inside and outside flags
Definition: elementpointlistbase.hh:197
@ OUTSIDE
Definition: elementpointlistbase.hh:197
@ INSIDE
Definition: elementpointlistbase.hh:197
const LocalCoordinateType & localPoint(size_t i) const
obtain local coordinates of i-th integration point
Definition: elementpointlistbase.hh:263
IntegrationTraits::CoordinateType CoordinateType
Definition: elementpointlistbase.hh:211
bool twisted() const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:318
size_t localCachingPoint(const size_t quadraturePoint) const
Definition: elementpointlistbase.hh:312
GridPartImp GridPartType
type of the grid partition
Definition: elementpointlistbase.hh:194
int twistId() const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:321
ElementPointListBase(const GeometryType &elementGeo, const int localFaceIndex, const QuadratureKeyType &quadKey)
constructor
Definition: elementpointlistbase.hh:238
static const int dimension
dimension of the grid
Definition: elementpointlistbase.hh:206
size_t cachingPoint(const size_t quadraturePoint) const
Definition: elementpointlistbase.hh:307
const IntegrationPointListType & quadImp() const
obtain the actual implementation of the quadrature
Definition: elementpointlistbase.hh:338
static GeometryType getFaceGeometry(const GeometryType &elementGeo, const int face)
Definition: elementpointlistbase.hh:344
size_t nop() const
obtain the number of integration points
Definition: elementpointlistbase.hh:248
int localFaceIndex() const
Definition: elementpointlistbase.hh:326
IntegrationTraits::IntegrationPointListType IntegrationPointListType
type of the integration point list
Definition: elementpointlistbase.hh:209
GeometryType geometry() const
obtain GeometryType for this integration point list
Definition: elementpointlistbase.hh:284
int cachingPointStart() const
Definition: elementpointlistbase.hh:324
IntegrationPointListType ::QuadratureKeyType QuadratureKeyType
Definition: elementpointlistbase.hh:214
int nCachingPoints() const
Definition: elementpointlistbase.hh:323
static const int codimension
codimension of the element integration point list
Definition: elementpointlistbase.hh:200
size_t id() const
obtain the identifier of the integration point list
Definition: elementpointlistbase.hh:270
GridPartType::ctype RealType
coordinate type
Definition: elementpointlistbase.hh:203
int order() const
obtain order of the integration point list
Definition: elementpointlistbase.hh:277
GeometryType elementGeometry() const
obtain GeometryType of the corresponding codim-0 the integration point list belongs to
Definition: elementpointlistbase.hh:302
Side
inside and outside flags
Definition: elementpointlistbase.hh:31
ElementPointListBase(const IntegrationPointListType &ipList)
constructor
Definition: elementpointlistbase.hh:64
GridPartType::ctype RealType
coordinate type
Definition: elementpointlistbase.hh:37
int nCachingPoints() const
Definition: elementpointlistbase.hh:165
IntegrationPointListType ::QuadratureKeyType QuadratureKeyType
Definition: elementpointlistbase.hh:48
size_t cachingPoint(const size_t quadraturePoint) const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:143
int twistId() const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:158
GeometryType elementGeometry() const
obtain GeometryType of the corresponding codim-0 the integration point list belongs to
Definition: elementpointlistbase.hh:137
GeometryType type() const
Definition: elementpointlistbase.hh:112
bool twisted() const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:155
GeometryType geometry() const
Definition: elementpointlistbase.hh:105
IntegrationPointListType quad_
Definition: elementpointlistbase.hh:181
int order() const
obtain order of the integration point list
Definition: elementpointlistbase.hh:98
ElementPointListBase(const GeometryType &geometry, const QuadratureKeyType &quadKey)
constructor
Definition: elementpointlistbase.hh:55
const IntegrationPointListType & quadImp() const
obtain the actual implementation of the quadrature
Definition: elementpointlistbase.hh:175
IntegrationTraits::IntegrationPointListType IntegrationPointListType
type of the integration point list
Definition: elementpointlistbase.hh:43
const LocalCoordinateType & localPoint(size_t i) const
obtain local coordinates of i-th integration point
Definition: elementpointlistbase.hh:84
IntegrationPointListType::CoordinateType LocalCoordinateType
Definition: elementpointlistbase.hh:46
int localFaceIndex() const
Definition: elementpointlistbase.hh:160
size_t id() const
obtain the identifier of the integration point list
Definition: elementpointlistbase.hh:91
size_t localCachingPoint(const size_t quadraturePoint) const
convenience implementation for Dune::Fem::CachingInterface
Definition: elementpointlistbase.hh:149
GeometryType geometryType() const
Definition: elementpointlistbase.hh:119
GridPartImp GridPartType
type of the grid partition
Definition: elementpointlistbase.hh:28
int cachingPointStart() const
Definition: elementpointlistbase.hh:166
IntegrationTraits::CoordinateType CoordinateType
Definition: elementpointlistbase.hh:45
size_t nop() const
obtain the number of integration points
Definition: elementpointlistbase.hh:69