Triangle.cc 6 KB
Newer Older
1
2
3
4
5
#include "Triangle.h"
#include "DOFAdmin.h"
#include "Mesh.h"
#include "CoarseningManager.h"
#include "FixVec.h"
6
#include "ElementDofIterator.h"
7
8
9

namespace AMDiS {

Thomas Witkowski's avatar
Thomas Witkowski committed
10
11
12
  const int Triangle::vertexOfEdge[3][2] = {{1, 2}, {2, 0}, {0, 1}};
  const int Triangle::sideOfChild[2][3] = {{-1, 2, 0}, {2, -1, 1}};
  const int Triangle::vertexOfParent[2][3] = {{2, 0, -1}, {1, 2, -1}};
13
14
15
16

  bool Triangle::hasSide(Element* sideElem) const
  {
    FUNCNAME("Triangle::hasSide");
17
    TEST_EXIT_DBG(sideElem->isLine())("called for sideElem-type != Line\n");
18
19
20
21
    ERROR_EXIT("not yet\n");
    return false;
  }

22
23
  int Triangle::getVertexOfPosition(GeoIndex position, int positionIndex,
				    int vertexIndex) const 
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
  {
    FUNCNAME("Triangle::getVertexOfPosition");
    switch(position) {
    case VERTEX:
      return positionIndex;
      break;
    case EDGE:
      return vertexOfEdge[positionIndex][vertexIndex];
      break;
    case CENTER:
      return vertexIndex;
      break;
    default:
      ERROR_EXIT("invalid position\n");
      return 0;
    }
  }

42
43
  const FixVec<int, WORLD>& Triangle::sortFaceIndices(int face, 
						      FixVec<int,WORLD> *vec) const
44
45
46
47
48
  {
    static MatrixOfFixVecs<FixVec<int,WORLD> > *sorted_2d = NULL;

    int no = 0;
    FixVec<int,WORLD> *val = NULL;
49
    const int *vof = vertexOfEdge[face];
50

51
52
    if (NULL == sorted_2d) {
      sorted_2d = new MatrixOfFixVecs<FixVec<int,WORLD> >(2, 3, 2, NO_INIT);
53

54
55
56
57
58
59
      (*sorted_2d)[1][0][1] = (*sorted_2d)[1][1][0] =
	(*sorted_2d)[2][0][0] = (*sorted_2d)[2][1][1] = 0;
      (*sorted_2d)[0][0][0] = (*sorted_2d)[0][1][1] =
	(*sorted_2d)[2][0][1] = (*sorted_2d)[2][1][0] = 1;
      (*sorted_2d)[0][0][1] = (*sorted_2d)[0][1][0] =
	(*sorted_2d)[1][0][0] = (*sorted_2d)[1][1][1] = 2;
60
61
62
63
64
65
66
67
68
69
    }

    if (dof[vof[0]][0] < dof[vof[1]][0])
      no = 0;
    else
      no = 1;

    if (vec) {
      val = vec;
      *val = (*sorted_2d)[face][no];
70
    } else {
71
      val = &((*sorted_2d)[face][no]);
72
    }
73

74
    return *(const_cast<const FixVec<int,WORLD>* >(val));
75
76
  }

77

78
  void Triangle::getVertexDofs(FiniteElemSpace* feSpace, 
79
			       BoundaryObject bound,
80
			       DofContainer& dofs) const
81
82
83
  {
    FUNCNAME("Triangle::getVertexDofs()");

84
85
86
87
88
89
90
    if (bound.subObj == VERTEX) {
      dofs.push_back(dof[bound.ithObj]);
      return;
    }

    TEST_EXIT_DBG(bound.subObj == EDGE)("This should not happen!\n");

91
92
93
    BoundaryObject nextBound = bound;

    switch (bound.ithObj) {
94
95
    case 0: 
      {
96
97
98
99
100
101
102
103
104
105
106
107
108
109
	if (child[1] && child[1]->getFirstChild()) {
 	  if (bound.reverseMode) {
 	    nextBound.ithObj = 1;
 	    child[1]->getSecondChild()->getVertexDofs(feSpace, nextBound, dofs);
 	    dofs.push_back(child[1]->getFirstChild()->getDOF(2));	    
 	    nextBound.ithObj = 0;
 	    child[1]->getFirstChild()->getVertexDofs(feSpace, nextBound, dofs);	    
 	  } else {
	    nextBound.ithObj = 0;
	    child[1]->getFirstChild()->getVertexDofs(feSpace, nextBound, dofs);	    
	    dofs.push_back(child[1]->getFirstChild()->getDOF(2));	    
	    nextBound.ithObj = 1;
	    child[1]->getSecondChild()->getVertexDofs(feSpace, nextBound, dofs);
	  }
110
111
112
113
114
	}
      }
      break;
    case 1:
      {
115
116
117
118
119
120
121
122
123
124
125
126
127
128
	if (child[0] && child[0]->getFirstChild()) {
	  if (bound.reverseMode) {
	    nextBound.ithObj = 1;
	    child[0]->getSecondChild()->getVertexDofs(feSpace, nextBound, dofs);
	    dofs.push_back(child[0]->getFirstChild()->getDOF(2));
	    nextBound.ithObj = 0;
	    child[0]->getFirstChild()->getVertexDofs(feSpace, nextBound, dofs);
	  } else {
	    nextBound.ithObj = 0;
	    child[0]->getFirstChild()->getVertexDofs(feSpace, nextBound, dofs);
	    dofs.push_back(child[0]->getFirstChild()->getDOF(2));
	    nextBound.ithObj = 1;
	    child[0]->getSecondChild()->getVertexDofs(feSpace, nextBound, dofs);
	  }
129
130
131
132
133
	}
      }
      break;
    case 2:      
      if (child[0]) {
134
135
136
137
138
139
140
141
142
143
144
145
146
	if (bound.reverseMode) {
	  nextBound.ithObj = 1;
	  child[1]->getVertexDofs(feSpace, nextBound, dofs);
	  dofs.push_back(child[0]->getDOF(2));	  
	  nextBound.ithObj = 0;
	  child[0]->getVertexDofs(feSpace, nextBound, dofs);
	} else {
	  nextBound.ithObj = 0;
	  child[0]->getVertexDofs(feSpace, nextBound, dofs);	  
	  dofs.push_back(child[0]->getDOF(2));	  
	  nextBound.ithObj = 1;
	  child[1]->getVertexDofs(feSpace, nextBound, dofs);
	}
147
148
149
150
151
152
153
154
      }
      break;      
    default:
      ERROR_EXIT("Should never happen!\n");
    }
  }


155
  void Triangle::getNonVertexDofs(FiniteElemSpace* feSpace, 
156
				  BoundaryObject bound,
157
158
159
160
				  DofContainer& dofs) const
  {
    FUNCNAME("Triange::getNonVertexDofs()");

161
162
163
164
165
    if (bound.subObj == VERTEX)
      return;

    TEST_EXIT_DBG(bound.subObj == EDGE)("This should not happen!\n");

166
    bool addThisEdge = false;
167
    BoundaryObject nextBound = bound;
168

169
    switch (bound.ithObj) {
170
    case 0:
171
172
173
174
      if (child[1]) {
	nextBound.ithObj = 2;
	child[1]->getNonVertexDofs(feSpace, nextBound, dofs);
      } else {
175
	addThisEdge = true;
176
      }
177
178
179

      break;
    case 1:
180
181
182
183
      if (child[0]) {
	nextBound.ithObj = 2;
	child[0]->getNonVertexDofs(feSpace, nextBound, dofs);
      } else {
184
	addThisEdge = true;
185
      }
186
187
188
189

      break;
    case 2:
      if (child[0]) {
190
191
192
193
194
195
196
197
198
199
200
	if (bound.reverseMode) {
	  nextBound.ithObj = 1;
	  child[1]->getNonVertexDofs(feSpace, nextBound, dofs);
	  nextBound.ithObj = 0;
	  child[0]->getNonVertexDofs(feSpace, nextBound, dofs);
	} else {
	  nextBound.ithObj = 0;
	  child[0]->getNonVertexDofs(feSpace, nextBound, dofs);
	  nextBound.ithObj = 1;
	  child[1]->getNonVertexDofs(feSpace, nextBound, dofs);
	}
201
202
203
204
205
206
207
208
209
210
      } else {
	addThisEdge = true;
      } 

      break;
    default:
      ERROR_EXIT("Should never happen!\n");
    }

    if (addThisEdge) {
211
      DofContainer addDofs;
212
213
214
215
      ElementDofIterator elDofIter(feSpace, true);
      elDofIter.reset(this);
      do {
	if (elDofIter.getCurrentPos() == 1 && 
216
	    elDofIter.getCurrentElementPos() == bound.ithObj)
217
	  addDofs.push_back(elDofIter.getDofPtr());	
218
      } while(elDofIter.next());      
219
220
221
222
223
224
225
226


      if (bound.reverseMode)
 	for (int i = addDofs.size() - 1; i >= 0; i--)
 	  dofs.push_back(addDofs[i]);
      else
 	for (unsigned int i = 0; i < addDofs.size(); i++) 
 	  dofs.push_back(addDofs[i]);      
227
228
229
230
    }
  }


231
}