MeshStructure.cc 5.68 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
#include "MeshStructure.h"
#include "MeshStructure_ED.h"
#include "PartitionElementData.h"
#include "Mesh.h"
#include "Element.h"
#include "Traverse.h"
#include "ElInfo.h"
#include "RefinementManager.h"

namespace AMDiS {

12
  const int MeshStructure::unsignedLongSize = sizeof(unsigned long int) * 8;
13

14
15
  void MeshStructure::insertElement(bool isLeaf) 
  {
16
    // overflow? -> next index
17
18
19
20
    if (pos >= unsignedLongSize) {
      code.push_back(currentCode);
      pos = 0;
      currentCode = 0;
21
22
23
    }

    // insert element in binary code
24
    if (!isLeaf) {
25
      unsigned long int one = 1;
26
      currentCode += (one << pos);
27
28
    } 

29
30
    pos++;
    nElements++;
31
32
33
34
  }

  void MeshStructure::clear()
  {
35
36
37
38
39
    currentCode = 0;
    code.resize(0);
    pos = 0;
    nElements = 0;
    currentElement = 0;
40
41
42
43
44
  }

  void MeshStructure::init(Mesh *mesh) 
  {
    clear();
45

46
47
    TraverseStack stack;
    ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
48
    while (elInfo) {
49
50
51
52
53
54
55
      insertElement(elInfo->getElement()->isLeaf());
      elInfo = stack.traverseNext(elInfo);
    }
  
    commit();
  }

56
  void MeshStructure::init(Element *el, int ithSide, int elType, bool reverseOrder)
57
58
59
60
61
  {
    FUNCNAME("MeshStructure::init()");

    clear();

62
    addAlongSide(el, ithSide, elType, reverseOrder);
63
64
65
66

    commit();    
  }

67
68
  void MeshStructure::addAlongSide(Element *el, int ithSide, int elType, 
				   bool reverseOrder)
69
70
71
72
73
74
75
  {
    insertElement(el->isLeaf());
    
    if (!el->isLeaf()) {
      int s1 = el->getSideOfChild(0, ithSide, elType);
      int s2 = el->getSideOfChild(1, ithSide, elType);

76
77
78
79
80
81
82
83
84
85
86
      if (!reverseOrder) {
	if (s1 != -1) 
	  addAlongSide(el->getFirstChild(), s1, el->getChildType(elType), reverseOrder);
	if (s2 != -1)
	  addAlongSide(el->getSecondChild(), s2, el->getChildType(elType), reverseOrder);
      } else {
	if (s2 != -1)
	  addAlongSide(el->getSecondChild(), s2, el->getChildType(elType), reverseOrder);
	if (s1 != -1) 
	  addAlongSide(el->getFirstChild(), s1, el->getChildType(elType), reverseOrder);
      }
87
88
89
    }
  }

90
91
92
93
94
95
  void MeshStructure::reset() 
  {
    currentIndex = 0;
    currentCode = code[0];
    pos = 0;
    currentElement = 0;
96
97
98
99
  }

  bool MeshStructure::nextElement(MeshStructure *insert)
  {
100
101
102
    FUNCNAME("MeshStructure::nextElement()");

    if (insert)
103
104
      insert->insertElement(isLeafElement());

105
106
    pos++;
    currentElement++;
107

108
109
    if (currentElement >= nElements) 
      return false;
110

111
112
113
114
115
116
    if (pos >= unsignedLongSize) {
      currentIndex++;
      TEST_EXIT_DBG(currentIndex < static_cast<int>(code.size()))
	("End of structure reached!\n");
      pos = 0;
      currentCode = code[currentIndex];
117
    } else {
118
      currentCode >>= 1;
119
    }
120

121
122
123
124
125
    return true;
  }

