From 1f9daa76bc54b81ef8b05ac57b94134bf2bf6999 Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Wed, 12 Jan 2011 09:18:19 +0000
Subject: [PATCH] Added some parallel debuging test for periodic boundary
conditions.
---
AMDiS/src/AdaptInfo.h | 4 +-
AMDiS/src/Boundary.h | 2 +-
AMDiS/src/CoarseningManager.cc | 4 +-
AMDiS/src/CoarseningManager1d.cc | 2 +-
AMDiS/src/DOFAdmin.cc | 6 +-
AMDiS/src/DOFIterator.h | 2 +-
AMDiS/src/Error.hh | 2 +-
AMDiS/src/Global.h | 11 -
AMDiS/src/MathFunctions.h | 14 +-
AMDiS/src/Mesh.cc | 45 ++--
AMDiS/src/Mesh.h | 28 +++
AMDiS/src/Operator.cc | 2 +-
AMDiS/src/Parameters.cc | 2 +-
AMDiS/src/RefinementManager1d.cc | 2 +-
AMDiS/src/RefinementManager2d.cc | 4 +-
AMDiS/src/RefinementManager3d.cc | 7 +-
AMDiS/src/ResidualEstimator.cc | 4 +-
AMDiS/src/Tetrahedron.h | 2 +-
AMDiS/src/Triangle.h | 2 +-
AMDiS/src/VertexVector.cc | 8 +-
AMDiS/src/io/MacroInfo.cc | 8 +-
AMDiS/src/io/MacroReader.cc | 12 +-
AMDiS/src/io/PovrayWriter.cc | 24 +-
AMDiS/src/parallel/ElementObjectData.cc | 130 +++++++----
AMDiS/src/parallel/ElementObjectData.h | 22 +-
AMDiS/src/parallel/InteriorBoundary.cc | 7 +-
AMDiS/src/parallel/MeshDistributor.cc | 267 ++++++++++++----------
AMDiS/src/parallel/ParMetisPartitioner.cc | 8 +-
AMDiS/src/parallel/ParallelDebug.cc | 96 +++++++-
AMDiS/src/parallel/ParallelDebug.h | 5 +
AMDiS/src/parallel/PetscSolver.cc | 11 +-
31 files changed, 479 insertions(+), 264 deletions(-)
diff --git a/AMDiS/src/AdaptInfo.h b/AMDiS/src/AdaptInfo.h
index 52b33ace..ac7b2f1c 100644
--- a/AMDiS/src/AdaptInfo.h
+++ b/AMDiS/src/AdaptInfo.h
@@ -367,7 +367,7 @@ namespace AMDiS {
/// Sets \ref timestepNumber.
inline void setTimestepNumber(int num)
{
- timestepNumber = AMDiS::min(nTimesteps,num);
+ timestepNumber = std::min(nTimesteps, num);
}
/// Returns \ref nTimesteps.
@@ -379,7 +379,7 @@ namespace AMDiS {
/// Sets \ref nTimesteps.
inline void setNumberOfTimesteps(int num)
{
- nTimesteps = AMDiS::max(0,num);
+ nTimesteps = std::max(0, num);
}
/// Increments \ref timestepNumber by 1;
diff --git a/AMDiS/src/Boundary.h b/AMDiS/src/Boundary.h
index 4cd2f773..40c7634b 100644
--- a/AMDiS/src/Boundary.h
+++ b/AMDiS/src/Boundary.h
@@ -43,4 +43,4 @@ namespace AMDiS {
}
-#endif // _Boundary_H_
+#endif
diff --git a/AMDiS/src/CoarseningManager.cc b/AMDiS/src/CoarseningManager.cc
index 14b26a68..6ba5f78d 100644
--- a/AMDiS/src/CoarseningManager.cc
+++ b/AMDiS/src/CoarseningManager.cc
@@ -45,7 +45,7 @@ namespace AMDiS {
if (el->getChild(0)) {
// interior node of the tree
- int mark = max(el->getChild(0)->getMark(), el->getChild(1)->getMark());
+ int mark = std::max(el->getChild(0)->getMark(), el->getChild(1)->getMark());
el->setMark(std::min(mark + 1, 0));
} else {
// leaf node of the tree
@@ -63,7 +63,7 @@ namespace AMDiS {
ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL);
while (elInfo) {
Element *el = elInfo->getElement();
- el->setMark(max(el->getMark(), 0));
+ el->setMark(std::max(el->getMark(), 0));
elInfo = stack.traverseNext(elInfo);
}
}
diff --git a/AMDiS/src/CoarseningManager1d.cc b/AMDiS/src/CoarseningManager1d.cc
index 37a62db0..8d4972d1 100644
--- a/AMDiS/src/CoarseningManager1d.cc
+++ b/AMDiS/src/CoarseningManager1d.cc
@@ -57,7 +57,7 @@ namespace AMDiS {
int mark0 = coarsenRecursive(dynamic_cast<Line*>(const_cast<Element*>(parent->getChild(0))));
int mark1 = coarsenRecursive(dynamic_cast<Line*>(const_cast<Element*>( parent->getChild(1))));
- mark = max(mark0, mark1);
+ mark = std::max(mark0, mark1);
child[0] = dynamic_cast<Line*>(const_cast<Element*>( parent->getChild(0)));
child[1] = dynamic_cast<Line*>(const_cast<Element*>( parent->getChild(1)));
diff --git a/AMDiS/src/DOFAdmin.cc b/AMDiS/src/DOFAdmin.cc
index dfc07484..67c604db 100644
--- a/AMDiS/src/DOFAdmin.cc
+++ b/AMDiS/src/DOFAdmin.cc
@@ -167,7 +167,7 @@ namespace AMDiS {
usedCount++;
if (holeCount > 0)
holeCount--;
- sizeUsed = max(sizeUsed, dof + 1);
+ sizeUsed = std::max(sizeUsed, dof + 1);
return dof;
}
@@ -182,7 +182,7 @@ namespace AMDiS {
if (old > minsize)
return;
- int newval = max(minsize, static_cast<int>((dofFree.size() + sizeIncrement)));
+ int newval = std::max(minsize, static_cast<int>((dofFree.size() + sizeIncrement)));
size = newval;
@@ -324,7 +324,7 @@ namespace AMDiS {
for (std::list<DOFContainer*>::iterator dc = dofContainerList.begin();
dc != dofContainerList.end(); ++dc)
- (*dc)->compressDofContainer(n, newDofIndex);
+ (*dc)->compressDofContainer(n, newDofIndex);
}
diff --git a/AMDiS/src/DOFIterator.h b/AMDiS/src/DOFIterator.h
index ea97dcef..815fe658 100644
--- a/AMDiS/src/DOFIterator.h
+++ b/AMDiS/src/DOFIterator.h
@@ -173,7 +173,7 @@ namespace AMDiS {
}
/// Dereferntiation of the \ref dofFreeIterator
- virtual bool isDOFFree()
+ virtual bool isDofFree()
{
return *dofFreeIterator;
}
diff --git a/AMDiS/src/Error.hh b/AMDiS/src/Error.hh
index 7fcbaecd..ae26ec46 100644
--- a/AMDiS/src/Error.hh
+++ b/AMDiS/src/Error.hh
@@ -77,7 +77,7 @@ namespace AMDiS {
for (int i = 0; i < nPoints; i++) {
err = u_vec[i] > uh_vec[i] ? u_vec[i] - uh_vec[i] : uh_vec[i] - u_vec[i];
- maxErr = max(maxErr, err);
+ maxErr = std::max(maxErr, err);
}
elInfo = stack.traverseNext(elInfo);
diff --git a/AMDiS/src/Global.h b/AMDiS/src/Global.h
index ebefbcb1..64e17af7 100644
--- a/AMDiS/src/Global.h
+++ b/AMDiS/src/Global.h
@@ -117,17 +117,6 @@ namespace AMDiS {
// ===== some simple template functions ======================================
- /// template<typename T> T max(T a,T b ) {return ((a) > (b) ? (a) : (b));}
- template<typename T,typename S> inline T max(T a,S b )
- {
- return ((a) > (b) ? (a) : (b));
- }
-
- template<typename T> inline T min(T a,T b)
- {
- return ((a) < (b)) ? (a) : (b);
- }
-
template<typename T> inline T abs(T a)
{
return ((a) >= 0 ? (a) : -(a));
diff --git a/AMDiS/src/MathFunctions.h b/AMDiS/src/MathFunctions.h
index 26063236..aa354f5e 100644
--- a/AMDiS/src/MathFunctions.h
+++ b/AMDiS/src/MathFunctions.h
@@ -42,15 +42,17 @@ namespace AMDiS {
return 0;
}
- inline double Phi1ToR(double p1, double eps) {
- double x = max(-1.0 + numeric_limits< double >::epsilon(),
- min(1.0 - numeric_limits< double >::epsilon(), p1));
+ inline double Phi1ToR(double p1, double eps)
+ {
+ double x = std::max(-1.0 + numeric_limits<double>::epsilon(),
+ std::min(1.0 - numeric_limits<double>::epsilon(), p1));
return eps / 3.0 * log((1 + x) / (1 - x)) * 0.5;
}
- inline double Phi2ToR(double p2, double eps) {
- double x = max(-1.0 + numeric_limits< double >::epsilon(),
- min(1.0 - numeric_limits< double >::epsilon(), 1 + 2 * p2));
+ inline double Phi2ToR(double p2, double eps)
+ {
+ double x = std::max(-1.0 + numeric_limits<double>::epsilon(),
+ std::min(1.0 - numeric_limits<double>::epsilon(), 1 + 2 * p2));
return eps / 3.0 * log( (1 + x) / (1 - x) );
}
diff --git a/AMDiS/src/Mesh.cc b/AMDiS/src/Mesh.cc
index cd2a656a..f140ca49 100644
--- a/AMDiS/src/Mesh.cc
+++ b/AMDiS/src/Mesh.cc
@@ -102,9 +102,12 @@ namespace AMDiS {
changeIndex(0),
final_lambda(dimension, DEFAULT_VALUE, 0.0)
{
-
FUNCNAME("Mesh::Mesh()");
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+ nParallelPreRefinements = 0;
+#endif
+
// set default element prototype
switch(dim) {
case 1:
@@ -1248,42 +1251,41 @@ namespace AMDiS {
// nMacroElements * 2^gr >= nProcs * 10
// => gr = log_2(nProcs * 10 / nMacroElements)
- double tmp = 10.0 * MPI::COMM_WORLD.Get_size() / nMacroElements;
- int nrRefine = ceil(log(tmp) / log(2));
- if (nrRefine < 0)
- nrRefine = 0;
+ double scale = 10.0 * MPI::COMM_WORLD.Get_size() / nMacroElements;
+ nParallelPreRefinements = static_cast<int>(std::max(0.0, ceil(log(scale) / log(2))));
if (dim == 3) {
- int newElType = (elType + nrRefine) % 3;
+ int newElType = (elType + nParallelPreRefinements) % 3;
switch (newElType) {
case 1:
- if (nrRefine > 0)
- nrRefine--;
+ if (nParallelPreRefinements > 0)
+ nParallelPreRefinements--;
else
- nrRefine = 2;
+ nParallelPreRefinements = 2;
break;
case 2:
- nrRefine++;
+ nParallelPreRefinements++;
break;
}
- TEST_EXIT((elType + nrRefine) % 3 == 0)("This should not happen!\n");
+ TEST_EXIT((elType + nParallelPreRefinements) % 3 == 0)
+ ("This should not happen!\n");
}
// === Check if number of pre refinements is set ini init file. ===
- int nPreRefine = 0;
- GET_PARAMETER(0, "parallel->pre refine", "%d", &nPreRefine);
- if (nPreRefine != 0) {
+ int tmp = 0;
+ GET_PARAMETER(0, "parallel->pre refine", "%d", &tmp);
+ if (tmp != 0) {
MSG("Calculated %d pre refines to be useful, but %d are set in init file!\n",
- nrRefine, nPreRefine);
- nrRefine = nPreRefine;
+ nParallelPreRefinements, tmp);
+ nParallelPreRefinements = tmp;
}
// === If we do not need to refine the mesh, return back. ===
- if (nrRefine == 0)
+ if (nParallelPreRefinements == 0)
return;
@@ -1323,7 +1325,7 @@ namespace AMDiS {
else
refManager = new RefinementManager3d();
- refManager->globalRefine(&testMesh, nrRefine);
+ refManager->globalRefine(&testMesh, nParallelPreRefinements);
delete refManager;
Lagrange* basFcts = Lagrange::getLagrange(dim, 1);
@@ -1348,21 +1350,20 @@ namespace AMDiS {
int globalRefinements = 0;
GET_PARAMETER(0, name + "->global refinements", "%d", &globalRefinements);
- if (globalRefinements < nrRefine)
+ if (globalRefinements < nParallelPreRefinements)
globalRefinements = 0;
else
- globalRefinements -= nrRefine;
+ globalRefinements -= nParallelPreRefinements;
std::stringstream oss;
oss << globalRefinements;
-
ADD_PARAMETER(0, name + "->global refinements", oss.str().c_str());
// === Print a note to the screen that another mesh file will be used. ===
MSG("The macro mesh file \"%s\" was refined %d times and stored to file \"%s\".\n",
- macroFilename.c_str(), nrRefine, newMacroFilename.str().c_str());
+ macroFilename.c_str(), nParallelPreRefinements, newMacroFilename.str().c_str());
macroFilename = newMacroFilename.str();
if (periodicFilename != "")
diff --git a/AMDiS/src/Mesh.h b/AMDiS/src/Mesh.h
index 93c9b8e7..0b0a0bea 100644
--- a/AMDiS/src/Mesh.h
+++ b/AMDiS/src/Mesh.h
@@ -608,6 +608,21 @@ namespace AMDiS {
///
void deleteMeshStructure();
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+ /// In parallel computations the level of all macro elements is equal to the number
+ /// of global pre refinements, \ref nParallelPreRefinements.
+ inline int getMacroElementLevel()
+ {
+ return nParallelPreRefinements;
+ }
+#else
+ /// In sequentiel computations the level of all macro elements is always 0.
+ inline int getMacroElementLevel()
+ {
+ return 0;
+ }
+#endif
+
public:
///
static const Flag FILL_NOTHING;
@@ -702,6 +717,12 @@ namespace AMDiS {
void checkParallelMacroFile(std::string ¯oFilename,
std::string &periodicFilename,
int check);
+
+ /// Returns \ref nParallelPreRefinements
+ int getParallelPreRefinements() const
+ {
+ return nParallelPreRefinements;
+ }
#endif
protected:
@@ -855,6 +876,13 @@ namespace AMDiS {
*/
long changeIndex;
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+ /// In parallel computations the mesh may be globally prerefined to achieve a fine
+ /// enought starting mesh for the given number of ranks. The value of the variable
+ /// will be defined in function \ref checkParallelMacroFile.
+ int nParallelPreRefinements;
+#endif
+
protected:
/// for findElement-Fcts
DimVec<double> final_lambda;
diff --git a/AMDiS/src/Operator.cc b/AMDiS/src/Operator.cc
index 6ef9f206..57b8d76d 100644
--- a/AMDiS/src/Operator.cc
+++ b/AMDiS/src/Operator.cc
@@ -53,7 +53,7 @@ namespace AMDiS {
for (int i = 0; i < static_cast<int>(terms->size()); i++) {
OperatorTerm *term = (*terms)[i];
- maxTermDegree = max(maxTermDegree, term->degree);
+ maxTermDegree = std::max(maxTermDegree, term->degree);
}
return psiDegree + phiDegree - order + maxTermDegree;
diff --git a/AMDiS/src/Parameters.cc b/AMDiS/src/Parameters.cc
index c04502a7..9559167c 100644
--- a/AMDiS/src/Parameters.cc
+++ b/AMDiS/src/Parameters.cc
@@ -50,7 +50,7 @@ namespace AMDiS {
if (Parameters::singlett->paramInfo) {
if (Parameters::singlett->paramInfo > 1)
- info = max(info, Parameters::singlett->paramInfo - 1);
+ info = std::max(info, Parameters::singlett->paramInfo - 1);
} else {
info = 0;
}
diff --git a/AMDiS/src/RefinementManager1d.cc b/AMDiS/src/RefinementManager1d.cc
index d26cf503..6bc99dad 100644
--- a/AMDiS/src/RefinementManager1d.cc
+++ b/AMDiS/src/RefinementManager1d.cc
@@ -41,7 +41,7 @@ namespace AMDiS {
child[0] = dynamic_cast<Line*>(mesh->createNewElement(el));
child[1] = dynamic_cast<Line*>(mesh->createNewElement(el));
- int mark = max(0, el->getMark() - 1);
+ int mark = std::max(0, el->getMark() - 1);
child[0]->setMark(mark);
child[1]->setMark(mark);
el->setMark(0);
diff --git a/AMDiS/src/RefinementManager2d.cc b/AMDiS/src/RefinementManager2d.cc
index d730ca9f..601a67cf 100644
--- a/AMDiS/src/RefinementManager2d.cc
+++ b/AMDiS/src/RefinementManager2d.cc
@@ -257,7 +257,7 @@ namespace AMDiS {
child[0] = dynamic_cast<Triangle*>(mesh->createNewElement(el));
child[1] = dynamic_cast<Triangle*>(mesh->createNewElement(el));
- int newMark = max(0, el->getMark() - 1);
+ int newMark = std::max(0, el->getMark() - 1);
child[0]->setMark(newMark);
child[1]->setMark(newMark);
el->setMark(0);
@@ -334,7 +334,7 @@ namespace AMDiS {
int opp_vertex = (*elInfo)->getOppVertex(2);
ElInfo *neigh_info = stack->traverseNeighbour2d(*elInfo, 2);
- neigh_info->getElement()->setMark(max(neigh_info->getElement()->getMark(), 1));
+ neigh_info->getElement()->setMark(std::max(neigh_info->getElement()->getMark(), 1));
neigh_info = refineFunction(neigh_info);
// === Now go back to the original element and refine the patch. ===
diff --git a/AMDiS/src/RefinementManager3d.cc b/AMDiS/src/RefinementManager3d.cc
index b5f2ad0e..e3cde25c 100644
--- a/AMDiS/src/RefinementManager3d.cc
+++ b/AMDiS/src/RefinementManager3d.cc
@@ -25,6 +25,7 @@
#include "DOFVector.h"
#include "PeriodicBC.h"
#include "VertexVector.h"
+#include "Debug.h"
namespace AMDiS {
@@ -41,7 +42,7 @@ namespace AMDiS {
child[0] = dynamic_cast<Tetrahedron*>(mesh->createNewElement(el));
child[1] = dynamic_cast<Tetrahedron*>(mesh->createNewElement(el));
- int mark = max(0, el->getMark() - 1);
+ int mark = std::max(0, el->getMark() - 1);
child[0]->setMark(mark);
child[1]->setMark(mark);
el->setMark(0);
@@ -249,6 +250,8 @@ namespace AMDiS {
int RefinementManager3d::newCoordsFct(ElInfo *el_info)
{
+ FUNCNAME("RefinementManager3d::newCoordsFct()");
+
Element *el = el_info->getElement();
DegreeOfFreedom *edge[2];
ElInfo *elinfo = el_info;
@@ -474,7 +477,7 @@ namespace AMDiS {
// Only 0 can be a compatible commen refinement edge. Thus, neigh has not
// a compatible refinement edge. Refine it first.
- neigh->setMark(max(neigh->getMark(), 1));
+ neigh->setMark(std::max(neigh->getMark(), 1));
neighInfo = refineFunction(neighInfo);
diff --git a/AMDiS/src/ResidualEstimator.cc b/AMDiS/src/ResidualEstimator.cc
index 52f05c6b..2b9e8fbc 100644
--- a/AMDiS/src/ResidualEstimator.cc
+++ b/AMDiS/src/ResidualEstimator.cc
@@ -256,7 +256,7 @@ namespace AMDiS {
el->setEstimation(est_el, row);
el->setMark(0);
est_sum += est_el;
- est_max = max(est_max, est_el);
+ est_max = std::max(est_max, est_el);
}
@@ -295,7 +295,7 @@ namespace AMDiS {
}
double v = C3 * det * result;
est_t_sum += v;
- est_t_max = max(est_t_max, v);
+ est_t_max = std::max(est_t_max, v);
}
}
}
diff --git a/AMDiS/src/Tetrahedron.h b/AMDiS/src/Tetrahedron.h
index 1a3f274f..0f77b8bd 100644
--- a/AMDiS/src/Tetrahedron.h
+++ b/AMDiS/src/Tetrahedron.h
@@ -274,7 +274,7 @@ namespace AMDiS {
DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
- DofEdge edge = std::make_pair(min(dof0, dof1), max(dof0, dof1));
+ DofEdge edge = std::make_pair(std::min(dof0, dof1), std::max(dof0, dof1));
return edge;
}
diff --git a/AMDiS/src/Triangle.h b/AMDiS/src/Triangle.h
index ef141a01..f6a7f045 100644
--- a/AMDiS/src/Triangle.h
+++ b/AMDiS/src/Triangle.h
@@ -178,7 +178,7 @@ namespace AMDiS {
DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
- DofEdge edge = std::make_pair(min(dof0, dof1), max(dof0, dof1));
+ DofEdge edge = std::make_pair(std::min(dof0, dof1), std::max(dof0, dof1));
return edge;
}
diff --git a/AMDiS/src/VertexVector.cc b/AMDiS/src/VertexVector.cc
index 3544a529..886baabd 100644
--- a/AMDiS/src/VertexVector.cc
+++ b/AMDiS/src/VertexVector.cc
@@ -26,29 +26,33 @@ namespace AMDiS {
const_cast<DOFAdmin*>(admin)->addDOFContainer(this);
}
+
VertexVector::~VertexVector()
{
const_cast<DOFAdmin*>(admin)->removeDOFIndexed(this);
const_cast<DOFAdmin*>(admin)->removeDOFContainer(this);
}
+
void VertexVector::set(DegreeOfFreedom val)
{
DOFIteratorBase it(const_cast<DOFAdmin*>(admin), USED_DOFS);
for (it.reset(); !it.end(); ++it)
- if (!it.isDOFFree())
+ if (!it.isDofFree())
operator[](it.getDOFIndex()) = val;
}
+
void VertexVector::print()
{
DOFIteratorBase it(const_cast<DOFAdmin*>(admin), USED_DOFS);
for (it.reset(); !it.end(); ++it)
- if (!it.isDOFFree())
+ if (!it.isDofFree())
std::cout << "DOF " << it.getDOFIndex() << " -> "
<< operator[](it.getDOFIndex()) << "\n";
}
+
void VertexVector::changeDofIndices(std::map<int, int>& dofIndexMap)
{
FUNCNAME("VertexVector::changedofIndices()");
diff --git a/AMDiS/src/io/MacroInfo.cc b/AMDiS/src/io/MacroInfo.cc
index c0ae8aea..a1a853ae 100644
--- a/AMDiS/src/io/MacroInfo.cc
+++ b/AMDiS/src/io/MacroInfo.cc
@@ -525,19 +525,19 @@ namespace AMDiS {
int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j + 2) % 3];
- bound[j1] = max(bound[j1], (*melIt)->getBoundary(j));
- bound[j2] = max(bound[j2], (*melIt)->getBoundary(j));
+ bound[j1] = std::max(bound[j1], (*melIt)->getBoundary(j));
+ bound[j2] = std::max(bound[j2], (*melIt)->getBoundary(j));
} else if ((*melIt)->getBoundary(j) <= NEUMANN) {
int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j + 2) % 3];
if (bound[j1] != INTERIOR)
- bound[j1] = max(bound[j1], (*melIt)->getBoundary(j));
+ bound[j1] = std::max(bound[j1], (*melIt)->getBoundary(j));
else
bound[j1] = (*melIt)->getBoundary(j);
if (bound[j2] != INTERIOR)
- bound[j2] = max(bound[j2], (*melIt)->getBoundary(j));
+ bound[j2] = std::max(bound[j2], (*melIt)->getBoundary(j));
else
bound[j2] = (*melIt)->getBoundary(j);
}
diff --git a/AMDiS/src/io/MacroReader.cc b/AMDiS/src/io/MacroReader.cc
index 31594995..484c6467 100644
--- a/AMDiS/src/io/MacroReader.cc
+++ b/AMDiS/src/io/MacroReader.cc
@@ -495,23 +495,19 @@ namespace AMDiS {
periodic[it->elementSide] = true;
}
- for (k = 0; k < mesh->getGeo(EDGE); k++) {
- /*********************************************************************/
- /* check for not counted edges */
- /*********************************************************************/
+ for (k = 0; k < mesh->getGeo(EDGE); k++) {
+ // === Check for not counted edges. ===
n_neigh = 1;
if (newEdge(mesh, (*(mel + i)), k, &n_neigh)) {
mesh->incrementNumberOfEdges(1);
- max_n_neigh = max(max_n_neigh, n_neigh);
+ max_n_neigh = std::max(max_n_neigh, n_neigh);
}
}
for (k = 0; k < mesh->getGeo(NEIGH); k++) {
neigh = (*(mel + i))->getNeighbour(k);
- /*********************************************************************/
- /* face is counted and dof is added by the element with bigger index */
- /*********************************************************************/
+ // === Face is counted and dof is added by the element with bigger index. ===
if (neigh && (neigh->getIndex() > (*(mel + i))->getIndex()))
continue;
diff --git a/AMDiS/src/io/PovrayWriter.cc b/AMDiS/src/io/PovrayWriter.cc
index 32b94e8d..1c1d52c5 100644
--- a/AMDiS/src/io/PovrayWriter.cc
+++ b/AMDiS/src/io/PovrayWriter.cc
@@ -59,8 +59,8 @@ namespace AMDiS {
min_value = max_value = (*values)[0];
for (int index = 1; index < values->getSize(); index++) {
- min_value = min(min_value, (*values)[index]);
- max_value = max(max_value, (*values)[index]);
+ min_value = std::min(min_value, (*values)[index]);
+ max_value = std::max(max_value, (*values)[index]);
}
/* map DOFs to values */
@@ -126,13 +126,13 @@ namespace AMDiS {
for (it.reset(); !it.end(); ++it) {
// initialize mit and max values with coordinates of first vertex
std::list<VertexInfo>::iterator it2 = it->begin();
- bBox->minx = min(bBox->minx, it2->coords[0]);
- bBox->maxx = max(bBox->maxx, it2->coords[0]);
- bBox->miny = min(bBox->miny, it2->coords[1]);
- bBox->maxy = max(bBox->maxy, it2->coords[1]);
+ bBox->minx = std::min(bBox->minx, it2->coords[0]);
+ bBox->maxx = std::max(bBox->maxx, it2->coords[0]);
+ bBox->miny = std::min(bBox->miny, it2->coords[1]);
+ bBox->maxy = std::max(bBox->maxy, it2->coords[1]);
if (dow == 3) {
- bBox->minz = min(bBox->minz, it2->coords[2]);
- bBox->maxz = max(bBox->maxz, it2->coords[2]);
+ bBox->minz = std::min(bBox->minz, it2->coords[2]);
+ bBox->maxz = std::max(bBox->maxz, it2->coords[2]);
}
}
@@ -149,8 +149,8 @@ namespace AMDiS {
// determines the color for a given value/weight
string getColorString(double &value)
{
- value = max(0.0, value);
- value = min(value, 1.0);
+ value = std::max(0.0, value);
+ value = std::min(value, 1.0);
// rot: 1,0,0
// blau: 0,0,1
@@ -324,8 +324,8 @@ namespace AMDiS {
for (std::list<VertexInfo>::iterator it2 = it->begin(); it2 != it->end(); ++it2) {
// test: use y-coordinate to compute color
double redValue = it2->coords[1];
- redValue = max(0., redValue);
- redValue = min(redValue, 1.);
+ redValue = std::max(0., redValue);
+ redValue = std::min(redValue, 1.);
out << "\n\t\ttexture{ pigment{ rgb"<< getColorString(redValue) <<" } }" << ",";
}
diff --git a/AMDiS/src/parallel/ElementObjectData.cc b/AMDiS/src/parallel/ElementObjectData.cc
index b005e1fc..a97b17e5 100644
--- a/AMDiS/src/parallel/ElementObjectData.cc
+++ b/AMDiS/src/parallel/ElementObjectData.cc
@@ -44,6 +44,8 @@ namespace AMDiS {
if (mesh->isPeriodicAssociation(elInfo->getBoundary(EDGE, i))) {
// The current element's i-th edge is periodic.
Element *neigh = elInfo->getNeighbour(i);
+
+ TEST_EXIT_DBG(neigh)("Should not happen!\n");
DofEdge edge0 = el->getEdge(i);
DofEdge edge1 = neigh->getEdge(elInfo->getOppVertex(i));
@@ -72,6 +74,8 @@ namespace AMDiS {
if (mesh->isPeriodicAssociation(elInfo->getBoundary(FACE, i))) {
// The current element's i-th face is periodic.
Element *neigh = elInfo->getNeighbour(i);
+
+ TEST_EXIT_DBG(neigh)("Should not happen!\n");
DofFace face0 = el->getFace(i);
DofFace face1 = neigh->getFace(elInfo->getOppVertex(i));
@@ -132,21 +136,8 @@ namespace AMDiS {
if (periodicVertices.size() == 0)
return;
-
-
- // === Search for an unsed boundary index. ===
+
- BoundaryType newPeriodicBoundaryType = 0;
- for (map<BoundaryType, VertexVector*>::iterator it = mesh->getPeriodicAssociations().begin();
- it != mesh->getPeriodicAssociations().end(); ++it)
- newPeriodicBoundaryType = min(newPeriodicBoundaryType, it->first);
-
- TEST_EXIT_DBG(newPeriodicBoundaryType < 0)("Should not happen!\n");
- newPeriodicBoundaryType--;
- mesh->getPeriodicAssociations()[newPeriodicBoundaryType] =
- new VertexVector(feSpace->getAdmin(), "");
-
-
// === Get all vertex DOFs that have multiple periodic associations. ===
// We group all vertices together, that have either two or three periodic
@@ -223,10 +214,12 @@ namespace AMDiS {
("Should not happen!\n");
TEST_EXIT_DBG(periodicVertices.count(make_pair(dof3, dof0)) == 0)
("Should not happen!\n");
-
- periodicVertices[make_pair(dof0, dof3)] = newPeriodicBoundaryType;
- periodicVertices[make_pair(dof3, dof0)] = newPeriodicBoundaryType;
-
+
+ periodicVertices[make_pair(dof0, dof3)] =
+ provideConnectedPeriodicBoundary(type0, type1);
+ periodicVertices[make_pair(dof3, dof0)] =
+ provideConnectedPeriodicBoundary(type0, type1);
+
for (unsigned int j = i + 1; j < multPeriodicDof2.size(); j++)
if (multPeriodicDof2[j] == dof3)
multPeriodicDof2[j] = -1;
@@ -235,6 +228,7 @@ namespace AMDiS {
if (multPeriodicDof3.size() > 0) {
int nMultPeriodicDofs = multPeriodicDof3.size();
+
for (int i = 0; i < nMultPeriodicDofs; i++) {
for (int j = i + 1; j < nMultPeriodicDofs; j++) {
pair<DegreeOfFreedom, DegreeOfFreedom> perDofs0 =
@@ -243,20 +237,15 @@ namespace AMDiS {
make_pair(multPeriodicDof3[j], multPeriodicDof3[i]);
if (periodicVertices.count(perDofs0) == 0) {
- TEST_EXIT_DBG(periodicVertices.count(perDofs1) == 0)
- ("Should not happen!\n");
-
- periodicVertices[perDofs0] = newPeriodicBoundaryType;
- periodicVertices[perDofs1] = newPeriodicBoundaryType;
- newPeriodicBoundaryType--;
- mesh->getPeriodicAssociations()[newPeriodicBoundaryType] =
- new VertexVector(feSpace->getAdmin(), "");
+ BoundaryType b = getNewBoundaryType();
+ periodicVertices[perDofs0] = b;
+ periodicVertices[perDofs1] = b;
}
}
}
}
-
-
+
+
// === Get all edges that have multiple periodic associations (3D only!). ===
for (map<DofEdge, std::set<DofEdge> >::iterator it = periodicEdgeAssoc.begin();
@@ -265,16 +254,24 @@ namespace AMDiS {
TEST_EXIT_DBG(mesh->getDim() == 3)("Should not happen!\n");
TEST_EXIT_DBG(it->second.size() == 2)("Should not happen!\n");
- pair<DofEdge, DofEdge> perEdge0 =
- make_pair(*(it->second.begin()), *(++(it->second.begin())));
- pair<DofEdge, DofEdge> perEdge1 =
- make_pair(perEdge0.second, perEdge0.first);
+ DofEdge edge0 = it->first;
+ DofEdge edge1 = *(it->second.begin());
+ DofEdge edge2 = *(++(it->second.begin()));
- periodicEdges[perEdge0] = newPeriodicBoundaryType;
- periodicEdges[perEdge1] = newPeriodicBoundaryType;
- newPeriodicBoundaryType--;
- mesh->getPeriodicAssociations()[newPeriodicBoundaryType] =
- new VertexVector(feSpace->getAdmin(), "");
+ pair<DofEdge, DofEdge> perEdge0 = make_pair(edge1, edge2);
+ pair<DofEdge, DofEdge> perEdge1 = make_pair(edge2, edge1);
+
+ TEST_EXIT_DBG(periodicEdges.count(make_pair(edge0, edge1)) == 1)
+ ("Should not happen!\n");
+ TEST_EXIT_DBG(periodicEdges.count(make_pair(edge1, edge0)) == 1)
+ ("Should not happen!\n");
+
+ BoundaryType type0 = periodicEdges[make_pair(edge0, edge1)];
+ BoundaryType type1 = periodicEdges[make_pair(edge0, edge2)];
+ BoundaryType type2 = provideConnectedPeriodicBoundary(type0, type1);
+
+ periodicEdges[perEdge0] = type2;
+ periodicEdges[perEdge1] = type2;
}
}
@@ -314,6 +311,59 @@ namespace AMDiS {
}
+ BoundaryType ElementObjects::getNewBoundaryType()
+ {
+ FUNCNAME("ElementObjects::getNewBoundaryType()");
+
+ BoundaryType newPeriodicBoundaryType = 0;
+ for (map<BoundaryType, VertexVector*>::iterator it = mesh->getPeriodicAssociations().begin();
+ it != mesh->getPeriodicAssociations().end(); ++it)
+ newPeriodicBoundaryType = std::min(newPeriodicBoundaryType, it->first);
+
+ TEST_EXIT_DBG(newPeriodicBoundaryType < 0)("Should not happen!\n");
+ newPeriodicBoundaryType--;
+
+ mesh->getPeriodicAssociations()[newPeriodicBoundaryType] =
+ new VertexVector(feSpace->getAdmin(), "");
+
+ return newPeriodicBoundaryType;
+ }
+
+
+ BoundaryType ElementObjects::provideConnectedPeriodicBoundary(BoundaryType b0,
+ BoundaryType b1)
+ {
+ FUNCNAME("ElementObjects::provideConnectedPeriodicBoundary()");
+
+ std::pair<BoundaryType, BoundaryType> bConn =
+ (b0 <= b1 ? std::make_pair(b0, b1) : std::make_pair(b1, b0));
+
+ if (bConnMap.count(bConn) == 0) {
+ BoundaryType newPeriodicBoundaryType = getNewBoundaryType();
+
+ VertexVector &vecB0 = mesh->getPeriodicAssociations(b0);
+ VertexVector &vecB1 = mesh->getPeriodicAssociations(b1);
+ VertexVector &vecC = mesh->getPeriodicAssociations(newPeriodicBoundaryType);
+
+ DOFIteratorBase it(const_cast<DOFAdmin*>(feSpace->getAdmin()), USED_DOFS);
+
+ for (it.reset(); !it.end(); ++it) {
+ if (!it.isDofFree()) {
+ TEST_EXIT_DBG(vecB1[vecB0[it.getDOFIndex()]] == vecB0[vecB1[it.getDOFIndex()]])
+ ("Should not happen!\n");
+
+ vecC[it.getDOFIndex()] = vecB1[vecB0[it.getDOFIndex()]];
+ }
+ }
+
+ bConnMap[bConn] = newPeriodicBoundaryType;
+ }
+
+
+ return bConnMap[bConn];
+ }
+
+
void ElementObjects::createRankData(map<int, int>& macroElementRankMap)
{
FUNCNAME("ElementObjects::createRankData()");
@@ -329,7 +379,7 @@ namespace AMDiS {
if (it2->elIndex > vertexInRank[it->first][elementInRank].elIndex)
vertexInRank[it->first][elementInRank] = *it2;
- vertexOwner[it->first] = max(vertexOwner[it->first], elementInRank);
+ vertexOwner[it->first] = std::max(vertexOwner[it->first], elementInRank);
}
}
@@ -345,7 +395,7 @@ namespace AMDiS {
if (it2->elIndex > edgeInRank[it->first][elementInRank].elIndex)
edgeInRank[it->first][elementInRank] = *it2;
- edgeOwner[it->first] = max(edgeOwner[it->first], elementInRank);
+ edgeOwner[it->first] = std::max(edgeOwner[it->first], elementInRank);
}
}
@@ -361,7 +411,7 @@ namespace AMDiS {
if (it2->elIndex > faceInRank[it->first][elementInRank].elIndex)
faceInRank[it->first][elementInRank] = *it2;
- faceOwner[it->first] = max(faceOwner[it->first], elementInRank);
+ faceOwner[it->first] = std::max(faceOwner[it->first], elementInRank);
}
}
}
diff --git a/AMDiS/src/parallel/ElementObjectData.h b/AMDiS/src/parallel/ElementObjectData.h
index bbe33bcc..716d4345 100644
--- a/AMDiS/src/parallel/ElementObjectData.h
+++ b/AMDiS/src/parallel/ElementObjectData.h
@@ -360,6 +360,21 @@ namespace AMDiS {
return faceLocalMap[data];
}
+ PerBoundMap<DegreeOfFreedom>::type& getPeriodicVertices()
+ {
+ return periodicVertices;
+ }
+
+ PerBoundMap<DofEdge>::type& getPeriodicEdges()
+ {
+ return periodicEdges;
+ }
+
+ PerBoundMap<DofFace>::type& getPeriodicFaces()
+ {
+ return periodicFaces;
+ }
+
/// Write the element database to disk.
void serialize(ostream &out);
@@ -400,6 +415,10 @@ namespace AMDiS {
faceLocalMap[elObj] = face;
}
+ BoundaryType provideConnectedPeriodicBoundary(BoundaryType b0,
+ BoundaryType b1);
+
+ BoundaryType getNewBoundaryType();
/// Some auxiliary function to write the element object database to disk.
void serialize(ostream &out, vector<ElementObjectData>& elVec);
@@ -480,7 +499,8 @@ namespace AMDiS {
/// iterators, i.e., if no iteration is running.
GeoIndex iterGeoPos;
- public:
+ std::map<std::pair<BoundaryType, BoundaryType>, BoundaryType> bConnMap;
+
// The following three data structures store periodic DOFs, edges and faces.
PerBoundMap<DegreeOfFreedom>::type periodicVertices;
PerBoundMap<DofEdge>::type periodicEdges;
diff --git a/AMDiS/src/parallel/InteriorBoundary.cc b/AMDiS/src/parallel/InteriorBoundary.cc
index 612d82d7..539a1519 100644
--- a/AMDiS/src/parallel/InteriorBoundary.cc
+++ b/AMDiS/src/parallel/InteriorBoundary.cc
@@ -57,12 +57,17 @@ namespace AMDiS {
TEST_EXIT_DBG(localDofs0[el0_v1] == localDofs1[el1_v0] ||
localDofs0[el0_v1] == localDofs1[el1_v1])
("This should not happen!\n");
-
+
if (localDofs0[el0_v0] != localDofs1[el1_v0])
otherMode = true;
} else {
if (mesh->associated(localDofs0[el0_v0], localDofs1[el1_v0]) == false)
otherMode = true;
+
+ if ((elIndex == 28 && otherBound.elIndex == 41) ||
+ (elIndex == 41 && otherBound.elIndex == 28)) {
+ MSG("HERE TEST B (%d %d): %d\n", el0_v0, el1_v0, otherMode);
+ }
}
}
diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc
index f3490093..cf1b28a7 100644
--- a/AMDiS/src/parallel/MeshDistributor.cc
+++ b/AMDiS/src/parallel/MeshDistributor.cc
@@ -47,6 +47,7 @@ namespace AMDiS {
using boost::lexical_cast;
using namespace boost::filesystem;
+ using namespace std;
inline bool cmpDofsByValue(const DegreeOfFreedom* dof1, const DegreeOfFreedom* dof2)
{
@@ -54,7 +55,7 @@ namespace AMDiS {
}
- MeshDistributor::MeshDistributor(std::string str)
+ MeshDistributor::MeshDistributor(string str)
: probStat(0),
name(str),
feSpace(NULL),
@@ -118,15 +119,15 @@ namespace AMDiS {
macroElIndexMap.clear();
macroElIndexTypeMap.clear();
- std::map<int, bool>& elementInRank = partitioner->getElementInRank();
- for (std::vector<MacroElement*>::iterator it = allMacroElements.begin();
+ map<int, bool>& elementInRank = partitioner->getElementInRank();
+ for (vector<MacroElement*>::iterator it = allMacroElements.begin();
it != allMacroElements.end(); ++it) {
elementInRank[(*it)->getIndex()] = false;
macroElIndexMap[(*it)->getIndex()] = (*it)->getElement();
macroElIndexTypeMap[(*it)->getIndex()] = (*it)->getElType();
}
- for (std::deque<MacroElement*>::iterator it = mesh->getMacroElements().begin();
+ for (deque<MacroElement*>::iterator it = mesh->getMacroElements().begin();
it != mesh->getMacroElements().end(); ++it)
elementInRank[(*it)->getIndex()] = true;
@@ -183,6 +184,33 @@ namespace AMDiS {
ParallelDebug::printBoundaryInfo(*this);
#endif
+ for (deque<MacroElement*>::iterator it = mesh->firstMacroElement();
+ it != mesh->endOfMacroElements(); ++it) {
+ for (int i = 0; i < mesh->getGeo(NEIGH); i++) {
+ if ((*it)->getNeighbour(i) &&
+ mesh->isPeriodicAssociation((*it)->getBoundary(i))) {
+ (*it)->getNeighbour(i)->setNeighbour((*it)->getOppVertex(i), NULL);
+ (*it)->setNeighbour(i, NULL);
+ (*it)->setBoundary(i, 0);
+
+ macroElementNeighbours[(*it)->getIndex()][i] = -1;
+ }
+ }
+ }
+
+ for (vector<MacroElement*>::iterator it = allMacroElements.begin();
+ it != allMacroElements.end(); ++it) {
+ for (int i = 0; i < mesh->getGeo(NEIGH); i++) {
+ if ((*it)->getNeighbour(i) &&
+ mesh->isPeriodicAssociation((*it)->getBoundary(i))) {
+ (*it)->getNeighbour(i)->setNeighbour((*it)->getOppVertex(i), NULL);
+ (*it)->setNeighbour(i, NULL);
+ (*it)->setBoundary(i, 0);
+
+ macroElementNeighbours[(*it)->getIndex()][i] = -1;
+ }
+ }
+ }
// === Remove all macro elements that are not part of the rank partition. ===
@@ -195,7 +223,7 @@ namespace AMDiS {
// We have to remove the VertexVectors, which contain periodic assoiciations,
// because they are not valid anymore after some macro elements have been removed
// and the corresponding DOFs were deleted.
- for (std::map<BoundaryType, VertexVector*>::iterator it = mesh->getPeriodicAssociations().begin();
+ for (map<BoundaryType, VertexVector*>::iterator it = mesh->getPeriodicAssociations().begin();
it != mesh->getPeriodicAssociations().end(); ++it)
const_cast<DOFAdmin&>(mesh->getDofAdmin(0)).removeDOFContainer(dynamic_cast<DOFContainer*>(it->second));
@@ -221,17 +249,14 @@ namespace AMDiS {
// === Create periodic dof mapping, if there are periodic boundaries. ===
- createPeriodicMap();
+ createPeriodicMap();
-#if (DEBUG != 0)
- ParallelDebug::testPeriodicBoundary(*this);
-#endif
// === Global refinements. ===
-
+
int globalRefinement = 0;
GET_PARAMETER(0, mesh->getName() + "->global refinements", "%d", &globalRefinement);
-
+
if (globalRefinement > 0) {
refineManager->globalRefine(mesh, globalRefinement);
@@ -260,7 +285,7 @@ namespace AMDiS {
FUNCNAME("MeshDistributor::addProblemVec()");
if (feSpace != NULL) {
- std::vector<FiniteElemSpace*> feSpaces = probVec->getFeSpaces();
+ vector<FiniteElemSpace*> feSpaces = probVec->getFeSpaces();
for (unsigned int i = 0; i < feSpaces.size(); i++) {
TEST_EXIT(feSpace == feSpaces[i])
("Parallelizaton is not supported for multiple FE spaces!\n");
@@ -294,7 +319,7 @@ namespace AMDiS {
GET_PARAMETER(0, probVec->getName() + "->output->write serialization", "%d",
&writeSerialization);
if (writeSerialization && !writeSerializationFile) {
- std::string filename = "";
+ string filename = "";
GET_PARAMETER(0, name + "->output->serialization filename", &filename);
TEST_EXIT(filename != "")
@@ -312,11 +337,11 @@ namespace AMDiS {
GET_PARAMETER(0, probVec->getName() + "->input->read serialization", "%d",
&readSerialization);
if (readSerialization) {
- std::string filename = "";
+ string filename = "";
GET_PARAMETER(0, probVec->getName() + "->input->serialization filename", &filename);
- filename += ".p" + lexical_cast<std::string>(mpiRank);
+ filename += ".p" + lexical_cast<string>(mpiRank);
MSG("Start deserialization with %s\n", filename.c_str());
- std::ifstream in(filename.c_str());
+ ifstream in(filename.c_str());
TEST_EXIT(!in.fail())("Could not open deserialization file: %s\n",
filename.c_str());
@@ -331,7 +356,7 @@ namespace AMDiS {
TEST_EXIT(filename != "")
("No filename defined for parallel deserialization file!\n");
- std::string rankFilename = filename + ".p" + lexical_cast<std::string>(mpiRank);
+ string rankFilename = filename + ".p" + lexical_cast<string>(mpiRank);
in.open(rankFilename.c_str());
TEST_EXIT(!in.fail())("Could not open parallel deserialization file: %s\n",
@@ -378,11 +403,11 @@ namespace AMDiS {
void MeshDistributor::synchVector(DOFVector<double> &vec)
{
- StdMpi<std::vector<double> > stdMpi(mpiComm);
+ StdMpi<vector<double> > stdMpi(mpiComm);
for (RankToDofContainer::iterator sendIt = sendDofs.begin();
sendIt != sendDofs.end(); ++sendIt) {
- std::vector<double> dofs;
+ vector<double> dofs;
int nSendDofs = sendIt->second.size();
dofs.reserve(nSendDofs);
@@ -408,11 +433,11 @@ namespace AMDiS {
void MeshDistributor::synchVector(SystemVector &vec)
{
int nComponents = vec.getSize();
- StdMpi<std::vector<double> > stdMpi(mpiComm);
+ StdMpi<vector<double> > stdMpi(mpiComm);
for (RankToDofContainer::iterator sendIt = sendDofs.begin();
sendIt != sendDofs.end(); ++sendIt) {
- std::vector<double> dofs;
+ vector<double> dofs;
int nSendDofs = sendIt->second.size();
dofs.reserve(nComponents * nSendDofs);
@@ -517,10 +542,6 @@ namespace AMDiS {
{
FUNCNAME("MeshDistributor::checkMeshChange()");
-#if (DEBUG != 0)
- debug::writeMesh(feSpace, -1, debugOutputDir + "before_check_mesh");
-#endif
-
double first = MPI::Wtime();
// === If mesh has not been changed on all ranks, return. ===
@@ -616,11 +637,11 @@ namespace AMDiS {
// === Create mesh structure codes for all ranks boundary elements. ===
- std::map<int, MeshCodeVec> sendCodes;
+ map<int, MeshCodeVec> sendCodes;
for (RankToBoundMap::iterator it = allBound.begin(); it != allBound.end(); ++it) {
- for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
+ for (vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
MeshStructure elCode;
elCode.init(boundIt->rankObj);
@@ -642,7 +663,7 @@ namespace AMDiS {
MeshCodeVec &recvCodes = stdMpi.getRecvData()[it->first];
int i = 0;
- for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
+ for (vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
MeshStructure elCode;
@@ -748,8 +769,8 @@ namespace AMDiS {
Element *child0 = el->getFirstChild();
Element *child1 = el->getSecondChild();
if (reverseMode) {
- std::swap(s0, s1);
- std::swap(child0, child1);
+ swap(s0, s1);
+ swap(child0, child1);
}
// === We know that the edge/face is contained in only one of the children. ===
@@ -832,8 +853,8 @@ namespace AMDiS {
Element *child0 = el->getFirstChild();
Element *child1 = el->getSecondChild();
if (reverseMode) {
- std::swap(s0, s1);
- std::swap(child0, child1);
+ swap(s0, s1);
+ swap(child0, child1);
}
@@ -889,7 +910,7 @@ namespace AMDiS {
}
- void MeshDistributor::serialize(std::ostream &out, DofContainer &data)
+ void MeshDistributor::serialize(ostream &out, DofContainer &data)
{
int vecSize = data.size();
SerUtil::serialize(out, vecSize);
@@ -900,8 +921,8 @@ namespace AMDiS {
}
- void MeshDistributor::deserialize(std::istream &in, DofContainer &data,
- std::map<int, const DegreeOfFreedom*> &dofMap)
+ void MeshDistributor::deserialize(istream &in, DofContainer &data,
+ map<int, const DegreeOfFreedom*> &dofMap)
{
FUNCNAME("MeshDistributor::deserialize()");
@@ -921,7 +942,7 @@ namespace AMDiS {
}
- void MeshDistributor::serialize(std::ostream &out, RankToDofContainer &data)
+ void MeshDistributor::serialize(ostream &out, RankToDofContainer &data)
{
int mapSize = data.size();
SerUtil::serialize(out, mapSize);
@@ -933,8 +954,8 @@ namespace AMDiS {
}
- void MeshDistributor::deserialize(std::istream &in, RankToDofContainer &data,
- std::map<int, const DegreeOfFreedom*> &dofMap)
+ void MeshDistributor::deserialize(istream &in, RankToDofContainer &data,
+ map<int, const DegreeOfFreedom*> &dofMap)
{
data.clear();
@@ -954,13 +975,13 @@ namespace AMDiS {
elemWeights.clear();
- std::string filename = "";
+ string filename = "";
GET_PARAMETER(0, mesh->getName() + "->macro weights", &filename);
if (filename != "") {
MSG("Read macro weights from %s\n", filename.c_str());
- std::ifstream infile;
- infile.open(filename.c_str(), std::ifstream::in);
+ ifstream infile;
+ infile.open(filename.c_str(), ifstream::in);
while (!infile.eof()) {
int elNum, elWeight;
infile >> elNum;
@@ -993,13 +1014,13 @@ namespace AMDiS {
int nDofs = mesh->getDofAdmin(0).getUsedDofs();
int repartitioning = 0;
- std::vector<int> nDofsInRank(mpiSize);
+ vector<int> nDofsInRank(mpiSize);
mpiComm.Gather(&nDofs, 1, MPI_INT, &(nDofsInRank[0]), 1, MPI_INT, 0);
if (mpiRank == 0) {
int nOverallDofs = 0;
- int minDofs = std::numeric_limits<int>::max();
- int maxDofs = std::numeric_limits<int>::min();
+ int minDofs = numeric_limits<int>::max();
+ int maxDofs = numeric_limits<int>::min();
for (int i = 0; i < mpiSize; i++) {
nOverallDofs += nDofsInRank[i];
minDofs = std::min(minDofs, nDofsInRank[i]);
@@ -1027,7 +1048,7 @@ namespace AMDiS {
ParallelDebug::testDoubleDofs(mesh);
if (repartCounter == 0) {
- std::stringstream oss;
+ stringstream oss;
oss << debugOutputDir << "partitioning-" << repartCounter << ".vtu";
MSG("Write partitioning to %s\n", oss.str().c_str());
DOFVector<double> tmpa(feSpace, "tmp");
@@ -1063,8 +1084,8 @@ namespace AMDiS {
// === Create map that maps macro element indices to pointers to the ===
// === macro elements. ===
- std::map<int, MacroElement*> elIndexMap;
- for (std::vector<MacroElement*>::iterator it = allMacroElements.begin();
+ map<int, MacroElement*> elIndexMap;
+ for (vector<MacroElement*>::iterator it = allMacroElements.begin();
it != allMacroElements.end(); ++it)
elIndexMap[(*it)->getIndex()] = *it;
@@ -1073,9 +1094,9 @@ namespace AMDiS {
// === other ranks. ===
std::set<MacroElement*> newMacroEl;
- for (std::map<int, std::vector<int> >::iterator it = partitioner->getRecvElements().begin();
+ for (map<int, vector<int> >::iterator it = partitioner->getRecvElements().begin();
it != partitioner->getRecvElements().end(); ++it) {
- for (std::vector<int>::iterator elIt = it->second.begin();
+ for (vector<int>::iterator elIt = it->second.begin();
elIt != it->second.end(); ++elIt) {
TEST_EXIT_DBG(elIndexMap.count(*elIt) == 1)
("Could not find macro element %d\n", *elIt);
@@ -1089,7 +1110,7 @@ namespace AMDiS {
it != newMacroEl.end(); ++it)
allMacroEl.insert(*it);
- for (std::deque<MacroElement*>::iterator it = mesh->firstMacroElement();
+ for (deque<MacroElement*>::iterator it = mesh->firstMacroElement();
it != mesh->endOfMacroElements(); ++it)
allMacroEl.insert(*it);
@@ -1118,19 +1139,19 @@ namespace AMDiS {
// === Send and receive mesh structure codes. ===
- std::map<int, MeshCodeVec> sendCodes;
- std::map<int, std::vector<std::vector<double> > > sendValues;
+ map<int, MeshCodeVec> sendCodes;
+ map<int, vector<vector<double> > > sendValues;
- for (std::map<int, std::vector<int> >::iterator it = partitioner->getSendElements().begin();
+ for (map<int, vector<int> >::iterator it = partitioner->getSendElements().begin();
it != partitioner->getSendElements().end(); ++it) {
- for (std::vector<int>::iterator elIt = it->second.begin();
+ for (vector<int>::iterator elIt = it->second.begin();
elIt != it->second.end(); ++elIt) {
MeshStructure elCode;
elCode.init(mesh, *elIt);
sendCodes[it->first].push_back(elCode);
for (unsigned int i = 0; i < testVec.size(); i++) {
- std::vector<double> valVec;
+ vector<double> valVec;
elCode.getMeshStructureValues(mesh, *elIt, testVec[i], valVec);
sendValues[it->first].push_back(valVec);
}
@@ -1142,24 +1163,24 @@ namespace AMDiS {
stdMpi.recv(partitioner->getRecvElements());
stdMpi.startCommunication<uint64_t>(MPI_UNSIGNED_LONG);
- StdMpi<std::vector<std::vector<double> > > stdMpi2(mpiComm, true);
+ StdMpi<vector<vector<double> > > stdMpi2(mpiComm, true);
stdMpi2.send(sendValues);
stdMpi2.recv(partitioner->getRecvElements());
stdMpi2.startCommunication<double>(MPI_DOUBLE);
// === Adapte the new macro elements due to the received mesh structure codes. ===
- for (std::map<int, std::vector<int> >::iterator it = partitioner->getRecvElements().begin();
+ for (map<int, vector<int> >::iterator it = partitioner->getRecvElements().begin();
it != partitioner->getRecvElements().end(); ++it) {
MeshCodeVec &recvCodes = stdMpi.getRecvData()[it->first];
- std::vector<std::vector<double> > &recvValues = stdMpi2.getRecvData()[it->first];
+ vector<vector<double> > &recvValues = stdMpi2.getRecvData()[it->first];
int i = 0;
int j = 0;
TEST_EXIT_DBG(recvCodes.size() == it->second.size())
("Should not happen!\n");
- for (std::vector<int>::iterator elIt = it->second.begin();
+ for (vector<int>::iterator elIt = it->second.begin();
elIt != it->second.end(); ++elIt) {
recvCodes[i].fitMeshToStructure(mesh, refineManager, false, *elIt);
@@ -1202,9 +1223,9 @@ namespace AMDiS {
std::set<MacroElement*> deleteMacroElements;
- for (std::map<int, std::vector<int> >::iterator it = partitioner->getSendElements().begin();
+ for (map<int, vector<int> >::iterator it = partitioner->getSendElements().begin();
it != partitioner->getSendElements().end(); ++it) {
- for (std::vector<int>::iterator elIt = it->second.begin();
+ for (vector<int>::iterator elIt = it->second.begin();
elIt != it->second.end(); ++elIt) {
TEST_EXIT_DBG(elIndexMap.count(*elIt) == 1)
("Could not find macro element %d\n", *elIt);
@@ -1227,7 +1248,7 @@ namespace AMDiS {
mesh->dofCompress();
#if (DEBUG != 0)
- std::stringstream oss;
+ stringstream oss;
oss << debugOutputDir << "partitioning-" << repartCounter << ".vtu";
DOFVector<double> tmpa(feSpace, "tmp");
tmpa.set(mpiRank);
@@ -1331,7 +1352,7 @@ namespace AMDiS {
GeoIndex geoIndex = INDEX_OF_DIM(geoPos, mesh->getDim());
while (elObjects.iterate(geoIndex)) {
- std::map<int, ElementObjectData>& objData = elObjects.getIterateData();
+ map<int, ElementObjectData>& objData = elObjects.getIterateData();
if (objData.count(mpiRank) && objData.size() > 1) {
int owner = elObjects.getIterateOwner();
ElementObjectData& rankBoundEl = objData[mpiRank];
@@ -1350,13 +1371,13 @@ namespace AMDiS {
int edgeOfFace =
bound.rankObj.el->getEdgeOfFace(bound.rankObj.ithObj, edgeNo);
- bound.rankObj.excludedSubstructures.push_back(std::make_pair(EDGE, edgeOfFace));
+ bound.rankObj.excludedSubstructures.push_back(make_pair(EDGE, edgeOfFace));
}
}
if (owner == mpiRank) {
- for (std::map<int, ElementObjectData>::iterator it2 = objData.begin();
+ for (map<int, ElementObjectData>::iterator it2 = objData.begin();
it2 != objData.end(); ++it2) {
if (it2->first == mpiRank)
continue;
@@ -1399,8 +1420,8 @@ namespace AMDiS {
// === Create periodic boundary data structure. ===
- for (PerBoundMap<DegreeOfFreedom>::iterator it = elObjects.periodicVertices.begin();
- it != elObjects.periodicVertices.end(); ++it) {
+ for (PerBoundMap<DegreeOfFreedom>::iterator it = elObjects.getPeriodicVertices().begin();
+ it != elObjects.getPeriodicVertices().end(); ++it) {
if (elObjects.isInRank(it->first.first, mpiRank) == false)
continue;
@@ -1413,7 +1434,7 @@ namespace AMDiS {
ElementObjectData& perDofEl0 = elObjects.getElementsInRank(it->first.first)[mpiRank];
// MSG("DATA: %d %d %d\n", perDofEl0.elIndex, VERTEX, perDofEl0.ithObject);
- for (std::map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
+ for (map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
elIt != elObjects.getElementsInRank(it->first.second).end(); ++elIt) {
int otherElementRank = elIt->first;
@@ -1440,14 +1461,14 @@ namespace AMDiS {
}
- for (PerBoundMap<DofEdge>::iterator it = elObjects.periodicEdges.begin();
- it != elObjects.periodicEdges.end(); ++it) {
+ for (PerBoundMap<DofEdge>::iterator it = elObjects.getPeriodicEdges().begin();
+ it != elObjects.getPeriodicEdges().end(); ++it) {
if (elObjects.isInRank(it->first.first, mpiRank) == false)
continue;
ElementObjectData& perEdgeEl0 = elObjects.getElementsInRank(it->first.first)[mpiRank];
- for (std::map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
+ for (map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
elIt != elObjects.getElementsInRank(it->first.second).end(); ++elIt) {
int otherElementRank = elIt->first;
@@ -1479,8 +1500,8 @@ namespace AMDiS {
}
- for (PerBoundMap<DofFace>::iterator it = elObjects.periodicFaces.begin();
- it != elObjects.periodicFaces.end(); ++it) {
+ for (PerBoundMap<DofFace>::iterator it = elObjects.getPeriodicFaces().begin();
+ it != elObjects.getPeriodicFaces().end(); ++it) {
if (elObjects.isInRank(it->first.first, mpiRank) == false)
continue;
@@ -1491,7 +1512,7 @@ namespace AMDiS {
ElementObjectData& perEdgeEl0 = elObjects.getElementsInRank(it->first.first)[mpiRank];
- for (std::map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
+ for (map<int, ElementObjectData>::iterator elIt = elObjects.getElementsInRank(it->first.second).begin();
elIt != elObjects.getElementsInRank(it->first.second).end(); ++elIt) {
int otherElementRank = elIt->first;
@@ -1527,7 +1548,7 @@ namespace AMDiS {
// === bounds in the three boundary handling object. Eventually all the bound- ===
// === aries have to be in the same order on both ranks that share the bounday. ===
- StdMpi<std::vector<AtomicBoundary> > stdMpi(mpiComm);
+ StdMpi<vector<AtomicBoundary> > stdMpi(mpiComm);
stdMpi.send(myIntBoundary.boundary);
stdMpi.recv(otherIntBoundary.boundary);
stdMpi.startCommunication<int>(MPI_INT);
@@ -1649,7 +1670,7 @@ namespace AMDiS {
FUNCNAME("MeshDistributor::removeMacroElements()");
std::set<MacroElement*> macrosToRemove;
- for (std::deque<MacroElement*>::iterator it = mesh->firstMacroElement();
+ for (deque<MacroElement*>::iterator it = mesh->firstMacroElement();
it != mesh->endOfMacroElements(); ++it)
if (partitioner->getElementInRank()[(*it)->getIndex()] == false)
macrosToRemove.insert(*it);
@@ -1760,7 +1781,7 @@ namespace AMDiS {
#endif
int i = 0;
- StdMpi<std::vector<DegreeOfFreedom> > stdMpi(mpiComm, false);
+ StdMpi<vector<DegreeOfFreedom> > stdMpi(mpiComm, false);
for (RankToDofContainer::iterator sendIt = sendDofs.begin();
sendIt != sendDofs.end(); ++sendIt, i++) {
stdMpi.getSendData(sendIt->first).resize(0);
@@ -1804,7 +1825,7 @@ namespace AMDiS {
#if (DEBUG != 0)
- std::stringstream oss;
+ stringstream oss;
oss << debugOutputDir << "elementIndex-" << mpiRank << ".vtu";
debug::writeElementIndexMesh(mesh, oss.str());
@@ -1813,7 +1834,7 @@ namespace AMDiS {
ParallelDebug::testGlobalIndexByCoords(*this);
ParallelDebug::writeDebugFile(*this, debugOutputDir + "mpi-dbg", "dat");
-#if 0
+#if 1
MSG("------------- Debug information -------------\n");
MSG("nRankDofs = %d\n", nRankDofs);
MSG("nOverallDofs = %d\n", nOverallDofs);
@@ -1850,18 +1871,18 @@ namespace AMDiS {
periodicDof.clear();
periodicDofAssociations.clear();
- StdMpi<std::vector<int> > stdMpi(mpiComm, false);
+ StdMpi<vector<int> > stdMpi(mpiComm, false);
// === Each rank traverse its periodic boundaries and sends the dof indices ===
// === to the rank "on the other side" of the periodic boundary. ===
RankToDofContainer rankPeriodicDofs;
- std::map<int, std::vector<int> > rankToDofType;
+ map<int, vector<int> > rankToDofType;
for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
it != periodicBoundary.boundary.end(); ++it) {
- MSG("------------- WITH RANK %d ------------------\n", it->first);
+// MSG("------------- WITH RANK %d ------------------\n", it->first);
if (it->first == mpiRank) {
TEST_EXIT_DBG(it->second.size() % 2 == 0)("Should not happen!\n");
@@ -1874,9 +1895,9 @@ namespace AMDiS {
bound.neighObj.el->getVertexDofs(feSpace, bound.neighObj, dofs1);
bound.neighObj.el->getNonVertexDofs(feSpace, bound.neighObj, dofs1);
- MSG("BOUND-I %d-%d-%d WITH %d-%d-%d\n",
- bound.rankObj.elIndex, bound.rankObj.subObj, bound.rankObj.ithObj,
- bound.neighObj.elIndex, bound.neighObj.subObj, bound.neighObj.ithObj);
+// MSG("BOUND-I %d-%d-%d WITH %d-%d-%d\n",
+// bound.rankObj.elIndex, bound.rankObj.subObj, bound.rankObj.ithObj,
+// bound.neighObj.elIndex, bound.neighObj.subObj, bound.neighObj.ithObj);
TEST_EXIT_DBG(dofs0.size() == dofs1.size())("Should not happen!\n");
@@ -1890,7 +1911,7 @@ namespace AMDiS {
periodicDof[type][globalDof0] = globalDof1;
periodicDofAssociations[globalDof0].insert(type);
- MSG("SET(A TYPE %d) DOF %d -> %d\n", type, globalDof0, globalDof1);
+// MSG("SET(A TYPE %d) DOF %d -> %d\n", type, globalDof0, globalDof1);
}
}
}
@@ -1898,22 +1919,24 @@ namespace AMDiS {
} else {
// Create DOF indices on the boundary.
DofContainer& dofs = rankPeriodicDofs[it->first];
- for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
+ for (vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
int nDofs = dofs.size();
- MSG("BOUND-R (T = %d) %d-%d-%d WITH %d-%d-%d\n",
- boundIt->type,
- boundIt->rankObj.elIndex, boundIt->rankObj.subObj, boundIt->rankObj.ithObj,
- boundIt->neighObj.elIndex, boundIt->neighObj.subObj, boundIt->neighObj.ithObj);
+// MSG("BOUND-R (T = %d) %d-%d-%d WITH %d-%d-%d\n",
+// boundIt->type,
+// boundIt->rankObj.elIndex, boundIt->rankObj.subObj, boundIt->rankObj.ithObj,
+// boundIt->neighObj.elIndex, boundIt->neighObj.subObj, boundIt->neighObj.ithObj);
boundIt->rankObj.el->getVertexDofs(feSpace, boundIt->rankObj, dofs);
boundIt->rankObj.el->getNonVertexDofs(feSpace, boundIt->rankObj, dofs);
for (unsigned int i = 0; i < (dofs.size() - nDofs); i++) {
- MSG(" i = %d DOF = %d\n", nDofs + i, mapLocalGlobalDofs[*(dofs[nDofs + i])]);
+// WorldVector<double> c;
+// mesh->getDofIndexCoords(*(dofs[nDofs + i]), feSpace, c);
+// MSG(" SEND i = %d DOF = %d (%f %f %f)\n", nDofs + i, mapLocalGlobalDofs[*(dofs[nDofs + i])], c[0], c[1], c[2]);
rankToDofType[it->first].push_back(boundIt->type);
}
}
@@ -1921,7 +1944,7 @@ namespace AMDiS {
// Send the global indices to the rank on the other side.
stdMpi.getSendData(it->first).reserve(dofs.size());
for (unsigned int i = 0; i < dofs.size(); i++)
- stdMpi.getSendData(it->first).push_back(mapLocalGlobalDofs[*(dofs[i])]);
+ stdMpi.getSendData(it->first).push_back(mapLocalGlobalDofs[*(dofs[i])]);
// Receive from this rank the same number of dofs.
stdMpi.recv(it->first, dofs.size());
@@ -1932,8 +1955,6 @@ namespace AMDiS {
stdMpi.startCommunication<int>(MPI_INT);
- MSG("===============================\n");
-
// === The rank has received the dofs from the rank on the other side of ===
// === the boundary. Now it can use them to create the mapping between ===
// === the periodic dofs in this rank and the corresponding periodic ===
@@ -1943,10 +1964,10 @@ namespace AMDiS {
for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
it != periodicBoundary.boundary.end(); ++it) {
- MSG("------------- WITH RANK %d ------------------\n", it->first);
+// MSG("------------- WITH RANK %d ------------------\n", it->first);
DofContainer& dofs = rankPeriodicDofs[it->first];
- std::vector<int>& types = rankToDofType[it->first];
+ vector<int>& types = rankToDofType[it->first];
TEST_EXIT_DBG(dofs.size() == types.size())("Should not happen!\n");
@@ -1956,31 +1977,28 @@ namespace AMDiS {
int mapGlobalDofIndex = stdMpi.getRecvData(it->first)[i];
BoundaryType type = types[i];
+// WorldVector<double> c;
+// mesh->getDofIndexCoords(*(dofs[i]), feSpace, c);
+// MSG(" RECV i = %d DOF = %d (%f %f %f)\n", i, *(dofs[i]), c[0], c[1], c[2]);
+
+
// Check if this global dof with the corresponding boundary type was
// not added before by another periodic boundary from other rank.
if (periodicDofAssociations[globalDofIndex].count(type) == 0) {
- MSG("SET(B-%d TYPE %d) DOF %d -> %d\n", i, type, globalDofIndex, mapGlobalDofIndex);
+// MSG("SET(B-%d TYPE %d) DOF %d -> %d\n", i, type, globalDofIndex, mapGlobalDofIndex);
periodicDof[type][globalDofIndex] = mapGlobalDofIndex;
periodicDofAssociations[globalDofIndex].insert(type);
} else {
- MSG("ASSOC ALREADY SET FOR %d TYPE %d\n", i, type);
+ // MSG("ASSOC ALREADY SET FOR %d TYPE %d\n", i, type);
}
}
}
#if (DEBUG != 0)
- for (std::map<int, std::set<BoundaryType> >::iterator it = periodicDofAssociations.begin();
- it != periodicDofAssociations.end(); ++it) {
- WorldVector<double> c;
- mesh->getDofIndexCoords(it->first, feSpace, c);
- int nAssoc = it->second.size();
- TEST_EXIT_DBG(nAssoc == 1 || nAssoc == 3 || (mesh->getDim() == 3 && nAssoc == 7))
- ("Should not happen! DOF %d (%e %e %e) has %d periodic associations!\n",
- it->first, c[0], c[1], (mesh->getDim() == 2 ? 0.0 : c[2]), nAssoc);
- }
+ ParallelDebug::testPeriodicBoundary(*this);
#endif
}
@@ -1989,8 +2007,9 @@ namespace AMDiS {
{
FUNCNAME("MeshDistributor::createMacroElementInfo()");
- for (std::deque<MacroElement*>::iterator it = mesh->firstMacroElement();
+ for (deque<MacroElement*>::iterator it = mesh->firstMacroElement();
it != mesh->endOfMacroElements(); ++it) {
+
allMacroElements.push_back(*it);
macroElementNeighbours[(*it)->getIndex()].resize(mesh->getGeo(NEIGH));
@@ -2003,10 +2022,10 @@ namespace AMDiS {
void MeshDistributor::updateMacroElementInfo()
{
- std::deque<MacroElement*>& meshMacros = mesh->getMacroElements();
+ deque<MacroElement*>& meshMacros = mesh->getMacroElements();
for (unsigned int i = 0; i < allMacroElements.size(); i++) {
- for (std::deque<MacroElement*>::iterator it = meshMacros.begin();
+ for (deque<MacroElement*>::iterator it = meshMacros.begin();
it != meshMacros.end(); ++it) {
if ((*it)->getIndex() == allMacroElements[i]->getIndex() &&
*it != allMacroElements[i]) {
@@ -2018,7 +2037,7 @@ namespace AMDiS {
}
- void MeshDistributor::serialize(std::ostream &out)
+ void MeshDistributor::serialize(ostream &out)
{
FUNCNAME("MeshDistributor::serialize()");
@@ -2059,7 +2078,7 @@ namespace AMDiS {
}
- void MeshDistributor::deserialize(std::istream &in)
+ void MeshDistributor::deserialize(istream &in)
{
FUNCNAME("MeshDistributor::deserialize()");
@@ -2071,8 +2090,8 @@ namespace AMDiS {
// Create two maps: one from from element indices to the corresponding element
// pointers, and one map from Dof indices to the corresponding dof pointers.
- std::map<int, Element*> elIndexMap;
- std::map<int, const DegreeOfFreedom*> dofMap;
+ map<int, Element*> elIndexMap;
+ map<int, const DegreeOfFreedom*> dofMap;
ElementDofIterator elDofIter(feSpace);
TraverseStack stack;
ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
@@ -2114,7 +2133,7 @@ namespace AMDiS {
allMacroElements.resize(nSize);
for (int i = 0; i < nSize; i++) {
// We do not need the neighbour indices, but must still read them.
- std::vector<int> indices;
+ vector<int> indices;
allMacroElements[i] = new MacroElement(mesh->getDim());
allMacroElements[i]->writeNeighboursTo(&indices);
allMacroElements[i]->deserialize(in);
@@ -2128,7 +2147,7 @@ namespace AMDiS {
}
- void MeshDistributor::serialize(std::ostream &out, PeriodicDofMap &data)
+ void MeshDistributor::serialize(ostream &out, PeriodicDofMap &data)
{
int mapSize = data.size();
SerUtil::serialize(out, mapSize);
@@ -2143,12 +2162,12 @@ namespace AMDiS {
}
- void MeshDistributor::serialize(std::ostream &out, std::map<int, std::set<int> >& data)
+ void MeshDistributor::serialize(ostream &out, map<int, std::set<int> >& data)
{
int mapSize = data.size();
SerUtil::serialize(out, mapSize);
- for (std::map<int, std::set<int> >::iterator it = data.begin(); it != data.end(); ++it) {
+ for (map<int, std::set<int> >::iterator it = data.begin(); it != data.end(); ++it) {
int dof = it->first;
std::set<int> typeSet = it->second;
@@ -2158,7 +2177,7 @@ namespace AMDiS {
}
- void MeshDistributor::deserialize(std::istream &in, PeriodicDofMap &data)
+ void MeshDistributor::deserialize(istream &in, PeriodicDofMap &data)
{
data.clear();
@@ -2177,8 +2196,8 @@ namespace AMDiS {
}
- void MeshDistributor::deserialize(std::istream &in,
- std::map<int, std::set<int> >& data)
+ void MeshDistributor::deserialize(istream &in,
+ map<int, std::set<int> >& data)
{
data.clear();
@@ -2197,11 +2216,11 @@ namespace AMDiS {
}
- void MeshDistributor::writePartitioningMesh(std::string filename)
+ void MeshDistributor::writePartitioningMesh(string filename)
{
FUNCNAME("MeshDistributor::writePartitioningMesh()");
- std::map<int, double> vec;
+ map<int, double> vec;
TraverseStack stack;
ElInfo *elInfo = stack.traverseFirst(mesh, -1,
Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
diff --git a/AMDiS/src/parallel/ParMetisPartitioner.cc b/AMDiS/src/parallel/ParMetisPartitioner.cc
index 78d4d3f5..0c0181e6 100644
--- a/AMDiS/src/parallel/ParMetisPartitioner.cc
+++ b/AMDiS/src/parallel/ParMetisPartitioner.cc
@@ -288,7 +288,7 @@ namespace AMDiS {
if (elementInRank[index]) {
// get weight
float wgt = static_cast<float>(elemWeights[index]);
- maxWgt = max(wgt, maxWgt);
+ maxWgt = std::max(wgt, maxWgt);
// write float weight
TEST_EXIT_DBG(floatWgtsPos < floatWgts.size())("Should not happen!\n");
@@ -352,7 +352,7 @@ namespace AMDiS {
int kpartMax = 0;
for (int i = 0; i < nElements; i++)
if (wgts[i] < kpart)
- kpartMax = max(kpartMax, wgts[i]);
+ kpartMax = std::max(kpartMax, wgts[i]);
mpi::globalMax(kpartMax);
@@ -363,8 +363,8 @@ namespace AMDiS {
ssum = 0;
for (int i = 0; i < nElements; i++) {
- wgts[i] = min(wgts[i], kpartMax);
- wgts[i] = max(wgts[i], 1);
+ wgts[i] = std::min(wgts[i], kpartMax);
+ wgts[i] = std::max(wgts[i], 1);
smin = std::min(smin, wgts[i]);
smax = std::max(smax, wgts[i]);
diff --git a/AMDiS/src/parallel/ParallelDebug.cc b/AMDiS/src/parallel/ParallelDebug.cc
index f8cf20b6..4013113c 100644
--- a/AMDiS/src/parallel/ParallelDebug.cc
+++ b/AMDiS/src/parallel/ParallelDebug.cc
@@ -143,6 +143,24 @@ namespace AMDiS {
{
FUNCNAME("ParallelDebug::testPeriodicBoundary()");
+
+ // === 1. check: All periodic DOFs should have at least a correct number ===
+ // === of periodic associations. ===
+
+ for (map<int, std::set<BoundaryType> >::iterator it = pdb.periodicDofAssociations.begin();
+ it != pdb.periodicDofAssociations.end(); ++it) {
+ WorldVector<double> c;
+ pdb.mesh->getDofIndexCoords(it->first, pdb.feSpace, c);
+ int nAssoc = it->second.size();
+ TEST_EXIT_DBG(nAssoc == 1 || nAssoc == 3 || (pdb.mesh->getDim() == 3 && nAssoc == 7))
+ ("Should not happen! DOF %d (%e %e %e) has %d periodic associations!\n",
+ it->first, c[0], c[1], (pdb.mesh->getDim() == 2 ? 0.0 : c[2]), nAssoc);
+ }
+
+
+ // === 2. check: All periodic DOFs must be symmetric, i.e., if A is mapped ===
+ // === to B, then B must be mapped to A. ===
+
typedef MeshDistributor::PeriodicDofMap PeriodicDofMap;
StdMpi<PeriodicDofMap> stdMpi(pdb.mpiComm, true);
@@ -219,6 +237,80 @@ namespace AMDiS {
mpi::globalAdd(foundError);
TEST_EXIT(foundError == 0)("Error found on at least on rank!\n");
+
+
+ // === 3. check: On all edge and face periodic DOFs, at least on coordinate of ===
+ // === each periodic DOF pair must be equal (at least as long we consider ===
+ // === periodic boundaries only on rectangulars and boxes. ===
+
+ RankToCoords sendCoords;
+ std::map<int, std::vector<BoundaryType> > rankToDofType;
+
+ for (InteriorBoundary::RankToBoundMap::iterator it = pdb.periodicBoundary.boundary.begin();
+ it != pdb.periodicBoundary.boundary.end(); ++it) {
+ if (it->first == pdb.mpiRank)
+ continue;
+
+ for (vector<AtomicBoundary>::iterator boundIt = it->second.begin();
+ boundIt != it->second.end(); ++boundIt) {
+ if (boundIt->rankObj.subObj == VERTEX)
+ continue;
+
+ DofContainer dofs;
+ boundIt->rankObj.el->getVertexDofs(pdb.feSpace, boundIt->rankObj, dofs);
+ boundIt->rankObj.el->getNonVertexDofs(pdb.feSpace, boundIt->rankObj, dofs);
+
+ for (unsigned int i = 0; i < dofs.size(); i++) {
+ WorldVector<double> c;
+ pdb.mesh->getDofIndexCoords(*(dofs[i]), pdb.feSpace, c);
+ sendCoords[it->first].push_back(c);
+ rankToDofType[it->first].push_back(boundIt->type);
+ }
+ }
+ }
+
+ // Each rank must receive exactly the same number of coordinates as it sends
+ // to another rank.
+ RankToCoords recvCoords;
+ for (RankToCoords::iterator it = sendCoords.begin();
+ it != sendCoords.end(); ++it)
+ recvCoords[it->first].resize(it->second.size());
+
+
+ StdMpi<CoordsVec> stdMpiCoords(pdb.mpiComm, true);
+ stdMpiCoords.send(sendCoords);
+ stdMpiCoords.recv(recvCoords);
+ stdMpiCoords.startCommunication<double>(MPI_DOUBLE);
+
+
+ for (RankToCoords::iterator it = sendCoords.begin();
+ it != sendCoords.end(); ++it) {
+ for (unsigned int i = 0; i < it->second.size(); i++) {
+ WorldVector<double> &c0 = it->second[i];
+ WorldVector<double> &c1 = stdMpiCoords.getRecvData(it->first)[i];
+
+
+ int nEqual = 0;
+ for (int j = 0; j < pdb.mesh->getDim(); j++)
+ if (c0[j] == c1[j])
+ nEqual++;
+
+ if ((rankToDofType[it->first][i] >= -3 && nEqual < 2) ||
+ (rankToDofType[it->first][i] < -3 && nEqual == 0)) {
+ MSG("[DBG] %d-ith periodic DOF in boundary between ranks %d <-> %d is not correct!\n",
+ i, pdb.mpiRank, it->first);
+ MSG("[DBG] Coords on rank %d: %f %f %f\n",
+ pdb.mpiRank, c0[0], c0[1], (pdb.mesh->getDim() == 3 ? c0[2] : 0.0));
+ MSG("[DBG] Coords on rank %d: %f %f %f\n",
+ it->first, c1[0], c1[1], (pdb.mesh->getDim() == 3 ? c1[2] : 0.0));
+
+ foundError = 1;
+ }
+ }
+ }
+
+ mpi::globalAdd(foundError);
+ TEST_EXIT(foundError == 0)("Wrond periodic coordinates found on at least on rank!\n");
}
@@ -235,9 +327,6 @@ namespace AMDiS {
return;
}
- /// Defines a mapping type from rank numbers to sets of coordinates.
- typedef std::map<int, std::vector<WorldVector<double> > > RankToCoords;
-
/// Defines a mapping type from rank numbers to sets of DOFs.
typedef std::map<int, DofContainer> RankToDofContainer;
@@ -316,7 +405,6 @@ namespace AMDiS {
// === So we can check if also the coordinates of the communicated DOFs are ===
// === the same on both corresponding ranks. ===
- typedef std::vector<WorldVector<double> > CoordsVec;
StdMpi<CoordsVec> stdMpi(pdb.mpiComm, true);
stdMpi.send(sendCoords);
stdMpi.recv(recvCoords);
diff --git a/AMDiS/src/parallel/ParallelDebug.h b/AMDiS/src/parallel/ParallelDebug.h
index 101d2f21..1dc7e315 100644
--- a/AMDiS/src/parallel/ParallelDebug.h
+++ b/AMDiS/src/parallel/ParallelDebug.h
@@ -32,6 +32,11 @@ namespace AMDiS {
protected:
typedef MeshDistributor::RankToDofContainer RankToDofContainer;
+ typedef std::vector<WorldVector<double> > CoordsVec;
+
+ /// Defines a mapping type from rank numbers to sets of coordinates.
+ typedef std::map<int, CoordsVec> RankToCoords;
+
public:
/** \brief
* Tests the interior and the periodic boundaries on all ranks if their order
diff --git a/AMDiS/src/parallel/PetscSolver.cc b/AMDiS/src/parallel/PetscSolver.cc
index 202b77aa..652863ad 100644
--- a/AMDiS/src/parallel/PetscSolver.cc
+++ b/AMDiS/src/parallel/PetscSolver.cc
@@ -67,6 +67,8 @@ namespace AMDiS {
std::vector<double> values;
cols.reserve(300);
values.reserve(300);
+
+ std::vector<int> globalCols;
// === Traverse all rows of the dof matrix and insert row wise the values ===
// === to the petsc matrix. ===
@@ -123,6 +125,7 @@ namespace AMDiS {
// The col dof index is not periodic, simple add entry.
cols.push_back(colIndex);
values.push_back(value(*icursor));
+ globalCols.push_back(globalColDof);
}
}
}
@@ -136,7 +139,7 @@ namespace AMDiS {
std::set<int>& perAsc = meshDistributor->getPerDofAssociations(globalRowDof);
double scalFactor = 1.0 / (perAsc.size() + 1.0);
-
+
for (unsigned int i = 0; i < values.size(); i++)
values[i] *= scalFactor;
@@ -150,13 +153,15 @@ namespace AMDiS {
perCols.reserve(300);
std::vector<double> perValues;
perValues.reserve(300);
+
for (unsigned int i = 0; i < cols.size(); i++) {
int tmp = (cols[i] - dispAddCol) / dispMult;
- if (meshDistributor->isPeriodicDof(tmp, *perIt))
+ if (meshDistributor->isPeriodicDof(tmp, *perIt)) {
perCols.push_back((meshDistributor->getPeriodicMapping(*perIt, tmp) * dispMult) + dispAddCol);
- else
+ } else {
perCols.push_back(cols[i]);
+ }
perValues.push_back(values[i]);
}
--
GitLab