Liebe Gitlab-Nutzerin, lieber Gitlab-Nutzer,
es ist nun möglich sich mittels des ZIH-Logins/LDAP an unserem Dienst anzumelden. Die Konten der externen Nutzer:innen sind über den Reiter "Standard" erreichbar.
Die Administratoren


Dear Gitlab user,
it is now possible to log in to our service using the ZIH login/LDAP. The accounts of external users can be accessed via the "Standard" tab.
The administrators

Element.h 18.3 KB
Newer Older
1 2 3 4 5 6
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
7
// ==  TU Dresden                                                            ==
8
// ==                                                                        ==
9 10 11
// ==  Institut fr Wissenschaftliches Rechnen                               ==
// ==  Zellescher Weg 12-14                                                  ==
// ==  01069 Dresden                                                         ==
12 13 14 15
// ==  germany                                                               ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
16
// ==  https://gforge.zih.tu-dresden.de/projects/amdis/                      ==
17 18 19 20 21 22 23 24
// ==                                                                        ==
// ============================================================================

/** \file Element.h */

#ifndef AMDIS_ELEMENT_H
#define AMDIS_ELEMENT_H

25
#include "AMDiS_fwd.h"
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
#include "Global.h"
#include "RefinementManager.h"
#include "Serializable.h"
#include "ElementData.h"
#include "LeafData.h"

namespace AMDiS {

  template<typename T, GeoIndex d> class FixVec;

#define AMDIS_UNDEFINED  5

  /** \ingroup Triangulation 
   * \brief
   * Base class for Line, Triangle, Tetrahedron
   *
   * Elements in AMDiS are always simplices (a simplex is a Line in 1d, a 
   * Triangle in 2d and a Tetrahedron in 3d). 
   * We restrict ourselves here to simplicial meshes, for several reasons:
   * -# A simplex is one of the most simple geometric types and complex domains 
   *    may be approximated by a set of simplices quite easily.
   * -# Simplicial meshes allow local refinement without the need of 
   *    nonconforming meshes (hanging nodes), parametric elements, or mixture of
   *    element types (which is the case for quadrilateral meshes).
   * -# Polynomials of any degree are easily represented on a simplex using 
   *    local (barycentric) coordinates.
   *
   * A Line element and its refinement:
   *
   * <img src="line.png">
   *
   * A Triangle element and its refinement:
   *
   * <img src="triangle.png">
   *
   * A Tetrahedron element and its refinements:
   *
   * <img src="tetrahedron.png">
   */
  class Element : public Serializable
  {
  private:
68
    /// private standard constructor because an Element must know his Mesh
69
    Element() {}
70

71
  public:
72
    /// constructs an Element which belongs to Mesh
73 74
    Element(Mesh *);

75
    /// copy constructor
76 77
    Element(const Element& old);

78
    /// destructor
79 80
    virtual ~Element();

81
    ///
82 83
    void deleteElementDOFs();

84 85 86 87 88 89
    /** \brief
     * Clone this Element and return a reference to it. Because also the DOFs
     * are cloned, \ref Mesh::serializedDOfs must be used.
     */
    Element* cloneWithDOFs();

90 91 92 93
    /** \name getting methods
     * \{
     */

94
    /// Returns \ref child[0]
Thomas Witkowski's avatar
Thomas Witkowski committed
95 96
    inline Element* getFirstChild() const 
    {
97
      return child[0];
98
    }
99

100
    /// Returns \ref child[1]
Thomas Witkowski's avatar
Thomas Witkowski committed
101 102
    inline Element* getSecondChild() const 
    {
103
      return child[1];
104
    }
105

106
    /// Returns \ref child[i], i=0,1
Thomas Witkowski's avatar
Thomas Witkowski committed
107 108
    inline Element* getChild(int i) const 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
109
      TEST_EXIT_DBG(i == 0 || i == 1)("There is only child 0 or 1! (i = %d)\n", i);
110
      return child[i];
111
    }
112 113 114 115 116

