Tetrahedron.h 10.1 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
/******************************************************************************
 *
 * AMDiS - Adaptive multidimensional simulations
 *
 * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
 * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
 *
 * Authors: 
 * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
 *
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 *
 *
 * This file is part of AMDiS
 *
 * See also license.opensource.txt in the distribution.
 * 
 ******************************************************************************/
20
21


22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42

/** \file Tetrahedron.h */

#ifndef AMDIS_TETRAHEDRON_H
#define AMDIS_TETRAHEDRON_H

#include "Element.h"

namespace AMDiS {

  /** \ingroup Triangulation 
   * \brief
   * A Tetrahedron is a 3-dimensional Element. 
   *
   * A Tetrahedron and its refinements:
   *
   * <img src = "tetrahedron.png">
   */
  class Tetrahedron : public Element
  {
  public:
43
44
45
46
47
    
    /// calls base class contructor.
    Tetrahedron(Mesh* aMesh) 
      : Element(aMesh) 
    {}
48

49
    /// implements Element::clone
50
51
    inline Element *clone() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
52
      return new Tetrahedron(mesh); 
53
    }
54

55
    /// implements Element::getVertexOfEdge
56
57
    inline int getVertexOfEdge(int i, int j) const 
    {
58
      return vertexOfEdge[i][j];
59
    }
60

61
    /// implements Element::getVertexOfPosition
62
    int getVertexOfPosition(GeoIndex position,
63
64
			    int positionIndex,
			    int vertexIndex) const;
65
66


67
68
    virtual int getPositionOfVertex(int side, int vertex) const 
    {
69
      static const int positionOfVertex[4][4] = {{-1, 0, 1, 2},
70
71
72
					   {0, -1, 1, 2},
					   {0, 1, -1, 2},
					   {0, 1, 2, -1}};
73
      return positionOfVertex[side][vertex];
74
    }
75

76
    /// implements Element::getGeo
77
78
79
    inline int getGeo(GeoIndex i) const 
    {
      switch (i) {
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
      case VERTEX: case PARTS: case NEIGH:
	return 4;
	break;
      case EDGE:
	return 6;
      case FACE:
	return 4;
      case CENTER:
	return 1;
	break;
      case DIMEN:
	return 3;
	break;
      case BOUNDARY:
	return 14;
	break;
      case PROJECTION:
	return 10;
	break;
      default:
	ERROR_EXIT("invalid geo-index\n");
	return 0;
      }
103
    }
104

105
    /// implements Element::hasSide
106
107
    bool hasSide(Element* sideElem) const;

108
    /// implements Element::sortFaceIndices
109
    void sortFaceIndices(int face, FixVec<int,WORLD> &vec) const;
110

111
    /// Returns false because this element is a Tetrahedron.
112
113
    inline bool isLine() const 
    { 
114
      return false; 
115
    }
116

117
    /// Returns false because this element is a Tetrahedron.
118
119
    inline bool isTriangle() const 
    { 
120
      return false; 
121
    }
122
 
123
    /// Returns true because this element is a Tetrahedron.
124
125
    inline bool isTetrahedron() const 
    { 
126
      return true; 
127
    }
128
129
130
131
132
133

    /// Return the new element type of the children.
    inline int getChildType(int elType) const
    {
      return (elType + 1) % 3;
    }
134
 
135
136
    std::string getTypeName() const 
    { 
137
      return "Tetrahedron"; 
138
    }
139

140
    void getNodeDofs(const FiniteElemSpace* feSpace, 
141
		     BoundaryObject bound,
142
143
		     DofContainer& dofs,
		     bool baseDofPtr = false) const;
144

145
    void getNodeDofsAtFace(const FiniteElemSpace* feSpace, 
146
			   BoundaryObject bound,
147
148
			   DofContainer& dofs,
			   bool baseDofPtr) const;
149

150
    void getNodeDofsAtEdge(const FiniteElemSpace* feSpace, 
151
			   BoundaryObject bound,
152
153
			   DofContainer& dofs,
			   bool baseDofPtr) const;
154

155
    void getHigherOrderDofs(const FiniteElemSpace* feSpace,
156
			    BoundaryObject bound,
157
			    DofContainer& dofs,
158
			    bool baseDofPtr = false,
159
			    vector<GeoIndex>* dofGeoIndex = nullptr) const;
160

161
162
163
    void getSubBoundary(BoundaryObject bound, 
			vector<BoundaryObject> &subBound) const;

164
165
    void prepareNextBound(BoundaryObject &bound, int ithChild) const;

166
  public:
167
168
    /// nChildEdge[el_type][ichild][dir] gives local index of new edge on 
    /// child[ichild] part of face [2 + dir] on the parent.
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
    static const unsigned char nChildEdge[3][2][2];

