ElInfo.h 18.9 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  crystal growth group                                                  ==
// ==                                                                        ==
// ==  Stiftung caesar                                                       ==
// ==  Ludwig-Erhard-Allee 2                                                 ==
// ==  53175 Bonn                                                            ==
// ==  germany                                                               ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  http://www.caesar.de/cg/AMDiS                                         ==
// ==                                                                        ==
// ============================================================================

/** \file ElInfo.h */

#ifndef AMDIS_ELINFO_H
#define AMDIS_ELINFO_H

// ============================================================================
// ===== includes =============================================================
// ============================================================================
#include "Flag.h"
#include "Boundary.h"
#include "Global.h"
#include "FixVec.h"
#include "Element.h"

namespace AMDiS {

  // ============================================================================
  // ===== forward declarations =================================================
  // ============================================================================
  class MacroElement;
  class Mesh;
  class Element;
  class BasisFunction;
  class Projection;
  template<typename ReturnType, typename ArgumentType> class AbstractFunction;


  // ============================================================================
  // ===== class ElInfo =========================================================
  // ============================================================================

  /** \ingroup Traverse
   * \brief 
   * An ElInfo object holds informations wich are not stored in the corresponding
   * element. It is filled during mesh traversal by the traversal routines.
   * A fill flag determines which informations should be filled and which elements
   * should be visited. Since it is a
   * pure virtual base class for the dimension speciefic ElInfo classes, it must
   * not be instantiated directly.
   * \see ElInfo1d \see ElInfo2d \see ElInfo3d
   */

  class ElInfo
  {
    // ===== construtors, destructors =============================================
  protected: 
    /** \brief
     * Protected constructor. Avoids instatiation of the basis class
     */
    ElInfo();

    /** \brief
     * Protected constructor. Avoids instatiation of the basis class.
     * \param mesh pointer to the corresponding mesh.
     */
    ElInfo(Mesh *mesh);

  public:
    /** \brief
     * Virtual destructor because ElInfo is pure virtual.
     */
    virtual ~ElInfo();

    /** \brief
     * Assignement operator.
     * \param rhs right hand side.
     */
    ElInfo& operator=(const ElInfo& rhs) {
      mesh_ = rhs.mesh_;
      element_ = rhs.element_;
      parent_ = rhs.parent_;
      macroElement_ = rhs.macroElement_;
      fillFlag_ = rhs.fillFlag_;
93
      level = rhs.level;
94
      iChild = rhs.iChild;
95
96
97
98
      coord_ = rhs.coord_;
      boundary_ = rhs.boundary_;
      oppCoord_ = rhs.oppCoord_;
      neighbour_ = rhs.neighbour_;
99
      neighbourCoord_ = rhs.neighbourCoord_;
100
101
      oppVertex_ = rhs.oppVertex_;
      return *this;
102
    }
103
104
105
106
107
108
109
110
111
112
113
114

    // ===== getting-methods ======================================================
  public:
    /** \name getting methods
     * \{ 
     */

    /** \brief
     * Get ElInfo's \ref mesh_
     */
    inline Mesh* getMesh() const { 
      return mesh_; 
115
    }
116
117
118
119
120
121

    /** \brief
     * Get ElInfo's \ref macroElement_
     */
    inline MacroElement* getMacroElement() const { 
      return macroElement_; 
122
    }
123
124
125
126
127
128

    /** \brief
     * Get ElInfo's \ref element
     */
    inline Element* getElement() const { 
      return element_; 
129
    }
130
131
132
133
134
135

    /** \brief
     * Get ElInfo's \ref parent_
     */
    inline Element* getParent() const { 
      return parent_; 
136
    }
137
138
139
140
141
142

    /** \brief
     * Get ElInfo's \ref fillFlag_
     */
    inline Flag getFillFlag() const { 
      return fillFlag_; 
143
    }
144
145

    /** \brief
146
     * Get ElInfo's \ref level
147
     */
148
149
150
151
152
153
154
155
156
157
    inline int getLevel() const { 
      return level; 
    }

