Commit 32706a6b authored by Thomas Witkowski's avatar Thomas Witkowski

Solved edge interior boundary problem in domain parallelization.

parent 3119b7a7
......@@ -392,9 +392,8 @@ namespace AMDiS {
* element are put into the result list.
*/
virtual void getVertexDofs(FiniteElemSpace* feSpace,
int ith, int elType,
BoundaryObject bound,
DofContainer& dofs,
bool reverseMode,
bool parentVertices = false) const = 0;
/** \brief
......@@ -410,7 +409,7 @@ namespace AMDiS {
* \param[out] dofs All dofs are put to this dof list.
*/
virtual void getNonVertexDofs(FiniteElemSpace* feSpace,
int ith, int elType,
BoundaryObject bound,
DofContainer& dofs) const = 0;
......
......@@ -78,6 +78,9 @@ namespace AMDiS {
/// Defines type for a vector of DOF pointers.
typedef std::vector<const DegreeOfFreedom*> DofContainer;
/// Defines a type for global edge identification.
typedef std::pair<DegreeOfFreedom, DegreeOfFreedom> GlobalEdge;
/// returns the GeoIndex of d for dimension dim.
#define INDEX_OF_DIM(d, dim) (static_cast<GeoIndex>((d == dim) ? CENTER : d + 1))
......
......@@ -36,6 +36,22 @@ namespace AMDiS {
}
bool BoundaryObject::operator==(const BoundaryObject& other) const
{
return (other.elIndex == elIndex &&
other.subObj == subObj &&
other.ithObj == ithObj);
}
bool BoundaryObject::operator!=(const BoundaryObject& other) const
{
return (other.elIndex != elIndex ||
other.subObj != subObj ||
other.ithObj != ithObj);
}
AtomicBoundary& InteriorBoundary::getNewAtomic(int rank)
{
boundary[rank].resize(boundary[rank].size() + 1);
......
......@@ -25,17 +25,36 @@
#include <vector>
#include <map>
#include "MacroElement.h"
#include "AMDiS_fwd.h"
#include "MacroElement.h"
#include "Element.h"
namespace AMDiS {
/// Defines the geometrical objects that forms the boundary;
struct BoundaryObject {
BoundaryObject()
: reverseMode(false)
: reverseMode(false),
notIncludedSubStructures(0)
{}
BoundaryObject(Element *e, int eType, GeoIndex sObj, int iObj)
: el(e),
elIndex(e->getIndex()),
elType(eType),
subObj(sObj),
ithObj(iObj),
reverseMode(false),
notIncludedSubStructures(0)
{}
void setReverseMode(BoundaryObject &otherBound, FiniteElemSpace *feSpace);
bool operator==(const BoundaryObject& other) const;
bool operator!=(const BoundaryObject& other) const;
/// The macro element to which the boundary element corresponds to.
Element* el;
......@@ -64,7 +83,7 @@ namespace AMDiS {
bool reverseMode;
void setReverseMode(BoundaryObject &otherBound, FiniteElemSpace *feSpace);
std::vector<int> notIncludedSubStructures;
};
/** \brief
......
......@@ -23,6 +23,7 @@
#define AMDIS_LINE_H
#include "Element.h"
#include "InteriorBoundary.h"
namespace AMDiS {
......@@ -50,12 +51,12 @@ namespace AMDiS {
/// implements Element::getVertexOfPosition
virtual int getVertexOfPosition(GeoIndex position,
int positionIndex,
int vertexIndex) const;
int positionIndex,
int vertexIndex) const;
virtual int getPositionOfVertex(int side, int vertex) const
{
static int positionOfVertex[2][2] = {{0,-1},{-1,0}};
static int positionOfVertex[2][2] = {{0, -1}, {-1, 0}};
return positionOfVertex[side][vertex];
}
......@@ -153,21 +154,16 @@ namespace AMDiS {
return "Line";
}
void getVertexDofs(FiniteElemSpace* feSpace, int ithEdge, int elType,
DofContainer& dofs,
bool reverseMode,
void getVertexDofs(FiniteElemSpace*, BoundaryObject, DofContainer&,
bool parentVertices = false) const
{
FUNCNAME("Line::getVertexDofs()");
ERROR_EXIT("Not yet implemented!\n");
}
void getNonVertexDofs(FiniteElemSpace* feSpace, int ithEdge, int elType,
DofContainer& dofs) const
void getNonVertexDofs(FiniteElemSpace*, BoundaryObject, DofContainer&) const
{
FUNCNAME("Line::getNonVertexDofs()");
ERROR_EXIT("Not yet implemented!\n");
}
......
This diff is collapsed.
......@@ -262,7 +262,6 @@ namespace AMDiS {
out << "\t\t" << nVertices << "," << endl;
// For all DOFs of vertices, write the coordinates.
int counter = 0;
for (it.reset(); !it.end(); ++it) {
// for all vertex infos of this DOF
for (std::list<VertexInfo>::iterator it2 = it->begin(); it2 != it->end(); ++it2) {
......@@ -308,7 +307,6 @@ namespace AMDiS {
out << "\t\t" << nTextures << "," << endl;
// For all DOFs of vertices, write the coordinates.
int counter = 0;
for (it.reset(); !it.end(); ++it) {
// for all vertex infos of this DOF
for (std::list<VertexInfo>::iterator it2 = it->begin(); it2 != it->end(); ++it2) {
......
......@@ -25,6 +25,7 @@
#include <map>
#include "mpi.h"
#include "MeshStructure.h"
#include "InteriorBoundary.h"
namespace AMDiS {
......@@ -52,6 +53,21 @@ namespace AMDiS {
return data.size() * 2;
}
int intSizeOf(std::vector<AtomicBoundary> &data)
{
return data.size() * 3;
}
int intSizeOf(std::vector<BoundaryObject> &data)
{
return data.size() * 5;
}
int intSizeOf(std::vector<const DegreeOfFreedom*> &data)
{
return data.size();
}
void makeBuf(std::map<const DegreeOfFreedom*, DegreeOfFreedom> &data, int* buf)
{
int i = 0;
......@@ -138,22 +154,75 @@ namespace AMDiS {
}
}
void makeBuf(std::vector<AtomicBoundary> &data, int *buf)
{
for (unsigned int i = 0; i < data.size(); i++) {
buf[i * 3] = data[i].rankObj.elIndex;
buf[i * 3 + 1] = data[i].rankObj.subObj;
buf[i * 3 + 2] = data[i].rankObj.ithObj;
}
}
void makeFromBuf(std::vector<AtomicBoundary> &data, int *buf, int bufSize)
{
if (bufSize == 0)
return;
TEST_EXIT(bufSize % 3 == 0)("This should not happen!\n");
data.resize(bufSize / 3);
for (int i = 0; i < bufSize / 3; i++) {
data[i].rankObj.elIndex = buf[i * 3];
data[i].rankObj.subObj = static_cast<GeoIndex>(buf[i * 3 + 1]);
data[i].rankObj.ithObj = buf[i * 3 + 2];
}
}
void makeBuf(std::vector<BoundaryObject> &data, int *buf)
{
for (unsigned int i = 0; i < data.size(); i++) {
buf[i * 5] = data[i].elIndex;
buf[i * 5 + 1] = data[i].elType;
buf[i * 5 + 2] = data[i].subObj;
buf[i * 5 + 3] = data[i].ithObj;
buf[i * 5 + 4] = data[i].reverseMode;
}
}
void makeFromBuf(std::vector<BoundaryObject> &data, int *buf, int bufSize)
{
if (bufSize == 0)
return;
TEST_EXIT(bufSize % 5 == 0)("This should not happen!\n");
data.resize(bufSize / 5);
for (int i = 0; i < bufSize / 5; i++) {
data[i].elIndex = buf[i * 5];
data[i].elType = buf[i * 5 + 1];
data[i].subObj = static_cast<GeoIndex>(buf[i * 5 + 2]);
data[i].ithObj = buf[i * 5 + 3];
data[i].reverseMode = static_cast<bool>(buf[i * 5 + 4]);
}
}
template<typename SendT, typename RecvT>
class StdMpi
{
public:
StdMpi(MPI::Intracomm &comm) :
StdMpi(MPI::Intracomm &comm, bool b = true) :
mpiComm(comm),
commPrepared(false),
exchangeDataSize(true)
commPrepared(true),
exchangeDataSize(b)
{}
void prepareCommunication(bool b)
void clear()
{
sendData.clear();
recvData.clear();
exchangeDataSize = b;
sendDataSize.clear();
recvDataSize.clear();
commPrepared = true;
}
......@@ -180,6 +249,16 @@ namespace AMDiS {
}
}
std::map<int, SendT>& getSendData()
{
return sendData;
}
SendT& getSendData(int rank)
{
return sendData[rank];
}
void recv(int fromRank, int size = -1)
{
FUNCNAME("StdMpi::recv()");
......@@ -189,6 +268,16 @@ namespace AMDiS {
recvDataSize[fromRank] = size;
}
void recvFromAll(int size = -1)
{
FUNCNAME("StdMpi::recvFromAll()");
TEST_EXIT_DBG(commPrepared)("Communication is no prepared!\n");
for (int i = 0; i < mpiComm.Get_size(); i++)
recvDataSize[i] = size;
}
template<class T>
void recv(std::map<int, T> &fromRanks)
{
......@@ -198,7 +287,7 @@ namespace AMDiS {
for (typename std::map<int, T>::iterator it = fromRanks.begin();
it != fromRanks.end(); ++it)
recvDataSize[it->first] = -1;
recvDataSize[it->first] = (exchangeDataSize ? -1 : intSizeOf(it->second));
}
std::map<int, RecvT>& getRecvData()
......@@ -206,6 +295,11 @@ namespace AMDiS {
return recvData;
}
RecvT& getRecvData(int rank)
{
return recvData[rank];
}
void commDataSize()
{
MPI::Request request[sendData.size() + recvDataSize.size()];
......@@ -274,6 +368,7 @@ namespace AMDiS {
int i = 0;
for (std::map<int, int>::iterator recvIt = recvDataSize.begin();
recvIt != recvDataSize.end(); ++recvIt) {
makeFromBuf(recvData[recvIt->first], recvBuffers[i], recvIt->second);
delete [] recvBuffers[i];
i++;
......@@ -282,16 +377,6 @@ namespace AMDiS {
commPrepared = false;
}
void clear()
{
sendData.clear();
recvData.clear();
sendDataSize.clear();
recvDataSize.clear();
commPrepared = true;
}
protected:
///
MPI::Intracomm mpiComm;
......
......@@ -10,15 +10,19 @@ namespace AMDiS {
const unsigned char Tetrahedron::nChildEdge[3][2][2] = {{{5,4},{4,5}},
{{5,4},{5,4}},
{{5,4},{5,4}}};
const unsigned char Tetrahedron::nChildFace[3][2][2] = {{{1,2},{2,1}},
{{1,2},{1,2}},
{{1,2},{1,2}}};
const int Tetrahedron::childVertex[3][2][4] = {{{0,2,3,4},{1,3,2,4}},
{{0,2,3,4},{1,2,3,4}},
{{0,2,3,4},{1,2,3,4}}};
const int Tetrahedron::childEdge[3][2][6] = {{{1,2,0,5,5,4},{4,3,0,5,5,4}},
{{1,2,0,5,4,5},{3,4,0,5,4,5}},
{{1,2,0,5,4,5},{3,4,0,5,4,5}}};
const unsigned char Tetrahedron::adjacentChild[2][2] = {{0,1}, {1,0}};
const signed char Tetrahedron::childOrientation[3][2] = {{1,1},
......@@ -36,12 +40,12 @@ namespace AMDiS {
{1,2},
{1,3},
{2,3}};
const int Tetrahedron::vertexOfFace[4][3] = {{1,2,3},
{0,2,3},
{0,1,3},
{0,1,2}};
const int Tetrahedron::sideOfChild[3][2][4] = {{{-1, 3, 1, 2}, // type 0
{3, -1, 2, 1}},
{{-1, 3, 1, 2}, // type 1
......@@ -49,6 +53,13 @@ namespace AMDiS {
{{-1, 3, 1, 2}, // type 2
{3, -1, 1, 2}}};
const int Tetrahedron::edgeOfChild[3][2][6] = {{{2, 0, 1, -1, -1, 3}, // type 0
{2, -1, -1, 1, 0, 3}},
{{2, 0, 1, -1, -1, 3}, // type 1
{2, -1, -1, 0, 1, 3}},
{{2, 0, 1, -1, -1, 3}, // type 2
{2, -1, -1, 0, 1, 3}}};
const int Tetrahedron::vertexOfParent[3][2][4] = {{{0, 2, 3, -1}, // type 0
{1, 3, 2, -1}},
{{0, 2, 3, -1}, // type 1
......@@ -93,8 +104,8 @@ namespace AMDiS {
}
}
const FixVec<int, WORLD>&
Tetrahedron::sortFaceIndices(int face, FixVec<int,WORLD> *vec) const
const FixVec<int, WORLD>& Tetrahedron::sortFaceIndices(int face,
FixVec<int,WORLD> *vec) const
{
static MatrixOfFixVecs<FixVec<int,WORLD> > *sorted_3d = NULL;
FixVec<int,WORLD> *val = NULL;
......@@ -103,51 +114,51 @@ namespace AMDiS {
if (NULL == sorted_3d) {
sorted_3d = new MatrixOfFixVecs<FixVec<int,WORLD> >(3, 4, 7, NO_INIT);
(*sorted_3d)[0][0][0]=(*sorted_3d)[0][0][1]=
(*sorted_3d)[0][0][2]=(*sorted_3d)[1][0][0]=
(*sorted_3d)[1][0][1]=(*sorted_3d)[1][0][2]=
(*sorted_3d)[1][1][0]=(*sorted_3d)[1][2][1]=
(*sorted_3d)[1][3][0]=(*sorted_3d)[1][4][2]=
(*sorted_3d)[1][5][1]=(*sorted_3d)[1][6][2]=
(*sorted_3d)[2][0][0]=(*sorted_3d)[2][0][1]=
(*sorted_3d)[2][0][2]=(*sorted_3d)[2][1][0]=
(*sorted_3d)[2][2][1]=(*sorted_3d)[2][3][0]=
(*sorted_3d)[2][4][2]=(*sorted_3d)[2][5][1]=
(*sorted_3d)[2][6][2]=(*sorted_3d)[3][0][0]=
(*sorted_3d)[3][0][1]=(*sorted_3d)[3][0][2]=
(*sorted_3d)[3][1][0]=(*sorted_3d)[3][2][1]=
(*sorted_3d)[3][3][0]=(*sorted_3d)[3][4][2]=
(*sorted_3d)[3][5][1]=(*sorted_3d)[3][6][2]=0;
(*sorted_3d)[0][1][0]=(*sorted_3d)[0][2][1]=
(*sorted_3d)[0][3][0]=(*sorted_3d)[0][4][2]=
(*sorted_3d)[0][5][1]=(*sorted_3d)[0][6][2]=
(*sorted_3d)[2][1][2]=(*sorted_3d)[2][2][0]=
(*sorted_3d)[2][3][1]=(*sorted_3d)[2][4][1]=
(*sorted_3d)[2][5][2]=(*sorted_3d)[2][6][0]=
(*sorted_3d)[3][1][2]=(*sorted_3d)[3][2][0]=
(*sorted_3d)[3][3][1]=(*sorted_3d)[3][4][1]=
(*sorted_3d)[3][5][2]=(*sorted_3d)[3][6][0]=1;
(*sorted_3d)[0][1][2]=(*sorted_3d)[0][2][0]=
(*sorted_3d)[0][3][1]=(*sorted_3d)[0][4][1]=
(*sorted_3d)[0][5][2]=(*sorted_3d)[0][6][0]=
(*sorted_3d)[1][1][2]=(*sorted_3d)[1][2][0]=
(*sorted_3d)[1][3][1]=(*sorted_3d)[1][4][1]=
(*sorted_3d)[1][5][2]=(*sorted_3d)[1][6][0]=
(*sorted_3d)[3][1][1]=(*sorted_3d)[3][2][2]=
(*sorted_3d)[3][3][2]=(*sorted_3d)[3][4][0]=
(*sorted_3d)[3][5][0]=(*sorted_3d)[3][6][1]=2;
(*sorted_3d)[0][1][1]=(*sorted_3d)[0][2][2]=
(*sorted_3d)[0][3][2]=(*sorted_3d)[0][4][0]=
(*sorted_3d)[0][5][0]=(*sorted_3d)[0][6][1]=
(*sorted_3d)[1][1][1]=(*sorted_3d)[1][2][2]=
(*sorted_3d)[1][3][2]=(*sorted_3d)[1][4][0]=
(*sorted_3d)[1][5][0]=(*sorted_3d)[1][6][1]=
(*sorted_3d)[2][1][1]=(*sorted_3d)[2][2][2]=
(*sorted_3d)[2][3][2]=(*sorted_3d)[2][4][0]=
(*sorted_3d)[2][5][0]=(*sorted_3d)[2][6][1]=3;
(*sorted_3d)[0][0][0] = (*sorted_3d)[0][0][1] =
(*sorted_3d)[0][0][2] = (*sorted_3d)[1][0][0] =
(*sorted_3d)[1][0][1] = (*sorted_3d)[1][0][2] =
(*sorted_3d)[1][1][0] = (*sorted_3d)[1][2][1] =
(*sorted_3d)[1][3][0] = (*sorted_3d)[1][4][2] =
(*sorted_3d)[1][5][1] = (*sorted_3d)[1][6][2] =
(*sorted_3d)[2][0][0] = (*sorted_3d)[2][0][1] =
(*sorted_3d)[2][0][2] = (*sorted_3d)[2][1][0] =
(*sorted_3d)[2][2][1] = (*sorted_3d)[2][3][0] =
(*sorted_3d)[2][4][2] = (*sorted_3d)[2][5][1] =
(*sorted_3d)[2][6][2] = (*sorted_3d)[3][0][0] =
(*sorted_3d)[3][0][1] = (*sorted_3d)[3][0][2] =
(*sorted_3d)[3][1][0] = (*sorted_3d)[3][2][1] =
(*sorted_3d)[3][3][0] = (*sorted_3d)[3][4][2] =
(*sorted_3d)[3][5][1] = (*sorted_3d)[3][6][2] = 0;
(*sorted_3d)[0][1][0] = (*sorted_3d)[0][2][1] =
(*sorted_3d)[0][3][0] = (*sorted_3d)[0][4][2] =
(*sorted_3d)[0][5][1] = (*sorted_3d)[0][6][2] =
(*sorted_3d)[2][1][2] = (*sorted_3d)[2][2][0] =
(*sorted_3d)[2][3][1] = (*sorted_3d)[2][4][1] =
(*sorted_3d)[2][5][2] = (*sorted_3d)[2][6][0] =
(*sorted_3d)[3][1][2] = (*sorted_3d)[3][2][0] =
(*sorted_3d)[3][3][1] = (*sorted_3d)[3][4][1] =
(*sorted_3d)[3][5][2] = (*sorted_3d)[3][6][0] = 1;
(*sorted_3d)[0][1][2] = (*sorted_3d)[0][2][0] =
(*sorted_3d)[0][3][1] = (*sorted_3d)[0][4][1] =
(*sorted_3d)[0][5][2] = (*sorted_3d)[0][6][0] =
(*sorted_3d)[1][1][2] = (*sorted_3d)[1][2][0] =
(*sorted_3d)[1][3][1] = (*sorted_3d)[1][4][1] =
(*sorted_3d)[1][5][2] = (*sorted_3d)[1][6][0] =
(*sorted_3d)[3][1][1] = (*sorted_3d)[3][2][2] =
(*sorted_3d)[3][3][2] = (*sorted_3d)[3][4][0] =
(*sorted_3d)[3][5][0] = (*sorted_3d)[3][6][1] = 2;
(*sorted_3d)[0][1][1] = (*sorted_3d)[0][2][2] =
(*sorted_3d)[0][3][2] = (*sorted_3d)[0][4][0] =
(*sorted_3d)[0][5][0] = (*sorted_3d)[0][6][1] =
(*sorted_3d)[1][1][1] = (*sorted_3d)[1][2][2] =
(*sorted_3d)[1][3][2] = (*sorted_3d)[1][4][0] =
(*sorted_3d)[1][5][0] = (*sorted_3d)[1][6][1] =
(*sorted_3d)[2][1][1] = (*sorted_3d)[2][2][2] =
(*sorted_3d)[2][3][2] = (*sorted_3d)[2][4][0] =
(*sorted_3d)[2][5][0] = (*sorted_3d)[2][6][1] = 3;
}
int no = 0;
......@@ -159,7 +170,7 @@ namespace AMDiS {
no += 4;
if (!(no >= 1 && no <= 6))
ERROR_EXIT("can not sort face indices of element %d at face %d\n",
ERROR_EXIT("Cannot sort face indices of element %d at face %d\n",
getIndex(), face);
if (vec) {
......@@ -172,16 +183,37 @@ namespace AMDiS {
return(*(const_cast<const FixVec<int,WORLD>* >(val)));
}
void Tetrahedron::getVertexDofs(FiniteElemSpace* feSpace,
int ith, int elType,
BoundaryObject bound,
DofContainer& dofs,
bool reverseMode,
bool parentVertices) const
{
FUNCNAME("Tetrahedron::getVertexDofs()");
switch (bound.subObj) {
case EDGE:
// ERROR_EXIT("Not yet implemented: %d\n", bound.elIndex);
break;
case FACE:
getVertexDofsAtFace(feSpace, bound, dofs, parentVertices);
break;
default:
ERROR_EXIT("Should not happen!\n");
}
}
void Tetrahedron::getVertexDofsAtFace(FiniteElemSpace* feSpace,
BoundaryObject bound,
DofContainer& dofs,
bool parentVertices) const
{
if (parentVertices) {
switch (ith) {
ERROR_EXIT("Einbau notIncludedSubStructures!\n");
switch (bound.ithObj) {
case 0:
dofs.push_back(dof[1]); dofs.push_back(dof[2]); dofs.push_back(dof[3]);
break;
......@@ -199,91 +231,113 @@ namespace AMDiS {
}
}
switch (ith) {
if (!child[0])
return;
switch (bound.ithObj) {
case 0:
{
if (sideOfChild[elType][0][ith] != -1) {
std::cout << "ERROR 1 WITH elType = " << elType << "!\n";
exit(0);
}
if (child[1])
child[1]->getVertexDofs(feSpace,
sideOfChild[elType][1][ith],
(elType + 1) % 3,
dofs,
reverseMode);
BoundaryObject nextBound = bound;
prepareNextBound(nextBound, 1);
child[1]->getVertexDofs(feSpace, nextBound, dofs);
}
break;
case 1:
{