dune-grid  2.3beta2
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 >
24  class GridDefaultImplementation;
25 
26 
27 
28  namespace FacadeOptions
29  {
30 
33 
45  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp >
47  {
49  static const bool v = true;
50  };
51 
52  }
53 
54 
55 
56  //*****************************************************************************
57  //
58  // Geometry
59  // forwards the interface to the implementation
60  //
61  //*****************************************************************************
62 
100  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp >
101  class Geometry
102  {
103  #if DUNE_GRID_EXPERIMENTAL_GRID_EXTENSIONS
104  public:
105  #else
106  protected:
107  // give the GridDefaultImplementation class access to the realImp
108  friend class GridDefaultImplementation<
109  GridImp::dimension, GridImp::dimensionworld,
110  typename GridImp::ctype,
111  typename GridImp::GridFamily> ;
112  #endif
113  // type of underlying implementation, for internal use only
114  typedef GeometryImp< mydim, cdim, GridImp > Implementation;
115 
117  const Implementation &impl () const { return realGeometry; }
118 
119  public:
121  enum { dimension=GridImp::dimension };
123  enum { mydimension=mydim };
125  enum { coorddimension=cdim };
126 
128  enum { dimensionworld=GridImp::dimensionworld };
130  typedef typename GridImp::ctype ctype;
131 
133  typedef FieldVector<ctype, mydim> LocalCoordinate;
134 
136  typedef FieldVector< ctype, cdim > GlobalCoordinate;
137 
139  typedef typename Implementation::JacobianInverseTransposed JacobianInverseTransposed;
140 
143  typedef JacobianInverseTransposed Jacobian DUNE_DEPRECATED_MSG ( "type Geometry::Jacobian is deprecated, use Geometry::JacobianInverseTransposed instead." );
144 
146  typedef typename Implementation::JacobianTransposed JacobianTransposed;
147 
151  GeometryType type () const { return impl().type(); }
152 
154  bool affine() const { return impl().affine(); }
155 
162  int corners () const { return impl().corners(); }
163 
176  GlobalCoordinate corner ( int i ) const
177  {
178  return impl().corner( i );
179  }
180 
186  {
187  return impl().global( local );
188  }
189 
195  {
196  return impl().local( global );
197  }
198 
223  {
224  return impl().integrationElement( local );
225  }
226 
228  ctype volume () const
229  {
230  return impl().volume();
231  }
232 
244  {
245  return impl().center();
246  }
247 
257  const JacobianTransposed &
259  {
260  return impl().jacobianTransposed( local );
261  }
262 
284  {
285  return impl().jacobianInverseTransposed(local);
286  }
287 
289  explicit Geometry ( const Implementation &impl )
290  : realGeometry( impl )
291  {
292  deprecationWarning ( integral_constant< bool, storeReference >() );
293  }
294 
295  private:
297  const Geometry &operator= ( const Geometry &rhs );
298 
299  void DUNE_DEPRECATED_MSG( "This Dune::Geometry is still a reference to its implementation." )
300  deprecationWarning ( integral_constant< bool, true > ) {}
301 
302  void
303  deprecationWarning ( integral_constant< bool, false > ) {}
304 
305  protected:
307 
308  typename conditional< storeReference, const Implementation &, Implementation >::type realGeometry;
309  };
310 
311 
312 
313  //************************************************************************
314  // GEOMETRY Default Implementations
315  //*************************************************************************
316  //
317  // --GeometryDefault
318  //
320  template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>
322  {
323  public:
324  static const int mydimension = mydim;
325  static const int coorddimension = cdim;
326 
327  // save typing
328  typedef typename GridImp::ctype ctype;
329 
330  typedef FieldVector< ctype, mydim > LocalCoordinate;
331  typedef FieldVector< ctype, cdim > GlobalCoordinate;
332 
334  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
335 
337  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
338 
340  ctype volume () const
341  {
342  GeometryType type = asImp().type();
343 
344  // get corresponding reference element
345  const ReferenceElement< ctype , mydim > & refElement =
346  ReferenceElements< ctype, mydim >::general(type);
347 
348  LocalCoordinate localBaryCenter ( 0 );
349  // calculate local bary center
350  const int corners = refElement.size(0,0,mydim);
351  for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
352  localBaryCenter *= (ctype) (1.0/corners);
353 
354  // volume is volume of reference element times integrationElement
355  return refElement.volume() * asImp().integrationElement(localBaryCenter);
356  }
357 
360  {
361  GeometryType type = asImp().type();
362 
363  // get corresponding reference element
364  const ReferenceElement< ctype , mydim > & refElement =
365  ReferenceElements< ctype, mydim >::general(type);
366 
367  // center is (for now) the centroid of the reference element mapped to
368  // this geometry.
369  return asImp().global(refElement.position(0,0));
370  }
371 
372  private:
374  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
375  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
376  }; // end GeometryDefault
377 
378  template<int cdim, class GridImp, template<int,int,class> class GeometryImp>
379  class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
380  {
381  // my dimension
382  enum { mydim = 0 };
383 
384  public:
385  static const int mydimension = mydim;
386  static const int coorddimension = cdim;
387 
388  // save typing
389  typedef typename GridImp::ctype ctype;
390 
391  typedef FieldVector< ctype, mydim > LocalCoordinate;
392  typedef FieldVector< ctype, cdim > GlobalCoordinate;
393 
395  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
396 
398  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
399 
401  FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
402  {
403  return asImp().corner(0);
404  }
405 
407  FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
408  {
409  return FieldVector<ctype, mydim>();
410  }
411 
413  bool checkInside (const FieldVector<ctype, mydim>& ) const
414  {
415  return true;
416  }
417 
419  ctype volume () const
420  {
421  return 1.0;
422  }
423 
425  FieldVector<ctype, cdim> center () const
426  {
427  return asImp().corner(0);
428  }
429 
430  private:
431  // Barton-Nackman trick
432  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
433  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
434  }; // end GeometryDefault
435 
436 } // namespace Dune
437 
438 #endif // DUNE_GRID_GEOMETRY_HH