diff --git a/AMDiS/src/Element.h b/AMDiS/src/Element.h index a831cad07de58359b9ef5b155de0eafaae6fba78..9a304a1ef70425ba59e6a5fae8e7235435137711 100644 --- a/AMDiS/src/Element.h +++ b/AMDiS/src/Element.h @@ -370,6 +370,13 @@ namespace AMDiS { /// Returns whether Element has sideElem as one of its sides. virtual bool hasSide(Element *sideElem) const = 0; + /** \brief + * Returns for a given element type number the element type number of the children. + * For 1d and 2d this is always 0, because element type number are used in the + * 3d case only. + */ + virtual int getChildType(int elType) const = 0; + /** \brief * Traverses an edge/face of a given element (this includes also all children of * the element having the same edge/face). All vertex dofs alonge this edge/face diff --git a/AMDiS/src/Line.h b/AMDiS/src/Line.h index 6901c896d99163fdaa570df4663b2f2beaf67327..a32573f1898c1194c230bd9ebcd12adb1ddb1855 100644 --- a/AMDiS/src/Line.h +++ b/AMDiS/src/Line.h @@ -141,6 +141,12 @@ namespace AMDiS { { return false; } + + /// Element type number is not used in 1d, so return 0. + inline int getChildType(int) const + { + return 0; + } std::string getTypeName() const { diff --git a/AMDiS/src/MeshStructure.cc b/AMDiS/src/MeshStructure.cc index b79b10d1ea4178abbd018792bc8ef4b3481e2571..00ab6af411969c45cc02664eac04dbea59838a7e 100644 --- a/AMDiS/src/MeshStructure.cc +++ b/AMDiS/src/MeshStructure.cc @@ -54,6 +54,32 @@ namespace AMDiS { commit(); } + void MeshStructure::init(Element *el, int ithSide, int elType) + { + FUNCNAME("MeshStructure::init()"); + + clear(); + + addAlongSide(el, ithSide, elType); + + commit(); + } + + void MeshStructure::addAlongSide(Element *el, int ithSide, int elType) + { + insertElement(el->isLeaf()); + + if (!el->isLeaf()) { + int s1 = el->getSideOfChild(0, ithSide, elType); + int s2 = el->getSideOfChild(1, ithSide, elType); + + if (s1 != -1) + addAlongSide(el->getFirstChild(), s1, el->getChildType(elType)); + if (s2 != -1) + addAlongSide(el->getSecondChild(), s2, el->getChildType(elType)); + } + } + void MeshStructure::reset() { currentIndex = 0; diff --git a/AMDiS/src/MeshStructure.h b/AMDiS/src/MeshStructure.h index d5f312f8851cadf543016d47654e1bd64bb5c179..9e5c64879017c96ca07c91ef5c40f510dbcf264f 100644 --- a/AMDiS/src/MeshStructure.h +++ b/AMDiS/src/MeshStructure.h @@ -41,9 +41,11 @@ namespace AMDiS { void clear(); - /// Creates a mesh structure code from a Mesh object by traversing it in preorder. + /// Creates a mesh structure code from a mesh object by traversing it in preorder. void init(Mesh *mesh); + void init(Element *el, int ithSide, int elType); + void init(const std::vector<unsigned long int>& initCode, int n) { code = initCode; @@ -102,13 +104,13 @@ namespace AMDiS { bool cont = true; while (cont) { if (isLeafElement()) - MSG("0"); + std::cout << "0"; else - MSG("1"); + std::cout << "1"; cont = nextElement(); } - MSG("\n"); + std::cout << std::endl; } /// Returns the mesh structure code. @@ -131,6 +133,8 @@ namespace AMDiS { /// Insert a new element to the structure code. Is used by the init function. void insertElement(bool isLeaf); + void addAlongSide(Element *el, int ithSide, int elType); + /// Merges two mesh structure codes to one structure code. void merge(MeshStructure *structure1, MeshStructure *structure2, diff --git a/AMDiS/src/ParallelDomainBase.cc b/AMDiS/src/ParallelDomainBase.cc index 2fdde4b58c463846021fb827c393e8824412ab07..20026e4fd59ab8e66d0e6c3c3998aa8f20c6c532 100644 --- a/AMDiS/src/ParallelDomainBase.cc +++ b/AMDiS/src/ParallelDomainBase.cc @@ -18,6 +18,7 @@ #include "ElementFileWriter.h" #include "VertexVector.h" #include "StdMpi.h" +#include "MeshStructure.h" #include "petscksp.h" @@ -694,7 +695,7 @@ namespace AMDiS { stdMpi.send(sendIt->first, dofs); } - stdMpi.startCommunication(); + stdMpi.startCommunication<int>(); #endif #if 1 @@ -758,9 +759,56 @@ namespace AMDiS { void ParallelDomainBase::checkMeshChange() { + // === If mesh has not been changed, return. === + if (mesh->getChangeIndex() == lastMeshChangeIndex) return; + + // === Create mesh structure codes for all ranks boundary elements. === + + typedef std::vector<MeshStructure> MeshCodeVec; + std::map<int, MeshCodeVec > sendCodes; + + for (RankToBoundMap::iterator it = myIntBoundary.boundary.begin(); + it != myIntBoundary.boundary.end(); ++it) { + + for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin(); + boundIt != it->second.end(); ++boundIt) { + + MeshStructure elCode; + elCode.init(boundIt->rankObj.el, boundIt->rankObj.ithObj, + boundIt->rankObj.elType); + sendCodes[it->first].push_back(elCode); + } + } + + StdMpi<MeshCodeVec, MeshCodeVec> stdMpi(mpiComm); + stdMpi.prepareCommunication(true); + stdMpi.send(sendCodes); + + for (RankToBoundMap::iterator it = otherIntBoundary.boundary.begin(); + it != otherIntBoundary.boundary.end(); ++it) + stdMpi.recv(it->first); + + stdMpi.startCommunication<unsigned long int>(); + + +#if 0 + for (RankToBoundMap::iterator it = otherIntBoundary.boundary.begin(); + it != otherIntBoundary.boundary.end(); ++it) { + + for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin(); + boundIt != it->second.end(); ++boundIt) { + + MeshStructure elCode; + elCode.init(boundIt->rankObj.el, + boundIt->rankObj.ithObj, + boundIt->rankObj.elType); + } + } +#endif + std::cout << "MESH HAS BEEN CHANGED!" << std::endl; exit(0); @@ -1262,7 +1310,7 @@ namespace AMDiS { for (std::map<int, int>::iterator recvIt = recvNewDofs.begin(); recvIt != recvNewDofs.end(); ++recvIt, i++) stdMpi.recv(recvIt->first, recvIt->second * 2); - stdMpi.startCommunication(); + stdMpi.startCommunication<int>(); std::map<int, std::vector<std::pair<DegreeOfFreedom, DegreeOfFreedom> > > &dofMap = stdMpi.getRecvData(); // === Change dof indices at boundary from other ranks. === diff --git a/AMDiS/src/StdMpi.h b/AMDiS/src/StdMpi.h index 7365f557657711b636c5c581ed99cbb9fc5eb1b2..56ad1a0b5df431354995dbcdb5842955e77cd2d7 100644 --- a/AMDiS/src/StdMpi.h +++ b/AMDiS/src/StdMpi.h @@ -24,6 +24,7 @@ #include <map> #include "mpi.h" +#include "MeshStructure.h" namespace AMDiS { @@ -32,7 +33,16 @@ namespace AMDiS { return data.size() * 2; } - void makeIntBuf(std::map<const DegreeOfFreedom*, DegreeOfFreedom> &data, int* buf) + int intSizeOf(std::vector<MeshStructure> &data) + { + int s = 0; + for (int i = 0; i < data.size(); i++) + s += data[i].getCode().size() + 2; + + return s; + } + + void makeBuf(std::map<const DegreeOfFreedom*, DegreeOfFreedom> &data, int* buf) { int i = 0; for (std::map<const DegreeOfFreedom*, DegreeOfFreedom>::iterator it = data.begin(); @@ -42,8 +52,8 @@ namespace AMDiS { } } - void makeFromIntBuf(std::vector<std::pair<DegreeOfFreedom, DegreeOfFreedom> > &data, - int *buf, int bufSize) + void makeFromBuf(std::vector<std::pair<DegreeOfFreedom, DegreeOfFreedom> > &data, + int *buf, int bufSize) { if (bufSize == 0) return; @@ -58,7 +68,35 @@ namespace AMDiS { } while (i < bufSize); } + void makeBuf(std::vector<MeshStructure> &data, unsigned long int *buf) + { + int pos = 0; + for (int i = 0; i < data.size(); i++) { + buf[pos++] = data[i].getCode().size(); + buf[pos++] = data[i].getNumElements(); + for (int j = 0; j < data[i].getCode().size(); j++) + buf[pos++] = data[i].getCode()[j]; + } + } + void makeFromBuf(std::vector<MeshStructure> &data, unsigned long int *buf, int bufSize) + { + int pos = 0; + + while (pos < bufSize) { + int codeSize = buf[pos++]; + int nElements = buf[pos++]; + std::vector<unsigned long int> code; + code.reserve(codeSize); + for (int i = 0; i < codeSize; i++) + code[i] = buf[pos++]; + + MeshStructure meshCode; + meshCode.init(code, nElements); + + data.push_back(meshCode); + } + } template<typename SendT, typename RecvT> class StdMpi @@ -116,21 +154,48 @@ namespace AMDiS { return recvData; } + void commDataSize() + { + MPI::Request request[sendData.size() + recvDataSize.size()]; + std::vector<int> sendBuffers, recvBuffers; + int requestCounter = 0; + + int i = 0; + for (typename std::map<int, int>::iterator sendIt = sendDataSize.begin(); + sendIt != sendDataSize.end(); ++sendIt) { + sendBuffers.push_back(sendIt->second); + request[requestCounter++] = + mpiComm.Isend(&(sendBuffers[i]), 1, MPI_INT, sendIt->first, 0); + i++; + } + + for (std::map<int, int>::iterator recvIt = recvDataSize.begin(); + recvIt != recvDataSize.end(); ++recvIt) + request[requestCounter++] = + mpiComm.Irecv(&(recvIt->second), 1, MPI_INT, recvIt->first, 0); + + MPI::Request::Waitall(requestCounter, request); + } + + template<class T> void startCommunication() { FUNCNAME("StdMpi::startCommunication()"); TEST_EXIT_DBG(commPrepared)("Communication is not prepared!\n"); + if (exchangeDataSize) + commDataSize(); + MPI::Request request[sendData.size() + recvDataSize.size()]; int requestCounter = 0; - std::vector<int*> sendBuffers, recvBuffers; + std::vector<T*> sendBuffers, recvBuffers; for (typename std::map<int, SendT>::iterator sendIt = sendData.begin(); sendIt != sendData.end(); ++sendIt) { - int bufferSize = intSizeOf(sendIt->second); - int* buf = new int[bufferSize]; - makeIntBuf(sendIt->second, buf); + int bufferSize = sendDataSize[sendIt->first]; + T* buf = new T[bufferSize]; + makeBuf(sendIt->second, buf); request[requestCounter++] = mpiComm.Isend(buf, bufferSize, MPI_INT, sendIt->first, 0); @@ -140,7 +205,7 @@ namespace AMDiS { for (std::map<int, int>::iterator recvIt = recvDataSize.begin(); recvIt != recvDataSize.end(); ++recvIt) { - int* buf = new int[recvIt->second]; + T* buf = new T[recvIt->second]; request[requestCounter++] = mpiComm.Irecv(buf, recvIt->second, MPI_INT, recvIt->first, 0); @@ -157,7 +222,7 @@ namespace AMDiS { int i = 0; for (std::map<int, int>::iterator recvIt = recvDataSize.begin(); recvIt != recvDataSize.end(); ++recvIt) { - makeFromIntBuf(recvData[recvIt->first], recvBuffers[i], recvIt->second); + makeFromBuf(recvData[recvIt->first], recvBuffers[i], recvIt->second); delete [] recvBuffers[i]; i++; } diff --git a/AMDiS/src/Tetrahedron.h b/AMDiS/src/Tetrahedron.h index 329d2460a879da3528b07aa63cf910f6c36f4373..4de965e20cd7bf1035052ebad6d9634e5013bbca 100644 --- a/AMDiS/src/Tetrahedron.h +++ b/AMDiS/src/Tetrahedron.h @@ -123,6 +123,12 @@ namespace AMDiS { { return true; } + + /// Return the new element type of the children. + inline int getChildType(int elType) const + { + return (elType + 1) % 3; + } std::string getTypeName() const { diff --git a/AMDiS/src/Triangle.h b/AMDiS/src/Triangle.h index fe2f7bd922662a9a50e1e19d357c839796da4079..fa1053891847a83633f9ce3e98a4490b8f8bacc3 100644 --- a/AMDiS/src/Triangle.h +++ b/AMDiS/src/Triangle.h @@ -114,6 +114,12 @@ namespace AMDiS { return false; } + /// Element type number is not used in 2d, so return 0. + inline int getChildType(int) const + { + return 0; + } + /// implements Element::getSideOfChild() virtual int getSideOfChild(int child, int side, int) const {