dune-grid
2.2.0
|
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/exceptions.hh> 00012 #include <dune/common/typetraits.hh> 00013 00014 #include <dune/geometry/referenceelements.hh> 00015 #include <dune/geometry/genericgeometry/conversion.hh> 00016 00017 namespace Dune 00018 { 00019 00020 // External Forward Declarations 00021 // ----------------------------- 00022 00023 template< int dim, int dimworld, class ct, class GridFamily > 00024 class GridDefaultImplementation; 00025 00026 00027 00028 namespace FacadeOptions 00029 { 00030 00033 00045 template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp > 00046 struct StoreGeometryReference 00047 { 00049 static const bool v = true; 00050 }; 00051 00052 } 00053 00054 00055 00056 //***************************************************************************** 00057 // 00058 // Geometry 00059 // forwards the interface to the implementation 00060 // 00061 //***************************************************************************** 00062 00100 template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp > 00101 class Geometry 00102 { 00103 #if DUNE_GRID_EXPERIMENTAL_GRID_EXTENSIONS 00104 public: 00105 #else 00106 protected: 00107 // give the GridDefaultImplementation class access to the realImp 00108 friend class GridDefaultImplementation< 00109 GridImp::dimension, GridImp::dimensionworld, 00110 typename GridImp::ctype, 00111 typename GridImp::GridFamily> ; 00112 #endif 00113 // type of underlying implementation, for internal use only 00114 typedef GeometryImp< mydim, cdim, GridImp > Implementation; 00115 00116 #if 0 00117 00118 Implementation &impl () { return realGeometry; } 00119 #endif 00120 00121 const Implementation &impl () const { return realGeometry; } 00122 00123 public: 00125 enum { dimension=GridImp::dimension }; 00127 enum { mydimension=mydim }; 00129 enum { coorddimension=cdim }; 00130 00132 enum { dimensionworld=GridImp::dimensionworld }; 00134 typedef typename GridImp::ctype ctype; 00135 00137 typedef FieldVector<ctype, mydim> LocalCoordinate; 00138 00140 typedef FieldVector< ctype, cdim > GlobalCoordinate; 00141 00143 typedef FieldMatrix<ctype,cdim,mydim> Jacobian; 00144 //typedef typename Implementation :: Jacobian Jacobian; 00145 00147 typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed; 00148 //typedef typename Implementation :: JacobianTransposed JacobianTransposed; 00149 00153 GeometryType type () const { return impl().type(); } 00154 00156 bool affine() const { return impl().affine(); } 00157 00164 int corners () const { return impl().corners(); } 00165 00178 GlobalCoordinate corner ( int i ) const 00179 { 00180 return impl().corner( i ); 00181 } 00182 00187 GlobalCoordinate global (const LocalCoordinate& local) const 00188 { 00189 return impl().global( local ); 00190 } 00191 00196 LocalCoordinate local (const GlobalCoordinate& global) const 00197 { 00198 return impl().local( global ); 00199 } 00200 00224 ctype integrationElement (const LocalCoordinate& local) const 00225 { 00226 return impl().integrationElement( local ); 00227 } 00228 00230 ctype volume () const 00231 { 00232 return impl().volume(); 00233 } 00234 00245 GlobalCoordinate center () const 00246 { 00247 return impl().center(); 00248 } 00249 00259 const JacobianTransposed & 00260 jacobianTransposed ( const LocalCoordinate& local ) const 00261 { 00262 return impl().jacobianTransposed( local ); 00263 } 00264 00284 const Jacobian& jacobianInverseTransposed (const LocalCoordinate& local) const 00285 { 00286 return impl().jacobianInverseTransposed(local); 00287 } 00288 00290 explicit Geometry ( const Implementation &impl ) 00291 : realGeometry( impl ) 00292 { 00293 deprecationWarning ( integral_constant< bool, storeReference >() ); 00294 } 00295 00296 private: 00298 const Geometry &operator= ( const Geometry &rhs ); 00299 00300 void DUNE_DEPRECATED_MSG( "This Dune::Geometry is still a reference to its implementation." ) 00301 deprecationWarning ( integral_constant< bool, true > ) {} 00302 00303 void 00304 deprecationWarning ( integral_constant< bool, false > ) {} 00305 00306 protected: 00307 static const bool storeReference = FacadeOptions::StoreGeometryReference< mydim, cdim, GridImp, GeometryImp >::v; 00308 00309 typename SelectType< storeReference, const Implementation &, Implementation >::Type realGeometry; 00310 }; 00311 00312 00313 00314 //************************************************************************ 00315 // GEOMETRY Default Implementations 00316 //************************************************************************* 00317 // 00318 // --GeometryDefault 00319 // 00321 template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp> 00322 class GeometryDefaultImplementation 00323 { 00324 public: 00325 static const int mydimension = mydim; 00326 static const int coorddimension = cdim; 00327 00328 // save typing 00329 typedef typename GridImp::ctype ctype; 00330 00331 typedef FieldVector< ctype, mydim > LocalCoordinate; 00332 typedef FieldVector< ctype, cdim > GlobalCoordinate; 00333 00335 typedef FieldMatrix<ctype,cdim,mydim> Jacobian; 00336 00338 typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed; 00339 00341 ctype volume () const 00342 { 00343 GeometryType type = asImp().type(); 00344 00345 // get corresponding reference element 00346 const GenericReferenceElement< ctype , mydim > & refElement = 00347 GenericReferenceElements< ctype, mydim >::general(type); 00348 00349 LocalCoordinate localBaryCenter ( 0 ); 00350 // calculate local bary center 00351 const int corners = refElement.size(0,0,mydim); 00352 for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim); 00353 localBaryCenter *= (ctype) (1.0/corners); 00354 00355 // volume is volume of reference element times integrationElement 00356 return refElement.volume() * asImp().integrationElement(localBaryCenter); 00357 } 00358 00360 GlobalCoordinate center () const 00361 { 00362 GeometryType type = asImp().type(); 00363 00364 // get corresponding reference element 00365 const GenericReferenceElement< ctype , mydim > & refElement = 00366 GenericReferenceElements< ctype, mydim >::general(type); 00367 00368 // center is (for now) the centroid of the reference element mapped to 00369 // this geometry. 00370 return asImp().global(refElement.position(0,0)); 00371 } 00372 00373 private: 00375 GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);} 00376 const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);} 00377 }; // end GeometryDefault 00378 00379 template<int cdim, class GridImp, template<int,int,class> class GeometryImp> 00380 class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp> 00381 { 00382 // my dimension 00383 enum { mydim = 0 }; 00384 00385 public: 00386 static const int mydimension = mydim; 00387 static const int coorddimension = cdim; 00388 00389 // save typing 00390 typedef typename GridImp::ctype ctype; 00391 00392 typedef FieldVector< ctype, mydim > LocalCoordinate; 00393 typedef FieldVector< ctype, cdim > GlobalCoordinate; 00394 00396 typedef FieldMatrix<ctype,cdim,mydim> Jacobian; 00397 00399 typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed; 00400 00402 FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const 00403 { 00404 return asImp().corner(0); 00405 } 00406 00408 FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const 00409 { 00410 return FieldVector<ctype, mydim>(); 00411 } 00412 00414 bool checkInside (const FieldVector<ctype, mydim>& ) const 00415 { 00416 return true; 00417 } 00418 00420 ctype volume () const 00421 { 00422 return 1.0; 00423 } 00424 00426 FieldVector<ctype, cdim> center () const 00427 { 00428 return asImp().corner(0); 00429 } 00430 00431 private: 00432 // Barton-Nackman trick 00433 GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);} 00434 const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);} 00435 }; // end GeometryDefault 00436 00437 } // namespace Dune 00438 00439 #endif // DUNE_GRID_GEOMETRY_HH