dune-grid  2.5-git
common/geometry.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_GRID_GEOMETRY_HH
4 #define DUNE_GRID_GEOMETRY_HH
5 
10 #include <cassert>
11 
12 #include <dune/common/fmatrix.hh>
13 #include <dune/common/typetraits.hh>
14 
15 #include <dune/geometry/referenceelements.hh>
16 
17 namespace Dune
18 {
19 
20  // External Forward Declarations
21  // -----------------------------
22 
23  template< int dim, int dimworld, class ct, class GridFamily >
25 
26 
27  //*****************************************************************************
28  //
29  // Geometry
30  // forwards the interface to the implementation
31  //
32  //*****************************************************************************
33 
64  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp >
65  class Geometry
66  {
67  #if DUNE_GRID_EXPERIMENTAL_GRID_EXTENSIONS
68  public:
69  #else
70  protected:
71  // give the GridDefaultImplementation class access to the realImp
72  friend class GridDefaultImplementation<
73  GridImp::dimension, GridImp::dimensionworld,
74  typename GridImp::ctype,
75  typename GridImp::GridFamily> ;
76  #endif
77  // type of underlying implementation, for internal use only
78  typedef GeometryImp< mydim, cdim, GridImp > Implementation;
79 
81  const Implementation &impl () const { return realGeometry; }
82 
83  public:
85  enum { mydimension=mydim };
87  enum { coorddimension=cdim };
88 
90  typedef typename GridImp::ctype ctype;
91 
93  typedef FieldVector<ctype, mydim> LocalCoordinate;
94 
96  typedef FieldVector< ctype, cdim > GlobalCoordinate;
97 
107  typedef typename Implementation::JacobianInverseTransposed JacobianInverseTransposed;
108 
118  typedef typename Implementation::JacobianTransposed JacobianTransposed;
119 
123  GeometryType type () const { return impl().type(); }
124 
126  bool affine() const { return impl().affine(); }
127 
134  int corners () const { return impl().corners(); }
135 
148  GlobalCoordinate corner ( int i ) const
149  {
150  return impl().corner( i );
151  }
152 
157  GlobalCoordinate global (const LocalCoordinate& local) const
158  {
159  return impl().global( local );
160  }
161 
166  LocalCoordinate local (const GlobalCoordinate& global) const
167  {
168  return impl().local( global );
169  }
170 
194  ctype integrationElement (const LocalCoordinate& local) const
195  {
196  return impl().integrationElement( local );
197  }
198 
200  ctype volume () const
201  {
202  return impl().volume();
203  }
204 
215  GlobalCoordinate center () const
216  {
217  return impl().center();
218  }
219 
231  JacobianTransposed jacobianTransposed ( const LocalCoordinate& local ) const
232  {
233  return impl().jacobianTransposed( local );
234  }
235 
257  JacobianInverseTransposed jacobianInverseTransposed ( const LocalCoordinate &local ) const
258  {
259  return impl().jacobianInverseTransposed(local);
260  }
261  //===========================================================
265  //===========================================================
266 
268  explicit Geometry ( const Implementation &impl )
269  : realGeometry( impl )
270  {}
271 
273 
274  private:
276  const Geometry &operator= ( const Geometry &rhs );
277 
278  protected:
279 
281  };
282 
283 
284 
285  //************************************************************************
286  // GEOMETRY Default Implementations
287  //*************************************************************************
288  //
289  // --GeometryDefault
290  //
292  template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>
294  {
295  public:
296  static const int mydimension = mydim;
297  static const int coorddimension = cdim;
298 
299  // save typing
300  typedef typename GridImp::ctype ctype;
301 
302  typedef FieldVector< ctype, mydim > LocalCoordinate;
303  typedef FieldVector< ctype, cdim > GlobalCoordinate;
304 
306  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
307 
309  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
310 
312  ctype volume () const
313  {
314  GeometryType type = asImp().type();
315 
316  // get corresponding reference element
317  const ReferenceElement< ctype , mydim > & refElement =
318  ReferenceElements< ctype, mydim >::general(type);
319 
320  LocalCoordinate localBaryCenter ( 0 );
321  // calculate local bary center
322  const int corners = refElement.size(0,0,mydim);
323  for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
324  localBaryCenter *= (ctype) (1.0/corners);
325 
326  // volume is volume of reference element times integrationElement
327  return refElement.volume() * asImp().integrationElement(localBaryCenter);
328  }
329 
331  GlobalCoordinate center () const
332  {
333  GeometryType type = asImp().type();
334 
335  // get corresponding reference element
336  const ReferenceElement< ctype , mydim > & refElement =
337  ReferenceElements< ctype, mydim >::general(type);
338 
339  // center is (for now) the centroid of the reference element mapped to
340  // this geometry.
341  return asImp().global(refElement.position(0,0));
342  }
343 
344  private:
346  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
347  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
348  }; // end GeometryDefault
349 
350  template<int cdim, class GridImp, template<int,int,class> class GeometryImp>
351  class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
352  {
353  // my dimension
354  enum { mydim = 0 };
355 
356  public:
357  static const int mydimension = mydim;
358  static const int coorddimension = cdim;
359 
360  // save typing
361  typedef typename GridImp::ctype ctype;
362 
363  typedef FieldVector< ctype, mydim > LocalCoordinate;
364  typedef FieldVector< ctype, cdim > GlobalCoordinate;
365 
367  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
368 
370  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
371 
373  FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
374  {
375  return asImp().corner(0);
376  }
377 
379  FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
380  {
381  return FieldVector<ctype, mydim>();
382  }
383 
385  ctype volume () const
386  {
387  return 1.0;
388  }
389 
391  FieldVector<ctype, cdim> center () const
392  {
393  return asImp().corner(0);
394  }
395 
396  private:
397  // Barton-Nackman trick
398  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
399  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
400  }; // end GeometryDefault
401 
402 } // namespace Dune
403 
404 #endif // DUNE_GRID_GEOMETRY_HH
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:303
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:370
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:363
bool affine() const
Return true if the geometry mapping is affine and false otherwise.
Definition: common/geometry.hh:126
GridImp::ctype ctype
Definition: common/geometry.hh:300
Implementation::JacobianTransposed JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:118
FieldVector< ctype, cdim > center() const
return center of the geometry
Definition: common/geometry.hh:391
GeometryType
Type representing VTK&#39;s entity geometry types.
Definition: common.hh:178
GeometryType type() const
Return the type of the reference element. The type can be used to access the Dune::ReferenceElement.
Definition: common/geometry.hh:123
GlobalCoordinate corner(int i) const
Obtain a corner of the geometry.
Definition: common/geometry.hh:148
FieldVector< ctype, mydim > local(const FieldVector< ctype, cdim > &) const
return empty vector
Definition: common/geometry.hh:379
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:364
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:309
Definition: common/geometry.hh:24
GlobalCoordinate center() const
return center of geometry
Definition: common/geometry.hh:215
LocalCoordinate local(const GlobalCoordinate &global) const
Evaluate the inverse map .
Definition: common/geometry.hh:166
GeometryImp< mydim, cdim, GridImp > Implementation
Definition: common/geometry.hh:78
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:367
const Implementation & impl() const
return reference to the implementation
Definition: common/geometry.hh:81
Default implementation for class Geometry.
Definition: common/geometry.hh:293
int corners() const
Return the number of corners of the reference element.
Definition: common/geometry.hh:134
Implementation realGeometry
Definition: common/geometry.hh:280
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:306
Definition: common/geometry.hh:87
Include standard header files.
Definition: agrid.hh:59
FieldVector< ctype, cdim > GlobalCoordinate
type of the global coordinates
Definition: common/geometry.hh:96
GridImp::ctype ctype
Definition: common/geometry.hh:361
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the map .
Definition: common/geometry.hh:157
Wrapper class for geometries.
Definition: common/geometry.hh:65
ctype integrationElement(const LocalCoordinate &local) const
Return the factor appearing in the integral transformation formula.
Definition: common/geometry.hh:194
FieldVector< ctype, cdim > global(const FieldVector< ctype, mydim > &local) const
return the only coordinate
Definition: common/geometry.hh:373
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:302
ctype volume() const
return volume of the geometry
Definition: common/geometry.hh:385
Implementation::JacobianInverseTransposed JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:107
ctype volume() const
return volume of geometry
Definition: common/geometry.hh:200
GlobalCoordinate center() const
return center of the geometry
Definition: common/geometry.hh:331
ctype volume() const
return volume of the geometry
Definition: common/geometry.hh:312
GridImp::ctype ctype
define type used for coordinates in grid module
Definition: common/geometry.hh:90
Geometry(const Implementation &impl)
copy constructor from implementation
Definition: common/geometry.hh:268
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
Return the transposed of the Jacobian.
Definition: common/geometry.hh:231
FieldVector< ctype, mydim > LocalCoordinate
type of local coordinates
Definition: common/geometry.hh:93
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
Return inverse of transposed of Jacobian.
Definition: common/geometry.hh:257
Definition: common/geometry.hh:85