// // Software License for AMDiS // // Copyright (c) 2010 Dresden University of Technology // All rights reserved. // Authors: Simon Vey, Thomas Witkowski et al. // // This file is part of AMDiS // // See also license.opensource.txt in the distribution. #include "parallel/InteriorBoundary.h" #include "parallel/ElementObjectDatabase.h" #include "FiniteElemSpace.h" #include "BasisFunction.h" #include "Serializer.h" #include "VertexVector.h" namespace AMDiS { using namespace std; void InteriorBoundary::create(MPI::Intracomm &mpiComm, ElementObjectDatabase &elObjDb) { FUNCNAME("InteriorBoundary::clear()"); own.clear(); other.clear(); periodic.clear(); Mesh *mesh = elObjDb.getMesh(); TEST_EXIT_DBG(mesh)("Should not happen!\n"); int mpiRank = mpiComm.Get_rank(); // === Create interior boundary data structure. === for (int geoPos = 0; geoPos < mesh->getDim(); geoPos++) { GeoIndex geoIndex = INDEX_OF_DIM(geoPos, mesh->getDim()); while (elObjDb.iterate(geoIndex)) { map& objData = elObjDb.getIterateData(); if (!(objData.count(mpiRank) && objData.size() > 1)) continue; int owner = elObjDb.getIterateOwner(); ElementObjectData& rankBoundEl = objData[mpiRank]; AtomicBoundary bound; bound.maxLevel = elObjDb.getIterateMaxLevel(); bound.rankObj.el = elObjDb.getElementPtr(rankBoundEl.elIndex); bound.rankObj.elIndex = rankBoundEl.elIndex; bound.rankObj.elType = elObjDb.getElementType(rankBoundEl.elIndex); bound.rankObj.subObj = geoIndex; bound.rankObj.ithObj = rankBoundEl.ithObject; if (geoIndex == FACE) { for (int edgeNo = 0; edgeNo < 3; edgeNo++) { int edgeOfFace = bound.rankObj.el->getEdgeOfFace(bound.rankObj.ithObj, edgeNo); bound.rankObj.excludedSubstructures.push_back(make_pair(EDGE, edgeOfFace)); } } if (owner == mpiRank) { for (map::iterator it2 = objData.begin(); it2 != objData.end(); ++it2) { if (it2->first == mpiRank) continue; bound.neighObj.el = elObjDb.getElementPtr(it2->second.elIndex); bound.neighObj.elIndex = it2->second.elIndex; bound.neighObj.elType = elObjDb.getElementType(it2->second.elIndex); bound.neighObj.subObj = geoIndex; bound.neighObj.ithObj = it2->second.ithObject; bound.type = INTERIOR; AtomicBoundary& b = getNewOwn(it2->first); b = bound; if (geoIndex == EDGE) b.neighObj.reverseMode = elObjDb.getEdgeReverseMode(rankBoundEl, it2->second); if (geoIndex == FACE) b.neighObj.reverseMode = elObjDb.getFaceReverseMode(rankBoundEl, it2->second); } } else { TEST_EXIT_DBG(objData.count(owner) == 1) ("Should not happen!\n"); ElementObjectData& ownerBoundEl = objData[owner]; bound.neighObj.el = elObjDb.getElementPtr(ownerBoundEl.elIndex); bound.neighObj.elIndex = ownerBoundEl.elIndex; bound.neighObj.elType = -1; bound.neighObj.subObj = geoIndex; bound.neighObj.ithObj = ownerBoundEl.ithObject; bound.type = INTERIOR; AtomicBoundary& b = getNewOther(owner); b = bound; if (geoIndex == EDGE) b.rankObj.reverseMode = elObjDb.getEdgeReverseMode(rankBoundEl, ownerBoundEl); if (geoIndex == FACE) b.rankObj.reverseMode = elObjDb.getFaceReverseMode(rankBoundEl, ownerBoundEl); } } } // === Create periodic boundary data structure. === for (PerBoundMap::iterator it = elObjDb.getPeriodicVertices().begin(); it != elObjDb.getPeriodicVertices().end(); ++it) { if (elObjDb.isInRank(it->first.first, mpiRank) == false) continue; ElementObjectData& perDofEl0 = elObjDb.getElementsInRank(it->first.first)[mpiRank]; for (map::iterator elIt = elObjDb.getElementsInRank(it->first.second).begin(); elIt != elObjDb.getElementsInRank(it->first.second).end(); ++elIt) { int otherElementRank = elIt->first; ElementObjectData& perDofEl1 = elIt->second; AtomicBoundary bound; bound.rankObj.el = elObjDb.getElementPtr(perDofEl0.elIndex); bound.rankObj.elIndex = perDofEl0.elIndex; bound.rankObj.elType = elObjDb.getElementType(perDofEl0.elIndex); bound.rankObj.subObj = VERTEX; bound.rankObj.ithObj = perDofEl0.ithObject; bound.neighObj.el = elObjDb.getElementPtr(perDofEl1.elIndex); bound.neighObj.elIndex = perDofEl1.elIndex; bound.neighObj.elType = elObjDb.getElementType(perDofEl1.elIndex); bound.neighObj.subObj = VERTEX; bound.neighObj.ithObj = perDofEl1.ithObject; bound.type = it->second; AtomicBoundary& b = getNewPeriodic(otherElementRank); b = bound; } } for (PerBoundMap::iterator it = elObjDb.getPeriodicEdges().begin(); it != elObjDb.getPeriodicEdges().end(); ++it) { if (elObjDb.isInRank(it->first.first, mpiRank) == false) continue; ElementObjectData& perEdgeEl0 = elObjDb.getElementsInRank(it->first.first)[mpiRank]; for (map::iterator elIt = elObjDb.getElementsInRank(it->first.second).begin(); elIt != elObjDb.getElementsInRank(it->first.second).end(); ++elIt) { int otherElementRank = elIt->first; ElementObjectData& perEdgeEl1 = elIt->second; AtomicBoundary bound; bound.rankObj.el = elObjDb.getElementPtr(perEdgeEl0.elIndex); bound.rankObj.elIndex = perEdgeEl0.elIndex; bound.rankObj.elType = elObjDb.getElementType(perEdgeEl0.elIndex); bound.rankObj.subObj = EDGE; bound.rankObj.ithObj = perEdgeEl0.ithObject; bound.neighObj.el = elObjDb.getElementPtr(perEdgeEl1.elIndex); bound.neighObj.elIndex = perEdgeEl1.elIndex; bound.neighObj.elType = elObjDb.getElementType(perEdgeEl1.elIndex); bound.neighObj.subObj = EDGE; bound.neighObj.ithObj = perEdgeEl1.ithObject; bound.type = it->second; AtomicBoundary& b = getNewPeriodic(otherElementRank); b = bound; if (mpiRank > otherElementRank) b.neighObj.reverseMode = elObjDb.getEdgeReverseMode(perEdgeEl0, perEdgeEl1); else b.rankObj.reverseMode = elObjDb.getEdgeReverseMode(perEdgeEl0, perEdgeEl1); } } for (PerBoundMap::iterator it = elObjDb.getPeriodicFaces().begin(); it != elObjDb.getPeriodicFaces().end(); ++it) { if (elObjDb.isInRank(it->first.first, mpiRank) == false) continue; TEST_EXIT_DBG(elObjDb.getElements(it->first.first).size() == 1) ("Should not happen!\n"); TEST_EXIT_DBG(elObjDb.getElements(it->first.second).size() == 1) ("Should not happen!\n"); ElementObjectData& perFaceEl0 = elObjDb.getElementsInRank(it->first.first)[mpiRank]; for (map::iterator elIt = elObjDb.getElementsInRank(it->first.second).begin(); elIt != elObjDb.getElementsInRank(it->first.second).end(); ++elIt) { int otherElementRank = elIt->first; ElementObjectData& perFaceEl1 = elIt->second; AtomicBoundary bound; bound.rankObj.el = elObjDb.getElementPtr(perFaceEl0.elIndex); bound.rankObj.elIndex = perFaceEl0.elIndex; bound.rankObj.elType = elObjDb.getElementType(perFaceEl0.elIndex); bound.rankObj.subObj = FACE; bound.rankObj.ithObj = perFaceEl0.ithObject; bound.neighObj.el = elObjDb.getElementPtr(perFaceEl1.elIndex); bound.neighObj.elIndex = perFaceEl1.elIndex; bound.neighObj.elType = elObjDb.getElementType(perFaceEl1.elIndex); bound.neighObj.subObj = FACE; bound.neighObj.ithObj = perFaceEl1.ithObject; bound.type = it->second; AtomicBoundary& b = getNewPeriodic(otherElementRank); b = bound; if (mpiRank > otherElementRank) b.neighObj.reverseMode = elObjDb.getFaceReverseMode(perFaceEl0, perFaceEl1); else b.rankObj.reverseMode = elObjDb.getFaceReverseMode(perFaceEl0, perFaceEl1); } } // === Once we have this information, we must care about the order of the === // === atomic bounds in the three boundary handling object. Eventually === // === all the boundaries have to be in the same order on both ranks that === // === share the bounday. === StdMpi > stdMpi(mpiComm); stdMpi.send(own); stdMpi.recv(other); stdMpi.startCommunication(); // === The information about all neighbouring boundaries has been === // === received. So the rank tests if its own atomic boundaries are in === // === the same order. If not, the atomic boundaries are swaped to the === // === correct order. === for (RankToBoundMap::iterator rankIt = other.begin(); rankIt != other.end(); ++rankIt) { // === We have received from rank "rankIt->first" the ordered list of === // === element indices. Now, we have to sort the corresponding list in === // === this rank to get the same order. === for (unsigned int j = 0; j < rankIt->second.size(); j++) { // If the expected object is not at place, search for it. BoundaryObject &recvedBound = stdMpi.getRecvData()[rankIt->first][j].rankObj; if ((rankIt->second)[j].neighObj != recvedBound) { unsigned int k = j + 1; for (; k < rankIt->second.size(); k++) if ((rankIt->second)[k].neighObj == recvedBound) break; // The element must always be found, because the list is just in // another order. TEST_EXIT_DBG(k < rankIt->second.size())("Should never happen!\n"); // Swap the current with the found element. AtomicBoundary tmpBound = (rankIt->second)[k]; (rankIt->second)[k] = (rankIt->second)[j]; (rankIt->second)[j] = tmpBound; } } } // === Do the same for the periodic boundaries. === if (periodic.size() > 0) { stdMpi.clear(); RankToBoundMap sendBounds, recvBounds; for (RankToBoundMap::iterator rankIt = periodic.begin(); rankIt != periodic.end(); ++rankIt) { if (rankIt->first == mpiRank) continue; if (rankIt->first < mpiRank) sendBounds[rankIt->first] = rankIt->second; else recvBounds[rankIt->first] = rankIt->second; } stdMpi.send(sendBounds); stdMpi.recv(recvBounds); stdMpi.startCommunication(); for (RankToBoundMap::iterator rankIt = periodic.begin(); rankIt != periodic.end(); ++rankIt) { if (rankIt->first <= mpiRank) continue; for (unsigned int j = 0; j < rankIt->second.size(); j++) { BoundaryObject &recvRankObj = stdMpi.getRecvData()[rankIt->first][j].rankObj; BoundaryObject &recvNeighObj = stdMpi.getRecvData()[rankIt->first][j].neighObj; if (periodic[rankIt->first][j].neighObj != recvRankObj || periodic[rankIt->first][j].rankObj != recvNeighObj) { unsigned int k = j + 1; for (; k < rankIt->second.size(); k++) if (periodic[rankIt->first][k].neighObj == recvRankObj && periodic[rankIt->first][k].rankObj == recvNeighObj) break; // The element must always be found, because the list is just in // another order. TEST_EXIT_DBG(k < rankIt->second.size())("Should never happen!\n"); // Swap the current with the found element. AtomicBoundary tmpBound = (rankIt->second)[k]; (rankIt->second)[k] = (rankIt->second)[j]; (rankIt->second)[j] = tmpBound; } } } } // periodicBoundary.boundary.size() > 0 } void InteriorBoundary::serialize(ostream &out) { serialize(out, own); serialize(out, other); serialize(out, periodic); } void InteriorBoundary::serialize(ostream &out, RankToBoundMap& boundary) { FUNCNAME("InteriorBoundary::serialize()"); int mSize = boundary.size(); SerUtil::serialize(out, mSize); for (RankToBoundMap::iterator it = boundary.begin(); it != boundary.end(); ++it) { int rank = it->first; int boundSize = it->second.size(); SerUtil::serialize(out, rank); SerUtil::serialize(out, boundSize); for (int i = 0; i < boundSize; i++) { AtomicBoundary &bound = (it->second)[i]; SerUtil::serialize(out, bound.rankObj.elIndex); SerUtil::serialize(out, bound.rankObj.elType); SerUtil::serialize(out, bound.rankObj.subObj); SerUtil::serialize(out, bound.rankObj.ithObj); SerUtil::serialize(out, bound.rankObj.reverseMode); serializeExcludeList(out, bound.rankObj.excludedSubstructures); SerUtil::serialize(out, bound.neighObj.elIndex); SerUtil::serialize(out, bound.neighObj.elType); SerUtil::serialize(out, bound.neighObj.subObj); SerUtil::serialize(out, bound.neighObj.ithObj); SerUtil::serialize(out, bound.neighObj.reverseMode); serializeExcludeList(out, bound.neighObj.excludedSubstructures); SerUtil::serialize(out, bound.type); } } } void InteriorBoundary::deserialize(istream &in, Mesh *mesh) { map elIndexMap; mesh->getElementIndexMap(elIndexMap); deserialize(in, own, elIndexMap); deserialize(in, other, elIndexMap); deserialize(in, periodic, elIndexMap); } void InteriorBoundary::deserialize(istream &in, RankToBoundMap& boundary, map &elIndexMap) { FUNCNAME("InteriorBoundary::deserialize()"); int mSize = 0; SerUtil::deserialize(in, mSize); for (int i = 0; i < mSize; i++) { int rank = 0; int boundSize = 0; SerUtil::deserialize(in, rank); SerUtil::deserialize(in, boundSize); boundary[rank].resize(boundSize); for (int i = 0; i < boundSize; i++) { AtomicBoundary &bound = boundary[rank][i]; SerUtil::deserialize(in, bound.rankObj.elIndex); SerUtil::deserialize(in, bound.rankObj.elType); SerUtil::deserialize(in, bound.rankObj.subObj); SerUtil::deserialize(in, bound.rankObj.ithObj); SerUtil::deserialize(in, bound.rankObj.reverseMode); deserializeExcludeList(in, bound.rankObj.excludedSubstructures); SerUtil::deserialize(in, bound.neighObj.elIndex); SerUtil::deserialize(in, bound.neighObj.elType); SerUtil::deserialize(in, bound.neighObj.subObj); SerUtil::deserialize(in, bound.neighObj.ithObj); SerUtil::deserialize(in, bound.neighObj.reverseMode); deserializeExcludeList(in, bound.neighObj.excludedSubstructures); SerUtil::deserialize(in, bound.type); TEST_EXIT_DBG(elIndexMap.count(bound.rankObj.elIndex) == 1) ("Cannot find element with index %d for deserialization!\n", bound.rankObj.elIndex); TEST_EXIT_DBG(elIndexMap[bound.rankObj.elIndex]->getIndex() == bound.rankObj.elIndex)("Should not happen!\n"); bound.rankObj.el = elIndexMap[bound.rankObj.elIndex]; // For the case of periodic interior boundaries, a rank may have an // boundary with itself. In this case, also the pointer to the neighbour // object must be set correctly. if (elIndexMap.count(bound.neighObj.elIndex)) bound.neighObj.el = elIndexMap[bound.neighObj.elIndex]; else bound.neighObj.el = NULL; } } } AtomicBoundary& InteriorBoundary::getNewOwn(int rank) { int size = own[rank].size(); own[rank].resize(size + 1); return own[rank][size]; } AtomicBoundary& InteriorBoundary::getNewOther(int rank) { int size = other[rank].size(); other[rank].resize(size + 1); return other[rank][size]; } AtomicBoundary& InteriorBoundary::getNewPeriodic(int rank) { int size = periodic[rank].size(); periodic[rank].resize(size + 1); return periodic[rank][size]; } void InteriorBoundary::serializeExcludeList(ostream &out, ExcludeList &list) { int size = list.size(); SerUtil::serialize(out, size); for (int i = 0; i < size; i++) { SerUtil::serialize(out, list[i].first); SerUtil::serialize(out, list[i].second); } } void InteriorBoundary::deserializeExcludeList(istream &in, ExcludeList &list) { int size = 0; SerUtil::deserialize(in, size); list.resize(0); list.reserve(size); for (int i = 0; i < size; i++) { GeoIndex a; int b; SerUtil::deserialize(in, a); SerUtil::deserialize(in, b); list.push_back(make_pair(a, b)); } } }