    /** \brief
     * Get ElInfo's \ref displacement. Is used only during dual traverse.
     */
    inline int getDisplacement() const {
      return displacement;
    }
158

159
160
161
162
163
    /** \brief
     * Get ElInfo's \ref iChild
     */ 
    inline int getIChild() const {
      return iChild;
164
    }
165
166
167
168
169
170
171

    /** \brief
     * Get ElInfo's \ref coord_[i]. This is a WorldVector<double> filled with the
     * coordinates of the i-th vertex of element \ref el.
     */
    inline WorldVector<double>& getCoord(int i) { 
      return coord_[i]; 
172
    }
173
174
175
176
177
178
179

    /** \brief
     * Get ElInfo's \ref coord_[i]. This is a WorldVector<double> filled with the
     * coordinates of the i-th vertex of element \ref el.
     */
    inline const WorldVector<double>& getCoord(int i) const { 
      return coord_[i]; 
180
    }
181
182
183
184
185

    /** \brief
     * Get ElInfo's \ref coord_. This is a FixVec<WorldVector<double> > filled with the
     * coordinates of the all vertice of element \ref el.
     */
186
    inline FixVec<WorldVector<double>, VERTEX>& getCoords() { 
187
      return coord_; 
188
    }
189
190
191
192
193

    /** \brief
     * Get ElInfo's \ref coord_. This is a FixVec<WorldVector<double> > filled with the
     * coordinates of the all vertice of element \ref el.
     */
194
    inline const FixVec<WorldVector<double>, VERTEX>& getCoords() const { 
195
      return coord_; 
196
    }
197
198
199
200
201
202
203
204
205
206
207

    /** \brief
     * Get ElInfo's \ref oppCoord_[i]
     */
    inline WorldVector<double>& getOppCoord(int i) { 
      return oppCoord_[i]; 
    };

    /** \brief
     * Get ElInfo's \ref boundary_[i] 
     */
208
    inline BoundaryType getBoundary(int i) const { 
209
      return boundary_[i]; 
210
    }
211
212
213
214
215
216
217
218
219
220
221

    /** \brief
     * Get boundary type of i-th vertex/edge/face (pos).
     */
    BoundaryType getBoundary(GeoIndex pos, int i);

    /** \brief
     * Get ElInfo's \ref neighbour_[i]
     */
    inline Element* getNeighbour(int i) const { 
      return neighbour_[i]; 
222
    }
223

224
225
226
227
228
229
230
    /** \brief
     * Get ElInfo's \ref neighbourCoord_[i]
     */
    inline FixVec<WorldVector<double>, VERTEX> getNeighbourCoord(int i) const {
      return neighbourCoord_[i];
    }

231
232
233
234
235
    /** \brief
     * Get ElInfo's \ref oppVertex_[i] 
     */
    inline unsigned char getOppVertex(int i) const { 
      return oppVertex_[i]; 
236
    }
237
238
239

    virtual int getSideOfNeighbour(int i) { 
      return oppVertex_[i]; 
240
    }
241
242
243
244
245
246

    /** \brief
     * Get ElInfo's \ref det_
     */
    inline double getDet() const { 
      return det_; 
247
    }
248
249

    /** \brief
Thomas Witkowski's avatar
Thomas Witkowski committed
250
     * Returns \ref grdLambda
251
252
     */
    inline const DimVec<WorldVector<double> >& getGrdLambda() const { 
Thomas Witkowski's avatar
Thomas Witkowski committed
253
      return grdLambda; 
254
    }
255
256
257
258
259
260

    /** \brief
     * Returns \ref projection_[i]
     */
    inline Projection *getProjection(int i) const {
      return projection_[i];
261
    }
262
263
264
265
266
267

    /** \brief
     * Returns \ref parametric_
     */
    inline bool getParametric() { 
      return parametric_; 
268
269
    }

270
271
272
273
274
275
    /** \brief
     * Returns the element transformation matrix \ref subElemCoordsMat .
     * This value is set only during dual traverse.
     */
    inline DimMat<double> *getSubElemCoordsMat() const {
      return subElemCoordsMat;
276
    }
277
278
279
280
281
282
283
284
285
286
287
288
289
290

    /** \} */ 

    // ===== setting-methods ======================================================

    /** \name setting methods
     * \{ 
     */

    /** \brief
     * Set ElInfo's \ref mesh_
     */
    inline void setMesh(Mesh* aMesh) { 
      mesh_ = aMesh; 
291
    }
292
293
294
295
296
297

    /** \brief
     * Set ElInfo's \ref macroElement_
     */
    inline void setMacroElement(MacroElement* mel) { 
      macroElement_ = mel; 
298
    }
299
300
301
302
303
304

    /** \brief
     * Set ElInfo's \ref element
     */
    inline void setElement(Element* elem) { 
      element_ = elem; 
305
    }
306
307
308
309
310
311

    /** \brief
     * Set ElInfo's \ref parent_
     */
    inline void setParent(Element* elem) { 
      parent_ = elem; 
312
    }
313
314
315
316
317
318

    /** \brief
     * Set ElInfo's \ref fillFlag_
     */
    inline void setFillFlag(Flag flag) { 
      fillFlag_ = flag; 
319
    }
320
321
322
323
324
325

    /** \brief
     * Sets ElInfo's \ref coord_[i]. 
     */
    inline void setCoord(int i,WorldVector<double>& coord) { 
      coord_[i] = coord; 
326
    }
327
328
329
330
331
332

    /** \brief
     * Sets ElInfo's \ref coord. 
     */
    inline void setCoords(FixVec<WorldVector<double>,VERTEX >& coords) { 
      coord_ = coords; 
333
    }
334
335

    /** \brief
336
     * Set ElInfo's \ref level
337
338
     */
    inline void setLevel(int l) { 
339
340
341
342
343
344
345
346
347
      level = l; 
    }

    /** \brief
     * Set ElInfo's \ref displacement. Used only during dual traverse.
     */
    inline void setDisplacement(int d) {
      displacement = d;
    }
348
349
350
351
352
353

    /** \brief
     * Set ElInfo's \ref boundary_[i] 
     */
    inline void setBoundary(int i, BoundaryType t) { 
      boundary_[i] = newBound(boundary_[i], t);
354
    }
355
356
357
358
359
360

    /** \brief
     * Set \ref projection_[i] = p
     */
    inline void setProjection(int i, Projection *p) {
      projection_[i] = p;
361
    }
362
363
364
365
366
367

    /** \brief
     * Set \ref det_ = d
     */
    inline void setDet(double d) { 
      det_ = d; 
368
    }
369
370
371
372
373
374

    /** \brief
     * Set \ref parametric_ = param
     */
    inline void setParametric(bool param) { 
      parametric_ = param; 
375
376
    }

377
378
379
380
381
382
    /** \brief
     * Sets the element transformation matrix \ref subElemCoordsMat .
     * This value is used only during dual traverse.
     */
    inline void setSubElemCoordsMat(DimMat<double> *mat) {
      subElemCoordsMat = mat;
383
    }
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
  
    /** \} */

    // ===== other public methods =================================================
  

    /** \brief
     * Returns the absolute value of the determinant of the affine linear 
     * parametrization's Jacobian
     */
    virtual double calcDet() const;

    /** \brief
     * Used by non static method \ref calcDet(). Calculates the determinant
     * for a given vector of vertex coordinates.
     */
    static double calcDet(const FixVec<WorldVector<double>, VERTEX> &coords);

    /** \brief
     * Checks whether flag is set in ElInfo's \ref fillFlag_. If not, the program
     * exits.
     */
    void testFlag(const Flag& flag) const;

    /** \brief
409
     * Returns a pointer to a vector, which contains the world coordinates
410
411
412
413
414
     * of a point in barycentric coordinates lambda with respect to \ref element.
     * If world is not NULL the world coordinates are stored in this vector.
     * Otherwise the function itself provides memory for this vector. In this
     * case the vector is overwritten during the next call of coordToWorld.
     */
415
416
    virtual const WorldVector<double> *coordToWorld(const DimVec<double>& lambda,
						    WorldVector<double>* world) const;
417
418
419
  

    /** \brief
Thomas Witkowski's avatar
Thomas Witkowski committed
420
     * Fills ElInfo's \ref det_ and \ref grdLambda entries.
421
422
423
424
425
426
427
428
429
430
431
432
     */
    virtual void fillDetGrdLambda();

