Commit 9a7c2981 authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

Reengineering of parallel datastructures towards an efficient parallel preconditioner.

parent 2c48672b
This diff is collapsed.
This diff is collapsed.
This source diff could not be displayed because it is too large. You can view the blob instead.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
...@@ -91,12 +91,12 @@ namespace AMDiS { ...@@ -91,12 +91,12 @@ namespace AMDiS {
Element *el = smallElInfo->getElement(); Element *el = smallElInfo->getElement();
lastVecEl = lastMatEl = NULL; lastVecEl = lastMatEl = NULL;
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
if (smallElInfo == colElInfo) if (smallElInfo == colElInfo)
initElement(smallElInfo); initElement(smallElInfo);
else else
initElement(smallElInfo, largeElInfo); initElement(smallElInfo, largeElInfo);
}
if (el != lastMatEl || !operat->isOptimized()) { if (el != lastMatEl || !operat->isOptimized()) {
if (rememberElMat) if (rememberElMat)
......
...@@ -197,6 +197,9 @@ namespace AMDiS { ...@@ -197,6 +197,9 @@ namespace AMDiS {
/// ///
virtual int getEdgeOfFace(int face, int edge) const = 0; virtual int getEdgeOfFace(int face, int edge) const = 0;
///
virtual GlobalEdge getEdge(int localEdgeIndex) const = 0;
/// Returns the number of parts of type i in this element /// Returns the number of parts of type i in this element
virtual int getGeo(GeoIndex i) const = 0; virtual int getGeo(GeoIndex i) const = 0;
......
...@@ -84,10 +84,16 @@ namespace AMDiS { ...@@ -84,10 +84,16 @@ namespace AMDiS {
inline int getEdgeOfFace(int /* face */, int /*edge*/ ) const inline int getEdgeOfFace(int /* face */, int /*edge*/ ) const
{ {
ERROR_EXIT("This does not work in 2D!\n"); ERROR_EXIT("This does not work in 1D!\n");
return 0; return 0;
} }
GlobalEdge getEdge(int localEdgeIndex) const
{
ERROR_EXIT("This does not work in 1D!\n");
return GlobalEdge();
}
/// implements Element::sortFaceIndices /// implements Element::sortFaceIndices
const FixVec<int,WORLD>& sortFaceIndices(int face, const FixVec<int,WORLD>& sortFaceIndices(int face,
FixVec<int,WORLD> *vec) const; FixVec<int,WORLD> *vec) const;
......
...@@ -850,12 +850,11 @@ namespace AMDiS { ...@@ -850,12 +850,11 @@ namespace AMDiS {
void MacroInfo::fillBoundaryInfo(Mesh *mesh) void MacroInfo::fillBoundaryInfo(Mesh *mesh)
{ {
int i, j, k, nv = mesh->getNumberOfVertices(); int i, nv = mesh->getNumberOfVertices();
std::deque<MacroElement*>::iterator melIt; std::deque<MacroElement*>::iterator melIt;
BoundaryType *bound = new BoundaryType[nv]; std::vector<BoundaryType> bound(nv);
int dim = mesh->getDim();
switch(dim) { switch (mesh->getDim()) {
case 1: case 1:
break; break;
case 2: case 2:
...@@ -865,17 +864,17 @@ namespace AMDiS { ...@@ -865,17 +864,17 @@ namespace AMDiS {
for (i = 0, melIt = mesh->firstMacroElement(); for (i = 0, melIt = mesh->firstMacroElement();
melIt != mesh->endOfMacroElements(); melIt != mesh->endOfMacroElements();
++melIt, i++) { ++melIt, i++) {
for (j = 0; j < mesh->getGeo(NEIGH); j++) { for (int j = 0; j < mesh->getGeo(NEIGH); j++) {
if ((*melIt)->getBoundary(j) != INTERIOR) { if ((*melIt)->getBoundary(j) != INTERIOR) {
if ((*melIt)->getBoundary(j) >= DIRICHLET) { if ((*melIt)->getBoundary(j) >= DIRICHLET) {
int j1 = mel_vertex[i][(j+1)%3]; int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j+2)%3]; int j2 = mel_vertex[i][(j + 2) % 3];
bound[j1] = max(bound[j1], (*melIt)->getBoundary(j)); bound[j1] = max(bound[j1], (*melIt)->getBoundary(j));
bound[j2] = max(bound[j2], (*melIt)->getBoundary(j)); bound[j2] = max(bound[j2], (*melIt)->getBoundary(j));
} else if ((*melIt)->getBoundary(j) <= NEUMANN) { } else if ((*melIt)->getBoundary(j) <= NEUMANN) {
int j1 = mel_vertex[i][(j+1)%3]; int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j+2)%3]; int j2 = mel_vertex[i][(j + 2) % 3];
if (bound[j1] != INTERIOR) if (bound[j1] != INTERIOR)
bound[j1] = max(bound[j1], (*melIt)->getBoundary(j)); bound[j1] = max(bound[j1], (*melIt)->getBoundary(j));
...@@ -894,7 +893,7 @@ namespace AMDiS { ...@@ -894,7 +893,7 @@ namespace AMDiS {
for (i = 0, melIt = mesh->firstMacroElement(); for (i = 0, melIt = mesh->firstMacroElement();
melIt != mesh->endOfMacroElements(); melIt != mesh->endOfMacroElements();
++melIt, i++) ++melIt, i++)
for (j = 0; j < mesh->getGeo(VERTEX); j++) for (int j = 0; j < mesh->getGeo(VERTEX); j++)
(*melIt)->setBoundary(3 + j, bound[mel_vertex[i][j]]); (*melIt)->setBoundary(3 + j, bound[mel_vertex[i][j]]);
break; break;
...@@ -905,28 +904,26 @@ namespace AMDiS { ...@@ -905,28 +904,26 @@ namespace AMDiS {
for (i = 0, melIt = mesh->firstMacroElement(); for (i = 0, melIt = mesh->firstMacroElement();
melIt != mesh->endOfMacroElements(); melIt != mesh->endOfMacroElements();
++melIt, i++) { ++melIt, i++) {
for (j = 0; j < mesh->getGeo(NEIGH); j++) { for (int j = 0; j < mesh->getGeo(NEIGH); j++) {
for (k = 1; k < 4; k++) for (int k = 1; k < 4; k++)
bound[mel_vertex[i][(j+k)%4]] = bound[mel_vertex[i][(j + k) % 4]] =
((*melIt)->getBoundary(j) != INTERIOR) ? ((*melIt)->getBoundary(j) != INTERIOR) ?
newBound((*melIt)->getBoundary(j), newBound((*melIt)->getBoundary(j),
bound[mel_vertex[i][(j+k)%4]]) : bound[mel_vertex[i][(j + k) % 4]]) :
bound[mel_vertex[i][(j+k)%4]]; bound[mel_vertex[i][(j + k) % 4]];
} }
} }
for (i = 0, melIt = mesh->firstMacroElement(); for (i = 0, melIt = mesh->firstMacroElement();
melIt != mesh->endOfMacroElements(); melIt != mesh->endOfMacroElements();
++melIt, i++) ++melIt, i++)
for (j = 0; j < mesh->getGeo(VERTEX); j++) for (int j = 0; j < mesh->getGeo(VERTEX); j++)
(*melIt)->setBoundary(10 + j, bound[mel_vertex[i][j]]); (*melIt)->setBoundary(10 + j, bound[mel_vertex[i][j]]);
break; break;
default: default:
ERROR_EXIT("invalid dim\n"); ERROR_EXIT("invalid dim\n");
} }
delete [] bound;
} }
void MacroReader::computeNeighbours(Mesh *mesh) void MacroReader::computeNeighbours(Mesh *mesh)
...@@ -1518,8 +1515,6 @@ namespace AMDiS { ...@@ -1518,8 +1515,6 @@ namespace AMDiS {
// (wird fuer ALBERT-Routine write_macro benoetigt) // (wird fuer ALBERT-Routine write_macro benoetigt)
// ele wird nicht benoetigt (es kann NULL uebergeben werden) // ele wird nicht benoetigt (es kann NULL uebergeben werden)
void MacroReader::umbVkantMacro(Mesh *mesh, MacroElement* me, int ka, int *) void MacroReader::umbVkantMacro(Mesh *mesh, MacroElement* me, int ka, int *)
{ {
MacroElement* macr=new MacroElement(mesh->getDim()); MacroElement* macr=new MacroElement(mesh->getDim());
......
...@@ -1161,6 +1161,7 @@ namespace AMDiS { ...@@ -1161,6 +1161,7 @@ namespace AMDiS {
for (it = periodicAssociations.begin(); it != end; ++it) for (it = periodicAssociations.begin(); it != end; ++it)
if ((*(it->second))[dof1] == dof2) if ((*(it->second))[dof1] == dof2)
return true; return true;
return false; return false;
} }
...@@ -1186,6 +1187,7 @@ namespace AMDiS { ...@@ -1186,6 +1187,7 @@ namespace AMDiS {
} }
} }
} }
return false; return false;
} }
......
...@@ -542,6 +542,25 @@ namespace AMDiS { ...@@ -542,6 +542,25 @@ namespace AMDiS {
return periodicAssociations; return periodicAssociations;
} }
/// Returns the periodic association for a specific boundary type.
inline VertexVector& getPeriodicAssociations(BoundaryType b)
{
FUNCNAME("Mesh::getPeriodicAssociations()");
TEST_EXIT_DBG(periodicAssociations.count(b) == 1)
("There are no periodic assoications for boundary type %d!\n", b);
return (*(periodicAssociations[b]));
}
/// Returns whether the given boundary type is periodic, i.e., if there is
/// a periodic association for this boundary type.
inline bool isPeriodicAssociation(BoundaryType b)
{
return (periodicAssociations.count(b) == 1 ? true : false);
}
/// ///
bool associated(DegreeOfFreedom dof1, DegreeOfFreedom dof2); bool associated(DegreeOfFreedom dof1, DegreeOfFreedom dof2);
......
...@@ -266,6 +266,18 @@ namespace AMDiS { ...@@ -266,6 +266,18 @@ namespace AMDiS {
return edgeOfFace[face][edge]; return edgeOfFace[face][edge];
} }
GlobalEdge getEdge(int localEdgeIndex) const
{
FUNCNAME("Tetrahedron::getEdge()");
TEST_EXIT_DBG(localEdgeIndex >= 0 && localEdgeIndex < 6)("invalid edge\n");
DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
GlobalEdge edge = std::make_pair(min(dof0, dof1), max(dof0, dof1));
return edge;
}
protected: protected:
/// vertexOfEdge[i][j] is the local number of the j-vertex of edge i /// vertexOfEdge[i][j] is the local number of the j-vertex of edge i
......
...@@ -170,6 +170,18 @@ namespace AMDiS { ...@@ -170,6 +170,18 @@ namespace AMDiS {
return edge; return edge;
} }
GlobalEdge getEdge(int localEdgeIndex) const
{
FUNCNAME("Triangle::getEdge()");
TEST_EXIT_DBG(localEdgeIndex >= 0 && localEdgeIndex < 3)("invalid edge\n");
DegreeOfFreedom dof0 = dof[vertexOfEdge[localEdgeIndex][0]][0];
DegreeOfFreedom dof1 = dof[vertexOfEdge[localEdgeIndex][1]][0];
GlobalEdge edge = std::make_pair(min(dof0, dof1), max(dof0, dof1));
return edge;
}
std::string getTypeName() const std::string getTypeName() const
{ {
return "Triangle"; return "Triangle";
......
...@@ -4,6 +4,7 @@ ...@@ -4,6 +4,7 @@
#include "DOFVector.h" #include "DOFVector.h"
#include "Debug.h" #include "Debug.h"
#include "SystemVector.h" #include "SystemVector.h"
#include "VtkWriter.h"
#include "petscksp.h" #include "petscksp.h"
......
...@@ -67,6 +67,7 @@ namespace AMDiS { ...@@ -67,6 +67,7 @@ namespace AMDiS {
bool BoundaryObject::operator==(const BoundaryObject& other) const bool BoundaryObject::operator==(const BoundaryObject& other) const
{ {
return (other.elIndex == elIndex && return (other.elIndex == elIndex &&
// other.elType == elType &&
other.subObj == subObj && other.subObj == subObj &&
other.ithObj == ithObj); other.ithObj == ithObj);
} }
...@@ -75,11 +76,20 @@ namespace AMDiS { ...@@ -75,11 +76,20 @@ namespace AMDiS {
bool BoundaryObject::operator!=(const BoundaryObject& other) const bool BoundaryObject::operator!=(const BoundaryObject& other) const
{ {
return (other.elIndex != elIndex || return (other.elIndex != elIndex ||
// other.elType != elType ||
other.subObj != subObj || other.subObj != subObj ||
other.ithObj != ithObj); other.ithObj != ithObj);
} }
bool AtomicBoundary::operator==(const AtomicBoundary& other) const
{
return (rankObj == other.rankObj &&
neighObj == other.neighObj &&
type == other.type);
}
AtomicBoundary& InteriorBoundary::getNewAtomic(int rank) AtomicBoundary& InteriorBoundary::getNewAtomic(int rank)
{ {
boundary[rank].resize(boundary[rank].size() + 1); boundary[rank].resize(boundary[rank].size() + 1);
...@@ -87,6 +97,32 @@ namespace AMDiS { ...@@ -87,6 +97,32 @@ namespace AMDiS {
} }
bool InteriorBoundary::operator==(const InteriorBoundary& other) const
{
InteriorBoundary& other2 = const_cast<InteriorBoundary&>(other);
for (RankToBoundMap::const_iterator it = boundary.begin();
it != boundary.end(); ++it) {
if (other2.boundary.count(it->first) == 0)
return false;
if (other2.boundary[it->first].size() != it->second.size())
return false;
for (unsigned int i = 0; i < it->second.size(); i++) {
std::vector<AtomicBoundary>::iterator bIt =
find(other2.boundary[it->first].begin(),
other2.boundary[it->first].end(), it->second[i]);
if (bIt == other2.boundary[it->first].end())
return false;
}
}
return true;
}
void InteriorBoundary::serialize(std::ostream &out) void InteriorBoundary::serialize(std::ostream &out)
{ {
int mSize = boundary.size(); int mSize = boundary.size();
......
...@@ -38,7 +38,8 @@ namespace AMDiS { ...@@ -38,7 +38,8 @@ namespace AMDiS {
struct BoundaryObject { struct BoundaryObject {
BoundaryObject() BoundaryObject()
: reverseMode(false), : elType(0),
reverseMode(false),
excludedSubstructures(0) excludedSubstructures(0)
{} {}
...@@ -102,6 +103,12 @@ namespace AMDiS { ...@@ -102,6 +103,12 @@ namespace AMDiS {
* boundary goes through. * boundary goes through.
*/ */
struct AtomicBoundary { struct AtomicBoundary {
AtomicBoundary()
: type(INTERIOR)
{}
bool operator==(const AtomicBoundary& other) const;
/// The rank's part of the boundary. /// The rank's part of the boundary.
BoundaryObject rankObj; BoundaryObject rankObj;
...@@ -217,6 +224,10 @@ namespace AMDiS { ...@@ -217,6 +224,10 @@ namespace AMDiS {
/// Reads the state of an interior boundary from a file. /// Reads the state of an interior boundary from a file.
void deserialize(std::istream &in, std::map<int, Element*> &elIndexMap); void deserialize(std::istream &in, std::map<int, Element*> &elIndexMap);
/// Compares this interior boundaries with some other. The order of the
/// boundary elements within the object does not play a role.
bool operator==(const InteriorBoundary& other) const;
protected: protected:
void serializeExcludeList(std::ostream &out, ExcludeList &list); void serializeExcludeList(std::ostream &out, ExcludeList &list);
......
This diff is collapsed.
...@@ -42,6 +42,27 @@ namespace AMDiS { ...@@ -42,6 +42,27 @@ namespace AMDiS {
class ParMetisPartitioner; class ParMetisPartitioner;
struct ElementNeighbour {
ElementNeighbour()
: elIndex(-1),
ith(0),
boundaryIndex(INTERIOR)
{}
ElementNeighbour(int a, int b, BoundaryType c = INTERIOR)
: elIndex(a),
ith(b),
boundaryIndex(c)
{}
int elIndex;
int ith;
BoundaryType boundaryIndex;
};
class MeshDistributor class MeshDistributor
{ {
protected: protected:
...@@ -114,6 +135,12 @@ namespace AMDiS { ...@@ -114,6 +135,12 @@ namespace AMDiS {
return name; return name;
} }
/// Returns \ref feSpace.
inline const FiniteElemSpace* getFeSpace()
{
return feSpace;
}
/// Returns \ref nRankDOFs, the number of DOFs in the rank mesh. /// Returns \ref nRankDOFs, the number of DOFs in the rank mesh.
inline int getNumberRankDofs() inline int getNumberRankDofs()
{ {
...@@ -236,21 +263,7 @@ namespace AMDiS { ...@@ -236,21 +263,7 @@ namespace AMDiS {
*/ */
void createInteriorBoundaryInfo(); void createInteriorBoundaryInfo();
/** \brief void createBoundaryDataStructure();
* Deterimes the interior boundaries between ranks, that are based on the
* neighbourhood information of the macro elements. That means that in 2d the
* function search for all edge based interior boundaries and in 3d for all face
* based interior boundaries. This function cannot find boundaries of substructure
* elements, i.e. vertex boundaries in 2d and vertex and edge boundaries in 3d.
*/
void createMacroElementInteriorBoundaryInfo();
/** \brief
* Determines all interior boundaries between rank that consist of element's
* substructures. In 2d these may be only vertices. In 3d there may be also
* interior boundaries consisting of either a whole edge or of a single vertex.
*/
void createSubstructureInteriorBoundaryInfo();
/// Removes all macro elements from the mesh that are not part of ranks partition. /// Removes all macro elements from the mesh that are not part of ranks partition.
void removeMacroElements(); void removeMacroElements();
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment