common/geometry.hh

Go to the documentation of this file.
00001 #ifndef DUNE_GRID_GEOMETRY_HH
00002 #define DUNE_GRID_GEOMETRY_HH
00003 
00008 #include <cassert>
00009 
00010 #include <dune/common/fmatrix.hh>
00011 #include <dune/common/helpertemplates.hh>
00012 #include <dune/common/exceptions.hh>
00013 
00014 #include "referenceelements.hh"
00015 
00016 namespace Dune
00017 {
00018 
00019  // forward deklaration for volume implementation
00020  template<typename ctype, int dim> class ReferenceElement;
00021  template<typename ctype, int dim> class ReferenceElements;
00022   
00023 //*****************************************************************************
00024 //
00025 // Geometry
00026 // forwards the interface to the implementation
00027 //
00028 //*****************************************************************************
00029 
00064 template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>  
00065 class Geometry {
00066   // save typing
00067   typedef typename GridImp::ctype ct;
00068 protected:
00069   GeometryImp<mydim,cdim,GridImp> realGeometry;
00070   
00071   // type of underlying implementation, for internal use only 
00072   typedef GeometryImp<mydim,cdim,GridImp> ImplementationType;
00073 public:
00075   enum { dimension=GridImp::dimension  };
00077   enum { mydimension=mydim  };
00079   enum { coorddimension=cdim  };
00080 
00082   enum { dimensionworld=GridImp::dimensionworld  };
00084   typedef ct ctype;
00085   
00086 
00090   GeometryType type () const { return realGeometry.type(); };
00091 
00097   int corners () const { return realGeometry.corners(); };
00098 
00104   const FieldVector<ct, cdim>& operator[] (int i) const
00105     {
00106       return realGeometry[i];
00107     }
00108 
00113   FieldVector<ct, cdim> global (const FieldVector<ct, mydim>& local) const
00114     {
00115       return realGeometry.global(local);
00116     }
00117 
00122   FieldVector<ct, mydim> local (const FieldVector<ct, cdim>& global) const
00123     {
00124       return realGeometry.local(global);
00125     }
00126 
00128   bool checkInside (const FieldVector<ct, mydim>& local) const
00129     {
00130       return realGeometry.checkInside(local);
00131     }
00132 
00156   ct integrationElement (const FieldVector<ct, mydim>& local) const
00157   {
00158     return realGeometry.integrationElement(local);
00159   }
00160  
00162   ct volume () const
00163   {
00164     return realGeometry.volume();
00165   }
00166 
00183   const FieldMatrix<ct,mydim,mydim>& jacobianInverseTransposed (const FieldVector<ct, mydim>& local) const
00184   {
00185   return realGeometry.jacobianInverseTransposed(local);
00186   }
00187 
00188 public:    
00190   explicit Geometry(const GeometryImp<mydim,cdim,GridImp> & e) : realGeometry(e) {};
00191   
00192 protected:
00193   // give the GridDefaultImplementation class access to the realImp 
00194   friend class GridDefaultImplementation< 
00195             GridImp::dimension, GridImp::dimensionworld,
00196             typename GridImp::ctype,
00197             typename GridImp::GridFamily> ;
00198 
00200   GeometryImp<mydim,cdim,GridImp> & getRealImp() { return realGeometry; }
00202   const GeometryImp<mydim,cdim,GridImp> & getRealImp() const { return realGeometry; }
00203  
00204 protected:  
00206   Geometry(const Geometry& rhs) : realGeometry(rhs.realGeometry) {};
00208   Geometry & operator = (const Geometry& rhs) {
00209     realGeometry = rhs.realGeometry;
00210     return *this;
00211   };
00212 };
00213 
00214 
00215 //************************************************************************
00216 // GEOMETRY Default Implementations
00217 //*************************************************************************
00218 //
00219 // --GeometryDefault 
00220 // 
00222 template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>  
00223 class GeometryDefaultImplementation
00224 {
00225 public:
00226   // save typing
00227   typedef typename GridImp::ctype ctype;
00228  
00230   ctype volume () const 
00231   {
00232     GeometryType type = asImp().type();
00233 
00234     // get corresponding reference element
00235     const ReferenceElement< ctype , mydim > & refElement =
00236       ReferenceElements< ctype, mydim >::general(type);
00237 
00238     FieldVector<ctype,mydim> localBaryCenter (0.0);
00239     // calculate local bary center 
00240     const int corners = refElement.size(0,0,mydim);
00241     for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
00242     localBaryCenter *= (ctype) (1.0/corners);
00243 
00244     // volume is volume of reference element times integrationElement
00245     return refElement.volume() * asImp().integrationElement(localBaryCenter);
00246   }
00247   
00248 private:
00250   GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
00251   const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
00252 }; // end GeometryDefault 
00253 
00254 template<int cdim, class GridImp, template<int,int,class> class GeometryImp>  
00255 class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
00256 {
00257   // my dimension 
00258   enum { mydim = 0 };
00259 public:
00260   // save typing
00261   typedef typename GridImp::ctype ctype;
00262 
00264   FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
00265   {
00266     return asImp()[0];
00267   }
00268 
00270   FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
00271   {
00272     return FieldVector<ctype, mydim>();
00273   }
00274  
00276   bool checkInside (const FieldVector<ctype, mydim>& ) const
00277   {
00278     return true; 
00279   }
00280  
00282   ctype volume () const { return 1.0; }
00283   
00284 private:
00286   GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
00287   const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
00288 }; // end GeometryDefault 
00289 
00290 }
00291 #endif // DUNE_GRID_GEOMETRY_HH

Generated on 6 Nov 2008 with Doxygen (ver 1.5.6) [logfile].