dune-grid  2.2.0
mappings.hh
Go to the documentation of this file.
1 #ifndef DUNE_ALU3DGRIDMAPPINGS_HH
2 #define DUNE_ALU3DGRIDMAPPINGS_HH
3 
4 // System includes
5 #include <limits>
6 #include <cmath>
7 
8 // Dune includes
9 #include <dune/common/fvector.hh>
10 #include <dune/common/fmatrix.hh>
11 #include <dune/common/exceptions.hh>
12 
13 // Local includes
14 #include "alu3dinclude.hh"
15 
16 namespace Dune {
17 
18  static const alu3d_ctype ALUnumericEpsilon = 10.0 * std::numeric_limits< alu3d_ctype >::epsilon();
19 
20  template<int mydim, int coorddim, class GridImp>
21  class ALU3dGridGeometry;
22 
23  template< ALU3dGridElementType, class >
24  class ALU3dGrid;
25 
29  {
30  public:
31  typedef alu3d_ctype double_t[3];
32  typedef FieldVector<alu3d_ctype, 3> coord_t;
33  typedef FieldMatrix<alu3d_ctype, 3, 3> mat_t;
34  private:
35  static const double _epsilon ;
36 
37  // the internal mapping
38  alu3d_ctype a [8][3] ;
39  mat_t Df;
40  mat_t Dfi;
41  mat_t invTransposed_;
42  alu3d_ctype DetDf ;
43 
44  bool calcedDet_;
45  bool calcedLinear_;
46  bool calcedInv_;
47  bool affine_;
48 
49  void linear (const alu3d_ctype, const alu3d_ctype, const alu3d_ctype) ;
50  void linear (const coord_t&) ;
51  void inverse (const coord_t&) ;
52  public :
53  TrilinearMapping (const coord_t&, const coord_t&,
54  const coord_t&, const coord_t&,
55  const coord_t&, const coord_t&,
56  const coord_t&, const coord_t&);
57 
58  // only to call from geometry class
60 
62 
64  alu3d_ctype det (const coord_t&) ;
66  const mat_t& jacobianTransposed(const coord_t&);
67  void map2world (const coord_t&, coord_t&) const ;
68  void map2world (const alu3d_ctype , const alu3d_ctype , const alu3d_ctype ,
69  coord_t&) const ;
70  void world2map (const coord_t&, coord_t&) ;
71 
72  template <class vector_t>
73  void buildMapping(const vector_t&, const vector_t&,
74  const vector_t&, const vector_t&,
75  const vector_t&, const vector_t&,
76  const vector_t&, const vector_t&);
77 
78  // returns true if mapping is affine
79  inline bool affine () const { return affine_; }
80  };
81 
83  // NOTE: this class is different to the BilinearSurfaceMapping in
84  // ALUGrid, for example the reference elements differ
85  // here we have [0,1]^2 and in ALUGrid its [-1,1]^2
86  // also the point numbering is different
88  {
89  public:
90  // our coordinate types
91  typedef FieldVector<alu3d_ctype, 3> coord3_t;
92  typedef FieldVector<alu3d_ctype, 2> coord2_t;
93 
94  // type of coordinate vectors from elements
95  typedef alu3d_ctype double3_t[3];
96  protected:
97 
98  alu3d_ctype _n [3][3] ;
99 
100  static const double _epsilon ;
101 
102  bool _affine;
103 
104  public :
107 
110 
111  // returns true if mapping is affine
112  inline bool affine () const { return _affine ; }
113 
114  // calcuates normal
115  void normal(const coord2_t&, coord3_t&) const ;
116  void normal(const alu3d_ctype, const alu3d_ctype, coord3_t&)const;
117 
118  void negativeNormal(const coord2_t&, coord3_t&) const ;
119  void negativeNormal(const alu3d_ctype, const alu3d_ctype, coord3_t&)const;
120 
121  public:
122  // builds _b and _n, called from the constructors
123  // public because also used in faceutility
124  template <class vector_t>
125  void buildMapping (const vector_t & , const vector_t & ,
126  const vector_t & , const vector_t & );
127  protected:
128  // builds _b and _n, called from the constructors
129  // public because also used in faceutility
130  template <class vector_t>
131  void buildMapping (const vector_t & , const vector_t & ,
132  const vector_t & , const vector_t & ,
133  alu3d_ctype (&_b)[4][3] );
134  } ;
135 
136 
138  // NOTE: this class is different to the BilinearSurfaceMapping in
139  // ALUGrid, for example the reference elements differ
140  // here we have [0,1]^2 and in ALUGrid its [-1,1]^2
141  // also the point numbering is different
143  {
144  protected:
146 
147  using BaseType :: _n;
148  static const double _epsilon;
149 
150  // our coordinate types
151  typedef FieldVector<alu3d_ctype, 3> coord3_t;
152  typedef FieldVector<alu3d_ctype, 2> coord2_t;
153 
154  // type of coordinate vectors from elements
155  typedef alu3d_ctype double3_t[3];
156 
157  // type for helper matrices
158  typedef FieldMatrix<alu3d_ctype,3,3> mat3_t;
159 
160  // type for inverse matrices
161  typedef FieldMatrix<alu3d_ctype,2,3> matrix_t;
162 
163  // type for inverse matrices
164  typedef FieldMatrix<alu3d_ctype,3,2> inv_t;
165 
166  alu3d_ctype _b [4][3] ;
167 
168  mutable mat3_t Df,Dfi;
170  mutable matrix_t matrix_;
172 
173  mutable coord3_t normal_;
174  mutable coord3_t tmp_;
175 
176  mutable bool _calcedInv;
177  mutable bool _calcedTransposed;
178  mutable bool _calcedMatrix;
179 
180  public :
183 
185  BilinearSurfaceMapping (const coord3_t&, const coord3_t&,
186  const coord3_t&, const coord3_t&) ;
188  BilinearSurfaceMapping (const double3_t &, const double3_t &,
189  const double3_t &, const double3_t &) ;
192 
193  void inverse (const coord3_t&) const;
194  const inv_t& jacobianInverseTransposed(const coord2_t&) const ;
195 
196  const matrix_t& jacobianTransposed(const coord2_t&) const ;
197 
198  // calculates determinant of face mapping using the normal
199  alu3d_ctype det(const coord2_t&) const;
200 
201  // maps from local coordinates to global coordinates
202  void world2map(const coord3_t &, coord2_t & ) const;
203 
204  // maps form global coordinates to local (within reference element)
205  // coordinates
206  void map2world(const coord2_t&, coord3_t&) const ;
207  void map2world(const alu3d_ctype ,const alu3d_ctype , coord3_t&) const ;
208 
209  private:
210  void map2worldnormal(const alu3d_ctype, const alu3d_ctype, const alu3d_ctype , coord3_t&)const;
211  void map2worldlinear(const alu3d_ctype, const alu3d_ctype, const alu3d_ctype ) const;
212 
213  public:
214  // builds _b and _n, called from the constructors
215  // public because also used in faceutility
216  template <class vector_t>
217  void buildMapping (const vector_t & , const vector_t & ,
218  const vector_t & , const vector_t & );
219  } ;
220 
221 
222 
224  template< int cdim >
226  {
227  public:
229 
230  typedef FieldVector< ctype, cdim > world_t;
231  typedef FieldVector< ctype, 2 > map_t;
232 
233  typedef FieldMatrix< ctype, 2, cdim > matrix_t;
234  typedef FieldMatrix< ctype, cdim, 2 > inv_t;
235 
236  protected:
237  ctype _b [4][cdim];
238 
239  mutable ctype det_;
240  mutable matrix_t matrix_;
242 
243  mutable bool affine_;
244  mutable bool calcedMatrix_;
245  mutable bool calcedDet_;
246  mutable bool calcedInv_;
247 
248  public:
249  BilinearMapping ();
250  BilinearMapping ( const world_t &p0, const world_t &p1,
251  const world_t &p2, const world_t &p3 );
252  BilinearMapping ( const ctype (&p0)[ cdim ], const ctype (&p1)[ cdim ],
253  const ctype (&p2)[ cdim ], const ctype (&p3)[ cdim ] );
254 
255  bool affine () const;
256 
257  void world2map ( const world_t &, map_t & ) const;
258  void map2world ( const ctype x, const ctype y, world_t &w ) const;
259  void map2world ( const map_t &, world_t & ) const;
260 
261  ctype det ( const map_t & ) const;
262 
263  const matrix_t &jacobianTransposed ( const map_t & ) const;
264  const inv_t &jacobianInverseTransposed ( const map_t & ) const;
265 
266  template< class vector_t >
267  void buildMapping ( const vector_t &, const vector_t &,
268  const vector_t &, const vector_t & );
269 
270  protected:
271  static void multTransposedMatrix ( const matrix_t &, FieldMatrix< ctype, 2, 2 > & );
272  static void multMatrix ( const matrix_t &, const FieldMatrix< ctype, 2, 2 > &, inv_t & );
273 
274  void map2worldlinear ( const ctype, const ctype ) const;
275  void inverse ( const map_t & ) const;
276  };
277 
278 
279 
281  template< int cdim, int mydim >
283  {
284  public:
286 
287  typedef ctype double_t[ cdim ];
288 
289  typedef FieldVector< ctype, cdim > world_t;
290  typedef FieldVector< ctype, mydim > map_t;
291 
292  typedef FieldMatrix< ctype, mydim, cdim > matrix_t;
293  typedef FieldMatrix< ctype, cdim, mydim > inv_t;
294 
295  protected:
299 
301  mutable ctype _det;
302 
304  mutable bool _calcedInv;
305 
307  mutable bool _calcedDet;
308 
309  public:
311  LinearMapping ();
312 
314  LinearMapping (const LinearMapping &) ;
315 
316  // returns true if mapping is affine (which is always true)
317  inline bool affine () const { return true ; }
318 
319  // return reference to transposed jacobian
320  const matrix_t& jacobianTransposed(const map_t &) const ;
321 
322  // return reference to transposed jacobian inverse
323  const inv_t& jacobianInverseTransposed(const map_t &) const ;
324 
325  // calculates determinant of mapping
326  ctype det(const map_t&) const;
327 
328  // maps from local coordinates to global coordinates
329  void world2map(const world_t &, map_t &) const;
330 
331  // maps form global coordinates to local (within reference element)
332  // coordinates
333  void map2world(const map_t &, world_t &) const ;
334 
335  protected:
336  // calculate inverse
337  void inverse (const map_t&) const;
338 
339  // calculate inverse one codim one entity
340  void inverseCodimOne (const map_t&) const;
341 
342  // calculate determinant
343  void calculateDeterminant (const map_t&) const;
344 
345  void multTransposedMatrix(const matrix_t& matrix,
346  FieldMatrix<ctype, mydim, mydim>& result) const;
347 
348  void multMatrix ( const matrix_t& A,
349  const FieldMatrix< ctype, mydim, mydim> &B,
350  inv_t& ret ) const ;
351 
352  public:
353  // builds _b and _n, called from the constructors
354  // public because also used in faceutility
355  template <class vector_t>
356  void buildMapping (const vector_t & , const vector_t & ,
357  const vector_t & , const vector_t & );
358 
359  // builds _b and _n, called from the constructors
360  // public because also used in faceutility
361  template <class vector_t>
362  void buildMapping (const vector_t & , const vector_t & ,
363  const vector_t & );
364 
365  // builds _b and _n, called from the constructors
366  // public because also used in faceutility
367  template <class vector_t>
368  void buildMapping (const vector_t & , const vector_t & );
369 
370  template <class vector_t>
371  void buildMapping (const vector_t & );
372  } ;
373 
374 
376  //
377  // NonConforming Mappings
378  //
380 
381 
384  template< ALU3dGridElementType type, class Comm >
385  class NonConformingFaceMapping;
386 
388  template< class Comm >
389  struct NonConformingFaceMapping< tetra, Comm >
390  {
391  typedef FieldVector< alu3d_ctype, 3 > CoordinateType;
393 
394  NonConformingFaceMapping ( RefinementRuleType rule, int nChild )
395  : rule_( rule ), nChild_( nChild )
396  {}
397 
398  void child2parent ( const CoordinateType &childCoordinates,
399  CoordinateType &parentCoordinates) const;
400 
401  CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const;
402 
403  private:
404  void child2parentNosplit(const CoordinateType& childCoordinates,
405  CoordinateType& parentCoordinates) const;
406  void child2parentE01(const CoordinateType& childCoordinates,
407  CoordinateType& parentCoordinates) const;
408  void child2parentE12(const CoordinateType& childCoordinates,
409  CoordinateType& parentCoordinates) const;
410  void child2parentE20(const CoordinateType& childCoordinates,
411  CoordinateType& parentCoordinates) const;
412  void child2parentIso4(const CoordinateType& childCoordinates,
413  CoordinateType& parentCoordinates) const;
414 
415  RefinementRuleType rule_;
416  int nChild_;
417  };
418 
420  template< class Comm >
421  struct NonConformingFaceMapping< hexa, Comm >
422  {
423  typedef FieldVector< alu3d_ctype, 2 > CoordinateType;
425 
426  NonConformingFaceMapping ( RefinementRuleType rule, int nChild )
427  : rule_( rule ), nChild_( nChild )
428  {}
429 
430  void child2parent ( const CoordinateType &childCoordinates,
431  CoordinateType &parentCoordinates) const;
432 
433  CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const;
434 
435  private:
436  void child2parentNosplit(const CoordinateType& childCoordinates,
437  CoordinateType& parentCoordinates) const;
438  void child2parentIso4(const CoordinateType& childCoordinates,
439  CoordinateType& parentCoordinates) const;
440 
441  RefinementRuleType rule_;
442  int nChild_;
443  };
444 
445 } // end namespace Dune
446 
447 #if COMPILE_ALUGRID_INLINE
448  #include "mappings_imp.cc"
449 #endif
450 #endif