dune-grid  2.2.1
common/geometry.hh
Go to the documentation of this file.
1 #ifndef DUNE_GRID_GEOMETRY_HH
2 #define DUNE_GRID_GEOMETRY_HH
3 
8 #include <cassert>
9 
10 #include <dune/common/fmatrix.hh>
11 #include <dune/common/exceptions.hh>
12 #include <dune/common/typetraits.hh>
13 
14 #include <dune/geometry/referenceelements.hh>
15 #include <dune/geometry/genericgeometry/conversion.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 
116 #if 0
117 
118  Implementation &impl () { return realGeometry; }
119 #endif
120 
121  const Implementation &impl () const { return realGeometry; }
122 
123  public:
125  enum { dimension=GridImp::dimension };
127  enum { mydimension=mydim };
129  enum { coorddimension=cdim };
130 
132  enum { dimensionworld=GridImp::dimensionworld };
134  typedef typename GridImp::ctype ctype;
135 
137  typedef FieldVector<ctype, mydim> LocalCoordinate;
138 
140  typedef FieldVector< ctype, cdim > GlobalCoordinate;
141 
143  typedef FieldMatrix<ctype,cdim,mydim> Jacobian;
144  //typedef typename Implementation :: Jacobian Jacobian;
145 
147  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
148  //typedef typename Implementation :: JacobianTransposed JacobianTransposed;
149 
153  GeometryType type () const { return impl().type(); }
154 
156  bool affine() const { return impl().affine(); }
157 
164  int corners () const { return impl().corners(); }
165 
178  GlobalCoordinate corner ( int i ) const
179  {
180  return impl().corner( i );
181  }
182 
188  {
189  return impl().global( local );
190  }
191 
197  {
198  return impl().local( global );
199  }
200 
225  {
226  return impl().integrationElement( local );
227  }
228 
230  ctype volume () const
231  {
232  return impl().volume();
233  }
234 
246  {
247  return impl().center();
248  }
249 
259  const JacobianTransposed &
261  {
262  return impl().jacobianTransposed( local );
263  }
264 
285  {
286  return impl().jacobianInverseTransposed(local);
287  }
288 
290  explicit Geometry ( const Implementation &impl )
291  : realGeometry( impl )
292  {
293  deprecationWarning ( integral_constant< bool, storeReference >() );
294  }
295 
296  private:
298  const Geometry &operator= ( const Geometry &rhs );
299 
300  void DUNE_DEPRECATED_MSG( "This Dune::Geometry is still a reference to its implementation." )
301  deprecationWarning ( integral_constant< bool, true > ) {}
302 
303  void
304  deprecationWarning ( integral_constant< bool, false > ) {}
305 
306  protected:
308 
309  typename SelectType< storeReference, const Implementation &, Implementation >::Type realGeometry;
310  };
311 
312 
313 
314  //************************************************************************
315  // GEOMETRY Default Implementations
316  //*************************************************************************
317  //
318  // --GeometryDefault
319  //
321  template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>
323  {
324  public:
325  static const int mydimension = mydim;
326  static const int coorddimension = cdim;
327 
328  // save typing
329  typedef typename GridImp::ctype ctype;
330 
331  typedef FieldVector< ctype, mydim > LocalCoordinate;
332  typedef FieldVector< ctype, cdim > GlobalCoordinate;
333 
335  typedef FieldMatrix<ctype,cdim,mydim> Jacobian;
336 
338  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
339 
341  ctype volume () const
342  {
343  GeometryType type = asImp().type();
344 
345  // get corresponding reference element
346  const GenericReferenceElement< ctype , mydim > & refElement =
347  GenericReferenceElements< ctype, mydim >::general(type);
348 
349  LocalCoordinate localBaryCenter ( 0 );
350  // calculate local bary center
351  const int corners = refElement.size(0,0,mydim);
352  for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
353  localBaryCenter *= (ctype) (1.0/corners);
354 
355  // volume is volume of reference element times integrationElement
356  return refElement.volume() * asImp().integrationElement(localBaryCenter);
357  }
358 
361  {
362  GeometryType type = asImp().type();
363 
364  // get corresponding reference element
365  const GenericReferenceElement< ctype , mydim > & refElement =
366  GenericReferenceElements< ctype, mydim >::general(type);
367 
368  // center is (for now) the centroid of the reference element mapped to
369  // this geometry.
370  return asImp().global(refElement.position(0,0));
371  }
372 
373  private:
375  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
376  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
377  }; // end GeometryDefault
378 
379  template<int cdim, class GridImp, template<int,int,class> class GeometryImp>
380  class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
381  {
382  // my dimension
383  enum { mydim = 0 };
384 
385  public:
386  static const int mydimension = mydim;
387  static const int coorddimension = cdim;
388 
389  // save typing
390  typedef typename GridImp::ctype ctype;
391 
392  typedef FieldVector< ctype, mydim > LocalCoordinate;
393  typedef FieldVector< ctype, cdim > GlobalCoordinate;
394 
396  typedef FieldMatrix<ctype,cdim,mydim> Jacobian;
397 
399  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
400 
402  FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
403  {
404  return asImp().corner(0);
405  }
406 
408  FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
409  {
410  return FieldVector<ctype, mydim>();
411  }
412 
414  bool checkInside (const FieldVector<ctype, mydim>& ) const
415  {
416  return true;
417  }
418 
420  ctype volume () const
421  {
422  return 1.0;
423  }
424 
426  FieldVector<ctype, cdim> center () const
427  {
428  return asImp().corner(0);
429  }
430 
431  private:
432  // Barton-Nackman trick
433  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
434  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
435  }; // end GeometryDefault
436 
437 } // namespace Dune
438 
439 #endif // DUNE_GRID_GEOMETRY_HH