Liebe Gitlab-Nutzerin, lieber Gitlab-Nutzer,
es ist nun möglich sich mittels des ZIH-Logins/LDAP an unserem Dienst anzumelden. Die Konten der externen Nutzer:innen sind über den Reiter "Standard" erreichbar.
Die Administratoren


Dear Gitlab user,
it is now possible to log in to our service using the ZIH login/LDAP. The accounts of external users can be accessed via the "Standard" tab.
The administrators

DataCollector.cc 13.4 KB
Newer Older
1 2 3 4 5 6 7 8
#include "DataCollector.h"
#include "Traverse.h"
#include "DOFVector.h"
#include "SurfaceRegion_ED.h"
#include "ElementRegion_ED.h"
#include "Projection.h"

namespace AMDiS {
9

10 11 12 13 14 15 16 17 18
  DataCollector::DataCollector(const FiniteElemSpace *fe,
			       DOFVector<double> *val,
			       int l,
			       Flag flag,
			       bool (*writeElemFct)(ElInfo*))
    : values(val),
      level(l),
      traverseFlag(flag),
      feSpace(fe),
19 20 21 22 23 24 25
      nVertices(0),
      nElements(0),
      nInterpPoints(0),
      nConnection(0),
      elementDataCollected(false),
      valueDataCollected(false),
      periodicDataCollected(false),
26
      writeElem(writeElemFct)
27
  {
28
    FUNCNAME("DataCollector::DataCollector()");
29

30
    TEST_EXIT(feSpace)("No finite elem space defined!\n");
31
   
32
    // === get mesh  ===    
33
    mesh = feSpace->getMesh();
34 35

    // === get admin ===    
36
    localAdmin = const_cast<DOFAdmin*>(feSpace->getAdmin());
37
    // === create vertex info vector ===
38
    vertexInfos = new DOFVector<std::list<VertexInfo> >(feSpace, "vertex infos");
39

40 41 42
    interpPointInd = new DOFVector<int>(feSpace, "interpolation point indices");
    interpPointCoords = new DOFVector< std::list<WorldVector<double> > >(feSpace, "interpolation point coordinates");
    dofCoord = new DOFVector< std::list<WorldVector<double> > >(feSpace, "dof coords");
43

44
    dim = mesh->getDim();    
45
    nPreDofs = localAdmin->getNumberOfPreDOFs(VERTEX);
Thomas Witkowski's avatar
Thomas Witkowski committed
46
  } 
47

48

49 50
  DataCollector::~DataCollector()
  {
51 52 53 54
    delete vertexInfos;
    delete interpPointInd;
    delete interpPointCoords;
    delete dofCoord;
55 56
  }

57

Thomas Witkowski's avatar
Thomas Witkowski committed
58 59
  void DataCollector::fillAllData()
  {
60
    if (!elementDataCollected)
Thomas Witkowski's avatar
Thomas Witkowski committed
61
      startCollectingElementData();
62

63
    if (!periodicDataCollected)
Thomas Witkowski's avatar
Thomas Witkowski committed
64
      startCollectingPeriodicData();
65

66
    if (!valueDataCollected)
Thomas Witkowski's avatar
Thomas Witkowski committed
67 68 69
      startCollectingValueData();
  }

70

71
  void DataCollector::startCollectingElementData()
72
  {
73
    FUNCNAME("DataCollector::startCollectingElementData()");
74

75
    Flag flag = traverseFlag;
76 77 78 79 80 81
    flag |= 
      Mesh::FILL_NEIGH      | 
      Mesh::FILL_COORDS     |
      Mesh::FILL_OPP_COORDS |
      Mesh::FILL_BOUND;

82
    elements.resize(0, dim);
83 84 85 86

    TraverseStack stack;

    // Traverse elements to create continuous element indices
87
    ElInfo *elInfo = stack.traverseFirst(mesh, level, flag);
88
    while (elInfo) {
89 90
      if (!writeElem || writeElem(elInfo))
	outputIndices[elInfo->getElement()->getIndex()] = nElements++;
91 92 93 94
      elInfo = stack.traverseNext(elInfo);
    }

    // Traverse elements to create element information
95
    elInfo = stack.traverseFirst(mesh, level, flag);
96

97
    while (elInfo) {
98
      if (!writeElem || writeElem(elInfo))
99
	addElementData(elInfo);
100

101 102 103
      elInfo = stack.traverseNext(elInfo);
    }

104
    elementDataCollected = true;
105 106
  }

107

108
  void DataCollector::startCollectingValueData()
109
  {
110
    FUNCNAME("DataCollector::startCollectingValueData()");
111

112
    DOFVector<int>::Iterator intPointIt(interpPointInd, USED_DOFS);
113
    for (intPointIt.reset(); !intPointIt.end(); ++intPointIt)
114 115
      (*intPointIt) = -1;

116
    interpPoints.clear();
117

118 119 120
    basisFcts = const_cast<BasisFunction*>(feSpace->getBasisFcts());
    nBasisFcts = basisFcts->getNumber();
    localDOFs = new DegreeOfFreedom[nBasisFcts];
121
  
122 123 124 125
    TraverseStack stack;

    // Traverse elements to add value information and to mark
    // interpolation points.
126 127
    ElInfo *elInfo = stack.traverseFirst(mesh, level, 
					 traverseFlag | Mesh::FILL_COORDS);
128
    while (elInfo) {
129
      if (!writeElem || writeElem(elInfo))
Thomas Witkowski's avatar
Thomas Witkowski committed
130
	addValueData(elInfo);
131 132 133
      elInfo = stack.traverseNext(elInfo);
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
134 135
    // If there are interpolation points, add them to the corresponding 
    // data array.
136
    if (nInterpPoints > 0) {
Thomas Witkowski's avatar
Thomas Witkowski committed
137 138 139
      // Remove all interpolation marks and, instead, set to each
      // interpolation point its continous index starting from 0.
      int i = 0;
140 141
      for (intPointIt.reset(); !intPointIt.end(); ++intPointIt)
	if (*intPointIt == -3)
Thomas Witkowski's avatar
Thomas Witkowski committed
142 143 144
	  *intPointIt = i++;
      
      // Traverse elements to create interpolation values.
145
      elInfo = stack.traverseFirst(mesh, level, traverseFlag | Mesh::FILL_COORDS);
Thomas Witkowski's avatar
Thomas Witkowski committed
146
      while (elInfo) {
147
	if (!writeElem || writeElem(elInfo))
Thomas Witkowski's avatar
Thomas Witkowski committed
148 149 150
	  addInterpData(elInfo);
	elInfo = stack.traverseNext(elInfo);
      }
151
    }
152
   
153 154
    delete [] localDOFs;
    valueDataCollected = true;
155 156
  }

157

158
  void DataCollector::startCollectingPeriodicData()
159
  {    
160
    FUNCNAME("DataCollector::startCollectingPeriodicData()");
161

162
    periodicConnections.clear();
163 164
    
    TraverseStack stack;
165
    ElInfo *elInfo = stack.traverseFirst(mesh, level, traverseFlag);
166
    while (elInfo) {
167
      if (!writeElem || writeElem(elInfo)) {
168 169 170 171 172
	LeafDataPeriodic *ldp = dynamic_cast<LeafDataPeriodic*>
	  (elInfo->getElement()->
	   getElementData()->
	   getElementData(PERIODIC));
	
173
	if (ldp)
174
	  nConnection += dynamic_cast<LeafDataPeriodic*>(ldp)->getInfoList().size();
175

176
	periodicConnections.push_back(DimVec<bool>(dim, DEFAULT_VALUE, false));
177 178 179 180
      }
      elInfo = stack.traverseNext(elInfo);
    }   

181
    nConnection /= 2;
182

183
    periodicInfos.clear();
184

185
    Flag flag = traverseFlag;
186 187 188 189 190 191
    flag |= 
      Mesh::FILL_COORDS    |
      Mesh::FILL_OPP_COORDS|
      Mesh::FILL_NEIGH     | 
      Mesh::FILL_BOUND;

192
    elInfo = stack.traverseFirst(mesh, level, flag);
193
    while (elInfo) {
194
      if (!writeElem || writeElem(elInfo))
195 196 197 198
	addPeriodicData(elInfo);
      elInfo = stack.traverseNext(elInfo);
    }   

199
    periodicDataCollected = true;
200 201
  }

202

203
  void DataCollector::addElementData(ElInfo* elInfo)
204
  {
205
    FUNCNAME("DataCollector::addElementData()");
206

207
    const DegreeOfFreedom **dof = elInfo->getElement()->getDof();
208 209
    DegreeOfFreedom vertexDOF;
    WorldVector<double> vertexCoords;
210

211 212
    // MSG("ELEMENT: %d\n", elInfo->getElement()->getIndex());
    // MSG("DOFs: %d %d %d\n", dof[0][0], dof[1][0], dof[2][0]);
213

214
    // create ElementInfo
215
    ElementInfo elementInfo(dim);
216 217
    
    // read element region
218
    ElementData *ed = elInfo->getElement()->getElementData(ELEMENT_REGION);
219
    
Thomas Witkowski's avatar
Thomas Witkowski committed
220
    if (ed)
221
      elementInfo.elementRegion = dynamic_cast<ElementRegion_ED*>(ed)->getRegion();
Thomas Witkowski's avatar
Thomas Witkowski committed
222 223
    else 
      elementInfo.elementRegion = -1;    
224
   
225 226 227
    // read surface regions to element info
    ed = elInfo->getElement()->getElementData(SURFACE_REGION);
    elementInfo.surfaceRegions.set(-1);
228
    while (ed) {
229 230 231 232 233 234
      SurfaceRegion_ED *sr = dynamic_cast<SurfaceRegion_ED*>(ed);
      elementInfo.surfaceRegions[sr->getSide()] = sr->getRegion();
      ed = ed->getDecorated(SURFACE_REGION);
    }

    // for all vertices
235
    for (int i = 0; i < dim + 1; i++) {
236 237 238 239
      // get coords of this vertex
      vertexCoords = elInfo->getCoord(i);
      
      // get dof index of this vertex
240
      vertexDOF = dof[i][nPreDofs];
241

242
      // search for coords at this dof
Thomas Witkowski's avatar
Thomas Witkowski committed
243
      std::list<VertexInfo>::iterator it =
244
	  find((*vertexInfos)[vertexDOF].begin(), (*vertexInfos)[vertexDOF].end(),
245 246 247
	       vertexCoords);
      
      // coords not yet in list?
248
      if (it == (*vertexInfos)[vertexDOF].end()) {
249
	// create new vertex info and increment number of vertices
250
	VertexInfo newVertexInfo = {vertexCoords, nVertices};
251 252

	// add new vertex info to list
253
	(*vertexInfos)[vertexDOF].push_front(newVertexInfo);
254

255
	// set iterator to new vertex info
256
	it = (*vertexInfos)[vertexDOF].begin();	
257

258
	nVertices++;
259 260 261 262 263 264 265 266
      }
      
      // fill element info
      elementInfo.vertexInfo[i] = it;
      elementInfo.boundary[i] = elInfo->getBoundary(i);
      elementInfo.projection[i] = elInfo->getProjection(i);
      elementInfo.neighbour[i] = 
	elInfo->getNeighbour(i) ?
267
	outputIndices[elInfo->getNeighbour(i)->getIndex()] :
268 269 270
	-1;
    }
    
271
    if (dim == 3)
272
      elementInfo.type = elInfo->getType();
273 274

    // remember element info
275
    elements.push_back(elementInfo);
276 277
  }

278

279
  void DataCollector::addValueData(ElInfo *elInfo)
280
  {
281
    FUNCNAME("DataCollector::addValueData()");
282
    
283
    WorldVector<double> vertexCoords;
284
    basisFcts->getLocalIndices(elInfo->getElement(), localAdmin, localDOFs);
285 286

    // First, traverse all DOFs at the vertices of the element, determine
287 288 289
    // their coordinates and add them to the corresponding entry in dofCoords.
    for (int i = 0; i < mesh->getGeo(VERTEX); i++) {
      DegreeOfFreedom dofi = localDOFs[i];
Thomas Witkowski's avatar
Thomas Witkowski committed
290

291
      (*interpPointInd)[dofi] = -2; // mark as vertex
292 293
      
      // get coords of this vertex
294
      vertexCoords = elInfo->getCoord(i);
295

296
      // search for coords at this dof
Thomas Witkowski's avatar
Thomas Witkowski committed
297
      std::list<WorldVector<double> >::iterator it =
298 299
	  find((*dofCoord)[dofi].begin(),
	       (*dofCoord)[dofi].end(),
300
	       vertexCoords);
301 302
      
      // coords not yet in list?
303
      if (it == (*dofCoord)[dofi].end()) {
304
	// add new coords to list
305
	(*dofCoord)[dofi].push_back(vertexCoords);
306 307
      }
    }
308
   
309 310
    // Then, traverse all interpolation DOFs of the element, determine
    // their coordinates and add them to the corresponding entry in 
311 312 313
    // interpPointCoords.
    for (int i = mesh->getGeo(VERTEX); i < nBasisFcts; i++) {
      DegreeOfFreedom dofi = localDOFs[i];
314

315
      elInfo->coordToWorld(*basisFcts->getCoords(i), vertexCoords);
316
      
317
      if ((*interpPointInd)[dofi] == -1) {
318
	// mark as interpolation point
319
	(*interpPointInd)[dofi] = -3; 
320 321 322
	
	// search for interpolation point coordinates, and insert them to the 
	// dof-entry, if not contained in the list
Thomas Witkowski's avatar
Thomas Witkowski committed
323
	std::list<WorldVector<double> >::iterator it =
324
		find((*interpPointCoords)[dofi].begin(), 
325
		     (*interpPointCoords)[dofi].end(),
326
		     vertexCoords);
327
	
328
	if (it == (*interpPointCoords)[dofi].end()) {
329
	  (*interpPointCoords)[dofi].push_back(vertexCoords); 
330
	  nInterpPoints++;
331
	}
332
      }      
333
    }  
334 335
  }

336

337
  void DataCollector::addInterpData(ElInfo *elInfo)
338
  {
339
    FUNCNAME("DataCollector::addInterpData()");
340
    
341
    basisFcts->getLocalIndices(elInfo->getElement(), localAdmin, localDOFs);
342

343
    std::vector<int> elemInterpPoints(0);
344 345
    for (int i = mesh->getGeo(VERTEX); i < nBasisFcts; i++)
      elemInterpPoints.push_back((*interpPointInd)[localDOFs[i]]);
Thomas Witkowski's avatar
Thomas Witkowski committed
346
    
347
    interpPoints.push_back(elemInterpPoints); 
348 349
  }

350

351 352
  void DataCollector::addPeriodicData(ElInfo *elInfo) 
  {
353 354 355 356 357 358 359
    FUNCNAME("DataCollector::addPeriodicData");

    LeafDataPeriodic *ldp = dynamic_cast<LeafDataPeriodic*>
      (elInfo->getElement()->
       getElementData()->
       getElementData(PERIODIC));
    
360
    if (ldp) {
Thomas Witkowski's avatar
Thomas Witkowski committed
361
      std::list<LeafDataPeriodic::PeriodicInfo>::iterator it;
362
      
363 364 365 366
      for (it = dynamic_cast<LeafDataPeriodic*>(ldp)->getInfoList().begin();
	   it != dynamic_cast<LeafDataPeriodic*>(ldp)->getInfoList().end();
	   ++it) {

367 368
	int outputIndex = outputIndices[elInfo->getElement()->getIndex()];
	int neighIndex  = outputIndices[elInfo->
369 370 371
					 getNeighbour(it->elementSide)->
					 getIndex()];
	
372
	if (!periodicConnections[outputIndex][it->elementSide]) {
373
	  PeriodicInfo periodicInfo;
374
	  
375 376 377 378 379 380
	  periodicInfo.mode = it->periodicMode;
	  periodicInfo.type = it->type;
	  periodicInfo.outputIndex = outputIndex;
	  periodicInfo.neighIndex = neighIndex;
	  periodicInfo.vertexMap.clear();
	  
381 382
	  periodicConnections[outputIndex][it->elementSide] = true;
	  periodicConnections
383 384 385
	    [neighIndex][elInfo->getOppVertex(it->elementSide)] = true;
	    
  
386 387
	  for (int i = 0; i < dim; i++) {
	    int index1 = elInfo->getElement()->getVertexOfPosition(INDEX_OF_DIM(dim - 1, dim),
388 389
								   it->elementSide,
								   i);
390
	    int dof1 = elInfo->getElement()->getDof(index1, nPreDofs);
391
	    
392 393
	    for (int j = 0; j < dim; j++) {
	      int index2 = elInfo->getElement()->getVertexOfPosition(INDEX_OF_DIM(dim - 1, dim),
394 395
								     elInfo->getOppVertex(it->elementSide),
								     j);
396
	      int dof2 = elInfo->getNeighbour(it->elementSide)->getDof(index2, nPreDofs);
397
	      
398
	      if ((dof1 == dof2) || (mesh->associated(dof1, dof2))) {
399 400 401
		periodicInfo.vertexMap[index1] = index2;
		break;
	      } 
402 403
	    }
	  }
404
	  
405
	  periodicInfos.push_back(periodicInfo);
406
	}
407
      }
408 409 410
    }
  }

411

Thomas Witkowski's avatar
Thomas Witkowski committed
412
  std::list<ElementInfo>* DataCollector::getElementInfos()
413
  {        
414
    if (!elementDataCollected)
415 416
      startCollectingElementData();
      
417
    return &elements;
418 419
  }

420

Thomas Witkowski's avatar
Thomas Witkowski committed
421
  DOFVector< std::list<VertexInfo> >* DataCollector::getVertexInfos()
422
  {
423
    if (!elementDataCollected)
424 425
      startCollectingElementData();

426
    return vertexInfos;
427 428
  }

429

Thomas Witkowski's avatar
Thomas Witkowski committed
430
  std::list<PeriodicInfo>* DataCollector::getPeriodicInfos()
431
  {
432
    if (!periodicDataCollected)
433 434
      startCollectingPeriodicData();

435
    return &periodicInfos;
436 437
  }

438

439 440
  int DataCollector::getNumberVertices()
  {
441
    if (!elementDataCollected)
442 443
      startCollectingElementData();

444
    return nVertices;
445 446
  }

447

448 449
  int DataCollector::getNumberElements()
  {
450
    if (!elementDataCollected)
451 452
      startCollectingElementData();

453
    return nElements;
454 455
  }

456

457 458
  int DataCollector::getNumberInterpPoints()
  {
459
    if (!valueDataCollected)
460 461
      startCollectingValueData();

462
    return nInterpPoints;
463
  }
464

465 466 467
  
  int DataCollector::getNumberConnections()
  {
468
    if (!periodicDataCollected)
469 470
      startCollectingPeriodicData();

471
    return nConnection;
472 473 474 475 476
  }


  DOFVector<double>* DataCollector::getValues()
  {
477
    if (!valueDataCollected)
478 479
      startCollectingValueData();

480
    return values;
481 482
  }

483

Thomas Witkowski's avatar
Thomas Witkowski committed
484
  DOFVector< std::list<WorldVector<double> > >* DataCollector::getDofCoords()
485
  {
486
    if (!valueDataCollected)
487 488
      startCollectingValueData();

489
    return dofCoord;
490 491
  }

492

493 494
  DOFVector<int>* DataCollector::getInterpPointInd()
  {
495
    if (!valueDataCollected)
496 497
      startCollectingValueData();

498
    return interpPointInd;
499 500
  }

501

Thomas Witkowski's avatar
Thomas Witkowski committed
502
  DOFVector< std::list<WorldVector<double> > >* DataCollector::getInterpPointCoords()
503
  {
504
    if (!valueDataCollected)
505 506
      startCollectingValueData();

507
    return interpPointCoords;
508 509
  }

510

Thomas Witkowski's avatar
Thomas Witkowski committed
511
  std::vector< std::vector<int> >* DataCollector::getInterpPoints()
512
  {
513
    if (!valueDataCollected)
514 515
      startCollectingValueData();

516
    return &interpPoints;
517 518
  }
}