From 4873445c394074b5d1400a366ec5469808167f43 Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Tue, 20 Oct 2009 13:29:53 +0000
Subject: [PATCH] Fixed problem when using multiple periodic boundary
conditions in parallel domain decomposition.
---
AMDiS/src/Element.h | 29 ++++
AMDiS/src/ElementDofIterator.cc | 2 +-
AMDiS/src/ElementDofIterator.h | 4 +-
AMDiS/src/Global.h | 3 +
AMDiS/src/InteriorBoundary.cc | 28 ++--
AMDiS/src/InteriorBoundary.h | 14 +-
AMDiS/src/Line.h | 22 ++-
AMDiS/src/ParallelDomainBase.cc | 273 ++++++++++++--------------------
AMDiS/src/ParallelDomainBase.h | 28 ----
AMDiS/src/Tetrahedron.h | 28 ++--
AMDiS/src/Triangle.cc | 129 ++++++++++++++-
AMDiS/src/Triangle.h | 6 +
12 files changed, 322 insertions(+), 244 deletions(-)
diff --git a/AMDiS/src/Element.h b/AMDiS/src/Element.h
index bde0de6e..f5e07c3b 100644
--- a/AMDiS/src/Element.h
+++ b/AMDiS/src/Element.h
@@ -371,6 +371,35 @@ namespace AMDiS {
/// Returns whether Element has sideElem as one of its sides.
virtual bool hasSide(Element *sideElem) const = 0;
+ /** \brief
+ * Traverses an edge of a given element (this includes also all children of the
+ * element having the same edge). All vertex dofs alonge this edge are assembled
+ * and put together to a list.
+ *
+ * \param[in] feSpace FE space which is used to get the dofs.
+ * \param[in] ithEdge Defines the edge on which all the vertex dofs
+ * are assembled.
+ * \param[out] dofs List of dofs, where the result is stored.
+ * \param[in] parentVertices If true, also the two vertices of the parent
+ * element are put into the result list.
+ */
+ virtual void getVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs, bool parentVertices = 0) const = 0;
+
+ /** \brief
+ * Traverses an edge of a given element (this includes also all children of the
+ * element having the same edge). All non vertex dofs alonge this edge are
+ * assembled and put together to a list.
+ *
+ * \param[in] feSpace FE space which is used to get the dofs.
+ * \param[in] ithEdge Defines the edge on which all the non vertex
+ * dofs are assembled.
+ * \param[out] dofs All dofs are put to this dof list.
+ */
+ virtual void getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs) const = 0;
+
+
/** \} */
// ===== other public methods =================================================
diff --git a/AMDiS/src/ElementDofIterator.cc b/AMDiS/src/ElementDofIterator.cc
index 4e44a572..6fa9bba5 100644
--- a/AMDiS/src/ElementDofIterator.cc
+++ b/AMDiS/src/ElementDofIterator.cc
@@ -6,7 +6,7 @@
namespace AMDiS {
- void ElementDofIterator::reset(Element* el)
+ void ElementDofIterator::reset(const Element* el)
{
FUNCNAME("ElementDofIterator::reset()");
diff --git a/AMDiS/src/ElementDofIterator.h b/AMDiS/src/ElementDofIterator.h
index f4931a94..938a27a0 100644
--- a/AMDiS/src/ElementDofIterator.h
+++ b/AMDiS/src/ElementDofIterator.h
@@ -54,7 +54,7 @@ namespace AMDiS {
{}
/// Start a new traverse with the given element.
- void reset(Element* el);
+ void reset(const Element* el);
/// Go to next dof. Returns false, if there is dof anymore.
bool next();
@@ -116,7 +116,7 @@ namespace AMDiS {
int* orderPosition;
- Element* element;
+ const Element* element;
/// Current position (i.e., vertex, edge, face) of the traverse.
int pos;
diff --git a/AMDiS/src/Global.h b/AMDiS/src/Global.h
index c32a1504..76e9f459 100644
--- a/AMDiS/src/Global.h
+++ b/AMDiS/src/Global.h
@@ -74,6 +74,9 @@ namespace AMDiS {
/// datatype for degrees of freedom
typedef signed int DegreeOfFreedom;
+ /// Defines type for a vector of DOF pointers.
+ typedef std::vector<const DegreeOfFreedom*> DofContainer;
+
/// returns the GeoIndex of d for dimension dim.
#define INDEX_OF_DIM(d, dim) (static_cast<GeoIndex>((d == dim) ? CENTER : d + 1))
diff --git a/AMDiS/src/InteriorBoundary.cc b/AMDiS/src/InteriorBoundary.cc
index bf90070e..d751a1d0 100644
--- a/AMDiS/src/InteriorBoundary.cc
+++ b/AMDiS/src/InteriorBoundary.cc
@@ -21,13 +21,13 @@ namespace AMDiS {
for (int i = 0; i < boundSize; i++) {
AtomicBoundary &bound = (it->second)[i];
- SerUtil::serialize(out, bound.rankObject.elIndex);
- SerUtil::serialize(out, bound.rankObject.subObjAtBoundary);
- SerUtil::serialize(out, bound.rankObject.ithObjAtBoundary);
+ SerUtil::serialize(out, bound.rankObj.elIndex);
+ SerUtil::serialize(out, bound.rankObj.subObj);
+ SerUtil::serialize(out, bound.rankObj.ithObj);
- SerUtil::serialize(out, bound.neighbourObject.elIndex);
- SerUtil::serialize(out, bound.neighbourObject.subObjAtBoundary);
- SerUtil::serialize(out, bound.neighbourObject.ithObjAtBoundary);
+ SerUtil::serialize(out, bound.neighObj.elIndex);
+ SerUtil::serialize(out, bound.neighObj.subObj);
+ SerUtil::serialize(out, bound.neighObj.ithObj);
}
}
}
@@ -47,16 +47,16 @@ namespace AMDiS {
for (int i = 0; i < boundSize; i++) {
AtomicBoundary &bound = boundary[rank][i];
- SerUtil::deserialize(in, bound.rankObject.elIndex);
- SerUtil::deserialize(in, bound.rankObject.subObjAtBoundary);
- SerUtil::deserialize(in, bound.rankObject.ithObjAtBoundary);
+ SerUtil::deserialize(in, bound.rankObj.elIndex);
+ SerUtil::deserialize(in, bound.rankObj.subObj);
+ SerUtil::deserialize(in, bound.rankObj.ithObj);
- SerUtil::deserialize(in, bound.neighbourObject.elIndex);
- SerUtil::deserialize(in, bound.neighbourObject.subObjAtBoundary);
- SerUtil::deserialize(in, bound.neighbourObject.ithObjAtBoundary);
+ SerUtil::deserialize(in, bound.neighObj.elIndex);
+ SerUtil::deserialize(in, bound.neighObj.subObj);
+ SerUtil::deserialize(in, bound.neighObj.ithObj);
- bound.rankObject.el = elIndexMap[bound.rankObject.elIndex];
- bound.neighbourObject.el = NULL;
+ bound.rankObj.el = elIndexMap[bound.rankObj.elIndex];
+ bound.neighObj.el = NULL;
}
}
}
diff --git a/AMDiS/src/InteriorBoundary.h b/AMDiS/src/InteriorBoundary.h
index 79a48fa4..424f23ff 100644
--- a/AMDiS/src/InteriorBoundary.h
+++ b/AMDiS/src/InteriorBoundary.h
@@ -43,18 +43,18 @@ namespace AMDiS {
* Defines the geometrical object at the boundary. It must be "a part" of the
* macro element \ref el, i.e., either 1 (a vertex), 2 (an edge) or 3 (a face).
*/
- GeoIndex subObjAtBoundary;
+ GeoIndex subObj;
/** \brief
* Defines which of vertix, edge or face of the macro element is part of the
* boundary.
*
- * Example: If the macro element is a triangle, than \ref subObjAtBoundary may
- * be either 1 (vertex) or 2 (edge). Assume its the last one. So this variable
- * defines which of the three possible edges of the triangle is at the interior
+ * Example: If the macro element is a triangle, than \ref subObj may be either
+ * 1 (vertex) or 2 (edge). Assume its the last one. So this variable defines
+ * which of the three possible edges of the triangle is at the interior
* boundary.
*/
- int ithObjAtBoundary;
+ int ithObj;
};
/** \brief
@@ -63,10 +63,10 @@ namespace AMDiS {
*/
struct AtomicBoundary {
/// The rank's part of the boundary.
- BoundaryObject rankObject;
+ BoundaryObject rankObj;
/// The object on the other side of the boundary.
- BoundaryObject neighbourObject;
+ BoundaryObject neighObj;
};
/** \brief
diff --git a/AMDiS/src/Line.h b/AMDiS/src/Line.h
index cd6ef9d8..486c523a 100644
--- a/AMDiS/src/Line.h
+++ b/AMDiS/src/Line.h
@@ -124,19 +124,19 @@ namespace AMDiS {
return false;
}
- /// implements Element::isLine. Returns true because this element is a Line
+ /// Returns true because this element is a Line.
inline bool isLine() const
{
return true;
}
- /// implements Element::isTriangle. Returns false because this element is a Line
+ /// Returns false because this element is a Line.
inline bool isTriangle() const
{
return false;
}
- /// implements Element::isTetrahedron. Returns false because this element is a Line
+ /// Returns false because this element is a Line
inline bool isTetrahedron() const
{
return false;
@@ -147,6 +147,22 @@ namespace AMDiS {
return "Line";
}
+ void getVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs, bool parentVertices = 0) const
+ {
+ FUNCNAME("Line::getVertexDofs()");
+
+ ERROR_EXIT("Not yet implemented!\n");
+ }
+
+ void getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs) const
+ {
+ FUNCNAME("Line::getNonVertexDofs()");
+
+ ERROR_EXIT("Not yet implemented!\n");
+ }
+
protected:
/** \brief
* vertexOfEdge[i][j] is the local number of the j-th vertex of the i-th
diff --git a/AMDiS/src/ParallelDomainBase.cc b/AMDiS/src/ParallelDomainBase.cc
index cd1fceae..d41a4e6b 100644
--- a/AMDiS/src/ParallelDomainBase.cc
+++ b/AMDiS/src/ParallelDomainBase.cc
@@ -153,6 +153,7 @@ namespace AMDiS {
createPeriodicMap();
}
+ // exit(0);
#if (DEBUG != 0)
DbgTestCommonDofs(true);
@@ -982,16 +983,16 @@ namespace AMDiS {
// === Create information about the boundary between the two elements. ===
AtomicBoundary bound;
- bound.rankObject.el = element;
- bound.rankObject.elIndex = element->getIndex();
- bound.rankObject.subObjAtBoundary = EDGE;
- bound.rankObject.ithObjAtBoundary = i;
+ bound.rankObj.el = element;
+ bound.rankObj.elIndex = element->getIndex();
+ bound.rankObj.subObj = EDGE;
+ bound.rankObj.ithObj = i;
// Do not set a pointer to the element, because if will be deleted from
// mesh after partitioning and the pointer would become invalid.
- bound.neighbourObject.el = NULL;
- bound.neighbourObject.elIndex = elInfo->getNeighbour(i)->getIndex();
- bound.neighbourObject.subObjAtBoundary = EDGE;
- bound.neighbourObject.ithObjAtBoundary = elInfo->getSideOfNeighbour(i);
+ bound.neighObj.el = NULL;
+ bound.neighObj.elIndex = elInfo->getNeighbour(i)->getIndex();
+ bound.neighObj.subObj = EDGE;
+ bound.neighObj.ithObj = elInfo->getSideOfNeighbour(i);
// i == 2 => getSideOfNeighbour(i) == 2
TEST_EXIT_DBG(i != 2 || elInfo->getSideOfNeighbour(i) == 2)
@@ -1049,8 +1050,8 @@ namespace AMDiS {
int nSendInt = rankIt->second.size();
int* buffer = new int[nSendInt * 2];
for (int i = 0; i < nSendInt; i++) {
- buffer[i * 2] = (rankIt->second)[i].rankObject.elIndex;
- buffer[i * 2 + 1] = (rankIt->second)[i].rankObject.ithObjAtBoundary;
+ buffer[i * 2] = (rankIt->second)[i].rankObj.elIndex;
+ buffer[i * 2 + 1] = (rankIt->second)[i].rankObj.ithObj;
}
sendBuffers.push_back(buffer);
@@ -1086,12 +1087,12 @@ namespace AMDiS {
for (int j = 0; j < static_cast<int>(rankIt->second.size()); j++) {
// If the expected object is not at place, search for it.
- if ((rankIt->second)[j].neighbourObject.elIndex != recvBuffers[i][j * 2] ||
- (rankIt->second)[j].neighbourObject.ithObjAtBoundary != recvBuffers[i][j * 2 + 1]) {
+ if ((rankIt->second)[j].neighObj.elIndex != recvBuffers[i][j * 2] ||
+ (rankIt->second)[j].neighObj.ithObj != recvBuffers[i][j * 2 + 1]) {
int k = j + 1;
for (; k < static_cast<int>(rankIt->second.size()); k++)
- if ((rankIt->second)[k].neighbourObject.elIndex == recvBuffers[i][j * 2] &&
- (rankIt->second)[k].neighbourObject.ithObjAtBoundary == recvBuffers[i][j * 2 + 1])
+ if ((rankIt->second)[k].neighObj.elIndex == recvBuffers[i][j * 2] &&
+ (rankIt->second)[k].neighObj.ithObj == recvBuffers[i][j * 2 + 1])
break;
// The element must always be found, because the list is just in another order.
@@ -1130,8 +1131,8 @@ namespace AMDiS {
int nSendInt = rankIt->second.size();
int* buffer = new int[nSendInt * 2];
for (int i = 0; i < nSendInt; i++) {
- buffer[i * 2] = (rankIt->second)[i].rankObject.elIndex;
- buffer[i * 2 + 1] = (rankIt->second)[i].rankObject.ithObjAtBoundary;
+ buffer[i * 2] = (rankIt->second)[i].rankObj.elIndex;
+ buffer[i * 2 + 1] = (rankIt->second)[i].rankObj.ithObj;
}
sendBuffers.push_back(buffer);
@@ -1157,13 +1158,13 @@ namespace AMDiS {
if (rankIt->first > mpiRank) {
for (int j = 0; j < static_cast<int>(rankIt->second.size()); j++)
- if (periodicBoundary.boundary[rankIt->first][j].neighbourObject.elIndex != recvBuffers[i][j * 2] ||
- periodicBoundary.boundary[rankIt->first][j].neighbourObject.ithObjAtBoundary != recvBuffers[i][j * 2 + 1]) {
+ if (periodicBoundary.boundary[rankIt->first][j].neighObj.elIndex != recvBuffers[i][j * 2] ||
+ periodicBoundary.boundary[rankIt->first][j].neighObj.ithObj != recvBuffers[i][j * 2 + 1]) {
int k = j + 1;
for (; k < static_cast<int>(rankIt->second.size()); k++)
- if (periodicBoundary.boundary[rankIt->first][k].neighbourObject.elIndex == recvBuffers[i][j * 2] &&
- periodicBoundary.boundary[rankIt->first][k].neighbourObject.ithObjAtBoundary == recvBuffers[i][j * 2 + 1])
+ if (periodicBoundary.boundary[rankIt->first][k].neighObj.elIndex == recvBuffers[i][j * 2] &&
+ periodicBoundary.boundary[rankIt->first][k].neighObj.ithObj == recvBuffers[i][j * 2 + 1])
break;
// The element must always be found, because the list is just in another order.
@@ -1505,10 +1506,9 @@ namespace AMDiS {
for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
DofContainer dofs;
- addVertexDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs);
- addEdgeDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs);
+ boundIt->rankObj.el->getVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs);
+ boundIt->rankObj.el->getNonVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs);
+
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");
@@ -1528,10 +1528,8 @@ namespace AMDiS {
boundIt != it->second.end(); ++boundIt) {
DofContainer dofs;
- addEdgeDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs);
- addVertexDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs);
+ boundIt->rankObj.el->getNonVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs);
+ boundIt->rankObj.el->getVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs);
for (int i = static_cast<int>(dofs.size()) - 1; i >= 0; i--) {
TEST_EXIT_DBG(find(dofsToRecv.begin(), dofsToRecv.end(), dofs[i]) == dofsToRecv.end())
@@ -1666,121 +1664,7 @@ namespace AMDiS {
mapLocalToDofIndex[dofIt->second] = *(dofIt->first);
}
- void ParallelDomainBase::addVertexDofs(Element *el, int ithEdge, DofContainer& dofs,
- bool parentVertices)
- {
- FUNCNAME("ParallelDomainBase::addVertexDofs()");
-
- TEST_EXIT_DBG(mesh->getDim() == 2)("Not yet implemented for 3D!\n");
-
- if (parentVertices) {
- switch (ithEdge) {
- case 0:
- dofs.push_back(el->getDOF(1));
- break;
- case 1:
- dofs.push_back(el->getDOF(0));
- break;
- case 2:
- dofs.push_back(el->getDOF(0));
- break;
- default:
- ERROR_EXIT("Should never happen!\n");
- }
- }
-
- switch (ithEdge) {
- case 0:
- if (el->getSecondChild() && el->getSecondChild()->getFirstChild()) {
- addVertexDofs(el->getSecondChild()->getFirstChild(), 0, dofs);
- dofs.push_back(el->getSecondChild()->getFirstChild()->getDOF(2));
- addVertexDofs(el->getSecondChild()->getSecondChild(), 1, dofs);
- }
- break;
- case 1:
- if (el->getFirstChild() && el->getFirstChild()->getFirstChild()) {
- addVertexDofs(el->getFirstChild()->getFirstChild(), 0, dofs);
- dofs.push_back(el->getFirstChild()->getFirstChild()->getDOF(2));
- addVertexDofs(el->getFirstChild()->getSecondChild(), 1, dofs);
- }
- break;
- case 2:
- if (el->getFirstChild()) {
- addVertexDofs(el->getFirstChild(), 0, dofs);
- dofs.push_back(el->getFirstChild()->getDOF(2));
- addVertexDofs(el->getSecondChild(), 1, dofs);
- }
- break;
- default:
- ERROR_EXIT("Should never happen!\n");
- }
-
- if (parentVertices) {
- switch (ithEdge) {
- case 0:
- dofs.push_back(el->getDOF(2));
- break;
- case 1:
- dofs.push_back(el->getDOF(2));
- break;
- case 2:
- dofs.push_back(el->getDOF(1));
- break;
- default:
- ERROR_EXIT("Should never happen!\n");
- }
- }
- }
-
-
- void ParallelDomainBase::addEdgeDofs(Element *el, int ithEdge, DofContainer& dofs)
- {
- FUNCNAME("ParallelDomainBase::addEdgeDofs()");
-
- TEST_EXIT_DBG(mesh->getDim() == 2)("Not yet implemented for 3D!\n");
-
- bool addThisEdge = false;
-
- switch (ithEdge) {
- case 0:
- if (el->getSecondChild())
- addEdgeDofs(el->getSecondChild(), 2, dofs);
- else
- addThisEdge = true;
-
- break;
- case 1:
- if (el->getFirstChild())
- addEdgeDofs(el->getFirstChild(), 2, dofs);
- else
- addThisEdge = true;
-
- break;
- case 2:
- if (el->getFirstChild()) {
- addEdgeDofs(el->getFirstChild(), 0, dofs);
- addEdgeDofs(el->getSecondChild(), 1, dofs);
- } else {
- addThisEdge = true;
- }
-
- break;
- default:
- ERROR_EXIT("Should never happen!\n");
- }
-
- if (addThisEdge) {
- ElementDofIterator elDofIter(feSpace, true);
- elDofIter.reset(el);
- do {
- if (elDofIter.getCurrentPos() == 1 &&
- elDofIter.getCurrentElementPos() == ithEdge)
- dofs.push_back(elDofIter.getDofPtr());
- } while(elDofIter.next());
- }
- }
-
-
+
void ParallelDomainBase::createDofMemberInfo(DofToPartitions& partitionDofs,
DofContainer& rankOwnedDofs,
DofContainer& rankAllDofs,
@@ -1876,7 +1760,7 @@ namespace AMDiS {
// Clear all periodic dof mappings calculated before. We do it from scratch.
periodicDof.clear();
- MPI::Request request[periodicBoundary.boundary.size() * 2];
+ MPI::Request request[min(static_cast<int>(periodicBoundary.boundary.size() * 2), 4)];
int requestCounter = 0;
std::vector<int*> sendBuffers, recvBuffers;
@@ -1895,10 +1779,8 @@ namespace AMDiS {
for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
- addVertexDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs, true);
- addEdgeDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- dofs);
+ boundIt->rankObj.el->getVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs, true);
+ boundIt->rankObj.el->getNonVertexDofs(feSpace, boundIt->rankObj.ithObj, dofs);
}
// Send the global indices to the rank on the other side.
@@ -1922,14 +1804,16 @@ namespace AMDiS {
}
- MPI::Request::Waitall(requestCounter, request);
-
+ MPI::Request::Waitall(requestCounter, request);
// === The rank has received the dofs from the rank on the other side of ===
// === the boundary. Now it can use them to create the mapping between ===
// === the periodic dofs in this rank and the corresponding periodic ===
// === dofs from the other ranks. ===
+
+ std::map<DegreeOfFreedom, std::set<int> > dofFromRank;
+
int i = 0;
for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
it != periodicBoundary.boundary.end(); ++it) {
@@ -1939,39 +1823,74 @@ namespace AMDiS {
for (std::vector<AtomicBoundary>::iterator boundIt = it->second.begin();
boundIt != it->second.end(); ++boundIt) {
DofContainer tmpdofs;
- addEdgeDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- tmpdofs);
- addVertexDofs(boundIt->rankObject.el, boundIt->rankObject.ithObjAtBoundary,
- tmpdofs, true);
+ boundIt->rankObj.el->getNonVertexDofs(feSpace, boundIt->rankObj.ithObj, tmpdofs);
+ boundIt->rankObj.el->getVertexDofs(feSpace, boundIt->rankObj.ithObj, tmpdofs, true);
for (int j = static_cast<int>(tmpdofs.size()) - 1; j >= 0; j--)
dofs.push_back(tmpdofs[j]);
}
// Added the received dofs to the mapping.
- for (int j = 0; j < static_cast<int>(dofs.size()); j++)
- periodicDof[mapLocalGlobalDOFs[*(dofs[j])]].insert(recvBuffers[i][j]);
+ for (int j = 0; j < static_cast<int>(dofs.size()); j++) {
+ int globalDofIndex = mapLocalGlobalDOFs[*(dofs[j])];
+ periodicDof[globalDofIndex].insert(recvBuffers[i][j]);
+ dofFromRank[globalDofIndex].insert(it->first);
+ }
delete [] sendBuffers[i];
delete [] recvBuffers[i];
i++;
}
- // TODO search for 2 dof mappings.
-
- switch (mpiRank) {
- case 0:
- periodicDof[0].insert(3136);
- break;
- case 1:
- periodicDof[1024].insert(2080);
- break;
- case 2:
- periodicDof[2080].insert(1024);
- break;
- case 3:
- periodicDof[3136].insert(0);
- break;
- }
+ sendBuffers.clear();
+ recvBuffers.clear();
+
+ TEST_EXIT_DBG(mesh->getDim() == 2)
+ ("Periodic boundary corner problem must be generalized to 3d!\n");
+
+ requestCounter = 0;
+ for (std::map<DegreeOfFreedom, std::set<int> >::iterator it = dofFromRank.begin();
+ it != dofFromRank.end(); ++it) {
+ if (it->second.size() == 2) {
+ TEST_EXIT_DBG(periodicDof[it->first].size() == 2)("Missing periodic dof!\n");
+
+ int *sendbuf = new int[2];
+ sendbuf[0] = *(periodicDof[it->first].begin());
+ sendbuf[1] = *(++(periodicDof[it->first].begin()));
+
+ request[requestCounter++] =
+ mpiComm.Isend(sendbuf, 2, MPI_INT, *(it->second.begin()), 0);
+ request[requestCounter++] =
+ mpiComm.Isend(sendbuf, 2, MPI_INT, *(++(it->second.begin())), 0);
+
+ sendBuffers.push_back(sendbuf);
+
+ int *recvbuf1 = new int[2];
+ int *recvbuf2 = new int[2];
+
+ request[requestCounter++] =
+ mpiComm.Irecv(recvbuf1, 2, MPI_INT, *(it->second.begin()), 0);
+ request[requestCounter++] =
+ mpiComm.Irecv(recvbuf2, 2, MPI_INT, *(++(it->second.begin())), 0);
+
+ recvBuffers.push_back(recvbuf1);
+ recvBuffers.push_back(recvbuf2);
+ }
+ }
+
+ MPI::Request::Waitall(requestCounter, request);
+
+ i = 0;
+ for (std::map<DegreeOfFreedom, std::set<int> >::iterator it = dofFromRank.begin();
+ it != dofFromRank.end(); ++it) {
+ if (it->second.size() == 2) {
+ for (int k = 0; k < 2; k++)
+ for (int j = 0; j < 2; j++)
+ if (recvBuffers[i + k][j] != it->first)
+ periodicDof[it->first].insert(recvBuffers[i + k][j]);
+
+ i++;
+ }
+ }
}
@@ -2133,7 +2052,7 @@ namespace AMDiS {
int nSendInt = rankIt->second.size();
int* buffer = new int[nSendInt];
for (int i = 0; i < nSendInt; i++)
- buffer[i] = (rankIt->second)[i].rankObject.elIndex;
+ buffer[i] = (rankIt->second)[i].rankObj.elIndex;
sendBuffers.push_back(buffer);
@@ -2166,7 +2085,7 @@ namespace AMDiS {
for (int i = 0; i < static_cast<int>(rankIt->second.size()); i++) {
int elIndex1 = recvBuffers[bufCounter][i];
- int elIndex2 = otherIntBoundary.boundary[rankIt->first][i].neighbourObject.elIndex;
+ int elIndex2 = otherIntBoundary.boundary[rankIt->first][i].neighObj.elIndex;
TEST_EXIT(elIndex1 == elIndex2)("Wrong element index at interior boundary!\n");
}
@@ -2371,6 +2290,8 @@ namespace AMDiS {
void ParallelDomainBase::printMapPeriodic(int rank)
{
+ FUNCNAME("ParallelDomainBase::printMapPeriodic()");
+
if (rank == -1 || mpiRank == rank) {
std::cout << "====== DOF MAP PERIODIC ====== " << std::endl;
@@ -2382,12 +2303,14 @@ namespace AMDiS {
std::cout << *dofit << " ";
std::cout << std::endl;
- DegreeOfFreedom localdof;
+ DegreeOfFreedom localdof = -1;
for (DofMapping::iterator dofIt = mapLocalGlobalDOFs.begin();
dofIt != mapLocalGlobalDOFs.end(); ++dofIt)
if (dofIt->second == it->first)
localdof = dofIt->first;
+ TEST_EXIT(localdof != -1)("There is something wrong!\n");
+
WorldVector<double> coords;
mesh->getDofIndexCoords(localdof, feSpace, coords);
coords.print();
diff --git a/AMDiS/src/ParallelDomainBase.h b/AMDiS/src/ParallelDomainBase.h
index 84e0c198..06bfd814 100644
--- a/AMDiS/src/ParallelDomainBase.h
+++ b/AMDiS/src/ParallelDomainBase.h
@@ -50,9 +50,6 @@ namespace AMDiS {
public ProblemTimeInterface
{
private:
- /// Defines type for a vector of DOFs.
- typedef std::vector<const DegreeOfFreedom*> DofContainer;
-
/// Defines a mapping type from DOFs to rank numbers.
typedef std::map<const DegreeOfFreedom*, int> DofToRank;
@@ -238,31 +235,6 @@ namespace AMDiS {
DofIndexMap &rankOwnedDofsNewLocalIndex,
DofIndexMap &rankDofsNewGlobalIndex);
- /** \brief
- * Traverses an edge on a given element and all its child elements. All vertex
- * dofs alonge this edge are assembled and put together to a list.
- *
- * \param[in] el Element for the edge traverse.
- * \param[in] ithEdge Defines the edge on which all the vertex dofs
- * are assembled.
- * \param[out] dofs All dofs are put to this dof list.
- * \param[in] parentVertices If true, also the two vertices of the parent
- * element el are put into the result list.
- */
- void addVertexDofs(Element *el, int ithEdge, DofContainer& dofs,
- bool parentVertices = false);
-
- /** \brief
- * Traverses an edge on a given element and all its child elements. All non vertex
- * dofs alonge this edge are assembled and put together to a list.
- *
- * \param[in] el Element for the edge traverse.
- * \param[in] ithEdge Defines the edge on which all the non vertex
- * dofs are assembled.
- * \param[out] dofs All dofs are put to this dof list.
- */
- void addEdgeDofs(Element *el, int ithEdge, DofContainer& dofs);
-
/** \brief
* This function traverses the whole mesh, i.e. before it is really partitioned,
* and collects information about which DOF corresponds to which rank. Can only
diff --git a/AMDiS/src/Tetrahedron.h b/AMDiS/src/Tetrahedron.h
index a9e3195b..8a735399 100644
--- a/AMDiS/src/Tetrahedron.h
+++ b/AMDiS/src/Tetrahedron.h
@@ -105,25 +105,19 @@ namespace AMDiS {
const FixVec<int,WORLD>& sortFaceIndices(int face,
FixVec<int,WORLD> *vec) const;
- /// implements Element::isLine. Returns false because this element is a Tetrahedron
+ /// Returns false because this element is a Tetrahedron.
inline bool isLine() const
{
return false;
}
- /** \brief
- * implements Element::isTriangle. Returns false because this element is a
- * Tetrahedron
- */
+ /// Returns false because this element is a Tetrahedron.
inline bool isTriangle() const
{
return false;
}
- /** \brief
- * implements Element::isTetrahedron. Returns true because this element is a
- * Tetrahedron
- */
+ /// Returns true because this element is a Tetrahedron.
inline bool isTetrahedron() const
{
return true;
@@ -134,6 +128,22 @@ namespace AMDiS {
return "Tetrahedron";
}
+ void getVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs, bool parentVertices = 0) const
+ {
+ FUNCNAME("Line::getVertexDofs()");
+
+ ERROR_EXIT("Not yet implemented!\n");
+ }
+
+ void getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs) const
+ {
+ FUNCNAME("Line::getNonVertexDofs()");
+
+ ERROR_EXIT("Not yet implemented!\n");
+ }
+
public:
/** \brief
* nChildEdge[el_type][ichild][dir]
diff --git a/AMDiS/src/Triangle.cc b/AMDiS/src/Triangle.cc
index 1ef94d5d..5731f4b6 100644
--- a/AMDiS/src/Triangle.cc
+++ b/AMDiS/src/Triangle.cc
@@ -3,6 +3,7 @@
#include "Mesh.h"
#include "CoarseningManager.h"
#include "FixVec.h"
+#include "ElementDofIterator.h"
namespace AMDiS {
@@ -18,9 +19,8 @@ namespace AMDiS {
return false;
}
- int Triangle::getVertexOfPosition(GeoIndex position,
- int positionIndex,
- int vertexIndex) const
+ int Triangle::getVertexOfPosition(GeoIndex position, int positionIndex,
+ int vertexIndex) const
{
FUNCNAME("Triangle::getVertexOfPosition");
switch(position) {
@@ -39,8 +39,8 @@ namespace AMDiS {
}
}
- const FixVec<int, WORLD>&
- Triangle::sortFaceIndices(int face, FixVec<int,WORLD> *vec) const
+ const FixVec<int, WORLD>& Triangle::sortFaceIndices(int face,
+ FixVec<int,WORLD> *vec) const
{
static MatrixOfFixVecs<FixVec<int,WORLD> > *sorted_2d = NULL;
@@ -75,4 +75,123 @@ namespace AMDiS {
return *(const_cast<const FixVec<int,WORLD>* >(val));
}
+
+ void Triangle::getVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs, bool parentVertices) const
+ {
+ FUNCNAME("Triangle::getVertexDofs()");
+
+ if (parentVertices) {
+ switch (ithEdge) {
+ case 0:
+ dofs.push_back(dof[1]);
+ break;
+ case 1:
+ dofs.push_back(dof[0]);
+ break;
+ case 2:
+ dofs.push_back(dof[0]);
+ break;
+ default:
+ ERROR_EXIT("Should never happen!\n");
+ }
+ }
+
+ switch (ithEdge) {
+ case 0:
+ {
+ Element *child1 = child[1];
+ if (child1 && child1->getFirstChild()) {
+ child1->getFirstChild()->getVertexDofs(feSpace, 0, dofs);
+ dofs.push_back(child1->getFirstChild()->getDOF(2));
+ child1->getSecondChild()->getVertexDofs(feSpace, 1, dofs);
+ }
+ }
+ break;
+ case 1:
+ {
+ Element *child0 = child[0];
+ if (child0 && child0->getFirstChild()) {
+ child0->getFirstChild()->getVertexDofs(feSpace, 0, dofs);
+ dofs.push_back(child[0]->getFirstChild()->getDOF(2));
+ child0->getSecondChild()->getVertexDofs(feSpace, 1, dofs);
+ }
+ }
+ break;
+ case 2:
+ if (child[0]) {
+ child[0]->getVertexDofs(feSpace, 0, dofs);
+ dofs.push_back(child[0]->getDOF(2));
+ child[1]->getVertexDofs(feSpace, 1, dofs);
+ }
+ break;
+ default:
+ ERROR_EXIT("Should never happen!\n");
+ }
+
+ if (parentVertices) {
+ switch (ithEdge) {
+ case 0:
+ dofs.push_back(dof[2]);
+ break;
+ case 1:
+ dofs.push_back(dof[2]);
+ break;
+ case 2:
+ dofs.push_back(dof[1]);
+ break;
+ default:
+ ERROR_EXIT("Should never happen!\n");
+ }
+ }
+ }
+
+
+ void Triangle::getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs) const
+ {
+ FUNCNAME("Triange::getNonVertexDofs()");
+
+ bool addThisEdge = false;
+
+ switch (ithEdge) {
+ case 0:
+ if (child[1])
+ child[1]->getNonVertexDofs(feSpace, 2, dofs);
+ else
+ addThisEdge = true;
+
+ break;
+ case 1:
+ if (child[0])
+ child[0]->getNonVertexDofs(feSpace, 2, dofs);
+ else
+ addThisEdge = true;
+
+ break;
+ case 2:
+ if (child[0]) {
+ child[0]->getNonVertexDofs(feSpace, 0, dofs);
+ child[1]->getNonVertexDofs(feSpace, 1, dofs);
+ } else {
+ addThisEdge = true;
+ }
+
+ break;
+ default:
+ ERROR_EXIT("Should never happen!\n");
+ }
+
+ if (addThisEdge) {
+ ElementDofIterator elDofIter(feSpace, true);
+ elDofIter.reset(this);
+ do {
+ if (elDofIter.getCurrentPos() == 1 &&
+ elDofIter.getCurrentElementPos() == ithEdge)
+ dofs.push_back(elDofIter.getDofPtr());
+ } while(elDofIter.next());
+ }
+ }
+
+
}
diff --git a/AMDiS/src/Triangle.h b/AMDiS/src/Triangle.h
index 9a6c3c25..f149bd3d 100644
--- a/AMDiS/src/Triangle.h
+++ b/AMDiS/src/Triangle.h
@@ -150,6 +150,12 @@ namespace AMDiS {
return "Triangle";
}
+ void getVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs, bool parentVertices = 0) const;
+
+ void getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge,
+ DofContainer& dofs) const;
+
protected:
/** \brief
* vertexOfEdge[i][j] is the local number of the j-th vertex of the i-th
--
GitLab