    // ===== pure virtual functions. Must be overriden in sub-classes ============

    /** \brief
     * Returns a pointer to a vector, which contains the barycentric coordinates
     * with respect to \ref element of a point with world coordinates world.
     * The barycentric coordinates are stored in lambda. 
     * pure virtual => must be overriden in sub-class.
     */
    virtual const int worldToCoord(const WorldVector<double>& world, 
433
				   DimVec<double>* lambda) const = 0;
434
435
436
437
438
439
440
441
442
443
444
445

    /** \brief
     * Fills this ElInfo with macro element information of mel.
     * pure virtual => must be overriden in sub-class.
     */
    virtual void fillMacroInfo(const MacroElement *mel) = 0;

    /** \brief
     * Fills this ElInfo for the child ichild using hierarchy information and
     * parent data parentInfo.
     * pure virtual => must be overriden in sub-class.
     */
446
    virtual void fillElInfo(int ichild, const ElInfo *parentInfo) = 0;
447
448
449
450
451
452
453

    /** \brief
     * calculates the Jacobian of the barycentric coordinates on \element and stores
     * the matrix in grd_lam. The return value of the function is the absolute
     * value of the determinant of the affine linear paraetrization's Jacobian.
     * pure virtual => must be overriden in sub-class.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
454
    virtual double calcGrdLambda(DimVec<WorldVector<double> >& grd_lam) = 0;
455
456
457
458
459
460
461

    /** \brief
     * calculates a normal of the given side (1d,2d: edge, 3d: face) of \ref element.
     * Returns the absolute value of the determinant of the
     * transformation to the reference element.
     * pure virtual => must be overriden in sub-class.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
462
    virtual double getNormal(int side, WorldVector<double> &normal) = 0;
463
464
465
466
467
468
469
470
471
472


    /** \brief
     * calculates a normal of the element in dim of world = dim + 1.
     * Returns the absolute value of the determinant of the
     * transformation to the reference element.
     * pure virtual => must be overriden in sub-class.
     */
    virtual double getElementNormal(WorldVector<double> &elementNormal) const 
    {
473
474
      FUNCNAME("ElInfo::getElementNormal()");

475
476
477
      ERROR("virtual function not implemented in this sub-class ");
    
      return(0.0);
478
    }
479

480
    virtual void getRefSimplexCoords(DimMat<double> *coords) const = 0;
481

482
    virtual void getSubElementCoords(DimMat<double> *coords,
483
				     int iChild) const = 0;
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518

  protected:

    /** \brief
     *  Pointer to the current mesh
     */
    Mesh *mesh_;

    /** \brief 
     * Pointer to the current element
     */
    Element *element_;

    /** \brief 
     * \ref element is child of element parent_
     */
    Element *parent_;

    /** \brief 
     * \ref element is an element of the binary tree located at MacroElement 
     * macroElement_
     */
    MacroElement *macroElement_;

    /** \brief
     * Indicates wich elements will be called and wich information should be
     * present while mesh traversal.
     */
    Flag fillFlag_;

    /** \brief 
     * Level of the element. The level is zero for macro elements and the level
     * of the children is (level of the parent + 1). level_ is filled always by
     * the traversal routines.
     */
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
    unsigned char level;

    /** \brief
     * This value is set only during dual traverse, and only if this is the smaller
     * element of an element pair. Then the value describes the displacement of
     * this element to the first element in the refinement hierachy starting from
     * the bigger element of the element pair.
     *
     * Example, if the level difference of the two elements is 1:
     *
     *                 a
     *                / \
     *               b   c
     * 
     * The element b has displacement 0, the element c has displacement 1.
     *
     * Example, if the level difference of the two elements is 2:
     *
     *                 a
     *                / \
     *               /   \
     *              b     c
     *             / \   / \
     *            d   e f   g
     *
     * The displacements are as follows: d = 0, e = 1, f = 2, g = 3.
     */
    int displacement;
547

548
549
550
551
552
    /** \brief
     * This ElInfo is the iChild-th child of the parent element.
     */
    int iChild;

553
554
555
556
    /** \brief 
     * \ref coord_[i] is a WorldVector<double> storing the world coordinates of the
     * i-th vertex of element \ref element.
     */
557
    FixVec<WorldVector<double>, VERTEX> coord_;
558
559
560
561
562
563
564
565

