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