    /** \brief
     * nChildFace[el_type][ichild][dir]                                     
     * gives local index of sub-face on child[ichild] part of face [2+dir] on
     * the parent   
     */
    static const unsigned char nChildFace[3][2][2];

    /** \brief
     * childVertex[el_type][child][i] =
     * parent's local vertex index of new vertex i.
     * 4 stands for the newly generated vertex                            
     */
    static const int childVertex[3][2][4];
 
    /** \brief
     * childEdge[el_type][child][i] =
     * parent's local edge index of new edge i
     * new edge 2 is half of old edge 0,
     * new edges 4,5 are really new edges, and value is different:
     * childEdge[][][4,5] = index of same edge in other child
191
192
193
     *
     * Note: el_type, i.e., the first index of the array, defines the element type
     * of the parent elements.
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
     */
    static const int childEdge[3][2][6];

    /** \brief
     * adjacentChild[position][ichild]                                      
     * gives number of the adjacent child on a neighbour element:             
     *   position = 0  same position of element and neigh at refinement edge, 
     *   position = 1  different ...                                         
     */
    static const unsigned char adjacentChild[2][2];

    /** \brief
     * childOrientation[el_type][child] =
     *   +1 if orientation is not changed during refinement, 
     *   -1 if orientation is changed during refinement
     */
    static const signed char childOrientation[3][2];

212
    /// edgeOfDOFs[i][j]: gives the local index of edge with vertices i and j 
213
    static const unsigned char edgeOfDofs[4][4];
214

215
    ///
216
217
    static const int edgeOfFace[4][3];

218
    /// implements Element::getSideOfChild()
219
    int getSideOfChild(int child, int side, int elType = 0) const 
220
    {
221
      FUNCNAME_DBG("Tetrahedron::getSideOfChild()");
222
223
224

      TEST_EXIT_DBG(child == 0 || child == 1)("Child must be in (0,1)!\n");
      TEST_EXIT_DBG(side >= 0 && side <= 3)("Side must be between 0 and 3!\n");
225
226
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)
	("ElType must be between 0 and 2!\n");
227

228
      return sideOfChild[elType][child][side];
229
    }
230

231
232
    /// Returns for an edge number its local edge number on a child element. See
    /// \ref edgeOfChild for mor information.
233
234
    inline int getEdgeOfChild(int child, int edge, int elType) const
    {
235
      FUNCNAME_DBG("Tetrahedron::getEdgeOfChild()");
236
237
238

      TEST_EXIT_DBG(child == 0 || child == 1)("Child must be in (0,1)!\n");
      TEST_EXIT_DBG(edge >= 0 && edge <= 5)("Side must be between 0 and 3!\n");
239
240
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)
	("ElType must be between 0 and 2!\n");
241
242
243
244

      return edgeOfChild[elType][child][edge];
    }

245
246
    int getSubObjOfChild(int childnr, GeoIndex subObj, 
			 int ithObj, int elType) const
247
    {
248
      FUNCNAME_DBG("Tetrahedron::getSubObjOfChild()");
249
250
251
252
253
254
255
256
257

      TEST_EXIT_DBG(subObj == EDGE || subObj == FACE)("Not yet implemented!\n");

      if (subObj == FACE)
	return getSideOfChild(childnr, ithObj, elType);
      else
	return getEdgeOfChild(childnr, ithObj, elType);
    }