    /** \brief
     * Returns true if Element is a leaf element (\ref child[0] == NULL), returns
     * false otherwise.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
117 118
    inline const bool isLeaf() const 
    { 
119
      return (child[0] == NULL); 
120
    }
121

122
    /// Returns \ref dof[i][j] which is the j-th DOF of the i-th node of Element.
123
    const DegreeOfFreedom getDof(int i, int j) const 
Thomas Witkowski's avatar
Thomas Witkowski committed
124
    { 
125
      return dof[i][j];
126
    }
127

128
    /// Returns \ref dof[i] which is a pointer to the DOFs of the i-th node.
129
    const DegreeOfFreedom* getDof(int i) const 
Thomas Witkowski's avatar
Thomas Witkowski committed
130
    {
131
      return dof[i];
132
    }
133

134
    /// Returns a pointer to the DOFs of this Element
135
    const DegreeOfFreedom** getDof() const 
Thomas Witkowski's avatar
Thomas Witkowski committed
136
    {
137
      return const_cast<const DegreeOfFreedom**>(dof);
138
    }
139

140
    /// Returns \ref mesh of Element
Thomas Witkowski's avatar
Thomas Witkowski committed
141 142
    inline Mesh* getMesh() const 
    { 
143
      return mesh; 
144
    }
145 146 147 148 149 150 151

    /** \brief
     * Returns \ref elementData's error estimation, if Element is a leaf element
     * and has leaf data. 
     */
    inline double getEstimation(int row) const
    {
152
      if (isLeaf()) {
153
	TEST_EXIT_DBG(elementData)("leaf element without leaf data\n");
154
	ElementData *ld = elementData->getElementData(ESTIMATABLE);
155
	TEST_EXIT_DBG(ld)("leaf data not estimatable!\n");
156

157 158 159 160
	return dynamic_cast<LeafDataEstimatableInterface*>(ld)->getErrorEstimate(row);
      }	
      
      return 0.0;
161
    }
162 163 164 165 166

    /** \brief
     * Returns Element's coarsening error estimation, if Element is a leaf 
     * element and if it has leaf data and if this leaf data are coarsenable.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
167 168
    inline double getCoarseningEstimation(int row) 
    {
169
      if (isLeaf()) {
170
	TEST_EXIT_DBG(elementData)("leaf element without leaf data\n");
171
	ElementData *ld = elementData->getElementData(COARSENABLE);
172
	TEST_EXIT_DBG(ld)("element data not coarsenable!\n");
173

174
	return dynamic_cast<LeafDataCoarsenableInterface*>(ld)->getCoarseningErrorEstimate(row);
175
      }
176 177
      
      return 0.0;
178
    }
179

180
    /// Returns region of element if defined, -1 else.
181 182
    int getRegion() const;

183
    /// Returns local vertex number of the j-th vertex of the i-th edge
184 185 186 187 188 189 190
    virtual int getVertexOfEdge(int i, int j) const = 0; 

    /** \brief
     * Returns local vertex number of the vertexIndex-th vertex of the
     * positionIndex-th part of type position (vertex, edge, face)
     */
    virtual int getVertexOfPosition(GeoIndex position,
191 192
				    int positionIndex,
				    int vertexIndex) const = 0;
193

194
    ///
195 196
    virtual int getPositionOfVertex(int side, int vertex) const = 0;

197
    ///
198 199
    virtual int getEdgeOfFace(int face, int edge) const = 0;

200
    ///
201 202 203 204
    virtual DofEdge getEdge(int localEdgeIndex) const = 0;

    ///
    virtual DofFace getFace(int localFaceIndex) const = 0;
205

206
    /// Returns the number of parts of type i in this element
207 208
    virtual int getGeo(GeoIndex i) const = 0;

209
    /// Returns Element's \ref mark
Thomas Witkowski's avatar
Thomas Witkowski committed
210 211
    inline const signed char getMark() const 
    { 
212
      return mark;
213
    }
214

215
    /// Returns \ref newCoord[i]
Thomas Witkowski's avatar
Thomas Witkowski committed
216 217
    inline double getNewCoord(int i) const 
    {
218 219 220
	TEST_EXIT_DBG(newCoord)("newCoord = NULL\n");
	return (*newCoord)[i];
    }
221

222
    /// Returns Element's \ref index
Thomas Witkowski's avatar
Thomas Witkowski committed
223 224
    inline int getIndex() const 
    { 
225
      return index; 
226
    }
227

228
    /// Returns \ref newCoord
Thomas Witkowski's avatar
Thomas Witkowski committed
229 230
    inline WorldVector<double>* getNewCoord() const 
    { 
231
      return newCoord; 
232
    }
233 234 235 236 237 238 239

    /** \} */

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

240
    /// Sets \ref child[0]
Thomas Witkowski's avatar
Thomas Witkowski committed
241 242
    virtual void setFirstChild(Element *aChild) 
    {
243
      child[0] = aChild;
244
    }
245

246
    /// Sets \ref child[1]
Thomas Witkowski's avatar
Thomas Witkowski committed
247 248
    virtual void setSecondChild(Element *aChild) 
    {
249
      child[1] = aChild;
250
    }
251

252
    /// Sets \ref elementData of Element
Thomas Witkowski's avatar
Thomas Witkowski committed
253 254
    void setElementData(ElementData* ed) 
    {
255
      elementData = ed;
256
    }
257 258 259 260 261

