dune-grid  2.4
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:
87  enum { dimension=GridImp::dimension };
89  enum { mydimension=mydim };
91  enum { coorddimension=cdim };
92 
96  enum { dimensionworld=GridImp::dimensionworld };
98  typedef typename GridImp::ctype ctype;
99 
101  typedef FieldVector<ctype, mydim> LocalCoordinate;
102 
104  typedef FieldVector< ctype, cdim > GlobalCoordinate;
105 
115  typedef typename Implementation::JacobianInverseTransposed JacobianInverseTransposed;
116 
126  typedef typename Implementation::JacobianTransposed JacobianTransposed;
127 
131  GeometryType type () const { return impl().type(); }
132 
134  bool affine() const { return impl().affine(); }
135 
142  int corners () const { return impl().corners(); }
143 
156  GlobalCoordinate corner ( int i ) const
157  {
158  return impl().corner( i );
159  }
160 
165  GlobalCoordinate global (const LocalCoordinate& local) const
166  {
167  return impl().global( local );
168  }
169 
174  LocalCoordinate local (const GlobalCoordinate& global) const
175  {
176  return impl().local( global );
177  }
178 
202  ctype integrationElement (const LocalCoordinate& local) const
203  {
204  return impl().integrationElement( local );
205  }
206 
208  ctype volume () const
209  {
210  return impl().volume();
211  }
212 
223  GlobalCoordinate center () const
224  {
225  return impl().center();
226  }
227 
239  JacobianTransposed jacobianTransposed ( const LocalCoordinate& local ) const
240  {
241  return impl().jacobianTransposed( local );
242  }
243 
265  JacobianInverseTransposed jacobianInverseTransposed ( const LocalCoordinate &local ) const
266  {
267  return impl().jacobianInverseTransposed(local);
268  }
269  //===========================================================
273  //===========================================================
274 
276  explicit Geometry ( const Implementation &impl )
277  : realGeometry( impl )
278  {}
279 
281 
282  private:
284  const Geometry &operator= ( const Geometry &rhs );
285 
286  protected:
287 
289  };
290 
291 
292 
293  //************************************************************************
294  // GEOMETRY Default Implementations
295  //*************************************************************************
296  //
297  // --GeometryDefault
298  //
300  template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>
302  {
303  public:
304  static const int mydimension = mydim;
305  static const int coorddimension = cdim;
306 
307  // save typing
308  typedef typename GridImp::ctype ctype;
309 
310  typedef FieldVector< ctype, mydim > LocalCoordinate;
311  typedef FieldVector< ctype, cdim > GlobalCoordinate;
312 
314  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
315 
317  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
318 
320  ctype volume () const
321  {
322  GeometryType type = asImp().type();
323 
324  // get corresponding reference element
325  const ReferenceElement< ctype , mydim > & refElement =
326  ReferenceElements< ctype, mydim >::general(type);
327 
328  LocalCoordinate localBaryCenter ( 0 );
329  // calculate local bary center
330  const int corners = refElement.size(0,0,mydim);
331  for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
332  localBaryCenter *= (ctype) (1.0/corners);
333 
334  // volume is volume of reference element times integrationElement
335  return refElement.volume() * asImp().integrationElement(localBaryCenter);
336  }
337 
339  GlobalCoordinate center () const
340  {
341  GeometryType type = asImp().type();
342 
343  // get corresponding reference element
344  const ReferenceElement< ctype , mydim > & refElement =
345  ReferenceElements< ctype, mydim >::general(type);
346 
347  // center is (for now) the centroid of the reference element mapped to
348  // this geometry.
349  return asImp().global(refElement.position(0,0));
350  }
351 
352  private:
354  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
355  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
356  }; // end GeometryDefault
357 
358  template<int cdim, class GridImp, template<int,int,class> class GeometryImp>
359  class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
360  {
361  // my dimension
362  enum { mydim = 0 };
363 
364  public:
365  static const int mydimension = mydim;
366  static const int coorddimension = cdim;
367 
368  // save typing
369  typedef typename GridImp::ctype ctype;
370 
371  typedef FieldVector< ctype, mydim > LocalCoordinate;
372  typedef FieldVector< ctype, cdim > GlobalCoordinate;
373 
375  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
376 
378  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
379 
381  FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
382  {
383  return asImp().corner(0);
384  }
385 
387  FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
388  {
389  return FieldVector<ctype, mydim>();
390  }
391 
393  ctype volume () const
394  {
395  return 1.0;
396  }
397 
399  FieldVector<ctype, cdim> center () const
400  {
401  return asImp().corner(0);
402  }
403 
404  private:
405  // Barton-Nackman trick
406  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
407  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
408  }; // end GeometryDefault
409 
410 } // namespace Dune
411 
412 #endif // DUNE_GRID_GEOMETRY_HH
FieldVector< ctype, cdim > global(const FieldVector< ctype, mydim > &local) const
return the only coordinate
Definition: common/geometry.hh:381
Geometry(const Implementation &impl)
copy constructor from implementation
Definition: common/geometry.hh:276
FieldVector< ctype, mydim > LocalCoordinate
type of local coordinates
Definition: common/geometry.hh:101
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:378
Definition: common/geometry.hh:24
ctype volume() const
return volume of geometry
Definition: common/geometry.hh:208
Wrapper class for geometries.
Definition: common/geometry.hh:65
Definition: common/geometry.hh:87
FieldVector< ctype, cdim > center() const
return center of the geometry
Definition: common/geometry.hh:399
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:317
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:375
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
GeometryImp< mydim, cdim, GridImp > Implementation
Definition: common/geometry.hh:78
Include standard header files.
Definition: agrid.hh:59
bool affine() const
Return true if the geometry mapping is affine and false otherwise.
Definition: common/geometry.hh:134
Implementation realGeometry
Definition: common/geometry.hh:288
Definition: common/geometry.hh:96
GlobalCoordinate center() const
return center of geometry
Definition: common/geometry.hh:223
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
Return inverse of transposed of Jacobian.
Definition: common/geometry.hh:265
GridImp::ctype ctype
Definition: common/geometry.hh:369
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
Return the transposed of the Jacobian.
Definition: common/geometry.hh:239
GlobalCoordinate corner(int i) const
Obtain a corner of the geometry.
Definition: common/geometry.hh:156
Definition: common/geometry.hh:91
ctype integrationElement(const LocalCoordinate &local) const
Return the factor appearing in the integral transformation formula.
Definition: common/geometry.hh:202
Implementation::JacobianTransposed JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:126
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:371
Default implementation for class Geometry.
Definition: common/geometry.hh:301
FieldVector< ctype, mydim > local(const FieldVector< ctype, cdim > &) const
return empty vector
Definition: common/geometry.hh:387
GeometryType type() const
Return the name of the reference element. The type can be used to access the Dune::ReferenceElement.
Definition: common/geometry.hh:131
Implementation::JacobianInverseTransposed JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:115
const Implementation & impl() const
return reference to the implementation
Definition: common/geometry.hh:81
GridImp::ctype ctype
Definition: common/geometry.hh:308
GridImp::ctype ctype
define type used for coordinates in grid module
Definition: common/geometry.hh:98
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the map .
Definition: common/geometry.hh:165
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:314
ctype volume() const
return volume of the geometry
Definition: common/geometry.hh:320
int corners() const
Return the number of corners of the reference element.
Definition: common/geometry.hh:142
Definition: common/geometry.hh:89
FieldVector< ctype, cdim > GlobalCoordinate
type of the global coordinates
Definition: common/geometry.hh:104
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:372
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:310
static const int coorddimension
Definition: common/geometry.hh:305
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:311
static const int mydimension
Definition: common/geometry.hh:304
LocalCoordinate local(const GlobalCoordinate &global) const
Evaluate the inverse map .
Definition: common/geometry.hh:174
ctype volume() const
return volume of the geometry
Definition: common/geometry.hh:393
GlobalCoordinate center() const
return center of the geometry
Definition: common/geometry.hh:339