Commit 8dd34fa5 authored by Thomas Witkowski's avatar Thomas Witkowski

Repartitioning for mixed finite elements works now.

parent 5367e947
......@@ -108,8 +108,8 @@ namespace AMDiS {
{
FUNCNAME("DOFAdmin::freeDofIndex()");
TEST_EXIT_DBG(usedCount > 0)("no dofs in use\n");
TEST_EXIT_DBG(dof >= 0 && dof < size)("invalid dof index %d\n",dof);
TEST_EXIT_DBG(usedCount > 0)("No DOFs in use!\n");
TEST_EXIT_DBG(dof >= 0 && dof < size)("Invalid DOF index %d!\n", dof);
std::list<DOFIndexedBase*>::iterator di;
std::list<DOFIndexedBase*>::iterator end = dofIndexedList.end();
......@@ -277,8 +277,9 @@ namespace AMDiS {
for (it.reset(); !it.end(); ++it)
newDofIndex[it.getDOFIndex()] = 1;
// create a MONOTONE compress
int n = 0, last = 0;
for (int i = 0; i < size; i++) { /* create a MONOTONE compress */
for (int i = 0; i < size; i++) {
if (newDofIndex[i] == 1) {
newDofIndex[i] = n++;
last = i;
......
......@@ -631,12 +631,26 @@ namespace AMDiS {
void Element::getAllDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr)
bool baseDofPtr,
vector<GeoIndex>* dofGeoIndex)
{
getNodeDofs(feSpace, bound, dofs, baseDofPtr);
if (dofGeoIndex != NULL) {
// In the case dofGeoIndex should be filled, set all node DOFs to be
// vertex DOFs.
dofGeoIndex->resize(dofs.size());
for (unsigned int i = 0; i < dofs.size(); i++)
(*dofGeoIndex)[i] = VERTEX;
}
if (feSpace->getBasisFcts()->getDegree() > 1)
getHigherOrderDofs(feSpace, bound, dofs, baseDofPtr);
getHigherOrderDofs(feSpace, bound, dofs, baseDofPtr, dofGeoIndex);
if (dofGeoIndex) {
TEST_EXIT_DBG(dofs.size() == dofGeoIndex->size())
("Arrays do not fit together: %d %d\n", dofs.size(), dofGeoIndex->size());
}
}
}
......@@ -427,18 +427,24 @@ namespace AMDiS {
* dof* [\ref dof] of the element is inserted. If
* false, &(dof[.][n0]) is put to the result vector,
* with n0 beging the number of predofs.
* \param[out] dofGeoIndex Optional, the function can store to each DOF in
* the DofContainer dofs the geometric index, thus
* identifing the DOF to be a vertex, edge, face or
* center DOF.
*/
virtual void getHigherOrderDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr = false) const = 0;
bool baseDofPtr = false,
vector<GeoIndex>* dofGeoIndex = NULL) const = 0;
/// Combines \ref getNodeDofs and \ref getHigherOrderDofs to one function.
/// See parameter description there.
void getAllDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr = false);
bool baseDofPtr = false,
vector<GeoIndex>* dofGeoIndex = NULL);
/** \} */
......
......@@ -181,7 +181,7 @@ namespace AMDiS {
}
void getHigherOrderDofs(const FiniteElemSpace*, BoundaryObject,
DofContainer&, bool) const
DofContainer&, bool, vector<GeoIndex>*) const
{
FUNCNAME("Line::getHigherOrderDofs()");
ERROR_EXIT("Not yet implemented!\n");
......
......@@ -376,10 +376,13 @@ namespace AMDiS {
void Tetrahedron::getHigherOrderDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr) const
bool baseDofPtr,
vector<GeoIndex>* dofGeoIndex) const
{
FUNCNAME("Tetrahedron::getHigherOrderDofs()");
TEST_EXIT(dofGeoIndex != NULL)("Not yet implemented!\n");
switch (bound.subObj) {
case VERTEX:
return;
......@@ -403,14 +406,18 @@ namespace AMDiS {
if (bound.reverseMode) {
if (nextBound1.ithObj >= 0)
child[1]->getHigherOrderDofs(feSpace, nextBound1, dofs, baseDofPtr);
child[1]->getHigherOrderDofs(feSpace, nextBound1, dofs,
baseDofPtr, dofGeoIndex);
if (nextBound0.ithObj >= 0)
child[0]->getHigherOrderDofs(feSpace, nextBound0, dofs, baseDofPtr);
child[0]->getHigherOrderDofs(feSpace, nextBound0, dofs,
baseDofPtr, dofGeoIndex);
} else {
if (nextBound0.ithObj >= 0)
child[0]->getHigherOrderDofs(feSpace, nextBound0, dofs, baseDofPtr);
child[0]->getHigherOrderDofs(feSpace, nextBound0, dofs,
baseDofPtr, dofGeoIndex);
if (nextBound1.ithObj >= 0)
child[1]->getHigherOrderDofs(feSpace, nextBound1, dofs, baseDofPtr);
child[1]->getHigherOrderDofs(feSpace, nextBound1, dofs,
baseDofPtr, dofGeoIndex);
}
} else {
// Either the edge is not contained in further refined children, or
......@@ -451,12 +458,14 @@ namespace AMDiS {
if (childFace0 != -1) {
nextBound.ithObj = childFace0;
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs, baseDofPtr);
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
}
if (childFace1 != -1) {
nextBound.ithObj = childFace1;
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs, baseDofPtr);
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
}
} else {
ElementDofIterator elDofIter(feSpace, true);
......
......@@ -154,7 +154,8 @@ namespace AMDiS {
void getHigherOrderDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr = false) const;
bool baseDofPtr = false,
vector<GeoIndex>* dofGeoIndex = NULL) const;
void prepareNextBound(BoundaryObject &bound, int ithChild) const;
......
......@@ -115,16 +115,20 @@ namespace AMDiS {
if (bound.reverseMode) {
nextBound.ithObj = 1;
child[1]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 0;
child[1]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
} else {
nextBound.ithObj = 0;
child[1]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 1;
child[1]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
}
}
}
......@@ -136,16 +140,20 @@ namespace AMDiS {
if (bound.reverseMode) {
nextBound.ithObj = 1;
child[0]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 0;
child[0]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
} else {
nextBound.ithObj = 0;
child[0]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getFirstChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 1;
child[0]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getSecondChild()->getNodeDofs(feSpace, nextBound, dofs,
baseDofPtr);
}
}
}
......@@ -156,16 +164,16 @@ namespace AMDiS {
if (bound.reverseMode) {
nextBound.ithObj = 1;
child[1]->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getNodeDofs(feSpace, nextBound, dofs, baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 0;
child[0]->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getNodeDofs(feSpace, nextBound, dofs, baseDofPtr);
} else {
nextBound.ithObj = 0;
child[0]->getNodeDofs(feSpace, nextBound, dofs);
child[0]->getNodeDofs(feSpace, nextBound, dofs, baseDofPtr);
dofs.push_back(&(elDofs[2][n0]));
nextBound.ithObj = 1;
child[1]->getNodeDofs(feSpace, nextBound, dofs);
child[1]->getNodeDofs(feSpace, nextBound, dofs, baseDofPtr);
}
}
break;
......@@ -178,7 +186,8 @@ namespace AMDiS {
void Triangle::getHigherOrderDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr) const
bool baseDofPtr,
vector<GeoIndex>* dofGeoIndex) const
{
FUNCNAME("Triange::getHigherOrderDofs()");
......@@ -194,7 +203,8 @@ namespace AMDiS {
case 0:
if (child[1]) {
nextBound.ithObj = 2;
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
} else {
addThisEdge = true;
}
......@@ -203,7 +213,8 @@ namespace AMDiS {
case 1:
if (child[0]) {
nextBound.ithObj = 2;
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
} else {
addThisEdge = true;
}
......@@ -213,14 +224,18 @@ namespace AMDiS {
if (child[0]) {
if (bound.reverseMode) {
nextBound.ithObj = 1;
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
nextBound.ithObj = 0;
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
} else {
nextBound.ithObj = 0;
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[0]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
nextBound.ithObj = 1;
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs);
child[1]->getHigherOrderDofs(feSpace, nextBound, dofs,
baseDofPtr, dofGeoIndex);
}
} else {
addThisEdge = true;
......@@ -238,24 +253,31 @@ namespace AMDiS {
if (baseDofPtr) {
do {
if (elDofIter.getCurrentPos() == 1 &&
if (elDofIter.getPosIndex() == EDGE &&
elDofIter.getCurrentElementPos() == bound.ithObj)
addDofs.push_back(elDofIter.getBaseDof());
} while (elDofIter.nextStrict());
} else {
do {
if (elDofIter.getCurrentPos() == 1 &&
if (elDofIter.getPosIndex() == EDGE &&
elDofIter.getCurrentElementPos() == bound.ithObj)
addDofs.push_back(elDofIter.getDofPtr());
addDofs.push_back(elDofIter.getDofPtr());
} while (elDofIter.next());
}
if (bound.reverseMode)
for (int i = addDofs.size() - 1; i >= 0; i--)
if (bound.reverseMode) {
for (int i = addDofs.size() - 1; i >= 0; i--) {
dofs.push_back(addDofs[i]);
else
for (unsigned int i = 0; i < addDofs.size(); i++)
dofs.push_back(addDofs[i]);
if (dofGeoIndex != NULL)
dofGeoIndex->push_back(EDGE);
}
} else {
for (unsigned int i = 0; i < addDofs.size(); i++) {
dofs.push_back(addDofs[i]);
if (dofGeoIndex != NULL)
dofGeoIndex->push_back(EDGE);
}
}
}
}
......
......@@ -202,7 +202,8 @@ namespace AMDiS {
void getHigherOrderDofs(const FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool baseDofPtr = false) const;
bool baseDofPtr = false,
vector<GeoIndex>* dofGeoIndex = NULL) const;
protected:
/// vertexOfEdge[i][j] is the local number of the j-th vertex of the i-th
......
......@@ -1086,16 +1086,11 @@ namespace AMDiS {
{
FUNCNAME("MeshDistributor::repartitionMesh()");
TEST_EXIT(mesh->getNumberOfDOFAdmin() == 1)
("Only meshes with one DOFAdmin are supported!\n");
// === First we check if the rank with the maximum number of DOFs has at ===
// === least 20% more DOFs than the rank with the minimum number of DOFs. ===
// === In this case, the mesh will be repartition. ===
int repartitioning = 0;
vector<int> nDofsInRank(mpiSize);
int nDofs = mesh->getDofAdmin(0).getUsedDofs();
mpiComm.Gather(&nDofs, 1, MPI_INT, &(nDofsInRank[0]), 1, MPI_INT, 0);
......@@ -1121,7 +1116,6 @@ namespace AMDiS {
mpiComm.Bcast(&repartitioning, 1, MPI_INT, 0);
}
if (repartitioning == 0)
return;
......@@ -1199,6 +1193,8 @@ namespace AMDiS {
// === Add new macro elements to mesh. ===
TEST_EXIT(mesh->getGeo(FACE) || mesh->getGeo(CENTER))
("Mh, dass macht dann doch noch etwas Arbeit fr den Thomas!\n");
for (std::set<MacroElement*>::iterator it = newMacroEl.begin();
it != newMacroEl.end(); ++it) {
......@@ -1244,7 +1240,6 @@ namespace AMDiS {
}
}
StdMpi<MeshCodeVec> stdMpi(mpiComm, true);
stdMpi.send(sendCodes);
for (map<int, vector<int> >::iterator it = partitioner->getRecvElements().begin();
......@@ -1252,7 +1247,6 @@ namespace AMDiS {
stdMpi.recv(it->first);
stdMpi.startCommunication();
if (interchangeVectors.size() == 0)
WARNING("There are no interchange vectors defined!\n");
......@@ -1368,9 +1362,9 @@ namespace AMDiS {
MSG("Debug mode tests finished!\n");
#endif
// === In 3D we have to make some test, if the resulting mesh is valid. If ===
// === it is not valid, there is no possiblity yet to fix this problem, just ===
// === exit with an error message. ===
// === In 3D we have to make some test, if the resulting mesh is valid. ===
// === If it is not valid, there is no possiblity yet to fix this ===
// === problem, just exit with an error message. ===
check3dValidMesh();
......@@ -1516,9 +1510,11 @@ namespace AMDiS {
AtomicBoundary& b = rankIntBoundary.getNewAtomic(it2->first);
b = bound;
if (geoIndex == EDGE)
b.neighObj.reverseMode = elObjects.getEdgeReverseMode(rankBoundEl, it2->second);
b.neighObj.reverseMode =
elObjects.getEdgeReverseMode(rankBoundEl, it2->second);
if (geoIndex == FACE)
b.neighObj.reverseMode = elObjects.getFaceReverseMode(rankBoundEl, it2->second);
b.neighObj.reverseMode =
elObjects.getFaceReverseMode(rankBoundEl, it2->second);
}
} else {
......@@ -1538,9 +1534,11 @@ namespace AMDiS {
AtomicBoundary& b = otherIntBoundary.getNewAtomic(owner);
b = bound;
if (geoIndex == EDGE)
b.rankObj.reverseMode = elObjects.getEdgeReverseMode(rankBoundEl, ownerBoundEl);
b.rankObj.reverseMode =
elObjects.getEdgeReverseMode(rankBoundEl, ownerBoundEl);
if (geoIndex == FACE)
b.rankObj.reverseMode = elObjects.getFaceReverseMode(rankBoundEl, ownerBoundEl);
b.rankObj.reverseMode =
elObjects.getFaceReverseMode(rankBoundEl, ownerBoundEl);
}
}
}
......
......@@ -131,7 +131,7 @@ namespace AMDiS {
if (macrosProcessed.count(elIt->elIndex) == 1) {
TEST_EXIT_DBG(macroIndexMap.count(elIt->elIndex))
("Should not happen!\n");
Element *el0 = (*it)->getElement();
Element *el1 = macroIndexMap[elIt->elIndex]->getElement();
......@@ -141,17 +141,24 @@ namespace AMDiS {
BoundaryObject b1(el1, 0, EDGE, elIt->ithObject, false);
DofContainer dofs0, dofs1;
el0->getAllDofs(feSpace, b0, dofs0);
el1->getAllDofs(feSpace, b1, dofs1);
vector<GeoIndex> dofGeoIndex0, dofGeoIndex1;
el0->getAllDofs(feSpace, b0, dofs0, true, &dofGeoIndex0);
el1->getAllDofs(feSpace, b1, dofs1, true, &dofGeoIndex1);
#if (DEBUG != 0)
debug::testDofsByCoords(feSpace, dofs0, dofs1);
if (feSpaces.size() == 1)
debug::testDofsByCoords(feSpace, dofs0, dofs1);
else
TEST_EXIT_DBG(dofs0.size() == dofs1.size())("Should not happen!\n");
#endif
for (unsigned int i = 0; i < dofs0.size(); i++) {
mapDelDofs[dofs0[i]] = dofs1[i];
dofPosIndex[dofs0[i]] = EDGE;
dofPosIndex[dofs0[i]] = dofGeoIndex0[i];
TEST_EXIT_DBG(dofGeoIndex0[i] == dofGeoIndex1[i])
("Should not happen: %d %d\n",
dofGeoIndex0[i], dofGeoIndex1[i]);
}
break;
......@@ -182,10 +189,9 @@ namespace AMDiS {
BoundaryObject b0(el0, 0, FACE, i, reverseMode);
BoundaryObject b1(el1, 0, FACE, elIt->ithObject, false);
DofContainer dofs0, dofs1;
el0->getAllDofs(feSpace, b0, dofs0);
el1->getAllDofs(feSpace, b1, dofs1);
DofContainer dofs0, dofs1;
el0->getAllDofs(feSpace, b0, dofs0, true);
el1->getAllDofs(feSpace, b1, dofs1, true);
#if (DEBUG != 0)
debug::testDofsByCoords(feSpace, dofs0, dofs1);
......@@ -240,7 +246,8 @@ namespace AMDiS {
for (map<const DegreeOfFreedom*, const DegreeOfFreedom*>::iterator it = mapDelDofs.begin();
it != mapDelDofs.end(); ++it)
mesh->freeDof(const_cast<DegreeOfFreedom*>(it->first), VERTEX);
mesh->freeDof(const_cast<DegreeOfFreedom*>(it->first),
dofPosIndex[it->first]);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment