Tetrahedron.h 10 KB
Newer Older
1
2
3
4
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
5
// ==  http://www.amdis-fem.org                                              ==
6
7
// ==                                                                        ==
// ============================================================================
8
9
10
11
12
13
14
15
16
17
18
19
//
// Software License for AMDiS
//
// Copyright (c) 2010 Dresden University of Technology 
// All rights reserved.
// Authors: Simon Vey, Thomas Witkowski et al.
//
// 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

/** \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:
41
42
43
44
45
    
    /// calls base class contructor.
    Tetrahedron(Mesh* aMesh) 
      : Element(aMesh) 
    {}
46

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

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

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


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

74
    /// implements Element::getGeo
75
76
77
    inline int getGeo(GeoIndex i) const 
    {
      switch (i) {
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
      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;
      }
101
    }
102

103
    /// implements Element::hasSide
104
105
    bool hasSide(Element* sideElem) const;

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

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

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

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

139
    void getNodeDofs(const FiniteElemSpace* feSpace, 
140
141
		     BoundaryObject bound,
		     DofContainer& dofs) const;
142

143
    void getNodeDofsAtFace(const FiniteElemSpace* feSpace, 
144
145
			   BoundaryObject bound,
			   DofContainer& dofs) const;
146

147
    void getNodeDofsAtEdge(const FiniteElemSpace* feSpace, 
148
149
			   BoundaryObject bound,
			   DofContainer& dofs) const;
150

151
    void getHigherOrderDofs(const FiniteElemSpace* feSpace,
152
153
			    BoundaryObject bound,
			    DofContainer& dofs) const;
154

155
156
    void prepareNextBound(BoundaryObject &bound, int ithChild) const;

157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
  public:
    /** \brief
     * nChildEdge[el_type][ichild][dir]                                      
     * gives local index of new edge on child[ichild] part of face [2+dir] on
     * the parent  
     */
    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
185
186
187
     *
     * Note: el_type, i.e., the first index of the array, defines the element type
     * of the parent elements.
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
     */
    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];

206
    /// edgeOfDOFs[i][j]: gives the local index of edge with vertices i and j 
207
    static const unsigned char edgeOfDofs[4][4];
208

209
    ///
210
211
    static const int edgeOfFace[4][3];

212
    /// implements Element::getSideOfChild()
213
    int getSideOfChild(int child, int side, int elType = 0) const 
214
    {
215
      FUNCNAME("Tetrahedron::getSideOfChild()");
216
217
218
219
220

      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");
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)("ElType must be between 0 and 2!\n");

221
      return sideOfChild[elType][child][side];
222
    }
223

224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
    /** \brief
     * Returns for an edge number its local edge number on a child element. See
     * \ref edgeOfChild for mor information.
     */
    inline int getEdgeOfChild(int child, int edge, int elType) const
    {
      FUNCNAME("Tetrahedron::getEdgeOfChild()");

      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");
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)("ElType must be between 0 and 2!\n");

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

    int getSubObjOfChild(int childnr, GeoIndex subObj, int ithObj, int elType) const
    {
      FUNCNAME("Tetrahedron::getSubObjOfChild()");

      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);
    }

251
    /// implements Element::getVertexOfParent()
252
    int getVertexOfParent(int child, int side, int elType = 0) const 
253
    {
254
      FUNCNAME("Tetrahedron::getVertexOfParent()");
255
256
257
258
259

      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");
      TEST_EXIT_DBG(elType >= 0 && elType <= 2)("ElType must be between 0 and 2!\n");

260
      return vertexOfParent[elType][child][side];
261
    }
262

263
264
    inline int getEdgeOfFace(int face, int edge) const 
    {
265
266
267
268
269
      FUNCNAME("Tetrahedron::getEdgeOfFace()");

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

270
      return edgeOfFace[face][edge];
271
    }
272
   
273
    DofEdge getEdge(int localEdgeIndex) const
274
275
    {
      FUNCNAME("Tetrahedron::getEdge()");
276
      TEST_EXIT_DBG(localEdgeIndex >= 0 && localEdgeIndex < 6)("Invalid edge!\n");
277
278
279

      DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
      DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
280
      DofEdge edge = std::make_pair(std::min(dof0, dof1), std::max(dof0, dof1));
281
282
283
284

      return edge;
    }

285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
    DofFace getFace(int localFaceIndex) const
    {
      FUNCNAME("Tetrahedron::getFace()");
      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);
    }

315
316
  protected:

317
    /// vertexOfEdge[i][j] is the local number of the j-vertex of edge i
318
319
    static const int vertexOfEdge[6][2]; 

320
    /// vertexOfFace[i][j] is the local number of the j-vertex of face i
321
322
    static const int vertexOfFace[4][3];

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

326
327
328
329
330
331
332
333
    /** \brief
     * 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.
     */
    static const int edgeOfChild[3][2][6];

Thomas Witkowski's avatar
Thomas Witkowski committed
334
    /// See \ref Element::getVertexOfParent for more information.
335
336
337
338
339
340
    static const int vertexOfParent[3][2][4];
  };

}

#endif // AMDIS_TETRAHEDRON_H