dune-grid
2.2.0
|
00001 #ifndef DUNE_ALU3DGRIDMAPPINGS_HH 00002 #define DUNE_ALU3DGRIDMAPPINGS_HH 00003 00004 // System includes 00005 #include <limits> 00006 #include <cmath> 00007 00008 // Dune includes 00009 #include <dune/common/fvector.hh> 00010 #include <dune/common/fmatrix.hh> 00011 #include <dune/common/exceptions.hh> 00012 00013 // Local includes 00014 #include "alu3dinclude.hh" 00015 00016 namespace Dune { 00017 00018 static const alu3d_ctype ALUnumericEpsilon = 10.0 * std::numeric_limits< alu3d_ctype >::epsilon(); 00019 00020 template<int mydim, int coorddim, class GridImp> 00021 class ALU3dGridGeometry; 00022 00023 template< ALU3dGridElementType, class > 00024 class ALU3dGrid; 00025 00028 class TrilinearMapping 00029 { 00030 public: 00031 typedef alu3d_ctype double_t[3]; 00032 typedef FieldVector<alu3d_ctype, 3> coord_t; 00033 typedef FieldMatrix<alu3d_ctype, 3, 3> mat_t; 00034 private: 00035 static const double _epsilon ; 00036 00037 // the internal mapping 00038 alu3d_ctype a [8][3] ; 00039 mat_t Df; 00040 mat_t Dfi; 00041 mat_t invTransposed_; 00042 alu3d_ctype DetDf ; 00043 00044 bool calcedDet_; 00045 bool calcedLinear_; 00046 bool calcedInv_; 00047 bool affine_; 00048 00049 void linear (const alu3d_ctype, const alu3d_ctype, const alu3d_ctype) ; 00050 void linear (const coord_t&) ; 00051 void inverse (const coord_t&) ; 00052 public : 00053 TrilinearMapping (const coord_t&, const coord_t&, 00054 const coord_t&, const coord_t&, 00055 const coord_t&, const coord_t&, 00056 const coord_t&, const coord_t&); 00057 00058 // only to call from geometry class 00059 TrilinearMapping () {} 00060 00061 TrilinearMapping (const TrilinearMapping &) ; 00062 00063 ~TrilinearMapping () {} 00064 alu3d_ctype det (const coord_t&) ; 00065 const mat_t& jacobianInverseTransposed(const coord_t&); 00066 const mat_t& jacobianTransposed(const coord_t&); 00067 void map2world (const coord_t&, coord_t&) const ; 00068 void map2world (const alu3d_ctype , const alu3d_ctype , const alu3d_ctype , 00069 coord_t&) const ; 00070 void world2map (const coord_t&, coord_t&) ; 00071 00072 template <class vector_t> 00073 void buildMapping(const vector_t&, const vector_t&, 00074 const vector_t&, const vector_t&, 00075 const vector_t&, const vector_t&, 00076 const vector_t&, const vector_t&); 00077 00078 // returns true if mapping is affine 00079 inline bool affine () const { return affine_; } 00080 }; 00081 00083 // NOTE: this class is different to the BilinearSurfaceMapping in 00084 // ALUGrid, for example the reference elements differ 00085 // here we have [0,1]^2 and in ALUGrid its [-1,1]^2 00086 // also the point numbering is different 00087 class SurfaceNormalCalculator 00088 { 00089 public: 00090 // our coordinate types 00091 typedef FieldVector<alu3d_ctype, 3> coord3_t; 00092 typedef FieldVector<alu3d_ctype, 2> coord2_t; 00093 00094 // type of coordinate vectors from elements 00095 typedef alu3d_ctype double3_t[3]; 00096 protected: 00097 00098 alu3d_ctype _n [3][3] ; 00099 00100 static const double _epsilon ; 00101 00102 bool _affine; 00103 00104 public : 00106 SurfaceNormalCalculator(); 00107 00108 SurfaceNormalCalculator (const SurfaceNormalCalculator &) ; 00109 ~SurfaceNormalCalculator () {} 00110 00111 // returns true if mapping is affine 00112 inline bool affine () const { return _affine ; } 00113 00114 // calcuates normal 00115 void normal(const coord2_t&, coord3_t&) const ; 00116 void normal(const alu3d_ctype, const alu3d_ctype, coord3_t&)const; 00117 00118 void negativeNormal(const coord2_t&, coord3_t&) const ; 00119 void negativeNormal(const alu3d_ctype, const alu3d_ctype, coord3_t&)const; 00120 00121 public: 00122 // builds _b and _n, called from the constructors 00123 // public because also used in faceutility 00124 template <class vector_t> 00125 void buildMapping (const vector_t & , const vector_t & , 00126 const vector_t & , const vector_t & ); 00127 protected: 00128 // builds _b and _n, called from the constructors 00129 // public because also used in faceutility 00130 template <class vector_t> 00131 void buildMapping (const vector_t & , const vector_t & , 00132 const vector_t & , const vector_t & , 00133 alu3d_ctype (&_b)[4][3] ); 00134 } ; 00135 00136 00138 // NOTE: this class is different to the BilinearSurfaceMapping in 00139 // ALUGrid, for example the reference elements differ 00140 // here we have [0,1]^2 and in ALUGrid its [-1,1]^2 00141 // also the point numbering is different 00142 class BilinearSurfaceMapping : public SurfaceNormalCalculator 00143 { 00144 protected: 00145 typedef SurfaceNormalCalculator BaseType; 00146 00147 using BaseType :: _n; 00148 static const double _epsilon; 00149 00150 // our coordinate types 00151 typedef FieldVector<alu3d_ctype, 3> coord3_t; 00152 typedef FieldVector<alu3d_ctype, 2> coord2_t; 00153 00154 // type of coordinate vectors from elements 00155 typedef alu3d_ctype double3_t[3]; 00156 00157 // type for helper matrices 00158 typedef FieldMatrix<alu3d_ctype,3,3> mat3_t; 00159 00160 // type for inverse matrices 00161 typedef FieldMatrix<alu3d_ctype,2,3> matrix_t; 00162 00163 // type for inverse matrices 00164 typedef FieldMatrix<alu3d_ctype,3,2> inv_t; 00165 00166 alu3d_ctype _b [4][3] ; 00167 00168 mutable mat3_t Df,Dfi; 00169 mutable inv_t invTransposed_; 00170 mutable matrix_t matrix_; 00171 mutable alu3d_ctype DetDf; 00172 00173 mutable coord3_t normal_; 00174 mutable coord3_t tmp_; 00175 00176 mutable bool _calcedInv; 00177 mutable bool _calcedTransposed; 00178 mutable bool _calcedMatrix; 00179 00180 public : 00182 BilinearSurfaceMapping (); 00183 00185 BilinearSurfaceMapping (const coord3_t&, const coord3_t&, 00186 const coord3_t&, const coord3_t&) ; 00188 BilinearSurfaceMapping (const double3_t &, const double3_t &, 00189 const double3_t &, const double3_t &) ; 00190 BilinearSurfaceMapping (const BilinearSurfaceMapping &) ; 00191 ~BilinearSurfaceMapping () {} 00192 00193 void inverse (const coord3_t&) const; 00194 const inv_t& jacobianInverseTransposed(const coord2_t&) const ; 00195 00196 const matrix_t& jacobianTransposed(const coord2_t&) const ; 00197 00198 // calculates determinant of face mapping using the normal 00199 alu3d_ctype det(const coord2_t&) const; 00200 00201 // maps from local coordinates to global coordinates 00202 void world2map(const coord3_t &, coord2_t & ) const; 00203 00204 // maps form global coordinates to local (within reference element) 00205 // coordinates 00206 void map2world(const coord2_t&, coord3_t&) const ; 00207 void map2world(const alu3d_ctype ,const alu3d_ctype , coord3_t&) const ; 00208 00209 private: 00210 void map2worldnormal(const alu3d_ctype, const alu3d_ctype, const alu3d_ctype , coord3_t&)const; 00211 void map2worldlinear(const alu3d_ctype, const alu3d_ctype, const alu3d_ctype ) const; 00212 00213 public: 00214 // builds _b and _n, called from the constructors 00215 // public because also used in faceutility 00216 template <class vector_t> 00217 void buildMapping (const vector_t & , const vector_t & , 00218 const vector_t & , const vector_t & ); 00219 } ; 00220 00221 00222 00224 template< int cdim > 00225 class BilinearMapping 00226 { 00227 public: 00228 typedef alu3d_ctype ctype; 00229 00230 typedef FieldVector< ctype, cdim > world_t; 00231 typedef FieldVector< ctype, 2 > map_t; 00232 00233 typedef FieldMatrix< ctype, 2, cdim > matrix_t; 00234 typedef FieldMatrix< ctype, cdim, 2 > inv_t; 00235 00236 protected: 00237 ctype _b [4][cdim]; 00238 00239 mutable ctype det_; 00240 mutable matrix_t matrix_; 00241 mutable inv_t invTransposed_; 00242 00243 mutable bool affine_; 00244 mutable bool calcedMatrix_; 00245 mutable bool calcedDet_; 00246 mutable bool calcedInv_; 00247 00248 public: 00249 BilinearMapping (); 00250 BilinearMapping ( const world_t &p0, const world_t &p1, 00251 const world_t &p2, const world_t &p3 ); 00252 BilinearMapping ( const ctype (&p0)[ cdim ], const ctype (&p1)[ cdim ], 00253 const ctype (&p2)[ cdim ], const ctype (&p3)[ cdim ] ); 00254 00255 bool affine () const; 00256 00257 void world2map ( const world_t &, map_t & ) const; 00258 void map2world ( const ctype x, const ctype y, world_t &w ) const; 00259 void map2world ( const map_t &, world_t & ) const; 00260 00261 ctype det ( const map_t & ) const; 00262 00263 const matrix_t &jacobianTransposed ( const map_t & ) const; 00264 const inv_t &jacobianInverseTransposed ( const map_t & ) const; 00265 00266 template< class vector_t > 00267 void buildMapping ( const vector_t &, const vector_t &, 00268 const vector_t &, const vector_t & ); 00269 00270 protected: 00271 static void multTransposedMatrix ( const matrix_t &, FieldMatrix< ctype, 2, 2 > & ); 00272 static void multMatrix ( const matrix_t &, const FieldMatrix< ctype, 2, 2 > &, inv_t & ); 00273 00274 void map2worldlinear ( const ctype, const ctype ) const; 00275 void inverse ( const map_t & ) const; 00276 }; 00277 00278 00279 00281 template< int cdim, int mydim > 00282 class LinearMapping 00283 { 00284 public: 00285 typedef alu3d_ctype ctype; 00286 00287 typedef ctype double_t[ cdim ]; 00288 00289 typedef FieldVector< ctype, cdim > world_t; 00290 typedef FieldVector< ctype, mydim > map_t; 00291 00292 typedef FieldMatrix< ctype, mydim, cdim > matrix_t; 00293 typedef FieldMatrix< ctype, cdim, mydim > inv_t; 00294 00295 protected: 00296 matrix_t _matrix; 00297 mutable inv_t _invTransposed; 00298 world_t _p0; 00299 00301 mutable ctype _det; 00302 00304 mutable bool _calcedInv; 00305 00307 mutable bool _calcedDet; 00308 00309 public: 00311 LinearMapping (); 00312 00314 LinearMapping (const LinearMapping &) ; 00315 00316 // returns true if mapping is affine (which is always true) 00317 inline bool affine () const { return true ; } 00318 00319 // return reference to transposed jacobian 00320 const matrix_t& jacobianTransposed(const map_t &) const ; 00321 00322 // return reference to transposed jacobian inverse 00323 const inv_t& jacobianInverseTransposed(const map_t &) const ; 00324 00325 // calculates determinant of mapping 00326 ctype det(const map_t&) const; 00327 00328 // maps from local coordinates to global coordinates 00329 void world2map(const world_t &, map_t &) const; 00330 00331 // maps form global coordinates to local (within reference element) 00332 // coordinates 00333 void map2world(const map_t &, world_t &) const ; 00334 00335 protected: 00336 // calculate inverse 00337 void inverse (const map_t&) const; 00338 00339 // calculate inverse one codim one entity 00340 void inverseCodimOne (const map_t&) const; 00341 00342 // calculate determinant 00343 void calculateDeterminant (const map_t&) const; 00344 00345 void multTransposedMatrix(const matrix_t& matrix, 00346 FieldMatrix<ctype, mydim, mydim>& result) const; 00347 00348 void multMatrix ( const matrix_t& A, 00349 const FieldMatrix< ctype, mydim, mydim> &B, 00350 inv_t& ret ) const ; 00351 00352 public: 00353 // builds _b and _n, called from the constructors 00354 // public because also used in faceutility 00355 template <class vector_t> 00356 void buildMapping (const vector_t & , const vector_t & , 00357 const vector_t & , const vector_t & ); 00358 00359 // builds _b and _n, called from the constructors 00360 // public because also used in faceutility 00361 template <class vector_t> 00362 void buildMapping (const vector_t & , const vector_t & , 00363 const vector_t & ); 00364 00365 // builds _b and _n, called from the constructors 00366 // public because also used in faceutility 00367 template <class vector_t> 00368 void buildMapping (const vector_t & , const vector_t & ); 00369 00370 template <class vector_t> 00371 void buildMapping (const vector_t & ); 00372 } ; 00373 00374 00376 // 00377 // NonConforming Mappings 00378 // 00380 00381 00384 template< ALU3dGridElementType type, class Comm > 00385 class NonConformingFaceMapping; 00386 00388 template< class Comm > 00389 struct NonConformingFaceMapping< tetra, Comm > 00390 { 00391 typedef FieldVector< alu3d_ctype, 3 > CoordinateType; 00392 typedef typename ALU3dImplTraits< tetra, Comm >::HfaceRuleType RefinementRuleType; 00393 00394 NonConformingFaceMapping ( RefinementRuleType rule, int nChild ) 00395 : rule_( rule ), nChild_( nChild ) 00396 {} 00397 00398 void child2parent ( const CoordinateType &childCoordinates, 00399 CoordinateType &parentCoordinates) const; 00400 00401 CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const; 00402 00403 private: 00404 void child2parentNosplit(const CoordinateType& childCoordinates, 00405 CoordinateType& parentCoordinates) const; 00406 void child2parentE01(const CoordinateType& childCoordinates, 00407 CoordinateType& parentCoordinates) const; 00408 void child2parentE12(const CoordinateType& childCoordinates, 00409 CoordinateType& parentCoordinates) const; 00410 void child2parentE20(const CoordinateType& childCoordinates, 00411 CoordinateType& parentCoordinates) const; 00412 void child2parentIso4(const CoordinateType& childCoordinates, 00413 CoordinateType& parentCoordinates) const; 00414 00415 RefinementRuleType rule_; 00416 int nChild_; 00417 }; 00418 00420 template< class Comm > 00421 struct NonConformingFaceMapping< hexa, Comm > 00422 { 00423 typedef FieldVector< alu3d_ctype, 2 > CoordinateType; 00424 typedef typename ALU3dImplTraits< hexa, Comm >::HfaceRuleType RefinementRuleType; 00425 00426 NonConformingFaceMapping ( RefinementRuleType rule, int nChild ) 00427 : rule_( rule ), nChild_( nChild ) 00428 {} 00429 00430 void child2parent ( const CoordinateType &childCoordinates, 00431 CoordinateType &parentCoordinates) const; 00432 00433 CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const; 00434 00435 private: 00436 void child2parentNosplit(const CoordinateType& childCoordinates, 00437 CoordinateType& parentCoordinates) const; 00438 void child2parentIso4(const CoordinateType& childCoordinates, 00439 CoordinateType& parentCoordinates) const; 00440 00441 RefinementRuleType rule_; 00442 int nChild_; 00443 }; 00444 00445 } // end namespace Dune 00446 00447 #if COMPILE_ALUGRID_INLINE 00448 #include "mappings_imp.cc" 00449 #endif 00450 #endif