dune-alugrid 2.8.0
Loading...
Searching...
No Matches
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
16namespace 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>
22
23 template<int dim, int dimworld, ALU3dGridElementType, class >
24 class ALU3dGrid;
25
29 {
30 public:
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
96 protected:
97
98 alu3d_ctype _n [3][3] ;
99
100 static const double _epsilon ;
101
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
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
167
168 mutable mat3_t Df,Dfi;
172
174 mutable coord3_t tmp_;
175
176 mutable bool _calcedInv;
177 mutable bool _calcedTransposed;
178 mutable bool _calcedMatrix;
179
180 public :
183
186 const coord3_t&, const coord3_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_;
242
243 mutable bool affine_;
244 mutable bool calcedMatrix_;
245 mutable bool calcedDet_;
246 mutable bool calcedInv_;
247
248 public:
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< int dim, int dimworld, ALU3dGridElementType type, class Comm >
386
388 template< int dim, int dimworld, class Comm >
389 class NonConformingFaceMapping< dim, dimworld, tetra, Comm >
390 {
391 public:
392 typedef FieldVector< alu3d_ctype, 3 > CoordinateType;
394
396 : rule_( rule ), nChild_( nChild )
397 {}
398
399 void child2parent ( const CoordinateType &childCoordinates,
400 CoordinateType &parentCoordinates) const;
401
402 CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const;
403
404 private:
405 void child2parentNosplit(const CoordinateType& childCoordinates,
406 CoordinateType& parentCoordinates) const;
407 void child2parentE01(const CoordinateType& childCoordinates,
408 CoordinateType& parentCoordinates) const;
409 void child2parentE12(const CoordinateType& childCoordinates,
410 CoordinateType& parentCoordinates) const;
411 void child2parentE20(const CoordinateType& childCoordinates,
412 CoordinateType& parentCoordinates) const;
413 void child2parentIso4(const CoordinateType& childCoordinates,
414 CoordinateType& parentCoordinates) const;
415
416 RefinementRuleType rule_;
417 int nChild_;
418 };
419
421 template< int dim, int dimworld, class Comm >
422 class NonConformingFaceMapping< dim, dimworld, hexa, Comm >
423 {
424 public:
425 typedef FieldVector< alu3d_ctype, 2 > CoordinateType;
427
429 : rule_( rule ), nChild_( nChild )
430 {}
431
432 void child2parent ( const CoordinateType &childCoordinates,
433 CoordinateType &parentCoordinates) const;
434
435 CoordinateType child2parent ( const FieldVector< alu3d_ctype, 2 > &childCoordinates ) const;
436
437 private:
438 void child2parentNosplit(const CoordinateType& childCoordinates,
439 CoordinateType& parentCoordinates) const;
440 void child2parentIso4(const CoordinateType& childCoordinates,
441 CoordinateType& parentCoordinates) const;
442
443 RefinementRuleType rule_;
444 int nChild_;
445 };
446
447} // end namespace Dune
448
449#if COMPILE_ALUGRID_INLINE
450 #include "mappings_imp.cc"
451#endif
452#endif
Definition: alu3dinclude.hh:63
@ hexa
Definition: topology.hh:12
@ tetra
Definition: topology.hh:12
double alu3d_ctype
Definition: alu3dinclude.hh:68
static const alu3d_ctype ALUnumericEpsilon
Definition: mappings.hh:18
Definition: alu3dinclude.hh:242
[ provides Dune::Grid ]
Definition: 3d/grid.hh:429
Definition: geometry.hh:632
Definition: mappings.hh:29
void map2world(const coord_t &, coord_t &) const
Definition: mappings_imp.cc:112
TrilinearMapping()
Definition: mappings.hh:59
void world2map(const coord_t &, coord_t &)
Definition: mappings_imp.cc:219
FieldMatrix< alu3d_ctype, 3, 3 > mat_t
Definition: mappings.hh:33
FieldVector< alu3d_ctype, 3 > coord_t
Definition: mappings.hh:32
const mat_t & jacobianTransposed(const coord_t &)
Definition: mappings_imp.cc:86
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 &)
bool affine() const
Definition: mappings.hh:79
alu3d_ctype det(const coord_t &)
Definition: mappings_imp.cc:167
alu3d_ctype double_t[3]
Definition: mappings.hh:31
~TrilinearMapping()
Definition: mappings.hh:63
const mat_t & jacobianInverseTransposed(const coord_t &)
Definition: mappings_imp.cc:93
A bilinear surface mapping.
Definition: mappings.hh:88
void buildMapping(const vector_t &, const vector_t &, const vector_t &, const vector_t &, alu3d_ctype(&_b)[4][3])
void negativeNormal(const coord2_t &, coord3_t &) const
Definition: mappings_imp.cc:348
void normal(const coord2_t &, coord3_t &) const
Definition: mappings_imp.cc:333
bool affine() const
Definition: mappings.hh:112
static const double _epsilon
Definition: mappings.hh:100
FieldVector< alu3d_ctype, 3 > coord3_t
Definition: mappings.hh:91
SurfaceNormalCalculator()
Constructor creating empty mapping with double , i.e. zero.
Definition: mappings_imp.cc:253
FieldVector< alu3d_ctype, 2 > coord2_t
Definition: mappings.hh:92
alu3d_ctype _n[3][3]
Definition: mappings.hh:98
void buildMapping(const vector_t &, const vector_t &, const vector_t &, const vector_t &)
~SurfaceNormalCalculator()
Definition: mappings.hh:109
bool _affine
Definition: mappings.hh:102
alu3d_ctype double3_t[3]
Definition: mappings.hh:95
A bilinear surface mapping.
Definition: mappings.hh:143
const matrix_t & jacobianTransposed(const coord2_t &) const
Definition: mappings_imp.cc:475
mat3_t Df
Definition: mappings.hh:168
FieldVector< alu3d_ctype, 3 > coord3_t
Definition: mappings.hh:151
BilinearSurfaceMapping()
Constructor creating empty mapping with double , i.e. zero.
Definition: mappings_imp.cc:369
matrix_t matrix_
Definition: mappings.hh:170
void buildMapping(const vector_t &, const vector_t &, const vector_t &, const vector_t &)
FieldMatrix< alu3d_ctype, 3, 3 > mat3_t
Definition: mappings.hh:158
inv_t invTransposed_
Definition: mappings.hh:169
mat3_t Dfi
Definition: mappings.hh:168
bool _calcedTransposed
Definition: mappings.hh:177
FieldVector< alu3d_ctype, 2 > coord2_t
Definition: mappings.hh:152
const inv_t & jacobianInverseTransposed(const coord2_t &) const
Definition: mappings_imp.cc:533
alu3d_ctype double3_t[3]
Definition: mappings.hh:155
~BilinearSurfaceMapping()
Definition: mappings.hh:191
void world2map(const coord3_t &, coord2_t &) const
Definition: mappings_imp.cc:560
static const double _epsilon
Definition: mappings.hh:148
coord3_t normal_
Definition: mappings.hh:173
FieldMatrix< alu3d_ctype, 2, 3 > matrix_t
Definition: mappings.hh:161
bool _calcedInv
Definition: mappings.hh:176
alu3d_ctype _b[4][3]
Definition: mappings.hh:166
alu3d_ctype det(const coord2_t &) const
Definition: mappings_imp.cc:498
coord3_t tmp_
Definition: mappings.hh:174
alu3d_ctype DetDf
Definition: mappings.hh:171
void map2world(const coord2_t &, coord3_t &) const
Definition: mappings_imp.cc:423
void inverse(const coord3_t &) const
Definition: mappings_imp.cc:507
bool _calcedMatrix
Definition: mappings.hh:178
FieldMatrix< alu3d_ctype, 3, 2 > inv_t
Definition: mappings.hh:164
SurfaceNormalCalculator BaseType
Definition: mappings.hh:145
A bilinear mapping.
Definition: mappings.hh:226
bool calcedMatrix_
Definition: mappings.hh:244
void buildMapping(const vector_t &, const vector_t &, const vector_t &, const vector_t &)
ctype _b[4][cdim]
Definition: mappings.hh:237
matrix_t matrix_
Definition: mappings.hh:240
FieldMatrix< ctype, 2, cdim > matrix_t
Definition: mappings.hh:233
BilinearMapping()
Definition: mappings_imp.cc:600
void map2world(const ctype x, const ctype y, world_t &w) const
Definition: mappings_imp.cc:645
bool affine() const
Definition: mappings_imp.cc:631
bool calcedDet_
Definition: mappings.hh:245
ctype det_
Definition: mappings.hh:239
bool affine_
Definition: mappings.hh:243
inv_t invTransposed_
Definition: mappings.hh:241
void world2map(const world_t &, map_t &) const
Definition: mappings_imp.cc:654
static void multMatrix(const matrix_t &, const FieldMatrix< ctype, 2, 2 > &, inv_t &)
Definition: mappings_imp.cc:790
static void multTransposedMatrix(const matrix_t &, FieldMatrix< ctype, 2, 2 > &)
Definition: mappings_imp.cc:774
FieldVector< ctype, 2 > map_t
Definition: mappings.hh:231
bool calcedInv_
Definition: mappings.hh:246
void inverse(const map_t &) const
Definition: mappings_imp.cc:826
FieldVector< ctype, cdim > world_t
Definition: mappings.hh:230
alu3d_ctype ctype
Definition: mappings.hh:228
const matrix_t & jacobianTransposed(const map_t &) const
Definition: mappings_imp.cc:733
void map2worldlinear(const ctype, const ctype) const
Definition: mappings_imp.cc:800
FieldMatrix< ctype, cdim, 2 > inv_t
Definition: mappings.hh:234
const inv_t & jacobianInverseTransposed(const map_t &) const
Definition: mappings_imp.cc:742
ctype det(const map_t &) const
Definition: mappings_imp.cc:686
A linear mapping.
Definition: mappings.hh:283
inv_t _invTransposed
storage for inverse of jacobian (transposed)
Definition: mappings.hh:297
void inverse(const map_t &) const
Definition: mappings_imp.cc:1054
void multMatrix(const matrix_t &A, const FieldMatrix< ctype, mydim, mydim > &B, inv_t &ret) const
Definition: mappings_imp.cc:1192
void inverseCodimOne(const map_t &) const
Definition: mappings_imp.cc:1119
FieldMatrix< ctype, mydim, cdim > matrix_t
Definition: mappings.hh:292
bool affine() const
Definition: mappings.hh:317
FieldVector< ctype, cdim > world_t
Definition: mappings.hh:289
void buildMapping(const vector_t &, const vector_t &)
alu3d_ctype ctype
Definition: mappings.hh:285
FieldMatrix< ctype, cdim, mydim > inv_t
Definition: mappings.hh:293
world_t _p0
Definition: mappings.hh:298
void calculateDeterminant(const map_t &) const
Definition: mappings_imp.cc:1084
void map2world(const map_t &, world_t &) const
Definition: mappings_imp.cc:1027
FieldVector< ctype, mydim > map_t
Definition: mappings.hh:290
ctype det(const map_t &) const
Definition: mappings_imp.cc:1257
void multTransposedMatrix(const matrix_t &matrix, FieldMatrix< ctype, mydim, mydim > &result) const
Definition: mappings_imp.cc:1173
LinearMapping()
Constructor creating empty mapping with double , i.e. zero.
Definition: mappings_imp.cc:853
void buildMapping(const vector_t &)
void buildMapping(const vector_t &, const vector_t &, const vector_t &, const vector_t &)
const matrix_t & jacobianTransposed(const map_t &) const
Definition: mappings_imp.cc:1271
const inv_t & jacobianInverseTransposed(const map_t &) const
Definition: mappings_imp.cc:1279
bool _calcedInv
true if inverse has been calculated
Definition: mappings.hh:304
void buildMapping(const vector_t &, const vector_t &, const vector_t &)
ctype double_t[cdim]
Definition: mappings.hh:287
bool _calcedDet
true if determinant has been calculated
Definition: mappings.hh:307
matrix_t _matrix
transformation matrix (transposed)
Definition: mappings.hh:296
ctype _det
P[0].
Definition: mappings.hh:301
void world2map(const world_t &, map_t &) const
Definition: mappings_imp.cc:1039
Definition: mappings.hh:385
FieldVector< alu3d_ctype, 3 > CoordinateType
Definition: mappings.hh:392
ALU3dImplTraits< tetra, Comm >::HfaceRuleType RefinementRuleType
Definition: mappings.hh:393
NonConformingFaceMapping(RefinementRuleType rule, int nChild)
Definition: mappings.hh:395
NonConformingFaceMapping(RefinementRuleType rule, int nChild)
Definition: mappings.hh:428
FieldVector< alu3d_ctype, 2 > CoordinateType
Definition: mappings.hh:425
ALU3dImplTraits< hexa, Comm >::HfaceRuleType RefinementRuleType
Definition: mappings.hh:426