Commit ed1ab023 authored by Thomas Witkowski's avatar Thomas Witkowski

Bugfix for assembling with 8 or more processors.

parent 4d620174
......@@ -44,7 +44,7 @@ available_tags=" CXX F77"
# ### BEGIN LIBTOOL CONFIG
# Libtool was configured on host deimos101:
# Libtool was configured on host p2d087:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
......@@ -6760,7 +6760,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac`
# End:
# ### BEGIN LIBTOOL TAG CONFIG: CXX
# Libtool was configured on host deimos101:
# Libtool was configured on host p2d087:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
......@@ -7065,7 +7065,7 @@ include_expsyms=""
# ### BEGIN LIBTOOL TAG CONFIG: F77
# Libtool was configured on host deimos101:
# Libtool was configured on host p2d087:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
......
......@@ -42,7 +42,7 @@ namespace AMDiS {
bool ElementDofIterator::next()
{
// First iteratore over the dofs of one element (vertex, edge, face).
// First iterate over the dofs of one element (vertex, edge, face).
dofPos++;
if (dofPos >= nDofs) {
......
......@@ -2982,24 +2982,17 @@ namespace AMDiS {
if (n < 1)
return;
const Element *el;
const DegreeOfFreedom *cd;
DegreeOfFreedom pd[19], cdi;
int node0, n0, i, typ, lr_set;
const DOFAdmin *admin;
el = list->getElement(0);
typ = list->getType(0);
admin = drv->getFESpace()->getAdmin();
DegreeOfFreedom pd[20], cdi;
const Element *el = list->getElement(0);
int typ = list->getType(0);
const DOFAdmin *admin = drv->getFESpace()->getAdmin();
basFct->getLocalIndices(el, admin, pd);
/****************************************************************************/
/* values on child[0] */
/****************************************************************************/
cd = basFct->getLocalIndices(el->getChild(0), admin, NULL);
const DegreeOfFreedom *cd = basFct->getLocalIndices(el->getChild(0), admin, NULL);
(*drv)[pd[0]] +=
(0.0625*(-(*drv)[cd[3]] + (*drv)[cd[12]] + (*drv)[cd[14]]
......@@ -3087,15 +3080,15 @@ namespace AMDiS {
/* adjust neighbour values */
/****************************************************************************/
node0 = drv->getFESpace()->getMesh()->getNode(FACE);
n0 = admin->getNumberOfPreDOFs(FACE);
int node0 = drv->getFESpace()->getMesh()->getNode(FACE);
int n0 = admin->getNumberOfPreDOFs(FACE);
for (i = 1; i < n; i++) {
for (int i = 1; i < n; i++) {
el = list->getElement(i);
typ = list->getType(i);
basFct->getLocalIndices(el, admin, pd);
lr_set = 0;
int lr_set = 0;
if (list->getNeighbourElement(i, 0) && list->getNeighbourNr(i, 0) < i)
lr_set = 1;
......
......@@ -320,14 +320,22 @@ namespace AMDiS {
if (isRankDof[*cursor]) {
r -= rstart * nComponents;
TEST_EXIT_DBG(r >= 0 && r < nRankRows)("Should not happen!\n");
#if (DEBUG != 0)
if (r < 0 || r >= nRankRows) {
std::cout << "ERROR in rank: " << mpiRank << std::endl;
std::cout << " Wrong r = " << r << " " << *cursor << " "
<< mapLocalGlobalDOFs[*cursor] << " "
<< nRankRows << std::endl;
ERROR_EXIT("Should not happen!\n");
}
#endif
for (icursor_type icursor = begin<nz>(cursor),
icend = end<nz>(cursor); icursor != icend; ++icursor) {
if (value(*icursor) != 0.0) {
int c = mapLocalGlobalDOFs[col(*icursor)] * nComponents + j;
if (c >= rstart * nComponents &&
c < rstart * nComponents + nRankRows)
d_nnz[r]++;
......@@ -409,10 +417,9 @@ namespace AMDiS {
it != recvDofs.end(); ++it) {
int nSend = sendMatrixEntry[it->first].size();
if (nSend > 0) {
if (nSend > 0)
request[requestCounter++] =
mpiComm.Isend(sendBuffers[i], nSend * 2, MPI_INT, it->first, 0);
}
i++;
}
......@@ -588,7 +595,8 @@ namespace AMDiS {
KSPCreate(PETSC_COMM_WORLD, &solver);
KSPSetOperators(solver, petscMatrix, petscMatrix, SAME_NONZERO_PATTERN);
KSPSetTolerances(solver, 1.e-8, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT);
KSPSetTolerances(solver, 0.0, 1e-8, PETSC_DEFAULT, PETSC_DEFAULT);
KSPSetType(solver, KSPBCGS);
KSPMonitorSet(solver, myKSPMonitor, PETSC_NULL, 0);
KSPSetFromOptions(solver);
......@@ -912,7 +920,7 @@ namespace AMDiS {
DofContainer rankAllDofs;
DofToRank boundaryDofs;
createDOFMemberInfo(partitionDOFs, rankDOFs, rankAllDofs, boundaryDofs);
createDOFMemberInfo(partitionDOFs, rankDOFs, rankAllDofs, boundaryDofs, vertexDof);
nRankDOFs = rankDOFs.size();
nOverallDOFs = partitionDOFs.size();
......@@ -966,7 +974,7 @@ namespace AMDiS {
// Stores to each rank a map from DOF pointers (which are all owned by the rank
// and lie on an interior boundary) to new global DOF indices.
std::map<int, std::map<const DegreeOfFreedom*, DegreeOfFreedom> > sendNewDofs;
std::map<int, DofIndexMap > sendNewDofs;
// Maps to each rank the number of new DOF indices this rank will receive from.
// All these DOFs are at an interior boundary on this rank, but are owned by
......@@ -1011,18 +1019,14 @@ namespace AMDiS {
int requestCounter = 0;
i = 0;
for (std::map<int, std::map<const DegreeOfFreedom*, DegreeOfFreedom> >::iterator
sendIt = sendNewDofs.begin();
sendIt != sendNewDofs.end();
++sendIt, i++) {
for (std::map<int, DofIndexMap >::iterator sendIt = sendNewDofs.begin();
sendIt != sendNewDofs.end(); ++sendIt, i++) {
int nSendDofs = sendIt->second.size() * 2;
sendBuffers[i] = new int[nSendDofs];
int c = 0;
for (std::map<const DegreeOfFreedom*, DegreeOfFreedom>::iterator
dofIt = sendIt->second.begin();
dofIt != sendIt->second.end();
++dofIt) {
for (DofIndexMap::iterator dofIt = sendIt->second.begin();
dofIt != sendIt->second.end(); ++dofIt) {
sendBuffers[i][c++] = *(dofIt->first);
sendBuffers[i][c++] = dofIt->second;
......@@ -1035,8 +1039,7 @@ namespace AMDiS {
i = 0;
for (std::map<int, int>::iterator recvIt = recvNewDofs.begin();
recvIt != recvNewDofs.end();
++recvIt, i++) {
recvIt != recvNewDofs.end(); ++recvIt, i++) {
int nRecvDOFs = recvIt->second * 2;
recvBuffers[i] = new int[nRecvDOFs];
......@@ -1142,9 +1145,12 @@ namespace AMDiS {
// === Traverse on interior boundaries and move all not ranked owned DOFs from ===
// === rankDOFs to boundaryDOFs. ===
RankToDofContainer oldSendDofs = sendDofs;
RankToDofContainer oldRecvDofs = recvDofs;
sendDofs.clear();
recvDofs.clear();
for (RankToBoundMap::iterator it = myIntBoundary.boundary.begin();
it != myIntBoundary.boundary.end(); ++it) {
......@@ -1154,6 +1160,11 @@ namespace AMDiS {
DofContainer dofs;
DofContainer &dofsToSend = sendDofs[it->first];
for (DofContainer::iterator iit = oldSendDofs[it->first].begin();
iit != oldSendDofs[it->first].end(); ++iit)
if (vertexDof[*iit])
dofsToSend.push_back(*iit);
switch (boundIt->rankObject.ithObjAtBoundary) {
case 0:
dofs.push_back(boundIt->rankObject.el->getDOF(1));
......@@ -1171,21 +1182,27 @@ namespace AMDiS {
ERROR_EXIT("Should never happen!\n");
}
/*
for (DofContainer::iterator dofIt = dofs.begin(); dofIt != dofs.end(); ++dofIt)
if (find(dofsToSend.begin(), dofsToSend.end(), *dofIt) == dofsToSend.end())
dofsToSend.push_back(*dofIt);
if (find(oldSendDofs[it->first].begin(), oldSendDofs[it->first].end(), *dofIt) != oldSendDofs[it->first].end())
dofsToSend.push_back(*dofIt);
*/
dofs.clear();
addAllVertexDOFs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
dofs);
addAllEdgeDOFs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
dofs);
for (int i = 0; i < static_cast<int>(dofs.size()); i++)
dofsToSend.push_back(dofs[i]);
for (int i = 0; i < static_cast<int>(dofs.size()); i++) {
TEST_EXIT_DBG(find(dofsToSend.begin(), dofsToSend.end(), dofs[i]) == dofsToSend.end())
("Should not happen!\n");
dofsToSend.push_back(dofs[i]);
}
}
}
for (RankToBoundMap::iterator it = otherIntBoundary.boundary.begin();
it != otherIntBoundary.boundary.end(); ++it) {
......@@ -1196,7 +1213,26 @@ namespace AMDiS {
DofContainer dofs;
DofContainer &dofsToRecv = recvDofs[it->first];
switch (boundIt->rankObject.ithObjAtBoundary) {
/*
for (DofContainer::iterator iit = dofsToRecv.begin(); iit != dofsToRecv.end(); ++iit) {
DofContainer::iterator eraseIt = find(rankDOFs.begin(), rankDOFs.end(), *iit);
if (eraseIt != rankDOFs.end())
rankDOFs.erase(eraseIt);
}
*/
for (DofContainer::iterator iit = oldRecvDofs[it->first].begin();
iit != oldRecvDofs[it->first].end(); ++iit)
if (vertexDof[*iit]) {
dofsToRecv.push_back(*iit);
DofContainer::iterator eraseIt = find(rankDOFs.begin(), rankDOFs.end(), *iit);
if (eraseIt != rankDOFs.end())
rankDOFs.erase(eraseIt);
}
/* switch (boundIt->rankObject.ithObjAtBoundary) {
case 0:
if (boundIt->neighbourObject.ithObjAtBoundary == 0) {
dofs.push_back(boundIt->rankObject.el->getDOF(2));
......@@ -1225,12 +1261,20 @@ namespace AMDiS {
for (DofContainer::iterator dofIt = dofs.begin(); dofIt != dofs.end(); ++dofIt) {
DofContainer::iterator eraseIt = find(rankDOFs.begin(), rankDOFs.end(), *dofIt);
if (eraseIt != rankDOFs.end())
rankDOFs.erase(eraseIt);
if (eraseIt != rankDOFs.end()) {
if (mpiRank == 4)
std::cout << "ERASE a: " << **eraseIt << std::endl;
rankDOFs.erase(eraseIt);
}
if (find(dofsToRecv.begin(), dofsToRecv.end(), *dofIt) == dofsToRecv.end())
dofsToRecv.push_back(*dofIt);
if (find(oldRecvDofs[it->first].begin(), oldRecvDofs[it->first].end(), *dofIt) != oldRecvDofs[it->first].end())
dofsToRecv.push_back(*dofIt);
}
*/
dofs.clear();
addAllEdgeDOFs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
dofs);
......@@ -1238,12 +1282,15 @@ namespace AMDiS {
dofs);
for (int i = static_cast<int>(dofs.size()) - 1; i >= 0; i--) {
TEST_EXIT_DBG(find(rankDOFs.begin(), rankDOFs.end(), dofs[i]) != rankDOFs.end())
("Should never happen!\n");
// TEST_EXIT_DBG(find(rankDOFs.begin(), rankDOFs.end(), dofs[i]) != rankDOFs.end())
// ("Should never happen!\n");
TEST_EXIT_DBG(find(dofsToRecv.begin(), dofsToRecv.end(), dofs[i]) == dofsToRecv.end())
("Should not happen!\n");
DofContainer::iterator eraseIt = find(rankDOFs.begin(), rankDOFs.end(), dofs[i]);
if (eraseIt != rankDOFs.end())
rankDOFs.erase(eraseIt);
rankDOFs.erase(eraseIt);
dofsToRecv.push_back(dofs[i]);
}
......@@ -1307,8 +1354,9 @@ namespace AMDiS {
sendBuffers[i] = new int[nSendDofs];
int c = 0;
for (DofContainer::iterator dofIt = sendIt->second.begin();
dofIt != sendIt->second.end(); ++dofIt)
dofIt != sendIt->second.end(); ++dofIt) {
sendBuffers[i][c++] = rankDofsNewGlobalIndex[*dofIt];
}
request[requestCounter++] =
mpiComm.Isend(sendBuffers[i], nSendDofs, MPI_INT, sendIt->first, 0);
......@@ -1336,7 +1384,7 @@ namespace AMDiS {
for (DofContainer::iterator dofIt = recvIt->second.begin();
dofIt != recvIt->second.end(); ++dofIt) {
rankDofsNewGlobalIndex[*dofIt] = recvBuffers[i][j];
rankDofsNewGlobalIndex[*dofIt] = recvBuffers[i][j];
isRankDof[rankDofsNewLocalIndex[*dofIt]] = false;
j++;
}
......@@ -1455,12 +1503,14 @@ namespace AMDiS {
void ParallelDomainBase::createDOFMemberInfo(DofToPartitions& partitionDofs,
DofContainer& rankOwnedDofs,
DofContainer& rankAllDofs,
DofToRank& boundaryDofs)
DofToRank& boundaryDofs,
DofToBool& vertexDof)
{
partitionDofs.clear();
rankOwnedDofs.clear();
rankAllDofs.clear();
boundaryDofs.clear();
vertexDof.clear();
// === Determine to each dof the set of partitions the dof belongs to. ===
......@@ -1474,6 +1524,11 @@ namespace AMDiS {
do {
// Determine to each dof the partition(s) it corresponds to.
partitionDofs[elDofIt.getDofPtr()].insert(partitionVec[element->getIndex()]);
if (elDofIt.getCurrentPos() == 0)
vertexDof[elDofIt.getDofPtr()] = true;
else
vertexDof[elDofIt.getDofPtr()] = false;
} while(elDofIt.next());
elInfo = stack.traverseNext(elInfo);
......@@ -1486,11 +1541,13 @@ namespace AMDiS {
for (DofToPartitions::iterator it = partitionDofs.begin();
it != partitionDofs.end(); ++it) {
bool isInRank = false;
// iterate over all partition the current DOF is part of.
for (std::set<int>::iterator itpart1 = it->second.begin();
itpart1 != it->second.end(); ++itpart1) {
if (*itpart1 == mpiRank) {
if (*itpart1 == mpiRank) {
rankAllDofs.push_back(it->first);
if (it->second.size() == 1) {
......@@ -1516,10 +1573,14 @@ namespace AMDiS {
boundaryDofs[it->first] = highestRank;
}
isInRank = true;
break;
}
}
if (!isInRank)
vertexDof.erase(it->first);
}
sort(rankAllDofs.begin(), rankAllDofs.end(), cmpDofsByValue);
......@@ -1657,7 +1718,14 @@ namespace AMDiS {
for (std::vector<const DegreeOfFreedom*>::iterator dofIt = it->second.begin();
dofIt != it->second.end(); ++dofIt) {
bool b = mesh->getDofIndexCoords(*dofIt, feSpace, coords);
TEST_EXIT(b)("Cannot find DOF in mesh!\n");
if (!b) {
std::cout << "ERROR (send): Cannot find DOF in mesh!" << std::endl;
std::cout << " mpi rank = " << mpiRank
<< " send to rank = " << it->first
<< " dof = " << *dofIt << " = " << **dofIt << std::endl;
ERROR_EXIT("Should not happen!\n");
}
sendCoords[it->first].push_back(coords);
}
......@@ -1668,7 +1736,15 @@ namespace AMDiS {
for (std::vector<const DegreeOfFreedom*>::iterator dofIt = it->second.begin();
dofIt != it->second.end(); ++dofIt) {
bool b = mesh->getDofIndexCoords(*dofIt, feSpace, coords);
TEST_EXIT(b)("Cannot find DOF in mesh!\n");
if (!b) {
std::cout << "ERROR (recv): Cannot find DOF in mesh!" << std::endl;
std::cout << " mpi rank = " << mpiRank
<< " send to rank = " << it->first
<< " dof = " << *dofIt << " = " << **dofIt << std::endl;
ERROR_EXIT("Should not happen!\n");
}
recvCoords[it->first].push_back(coords);
}
}
......
......@@ -64,8 +64,11 @@ namespace AMDiS {
/// Defines a mapping type from DOF indices to DOF indices.
typedef std::map<DegreeOfFreedom, DegreeOfFreedom> DofMapping;
/// Defines a mapping type from DOFs to boolean values.
typedef std::map<const DegreeOfFreedom*, bool> DofToBool;
/// Defines a mapping type from DOF indices to boolean values.
typedef std::map<DegreeOfFreedom, bool> DofToBool;
typedef std::map<DegreeOfFreedom, bool> DofIndexToBool;
/// Defines a mapping type from rank numbers to sets of coordinates.
typedef std::map<int, std::vector<WorldVector<double> > > RankToCoords;
......@@ -244,7 +247,8 @@ namespace AMDiS {
void createDOFMemberInfo(DofToPartitions& partitionDofs,
DofContainer& rankOwnedDofs,
DofContainer& rankAllDofs,
DofToRank& boundaryDofs);
DofToRank& boundaryDofs,
DofToBool& vertexDof);
void setDofMatrix(DOFMatrix* mat, int dispMult = 1,
int dispAddRow = 0, int dispAddCol = 0);
......@@ -403,7 +407,9 @@ namespace AMDiS {
* owned by the rank. Otherwise, its an interior boundary DOF that is owned by
* another rank.
*/
DofToBool isRankDof;
DofIndexToBool isRankDof;
DofToBool vertexDof;
int rstart;
......
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