    /** \brief
     * Sets \ref newCoord of Element. Needed by refinement, if Element has a
     * boundary edge on a curved boundary.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
262 263
    inline void setNewCoord(WorldVector<double>* coord) 
    {
264
      newCoord = coord;
265
    }
266

267
    /// Sets \ref mesh.
Thomas Witkowski's avatar
Thomas Witkowski committed
268 269
    inline void setMesh(Mesh *m) 
    {
270
      mesh = m;
271
    }
272

273
    /// Sets the pointer to the DOFs of the i-th node of Element
274
    inline void setDof(int pos, DegreeOfFreedom* p) 
Thomas Witkowski's avatar
Thomas Witkowski committed
275
    {
276
      dof[pos] = p;
277
    }
278 279 280 281 282 283 284

    /** \brief
     * Checks whether Element is a leaf element and whether it has leaf data.
     * If the checks don't fail, leaf data's error estimation is set to est.
     */
    inline void setEstimation(double est, int row)
    {
285
      if (isLeaf()) {
286
	TEST_EXIT_DBG(elementData)("leaf element without leaf data\n");
287
	ElementData *ld = elementData->getElementData(ESTIMATABLE);
288
	TEST_EXIT_DBG(ld)("leaf data not estimatable\n");
289 290 291

	dynamic_cast<LeafDataEstimatableInterface*>(ld)->
	  setErrorEstimate(row, est);
292
      } else {
293 294
	ERROR_EXIT("setEstimation only for leaf elements!\n");
      }
295
    }
296 297 298 299 300 301 302

    /** \brief
     * Sets Element's coarsening error estimation, if Element is a leaf element
     * and if it has leaf data and if this leaf data are coarsenable.
     */
    inline void setCoarseningEstimation(double est, int row)
    {
303
      if (isLeaf()) {
304
	TEST_EXIT_DBG(elementData)("leaf element without leaf data\n");
305
	ElementData *ld = elementData->getElementData(COARSENABLE);
306
	TEST_EXIT_DBG(ld)("leaf data not coarsenable\n");
307 308 309

	dynamic_cast<LeafDataCoarsenableInterface*>(ld)->
	  setCoarseningErrorEstimate(row, est);
310
      } else {
311 312
	ERROR_EXIT("setEstimation only for leaf elements!\n");
      }
313
    }
314

315
    /// Sets Elements \ref mark = mark + 1;
Thomas Witkowski's avatar
Thomas Witkowski committed
316 317
    inline void incrementMark() 
    {
318 319
      mark++;
    }
320

321
    /// Sets Elements \ref mark = mark - 1;
Thomas Witkowski's avatar
Thomas Witkowski committed
322 323
    inline void decrementMark() 
    {
324 325
      if (0 < mark) 
	mark--;
326
    }
327

328
    /// Sets Element's \ref mark
Thomas Witkowski's avatar
Thomas Witkowski committed
329 330
    inline void setMark(signed char m) 
    {
331
      mark = m;
332
    }
333 334 335 336 337 338 339 340 341 342 343 344 345

    /** \} */


    /** \name pure virtual methods 
     * \{ 
     */

    /** \brief
     * orient the vertices of edges/faces.
     * Used by Estimator for the jumps => same quadrature nodes from both sides!
     */
    virtual const FixVec<int,WORLD>& 
346
      sortFaceIndices(int face, FixVec<int, WORLD> *vec) const = 0;
347

348
    /// Returns a copy of itself. Needed by Mesh to create Elements by a prototype. 
349 350 351
    virtual Element *clone() = 0;