    /** \brief 
     * boundary_[i] is the BoundaryType of the i-th edge/face
     * for i=0,...,N_NEIGH - 1. In 3d 
     * (*boundary_)[N_FACES + i] is a pointer to the Boundary
     * object of the i-th edge, for i=0,..,N_EDGES - 1. It is
     * a pointer to NULL for an interior edge/face.
     */
566
    FixVec<BoundaryType, BOUNDARY> boundary_;
567
568
569
570
571
572
573
574
575
576

    /** \brief
     * Vector storing pointers to projections for each face, edge, vertex.
     */
    FixVec<Projection*, PROJECTION> projection_;

    /** \brief 
     * oppCoord_[i] coordinates of the i-th neighbour vertex opposite the
     * common edge/face.
     */
577
    FixVec<WorldVector<double>, NEIGH> oppCoord_;
578
579
580
581
582

    /** \brief 
     * neighbour_[i] pointer to the element at the edge/face with local index i.
     * It is a pointer to NULL for boundary edges/faces.
     */
583
584
585
586
587
588
589
    FixVec<Element*, NEIGH> neighbour_;

    /** \brief
     * neighbourCoord_[i][j] are the coordinate of the j-th vertex of the i-th
     * neighbour element with the common edge/face.
     */
    FixVec<FixVec<WorldVector<double>, VERTEX>, NEIGH> neighbourCoord_;
590
591
592
593
594
595

    /** \brief 
     * oppVertex_[i] is undefined if neighbour_[i] is a pointer to NULL. 
     * Otherwise it is the local index of the neighbour's vertex opposite the
     * common edge/face.
     */
596
    FixVec<unsigned char, NEIGH> oppVertex_;
597
598
599
600
601
602
603
604
605

    /** \brief
     * Elements determinant.
     */
    double det_;

    /** \brief
     * Gradient of lambda.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
606
    DimVec<WorldVector<double> > grdLambda;
607
608

    /** \brief
609
     * True, if this elInfo stores parametrized information. False, otherwise.
610
611
612
     */
    bool parametric_;

Thomas Witkowski's avatar
Thomas Witkowski committed
613
614
615
616
617
    /** \brief
     * Stores the world dimension.
     */
    int dimOfWorld;

618
619
620
621
622
623
624
625
    /** \brief
     * This is a transformation matrix used during dual traverse. It is set, if
     * the current element is the smaller element of an element pair in the traverse.
     * Then this matrix defines a mapping for points defined in barycentric 
     * coordinates on the larger element, to the barycentric coordinates of the smaller
     * element.
     */
    DimMat<double> *subElemCoordsMat;
626

627
628
629
630
631
632
633
634
635
    // ===== static public members ================================================
  public:
    /** \brief 
     * child_vertex[el_type][child][i] = father's local vertex index of new 
     * vertex i. 4 stands for the newly generated vertex .
     */
    static const int childVertex[3][2][4];

    /** \brief 
636
     * child_edge[el_type][child][i] = father's local edge index of new edge i.
637
638
639
640
641
642
643
     * new edge 2 is half of old edge 0, new edges 4,5 are really new edges, and
     * value is different: child_edge[][][4,5] = index of same edge in other
     * child.
     */
    static const int childEdge[3][2][6];

    /** \brief
644
645
     * Used to determine to which traverse an ElInfo belongs.
     * Stores the value for each thread, maximum value 16!.
646
     */
647
    static int traverseId[16];
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673

    /** \brief
     * Used to determine to which traverse an ElInfo belongs
     */
    static int traverseIdCounter;

    friend class ElInfo1d;
    friend class ElInfo2d;
    friend class ElInfo3d;
  };

}

#include "ElInfo1d.h"
#include "ElInfo2d.h"
#include "ElInfo3d.h"

#endif  // AMDIS_ELINFO_H