MeshStructure.cc 4.7 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#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"
#include "mpi.h"

namespace AMDiS {

  const int MeshStructure::unsignedLongSize_ = sizeof(unsigned long int) * 8;

  void MeshStructure::insertElement(bool isLeaf) {
    // overflow? -> next index
17
    if (pos_ >= unsignedLongSize_) {
18
19
20
21
22
23
      code_.push_back(currentCode_);
      pos_ = 0;
      currentCode_ = 0;
    }

    // insert element in binary code
24
    if (!isLeaf) {
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
      unsigned long int one = 1;
      currentCode_ += (one << pos_);
    } 

    pos_++;
    numElements_++;
  }

  void MeshStructure::clear()
  {
    currentCode_ = 0;
    code_.resize(0);
    pos_ = 0;
    numElements_ = 0;
    currentElement_ = 0;
  }

  void MeshStructure::init(Mesh *mesh) 
  {
    clear();
    TraverseStack stack;
    ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
47
    while (elInfo) {
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
      insertElement(elInfo->getElement()->isLeaf());
      elInfo = stack.traverseNext(elInfo);
    }
  
    commit();
  }

  void MeshStructure::reset() {
    currentIndex_ = 0;
    currentCode_ = code_[0];
    pos_ = 0;
    currentElement_ = 0;
  }

  bool MeshStructure::nextElement(MeshStructure *insert)
  {
    if(insert) {
      insert->insertElement(isLeafElement());
    }

    pos_++;
    currentElement_++;

    if(currentElement_ >= numElements_) return false;

    if(pos_ >= unsignedLongSize_) {
      currentIndex_++;
75
      TEST_EXIT_DBG(currentIndex_ < static_cast<int>(code_.size()))
76
77
78
79
80
81
82
83
84
85
86
	("end of structure reached\n");
      pos_ = 0;
      currentCode_ = code_[currentIndex_];
    } else {
      currentCode_ >>= 1;
    }
    return true;
  }

  bool MeshStructure::skipBranch(MeshStructure *insert)
  {
87
    if (isLeafElement()) {
88
89
90
91
      return nextElement(insert);
    } else {
      bool cont = nextElement(insert);
      cont = skipBranch(insert); // left branch
92
      TEST_EXIT_DBG(cont)("invalid structure\n");
93
94
95
96
97
98
99
100
101
102
103
104
105
106
      cont = skipBranch(insert); // righ branch
      return cont;
    }
  }

  void MeshStructure::merge(MeshStructure *structure1,
			    MeshStructure *structure2,
			    MeshStructure *result)
  {
    result->clear();
    structure1->reset();
    structure2->reset();

    bool cont = true;
107
    while (cont) {
108
      bool cont1, cont2;
109
      if (structure1->isLeafElement() == structure2->isLeafElement()) {
110
111
112
	cont1 = structure1->nextElement(result);
	cont2 = structure2->nextElement();
      } else {
113
	if (structure1->isLeafElement()) {
114
115
116
117
118
119
120
	  cont1 = structure1->nextElement();
	  cont2 = structure2->skipBranch(result);
	} else {
	  cont1 = structure1->skipBranch(result);
	  cont2 = structure2->nextElement();
	}
      }
121
      TEST_EXIT_DBG(cont1 == cont2)("structures don't match\n");
122
123
124
125
126
127
128
129
      cont = cont1;
    }

    result->commit();
  }

  void MeshStructure::fitMeshToStructure(Mesh *mesh,
					 RefinementManager *manager,
130
131
					 bool checkPartition,
					 bool debugMode) 
132
133
134
135
136
137
138
139
140
141
142
143
  {
    FUNCNAME("MeshStructure::fitMeshToStructure()");

    TraverseStack stack;
    ElInfo *elInfo;


    bool cont = true;

    // decorate leaf data
    reset();
    elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
144
    while (elInfo) {
145
      TEST_EXIT_DBG(cont)("unexpected structure code end!\n");
146
147
148

      Element *element = elInfo->getElement();

149
      if (isLeafElement()) {
150
	TEST_EXIT_DBG(element->isLeaf())("mesh finer than code\n");
151
152
      };

153
      if (element->isLeaf() && !isLeafElement()) {
154
155
156
157
158
	MeshStructure *structure = NEW MeshStructure();
	cont = skipBranch(structure);
	structure->commit();

	bool decorate = true;
159
	if (checkPartition) {
160
161
	  PartitionElementData *partitionData = dynamic_cast<PartitionElementData*>
	    (element->getElementData(PARTITION_ED));
162
	  TEST_EXIT_DBG(partitionData)("no partition element data\n");
163
	  PartitionStatus status = partitionData->getPartitionStatus();
164
	  if ((debugMode == false) && (status == OUT || status == UNDEFINED)) {
165
166
167
168
	    decorate = false;
	  }
	}

169
	if (decorate) {
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
	  MeshStructure_ED *elData = NEW MeshStructure_ED(element->getElementData());
	  elData->setStructure(structure);
	  element->setElementData(elData);
	} else {
	  DELETE structure;
	}
      } else {
	cont = nextElement();
      }

      elInfo = stack.traverseNext(elInfo);
    }

    // refine mesh
    bool finished;

    do {
      finished = true;
      elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL);
189
      while (elInfo) {
190
	Element *element = elInfo->getElement();
191
	if (element->getElementData(MESH_STRUCTURE) != NULL) {
192
193
194
195
196
197
198
199
	  element->setMark(1);
	  finished = false;
	} else {
	  element->setMark(0);
	}
	elInfo = stack.traverseNext(elInfo);
      }
      manager->refineMesh(mesh);
200
    } while (!finished);
201
202
203
  }

}