    /** \brief
Thomas Witkowski's avatar
Thomas Witkowski committed
352 353
     * Returns which side of child[childnr] corresponds to side sidenr of this 
     * Element. If the child has no corresponding side, the return value is negative. 
354 355 356
     */
    virtual int getSideOfChild(int childnr, int sidenr, int elType = 0) const = 0;

357 358 359 360 361 362 363 364 365 366 367 368
    /** \brief
     * Generalization of \ref getSideOfChild to arbitrary subObject. Thus, e.g., in 3d
     * we can ask for the local id of a verte, edge or face on the elements children.
     *
     * \param[in]  childnr    Either 0 or 1 for the left or right children.
     * \param[in]  subObj     Defines whether we ask for VERTEX, EDGE or FACE.
     * \param[in]  ithObj     Number of the object on the parent.
     * \param[in]  elType     Type of the element. Important only in 3D.
     */
    virtual int getSubObjOfChild(int childnr, GeoIndex subObj, int ithObj, 
				 int elType = 0) const = 0;

369 370 371 372 373 374 375
    /** \brief
     * Returns which vertex of elements parent corresponds to the vertexnr of
     * the element, if the element is the childnr-th child of the parent.
     * If the vertex is the ner vertex at the refinement edge, -1 is returned.
     */
    virtual int getVertexOfParent(int childnr, int vertexnr, int elType = 0) const = 0;

376
    /// Returns whether Element is a Line
377 378
    virtual bool isLine() const = 0;

379
    /// Returns whether Element is a Triangle
380 381
    virtual bool isTriangle() const = 0;

382
    /// Returns whether Element is a Tetrahedron
383 384
    virtual bool isTetrahedron() const = 0;

385
    /// Returns whether Element has sideElem as one of its sides.
386 387
    virtual bool hasSide(Element *sideElem) const = 0;

388 389 390 391 392 393 394
    /** \brief
     * Returns for a given element type number the element type number of the children.
     * For 1d and 2d this is always 0, because element type number are used in the 
     * 3d case only.
     */
    virtual int getChildType(int elType) const = 0;

395
    /** \brief
396 397 398
     * Traverses an edge/face of a given element (this includes also all children of 
     * the element having the same edge/face). All vertex dofs alonge this edge/face
     * are assembled and put together to a list. 
399 400
     *
     * \param[in]  feSpace         FE space which is used to get the dofs.
401 402
     * \param[in]  bound           Defines the edge/face of the element on which
     *                             all vertex dofs are assembled.
403 404
     * \param[out] dofs            List of dofs, where the result is stored.
     */
405
    virtual void getVertexDofs(FiniteElemSpace* feSpace, 
406
			       BoundaryObject bound,
407
			       DofContainer& dofs) const = 0;
408 409

    /** \brief
410 411 412
     * Traverses an edge/face of a given element (this includes also all children of 
     * the element having the same edge/face). All non vertex dofs alonge this edge/face 
     * are assembled and put together to a list.
413 414
     *
     * \param[in]  feSpace         FE space which is used to get the dofs.
415 416
     * \param[in]  bound           Defines the edge/face of the element on which
     *                             all non vertex dofs are assembled.
417 418
     * \param[out] dofs            All dofs are put to this dof list.
     */
419
    virtual void getNonVertexDofs(FiniteElemSpace* feSpace, 
420
				  BoundaryObject bound,
421 422 423
				  DofContainer& dofs) const = 0;


424 425 426 427
    /** \} */

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

428
    /// assignment operator
429
    Element& operator=(const Element& el);
430 431 432 433 434 435 436

    /** \brief
     * Checks whether the face with vertices dof[0],..,dof[DIM-1] is
     * part of mel's boundary. returns the opposite vertex if true, -1 else
     */
    int oppVertex(FixVec<DegreeOfFreedom*, DIMEN> pdof) const;

437
    /// Refines Element's leaf data
Thomas Witkowski's avatar
Thomas Witkowski committed
438 439
    inline void refineElementData(Element* child1, Element* child2, int elType = 0) 
    {
440 441 442
      if (elementData) {
	bool remove = elementData->refineElementData(this, child1, child2, elType);
	if (remove) {
443
	  ElementData *tmp = elementData->getDecorated();
Thomas Witkowski's avatar
Thomas Witkowski committed
444
	  delete elementData;
445 446 447
	  elementData = tmp;
	}
      }
448
    }
449

450
    /// Coarsens Element's leaf data
Thomas Witkowski's avatar
Thomas Witkowski committed
451 452
    inline void coarsenElementData(Element* child1, Element* child2, int elType = 0) 
    {
453 454
      ElementData *childData;
      childData = child1->getElementData();
455
      if (childData) {
456
	childData->coarsenElementData(this, child1, child2, elType);
Thomas Witkowski's avatar
Thomas Witkowski committed
457
	delete childData;
458 459 460
	child1->setElementData(NULL);
      }
      childData = child2->getElementData();
461
      if (childData) {
462
	childData->coarsenElementData(this, child2, child1, elType);
Thomas Witkowski's avatar
Thomas Witkowski committed
463
	delete childData;
464 465
	child2->setElementData(NULL);
      }
466
    }
467

468
    /// Returns pointer to \ref elementData
Thomas Witkowski's avatar
Thomas Witkowski committed
469 470
    inline ElementData* getElementData() const 
    {
471
      return elementData;
472
    }
473

474
    ///
Thomas Witkowski's avatar
Thomas Witkowski committed
475 476 477
    inline ElementData* getElementData(int typeID) const 
    {
      if (elementData)
478
	return elementData->getElementData(typeID);
Thomas Witkowski's avatar
Thomas Witkowski committed
479

480
      return NULL;
481
    }
482

483
    /// Deletes the \ref elementData with a specific typeID.
484
    bool deleteElementData(int typeID);
485 486 487 488 489 490 491 492