258
    /// implements Element::getVertexOfParent()
259
    int getVertexOfParent(int child, int side, int elType = 0) const 
260
    {
261
      FUNCNAME_DBG("Tetrahedron::getVertexOfParent()");
262
263
264

      TEST_EXIT_DBG(child == 0 || child == 1)("Child must be in (0,1)!\n");
      TEST_EXIT_DBG(side >= 0 && side <= 3)("Side must be between 0 and 3!\n");
265
266
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)
	("ElType must be between 0 and 2!\n");
267

268
      return vertexOfParent[elType][child][side];
269
    }
270

271
272
    inline int getEdgeOfFace(int face, int edge) const 
    {
273
      FUNCNAME_DBG("Tetrahedron::getEdgeOfFace()");
274
275
276
277

      TEST_EXIT_DBG(face >= 0 && face < 4)("Invalid face number!\n");
      TEST_EXIT_DBG(edge >= 0 && edge < 3)("Invalid edge number!\n");

278
      return edgeOfFace[face][edge];
279
    }
280
   
281
    DofEdge getEdge(int localEdgeIndex) const
282
    {
283
      FUNCNAME_DBG("Tetrahedron::getEdge()");
284
      TEST_EXIT_DBG(localEdgeIndex >= 0 && localEdgeIndex < 6)("Invalid edge!\n");
285
286
287

      DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
      DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
288
      DofEdge edge = std::make_pair(std::min(dof0, dof1), std::max(dof0, dof1));
289
290
291
292

      return edge;
    }

293
294
    DofFace getFace(int localFaceIndex) const
    {
295
      FUNCNAME_DBG("Tetrahedron::getFace()");
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
      TEST_EXIT_DBG(localFaceIndex >= 0 && localFaceIndex < 4)("Invalid face!\n");

      // Get the three DOFs of the face.
      DegreeOfFreedom dof0 = dof[vertexOfFace[localFaceIndex][0]][0];
      DegreeOfFreedom dof1 = dof[vertexOfFace[localFaceIndex][1]][0];
      DegreeOfFreedom dof2 = dof[vertexOfFace[localFaceIndex][2]][0];
      
      // Sort the three DOFs of the face with respect to their values.
      DegreeOfFreedom dofMin0, dofMin1, dofMin2;
      if (dof0 < dof1 && dof0 < dof2) {
	dofMin0 = dof0;
	dofMin1 = std::min(dof1, dof2);
	dofMin2 = std::max(dof1, dof2);
      } else if (dof1 < dof0 && dof1 < dof2) {
	dofMin0 = dof1;
	dofMin1 = std::min(dof0, dof2);
	dofMin2 = std::max(dof0, dof2);
      } else {
	TEST_EXIT_DBG(dof2 < dof0 && dof2 < dof1)("Should not happen!\n");
	dofMin0 = dof2;
	dofMin1 = std::min(dof0, dof1);
	dofMin2 = std::max(dof0, dof1);
      }
      
      return DofFace(dofMin0, dofMin1, dofMin2);
    }

323
324
  protected:

325
    /// vertexOfEdge[i][j] is the local number of the j-vertex of edge i
326
327
    static const int vertexOfEdge[6][2]; 

328
    /// vertexOfFace[i][j] is the local number of the j-vertex of face i
329
330
    static const int vertexOfFace[4][3];

Thomas Witkowski's avatar
Thomas Witkowski committed
331
    /// See \ref Element::getSideOfChild for more information.
332
    static const int sideOfChild[3][2][4];
333

334
335
336
337
    /// edgeOfChild[elType][i][j] is the local edge number of the j-th edge within
    /// the i-th children of an element of elType. If the value is -1, the edge is
    /// not included in the element's child. Note that the 0 edge is included in
    /// both children only by its half.
338
339
    static const int edgeOfChild[3][2][6];

Thomas Witkowski's avatar
Thomas Witkowski committed
340
    /// See \ref Element::getVertexOfParent for more information.
341
342
343
344
345
346
    static const int vertexOfParent[3][2][4];
  };

}

#endif // AMDIS_TETRAHEDRON_H