From eb0262c3b7ab8769db95917a1aa200cddc1edd84 Mon Sep 17 00:00:00 2001 From: Thomas Witkowski <thomas.witkowski@gmx.de> Date: Mon, 27 Apr 2009 06:45:01 +0000 Subject: [PATCH] Parallel AMDiS --- AMDiS/src/InteriorBoundary.cc | 7 ++ AMDiS/src/InteriorBoundary.h | 18 +-- AMDiS/src/ParallelDomainProblem.cc | 178 ++++++++++++++++++++++++----- AMDiS/src/ParallelDomainProblem.h | 9 ++ 4 files changed, 178 insertions(+), 34 deletions(-) diff --git a/AMDiS/src/InteriorBoundary.cc b/AMDiS/src/InteriorBoundary.cc index 999c3051..4da79698 100644 --- a/AMDiS/src/InteriorBoundary.cc +++ b/AMDiS/src/InteriorBoundary.cc @@ -1,4 +1,11 @@ #include "InteriorBoundary.h" namespace AMDiS { + + AtomicBoundary& InteriorBoundary::getNewAtomicBoundary(int rank) + { + boundary[rank].resize(boundary[rank].size() + 1); + return boundary[rank][boundary[rank].size() - 1]; + } + } diff --git a/AMDiS/src/InteriorBoundary.h b/AMDiS/src/InteriorBoundary.h index 43eef370..59f5065c 100644 --- a/AMDiS/src/InteriorBoundary.h +++ b/AMDiS/src/InteriorBoundary.h @@ -23,6 +23,7 @@ #define AMDIS_INTERIORBOUNDARY_H #include <vector> +#include <map> #include "MacroElement.h" @@ -32,7 +33,7 @@ namespace AMDiS { struct BoundaryObject { /// The macro element to which the boundary element corresponds to. - MacroElement& el; + Element* el; /** \brief * Defines the geometrical object at the boundary. It must be "a part" of the @@ -58,13 +59,10 @@ namespace AMDiS { */ struct AtomicBoundary { /// The rank's part of the boundary. - BoundaryObject& rankObject; + BoundaryObject rankObject; /// The object on the other side of the boundary. - BoundaryObject& neighbourObject; - - /// The rank number of the process at the other side of the boundary. - int neighbourRank; + BoundaryObject neighbourObject; }; /** \brief @@ -73,8 +71,12 @@ namespace AMDiS { */ class InteriorBoundary { public: - protected: - std::vector<AtomicBoundary> boundary; + InteriorBoundary() {} + + AtomicBoundary& getNewAtomicBoundary(int rank); + + public: + std::map<int, std::vector<AtomicBoundary> > boundary; }; } diff --git a/AMDiS/src/ParallelDomainProblem.cc b/AMDiS/src/ParallelDomainProblem.cc index 473bc60f..5b86db88 100644 --- a/AMDiS/src/ParallelDomainProblem.cc +++ b/AMDiS/src/ParallelDomainProblem.cc @@ -45,31 +45,26 @@ namespace AMDiS { // and now partition the mesh partitionMesh(adaptInfo); - std::map<int, std::set<int> > partitionDofs; + /// === Determine to each dof the set of partitions the dof belongs to. === + + std::map<const DegreeOfFreedom*, std::set<int> > partitionDofs; TraverseStack stack; ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); - int nLeaves = 0; while (elInfo) { Element *element = elInfo->getElement(); - // Hidde elements which are not part of ranks partition. - PartitionElementData *partitionData = - dynamic_cast<PartitionElementData*> - (element->getElementData(PARTITION_ED)); - if (partitionData->getPartitionStatus() != IN) { - } else { - } - // Determine to each dof the partition(s) it corresponds to. for (int i = 0; i < 3; i++) - partitionDofs[element->getDOF(i, 0)].insert(partitionVec[element->getIndex()]); + partitionDofs[element->getDOF(i)].insert(partitionVec[element->getIndex()]); elInfo = stack.traverseNext(elInfo); } - std::vector<int> rankDofs; - for (std::map<int, std::set<int> >::iterator it = partitionDofs.begin(); + /// === Determine the set of ranks dofs and the dofs ownership at the boundary. === + + std::vector<const DegreeOfFreedom*> rankDofs; + for (std::map<const DegreeOfFreedom*, std::set<int> >::iterator it = partitionDofs.begin(); it != partitionDofs.end(); ++it) { for (std::set<int>::iterator itpart1 = it->second.begin(); @@ -79,30 +74,64 @@ namespace AMDiS { if (it->second.size() == 1) { rankDofs.push_back(it->first); } else { + // This dof is at the ranks boundary. It is owned by the rank only if + // the rank number is the highest of all ranks containing this dof. + bool insert = true; + int highestRank = mpiRank; for (std::set<int>::iterator itpart2 = it->second.begin(); itpart2 != it->second.end(); ++itpart2) { - if (*itpart2 > mpiRank) { + if (*itpart2 > mpiRank) insert = false; - break; - } + + if (*itpart2 > highestRank) + highestRank = *itpart2; } - if (insert) { + + if (insert) rankDofs.push_back(it->first); - } + + boundaryDofs[it->first] = highestRank; } } } } - if (mpiRank == 1) { - std::cout << "RANKS dofs = "; - for (int i = 0; i < rankDofs.size(); i++) - std::cout << rankDofs[i] << " "; - std::cout << std::endl; + + // === Create interior boundary information === + + elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL | Mesh::FILL_NEIGH); + while (elInfo) { + Element *element = elInfo->getElement(); + + // Hidde elements which are not part of ranks partition. + PartitionElementData *partitionData = + dynamic_cast<PartitionElementData*>(element->getElementData(PARTITION_ED)); + if (partitionData->getPartitionStatus() == IN) { + for (int i = 0; i < 3; i++) { + if (!elInfo->getNeighbour(i)) + continue; + + PartitionElementData *neighbourPartitionData = + dynamic_cast<PartitionElementData*>(elInfo->getNeighbour(i)->getElementData(PARTITION_ED)); + if (neighbourPartitionData->getPartitionStatus() == OUT) { + AtomicBoundary& bound = interiorBoundary. + getNewAtomicBoundary(partitionVec[elInfo->getNeighbour(i)->getIndex()]); + bound.rankObject.el = element; + bound.rankObject.subObjAtBoundary = EDGE; + bound.rankObject.ithObjAtBoundary = i; + bound.neighbourObject.el = elInfo->getNeighbour(i); + bound.neighbourObject.subObjAtBoundary = EDGE; + bound.neighbourObject.ithObjAtBoundary = -1; + } + } + } + + elInfo = stack.traverseNext(elInfo); } + // === Remove all macro elements that are not part of the rank partition. === std::vector<MacroElement*> macrosToRemove; @@ -124,15 +153,15 @@ namespace AMDiS { int *gOrder = (int*)(malloc(sizeof(int) * rankDofs.size())); int *lOrder = (int*)(malloc(sizeof(int) * rankDofs.size())); - for (std::vector<int>::iterator it = rankDofs.begin(); + for (std::vector<const DegreeOfFreedom*>::iterator it = rankDofs.begin(); it != rankDofs.end(); ++it) { - gOrder[nRankDOFs++] = *it; + gOrder[nRankDOFs++] = (*it)[0]; } int rstart = 0; MPI_Scan(&nRankDOFs, &rstart, 1, MPI_INT, MPI_SUM, PETSC_COMM_WORLD); rstart -= nRankDOFs; - + for (int i = 0; i < nRankDOFs; i++) { lOrder[i] = rstart + i; } @@ -141,6 +170,103 @@ namespace AMDiS { free(gOrder); free(lOrder); + + /// === Create information which dof indices must be send and which must be received. === + + std::map<int, std::map<DegreeOfFreedom, DegreeOfFreedom> > sendNewDofs; + std::map<int, std::vector<DegreeOfFreedom> > recvNewDofs; + + for (std::map<const DegreeOfFreedom*, int>::iterator it = boundaryDofs.begin(); + it != boundaryDofs.end(); + ++it) { + if (it->second == mpiRank) { + int oldDofIndex = (it->first)[0]; + int newDofIndex = 0; + for (int i = 0; i < static_cast<int>(rankDofs.size()); i++) { + if (rankDofs[i] == it->first) { + newDofIndex = rstart + i; + break; + } + } + + for (std::set<int>::iterator itRanks = partitionDofs[it->first].begin(); + itRanks != partitionDofs[it->first].end(); + ++itRanks) { + if (*itRanks != mpiRank) { + sendNewDofs[*itRanks][oldDofIndex] = newDofIndex; + } + } + } else { + recvNewDofs[it->second].push_back((it->first)[0]); + } + } + + /// === Send and receive the dof indices at boundary. === + + std::vector<int*> sendBuffers(sendNewDofs.size()); + std::vector<int*> recvBuffers(recvNewDofs.size()); + + int i = 0; + for (std::map<int, std::map<DegreeOfFreedom, DegreeOfFreedom> >::iterator sendIt = sendNewDofs.begin(); + sendIt != sendNewDofs.end(); + ++sendIt, i++) { + sendBuffers[i] = new int[sendIt->second.size() * 2]; + int c = 0; + for (std::map<DegreeOfFreedom, DegreeOfFreedom>::iterator dofIt = sendIt->second.begin(); + dofIt != sendIt->second.end(); + ++dofIt, c += 2) { + sendBuffers[i][c] = dofIt->first; + sendBuffers[i][c + 1] = dofIt->second; + } + + mpiComm.Isend(sendBuffers[i], sendIt->second.size() * 2, MPI_INT, sendIt->first, 0); + } + + i = 0; + for (std::map<int, std::vector<DegreeOfFreedom> >::iterator recvIt = recvNewDofs.begin(); + recvIt != recvNewDofs.end(); + ++recvIt, i++) { + recvBuffers[i] = new int[recvIt->second.size() * 2]; + + mpiComm.Irecv(recvBuffers[i], recvIt->second.size() * 2, MPI_INT, recvIt->first, 0); + } + + + mpiComm.Barrier(); + + /// === Change dof indices at boundary from other ranks. === + + i = 0; + for (std::map<int, std::vector<DegreeOfFreedom> >::iterator recvIt = recvNewDofs.begin(); + recvIt != recvNewDofs.end(); + ++recvIt, i++) { + + for (int j = 0; j < static_cast<int>(recvIt->second.size()); j++) { + for (std::map<const DegreeOfFreedom*, int>::iterator dofIt = boundaryDofs.begin(); + dofIt != boundaryDofs.end(); + ++dofIt) { + if ((dofIt->first)[0] == recvBuffers[i][j * 2]) { + const_cast<DegreeOfFreedom*>(dofIt->first)[0] = recvBuffers[i][j * 2 + 1]; + break; + } + } + } + + delete [] recvBuffers[i]; + } + + i = 0; + for (std::map<int, std::map<DegreeOfFreedom, DegreeOfFreedom> >::iterator sendIt = sendNewDofs.begin(); + sendIt != sendNewDofs.end(); + ++sendIt, i++) { + delete [] sendBuffers[i]; + } + + /// === Change dof indices for rank partition. === + + for (int i = 0; i < static_cast<int>(rankDofs.size()); i++) { + const_cast<DegreeOfFreedom*>(rankDofs[i])[0] = rstart + i; + } } void ParallelDomainProblemBase::exitParallelization(AdaptInfo *adaptInfo) diff --git a/AMDiS/src/ParallelDomainProblem.h b/AMDiS/src/ParallelDomainProblem.h index ac369330..fe206322 100644 --- a/AMDiS/src/ParallelDomainProblem.h +++ b/AMDiS/src/ParallelDomainProblem.h @@ -29,6 +29,7 @@ #include "ProblemIterationInterface.h" #include "FiniteElemSpace.h" #include "AdaptInfo.h" +#include "InteriorBoundary.h" #include "petscao.h" #include "mpi.h" @@ -156,6 +157,14 @@ namespace AMDiS { /// Number of DOFs in the rank mesh. int nRankDOFs; + + /** \brief + * Defines the interioir boundaries of the domain that result from partitioning + * the whole mesh. + */ + InteriorBoundary interiorBoundary; + + std::map<const DegreeOfFreedom*, int> boundaryDofs; }; class ParallelDomainProblemScal : public ParallelDomainProblemBase -- GitLab