    /** \brief
     * Returns whether element is refined at side side
     * el1, el2 are the corresponding children. 
     * (not neccessarly the direct children!)
     * elementTyp is the type of this element (comes from ElInfo)
     */
    bool isRefinedAtSide(int side, Element *el1, Element *el2, 
493
			 unsigned char elementTyp = 255);
494

495
    /// Returns whether Element's \ref newCoord is set
Thomas Witkowski's avatar
Thomas Witkowski committed
496 497
    inline bool isNewCoordSet() const 
    { 
498
      return (newCoord != NULL);
499
    }
500

501
    /// Frees memory for \ref newCoord
502
    void eraseNewCoord();
503 504
 
    /// Serialize the element to a file.
505
    void serialize(std::ostream &out);
506

507
    /// Deserialize an element from a file.
508
    void deserialize(std::istream &in);
509

510
    ///
511 512
    int calcMemoryUsage();

513 514 515
    /// Sets Element's \ref dof pointer.
    void createNewDofPtrs();
 
516
  protected:
517
    /// Sets Element's \ref index. Used by friend class Mesh.
Thomas Witkowski's avatar
Thomas Witkowski committed
518 519
    inline void setIndex(int i) 
    {
520
      index = i;
521
    }
522

523
    /// Used by friend class Mesh while dofCompress
524
    void newDofFct1(const DOFAdmin*);
525

526
    /// Used by friend class Mesh while dofCompress
527 528 529 530 531 532 533
    void newDofFct2(const DOFAdmin*);

    /// Changes old dofs to negative new dofs
    void changeDofs1(const DOFAdmin* admin, int n0, int nd0, int nd, int pos);

    /// Changes negative new dofs to positive
    void changeDofs2(int n0, int nd0, int nd, int pos);
534 535 536 537 538 539

  protected:
    /** \brief
     * Pointers to the two children of interior elements of the tree. Pointers
     * to NULL for leaf elements.
     */
540
    Element *child[2];
541 542 543 544 545

    /** \brief
     * Vector of pointers to DOFs. These pointers must be available for elements
     * vertices (for the geometric description of the mesh). There my be pointers
     * for the edges, for faces and for the center of an element. They are 
546 547 548
     * ordered the following way: The first N_VERTICES entries correspond to the
     * DOFs at the vertices of the element. The next ones are those at the edges,
     * if present, then those at the faces, if present, and then those at the 
549 550
     * barycenter, if present.
     */
551
    DegreeOfFreedom **dof;
552 553 554 555 556 557

    /** \brief
     * Unique global index of the element. these indices are not strictly ordered
     * and may be larger than the number of elements in the binary tree (the list
     * of indices may have holes after coarsening).
     */
558
    int index;
559 560 561 562 563 564

    /** \brief
     * Marker for refinement and coarsening. if mark is positive for a leaf
     * element, this element is refined mark times. if mark is negative for
     * a leaf element, this element is coarsened -mark times.
     */
565
    signed char mark;
566 567 568 569 570
 
    /** \brief
     * If the element has a boundary edge on a curved boundary, this is a pointer
     * to the coordinates of the new vertex that is created due to the refinement
     * of the element, otherwise it is a NULL pointer. Thus coordinate 
571 572
     * information can be also produced by the traversal routines in the case of 
     * curved boundary.
573 574 575
     */
    WorldVector<double> *newCoord;

576
    /// Pointer to the Mesh this element belongs to
577
    Mesh* mesh;
578

579
    /// Pointer to Element's leaf data
580
    ElementData* elementData;
581

582 583 584 585 586 587
    /** \brief
     * This map is used for deletion of all DOFs of all elements of a mesh. Once
     * a DOF-vector (all DOFS at a node, edge, etc.) is deleted, its address is
     * added to this map to note not to delete it a second time.
     */
    static std::map<DegreeOfFreedom*, bool> deletedDOFs;
588 589 590

    friend class Mesh;
  };
591 592 593 594 595


  /// Writes the element hierarchie to a Graphviz dot-file. Using the dot-tool from
  /// Graphvis, this dot-file can be converted to a ps-file. Useful for debugging!
  void writeDotFile(Element *el, std::string filename, int maxLevels = -1);
596 597 598 599
}

#endif  // AMDIS_ELEMENT_H