  bool MeshStructure::skipBranch(MeshStructure *insert)
  {
126
127
    FUNCNAME("MeshStructure::skipBranch()");

128
    if (isLeafElement()) {
129
130
131
132
      return nextElement(insert);
    } else {
      bool cont = nextElement(insert);
      cont = skipBranch(insert); // left branch
133
      TEST_EXIT_DBG(cont)("Invalid structure!\n");
134
135
136
137
138
139
140
141
142
      cont = skipBranch(insert); // righ branch
      return cont;
    }
  }

  void MeshStructure::merge(MeshStructure *structure1,
			    MeshStructure *structure2,
			    MeshStructure *result)
  {
143
144
    FUNCNAME("MeshStructure::merge()");

145
146
147
148
149
    result->clear();
    structure1->reset();
    structure2->reset();

    bool cont = true;
150
    while (cont) {
151
      bool cont1, cont2;
152
      if (structure1->isLeafElement() == structure2->isLeafElement()) {
153
154
155
	cont1 = structure1->nextElement(result);
	cont2 = structure2->nextElement();
      } else {
156
	if (structure1->isLeafElement()) {
157
158
159
160
161
162
163
	  cont1 = structure1->nextElement();
	  cont2 = structure2->skipBranch(result);
	} else {
	  cont1 = structure1->skipBranch(result);
	  cont2 = structure2->nextElement();
	}
      }
164
      TEST_EXIT_DBG(cont1 == cont2)("Structures don't match!\n");
165
166
167
168
169
170
171
172
      cont = cont1;
    }

    result->commit();
  }

  void MeshStructure::fitMeshToStructure(Mesh *mesh,
					 RefinementManager *manager,
173
174
					 bool checkPartition,
					 bool debugMode) 
175
176
177
178
179
180
181
  {
    FUNCNAME("MeshStructure::fitMeshToStructure()");

    bool cont = true;

    // decorate leaf data
    reset();
182
183
    TraverseStack stack;
    ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
184
    while (elInfo) {
185
      TEST_EXIT(cont)("unexpected structure code end!\n");
186
187
188

      Element *element = elInfo->getElement();

189
190
      if (isLeafElement())
	TEST_EXIT_DBG(element->isLeaf())("mesh finer than code\n");
191

192
      if (element->isLeaf() && !isLeafElement()) {
Thomas Witkowski's avatar
Thomas Witkowski committed
193
	MeshStructure *structure = new MeshStructure();
194
195
196
197
	cont = skipBranch(structure);
	structure->commit();

	bool decorate = true;
198
	if (checkPartition) {
199
200
	  PartitionElementData *partitionData = dynamic_cast<PartitionElementData*>
	    (element->getElementData(PARTITION_ED));
201
	  TEST_EXIT_DBG(partitionData)("no partition element data\n");
202
	  PartitionStatus status = partitionData->getPartitionStatus();
203
	  if (debugMode == false && (status == OUT || status == UNDEFINED))
204
205
206
	    decorate = false;
	}

207
	if (decorate) {
Thomas Witkowski's avatar
Thomas Witkowski committed
208
	  MeshStructure_ED *elData = new MeshStructure_ED(element->getElementData());
209
210
211
	  elData->setStructure(structure);
	  element->setElementData(elData);
	} else {
Thomas Witkowski's avatar
Thomas Witkowski committed
212
	  delete structure;
213
214
215
216
217
218
219
220
221
222
223
224
225
226
	}
      } else {
	cont = nextElement();
      }

      elInfo = stack.traverseNext(elInfo);
    }

    // refine mesh
    bool finished;

    do {
      finished = true;
      elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL);
227
      while (elInfo) {
228
	Element *element = elInfo->getElement();
229
	if (element->getElementData(MESH_STRUCTURE) != NULL) {
230
231
232
233
234
235
236
237
	  element->setMark(1);
	  finished = false;
	} else {
	  element->setMark(0);
	}
	elInfo = stack.traverseNext(elInfo);
      }
      manager->refineMesh(mesh);
238
    } while (!finished);
239
240
241
  }

}