Skip to content
GitLab
Commit 08f87122 authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

First work on parallel domain decomposition for 3d.

parent bcaf907b
Hide whitespace changes
Inline Side-by-side
AMDiS/src/ElInfo.cc
View file @ 08f87122
......@@ -27,7 +27,7 @@ namespace AMDiS {
oppCoord_(mesh_->getDim(), NO_INIT),
neighbour_(mesh_->getDim(), NO_INIT),
neighbourCoord_(mesh_->getDim(), NO_INIT),
oppVertex_(mesh_->getDim(), NO_INIT),
oppVertex(mesh_->getDim(), NO_INIT),
grdLambda(mesh_->getDim(), NO_INIT),
refinementPath(0),
refinementPathLength(0)
......
AMDiS/src/ElInfo.h
View file @ 08f87122
......@@ -78,7 +78,7 @@ namespace AMDiS {
oppCoord_ = rhs.oppCoord_;
neighbour_ = rhs.neighbour_;
neighbourCoord_ = rhs.neighbourCoord_;
oppVertex_ = rhs.oppVertex_;
oppVertex = rhs.oppVertex;
return *this;
}
......@@ -192,15 +192,15 @@ namespace AMDiS {
return neighbourCoord_[i];
}
/// Get ElInfo's \ref oppVertex_[i]
/// Get ElInfo's \ref oppVertex[i]
inline unsigned char getOppVertex(int i) const
{
return oppVertex_[i];
return oppVertex[i];
}
virtual int getSideOfNeighbour(int i)
{
return oppVertex_[i];
return oppVertex[i];
}
/// Get ElInfo's \ref det_
......@@ -221,10 +221,10 @@ namespace AMDiS {
return projection_[i];
}
/// Returns \ref parametric_
/// Returns \ref parametric
inline bool getParametric()
{
return parametric_;
return parametric;
}
virtual void getSubElemCoordsMat(mtl::dense2D<double>& mat,
......@@ -305,10 +305,10 @@ namespace AMDiS {
det_ = d;
}
/// Set \ref parametric_ = param
/// Set \ref parametric = param
inline void setParametric(bool param)
{
parametric_ = param;
parametric = param;
}
/// Set \ref refinementPath
......@@ -481,11 +481,11 @@ namespace AMDiS {
FixVec<FixVec<WorldVector<double>, VERTEX>, NEIGH> neighbourCoord_;
/** \brief
* oppVertex_[i] is undefined if neighbour_[i] is a pointer to NULL.
* oppVertex[i] is undefined if neighbour_[i] is a pointer to NULL.
* Otherwise it is the local index of the neighbour's vertex opposite the
* common edge/face.
*/
FixVec<unsigned char, NEIGH> oppVertex_;
FixVec<unsigned char, NEIGH> oppVertex;
/// Elements determinant.
double det_;
......@@ -494,7 +494,7 @@ namespace AMDiS {
DimVec<WorldVector<double> > grdLambda;
/// True, if this elInfo stores parametrized information. False, otherwise.
bool parametric_;
bool parametric;
/// Stores the world dimension.
int dimOfWorld;
......
AMDiS/src/ElInfo1d.cc
View file @ 08f87122
......@@ -95,7 +95,7 @@ namespace AMDiS {
}
}
neighbour_[i] = nb;
oppVertex_[i] = nb ? i : -1;
oppVertex[i] = nb ? i : -1;
}
}
......@@ -240,11 +240,10 @@ namespace AMDiS {
const FixVec<WorldVector<double>, VERTEX> *old_coord = &(elInfoOld->coord_);
coord_[ichild] = (*old_coord)[ichild];
if (elem->isNewCoordSet()) {
if (elem->isNewCoordSet())
coord_[1 - ichild] = *(elem->getNewCoord());
} else {
else
coord_[1 - ichild] = ((*old_coord)[0] + (*old_coord)[1]) * 0.5;
}
}
if (fillFlag_.isSet(Mesh::FILL_NEIGH) || fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) {
......@@ -255,36 +254,32 @@ namespace AMDiS {
for (int i = 0; i < neighbours; i++) {
if (i != ichild) {
nb = const_cast<Element*>( elem->getChild(1-ichild));
if ( fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) {
nb = const_cast<Element*>(elem->getChild(1-ichild));
if (fillFlag_.isSet(Mesh::FILL_OPP_COORDS))
oppC = elInfoOld->coord_[i];
}
} else {
nb = const_cast<Element*>( elInfoOld->getNeighbour(i));
if (nb && fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) {
if (nb && fillFlag_.isSet(Mesh::FILL_OPP_COORDS))
oppC = elInfoOld->oppCoord_[i];
}
}
if (nb) {
while (nb->getChild(0)) { // make nb nearest element
if (fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) {
if (nb->isNewCoordSet()) {
if (nb->isNewCoordSet())
oppC = *(nb->getNewCoord());
} else {
else
oppC = (coord_[i] + oppC) * 0.5;
}
}
nb = const_cast<Element*>( nb->getChild(1-i));
}
if (fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) {
if (fillFlag_.isSet(Mesh::FILL_OPP_COORDS))
oppCoord_[i] = oppC;
}
}
neighbour_[i] = nb;
oppVertex_[i] = nb ? i : -1;
oppVertex[i] = nb ? i : -1;
}
}
......@@ -299,7 +294,8 @@ namespace AMDiS {
}
}
void ElInfo1d::getSubElemCoordsMat(mtl::dense2D<double>& mat, int basisFctsDegree) const
void ElInfo1d::getSubElemCoordsMat(mtl::dense2D<double>& mat,
int basisFctsDegree) const
{
switch (basisFctsDegree) {
case 1:
......
AMDiS/src/ElInfo2d.cc
View file @ 08f87122
......@@ -76,7 +76,7 @@ namespace AMDiS {
neighbour_[i] = macroNeighbour->getElement();
Element *nb = const_cast<Element*>(neighbour_[i]);
int edgeNo = oppVertex_[i] = mel->getOppVertex(i);
int edgeNo = oppVertex[i] = mel->getOppVertex(i);
if (nb->getFirstChild() && edgeNo != 2) {
......@@ -127,7 +127,7 @@ namespace AMDiS {
// In both cases the opp vertex number is 2, as one can see in the
// pictures above.
oppVertex_[i] = 2;
oppVertex[i] = 2;
if (fill_opp_coords) {
if (nb->isNewCoordSet()) {
......@@ -298,7 +298,7 @@ namespace AMDiS {
// cooresponding oppVertex.
neighbour_[2] = elInfoOld->neighbour_[1];
oppVertex_[2] = elInfoOld->oppVertex_[1];
oppVertex[2] = elInfoOld->oppVertex[1];
if (neighbour_[2] && fill_opp_coords) {
oppCoord_[2] = elInfoOld->oppCoord_[1];
......@@ -314,7 +314,7 @@ namespace AMDiS {
elem->getSecondChild()->getFirstChild()) {
neighbour_[1] = elem->getSecondChild()->getSecondChild();
oppVertex_[1] = 2;
oppVertex[1] = 2;
if (fill_opp_coords) {
if (elem->getSecondChild()->isNewCoordSet()) {
......@@ -330,7 +330,7 @@ namespace AMDiS {
}
} else {
neighbour_[1] = elem->getSecondChild();
oppVertex_[1] = 0;
oppVertex[1] = 0;
if (fill_opp_coords) {
oppCoord_[1] = elInfoOld->coord_[1];
......@@ -347,14 +347,14 @@ namespace AMDiS {
nb = elInfoOld->neighbour_[2];
if (nb) {
TEST(elInfoOld->oppVertex_[2] == 2)("invalid neighbour\n");
TEST(elInfoOld->oppVertex[2] == 2)("invalid neighbour\n");
TEST_EXIT_DBG(nb->getFirstChild())("missing first child?\n");
TEST_EXIT_DBG(nb->getSecondChild())("missing second child?\n");
nb = nb->getSecondChild();
if (nb->getFirstChild()) {
oppVertex_[0] = 2;
oppVertex[0] = 2;
if (fill_opp_coords) {
if (nb->isNewCoordSet()) {
......@@ -372,7 +372,7 @@ namespace AMDiS {
nb = nb->getFirstChild();
} else {
oppVertex_[0] = 1;
oppVertex[0] = 1;
if (fill_opp_coords) {
oppCoord_[0] = elInfoOld->oppCoord_[2];
......@@ -391,7 +391,7 @@ namespace AMDiS {
// cooresponding oppVertex.
neighbour_[2] = elInfoOld->neighbour_[0];
oppVertex_[2] = elInfoOld->oppVertex_[0];
oppVertex[2] = elInfoOld->oppVertex[0];
if (neighbour_[2] && fill_opp_coords) {
oppCoord_[2] = elInfoOld->oppCoord_[0];
......@@ -404,7 +404,7 @@ namespace AMDiS {
if (elem->getFirstChild()->getFirstChild()) {
neighbour_[0] = elem->getFirstChild()->getFirstChild();
oppVertex_[0] = 2;
oppVertex[0] = 2;
if (fill_opp_coords) {
if (elem->getFirstChild()->isNewCoordSet()) {
......@@ -420,7 +420,7 @@ namespace AMDiS {
}
} else {
neighbour_[0] = elem->getFirstChild();
oppVertex_[0] = 1;
oppVertex[0] = 1;
if (fill_opp_coords) {
oppCoord_[0] = elInfoOld->coord_[0];
......@@ -436,11 +436,11 @@ namespace AMDiS {
nb = elInfoOld->neighbour_[2];
if (nb) {
TEST(elInfoOld->oppVertex_[2] == 2)("invalid neighbour\n");
TEST(elInfoOld->oppVertex[2] == 2)("invalid neighbour\n");
TEST((nb = nb->getFirstChild()))("missing child?\n");
if (nb->getFirstChild()) {
oppVertex_[1] = 2;
oppVertex[1] = 2;
if (fill_opp_coords) {
if (nb->isNewCoordSet()) {
......@@ -458,7 +458,7 @@ namespace AMDiS {
nb = nb->getSecondChild();
} else {
oppVertex_[1] = 0;
oppVertex[1] = 0;
if (fill_opp_coords) {
oppCoord_[1] = elInfoOld->oppCoord_[2];
......@@ -476,11 +476,10 @@ namespace AMDiS {
if (fill_flag.isSet(Mesh::FILL_BOUND)) {
if (elInfoOld->getBoundary(2)) {
if (elInfoOld->getBoundary(2))
boundary_[5] = elInfoOld->getBoundary(2);
} else {
else
boundary_[5] = INTERIOR;
}
if (ichild == 0) {
boundary_[3] = elInfoOld->getBoundary(5);
......
AMDiS/src/ElInfo3d.cc
View file @ 08f87122
......@@ -32,9 +32,8 @@ namespace AMDiS {
if (fillFlag_.isSet(Mesh::FILL_COORDS) ||
fillFlag_.isSet(Mesh::FILL_DET) ||
fillFlag_.isSet(Mesh::FILL_GRD_LAMBDA)) {
for (int i = 0; i < vertices; i++) {
for (int i = 0; i < vertices; i++)
coord_[i] = mel->coord[i];
}
}
int neighbours = mesh_->getGeo(NEIGH);
......@@ -48,7 +47,7 @@ namespace AMDiS {
neighbour_[i] = const_cast<Element*>( mel->getNeighbour(i)->getElement());
nb = const_cast<Element*>( neighbour_[i]);
int k;
k = oppVertex_[i] = mel->getOppVertex(i);
k = oppVertex[i] = mel->getOppVertex(i);
if (nb->getChild(0) && (k < 2)) { /*make nb nearest element.*/
if (k == 1) {
......@@ -58,7 +57,7 @@ namespace AMDiS {
neighbour_[i] = const_cast<Element*>( nb->getChild(1));
nb = const_cast<Element*>( neighbour_[i]);
}
k = oppVertex_[i] = 3;
k = oppVertex[i] = 3;
if (fill_opp_coords.isAnySet()) {
/* always edge between vertices 0 and 1 is bisected! */
if (mnb->getElement()->isNewCoordSet())
......@@ -78,13 +77,11 @@ namespace AMDiS {
}
if (fillFlag_.isSet(Mesh::FILL_BOUND)) {
for (int i = 0; i < element_->getGeo(BOUNDARY); i++) {
boundary_[i] = mel->getBoundary(i);
}
for (int i = 0; i < element_->getGeo(PROJECTION); i++) {
for (int i = 0; i < element_->getGeo(BOUNDARY); i++)
boundary_[i] = mel->getBoundary(i);
for (int i = 0; i < element_->getGeo(PROJECTION); i++)
projection_[i] = mel->getProjection(i);
}
}
if (fillFlag_.isSet(Mesh::FILL_ORIENTATION)) {
......@@ -231,9 +228,8 @@ namespace AMDiS {
if (abs(det) < DBL_TOL) {
ERROR("det = %le; abort\n", det);
for (int i = 0; i <= dim; i++) {
for (int i = 0; i <= dim; i++)
(*lambda)[i] = 1.0 / dim;
}
return 0;
}
......@@ -276,7 +272,7 @@ namespace AMDiS {
for (int ineigh = 0; ineigh < neighbours; ineigh++) {
if ((nb = dynamic_cast<Tetrahedron*>(const_cast<Element*>(neighbour_[ineigh])))) {
int ov = oppVertex_[ineigh];
int ov = oppVertex[ineigh];
if (ov < 2 && nb->getFirstChild()) {
if (fill_opp_coords != Flag(0)) {
int k = -1;
......@@ -285,11 +281,9 @@ namespace AMDiS {
k = j;
if (k == -1) {
for (int j = 0; j < vertices; j++) {
if (mesh_->associated(element_->getDOF(j, 0), nb->getDOF(1 - ov, 0))) {
for (int j = 0; j < vertices; j++)
if (mesh_->associated(element_->getDOF(j, 0), nb->getDOF(1 - ov, 0)))
k = j;
}
}
}
TEST_EXIT_DBG(k >= 0)("neighbour dof not found\n");
......@@ -300,7 +294,7 @@ namespace AMDiS {
oppCoord_[ineigh][j] = (oppCoord_[ineigh][j] + coord_[k][j]) / 2;
}
neighbour_[ineigh] = dynamic_cast<Tetrahedron*>(const_cast<Element*>(nb->getChild(1-ov)));
oppVertex_[ineigh] = 3;
oppVertex[ineigh] = 3;
}
}
}
......@@ -389,16 +383,14 @@ namespace AMDiS {
fillFlag_.isSet(Mesh::FILL_DET) ||
fillFlag_.isSet(Mesh::FILL_GRD_LAMBDA)) {
for (int i = 0; i < 3; i++) {
for (int j = 0; j < dimOfWorld; j++) {
for (int j = 0; j < dimOfWorld; j++)
coord_[i][j] = elInfoOld->coord_[cv[i]][j];
}
}
if (el_old->getNewCoord()) {
coord_[3] = *(el_old->getNewCoord());
} else {
for (int j = 0; j < dimOfWorld; j++) {
for (int j = 0; j < dimOfWorld; j++)
coord_[3][j] = (elInfoOld->coord_[0][j] + elInfoOld->coord_[1][j]) / 2;
}
}
}
......@@ -423,21 +415,19 @@ namespace AMDiS {
oppCoord_[0]= *(nb->getNewCoord());
} else {
int k = cvg[ochild][1];
for (int j = 0; j < dimOfWorld; j++) {
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[0][j] = (elInfoOld->coord_[ochild][j] + elInfoOld->coord_[k][j]) / 2;
}
}
}
(*neigh_local)[0] = const_cast<Element*>( nb->getChild(1));
oppVertex_[0] = 3;
oppVertex[0] = 3;
} else {
if (fill_opp_coords.isAnySet()) {
for (int j = 0; j < dimOfWorld; j++) {
if (fill_opp_coords.isAnySet())
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[0][j] = elInfoOld->coord_[ochild][j];
}
}
(*neigh_local)[0] = nb;
oppVertex_[0] = 0;
oppVertex[0] = 0;
}
} else {
ERROR_EXIT("no other child");
......@@ -456,7 +446,7 @@ namespace AMDiS {
nbk = const_cast<Element*>( nb->getChild(k));
if (nbk->getDOF(0) == el_old->getDOF(ichild)) {
/* opp. vertex */
dof = const_cast<int*>(nb->getDOF(elInfoOld->oppVertex_[cv[i]]));
dof = const_cast<int*>(nb->getDOF(elInfoOld->oppVertex[cv[i]]));
if (dof == nbk->getDOF(1)) {
ov = 1;
......@@ -465,13 +455,13 @@ namespace AMDiS {
if (nbk->getNewCoord()) {
oppCoord_[i] = *(nbk->getNewCoord());
} else {
for (int j = 0; j < dimOfWorld; j++) {
oppCoord_[i][j] = (elInfoOld->oppCoord_[cv[i]][j] + elInfoOld->coord_[ichild][j]) / 2;
}
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[i][j] = (elInfoOld->oppCoord_[cv[i]][j] +
elInfoOld->coord_[ichild][j]) / 2;
}
}
(*neigh_local)[i] = nbk->getChild(0);
oppVertex_[i] = 3;
oppVertex[i] = 3;
break;
}
} else {
......@@ -482,13 +472,12 @@ namespace AMDiS {
ov = 2;
}
if (fill_opp_coords.isAnySet()) {
for (int j = 0; j < dimOfWorld; j++) {
if (fill_opp_coords.isAnySet())
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[i][j] = elInfoOld->oppCoord_[cv[i]][j];
}
}
(*neigh_local)[i] = nbk;
oppVertex_[i] = ov;
oppVertex[i] = ov;
break;
}
......@@ -501,7 +490,7 @@ namespace AMDiS {
if (nbk->getDOF(0) == el_old->getDOF(ichild) ||
mesh->associated(nbk->getDOF(0, 0), el_old->getDOF(ichild, 0))) {
/* opp. vertex */
dof = const_cast<int*>(nb->getDOF(elInfoOld->oppVertex_[cv[i]]));
dof = const_cast<int*>(nb->getDOF(elInfoOld->oppVertex[cv[i]]));
if (dof == nbk->getDOF(1) ||
mesh->associated(dof[0], nbk->getDOF(1, 0))) {
......@@ -511,13 +500,13 @@ namespace AMDiS {
if (nbk->getNewCoord()) {
oppCoord_[i] = *(nbk->getNewCoord());
} else {
for (int j = 0; j < dimOfWorld; j++) {
oppCoord_[i][j] = (elInfoOld->oppCoord_[cv[i]][j] + elInfoOld->coord_[ichild][j]) / 2;
}
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[i][j] = (elInfoOld->oppCoord_[cv[i]][j] +
elInfoOld->coord_[ichild][j]) / 2;
}
}
(*neigh_local)[i] = nbk->getChild(0);
oppVertex_[i] = 3;
oppVertex[i] = 3;
break;
}
} else {
......@@ -527,17 +516,17 @@ namespace AMDiS {
ov = 2;
}
if (fill_opp_coords.isAnySet()) {
for (int j = 0; j < dimOfWorld; j++) {
if (fill_opp_coords.isAnySet())
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[i][j] = elInfoOld->oppCoord_[cv[i]][j];
}
}
(*neigh_local)[i] = nbk;
oppVertex_[i] = ov;
oppVertex[i] = ov;
break;
}
} /* end for k */
TEST_EXIT_DBG(k < 2)("child not found with vertex\n");
}
} else {
......@@ -549,20 +538,17 @@ namespace AMDiS {
/*----- nb[3] is old neighbour neigh_old[ochild] ------------------------*/
if (((*neigh_local)[3] = (*neigh_old)[ochild])) {
oppVertex_[3] = elInfoOld->oppVertex_[ochild];
oppVertex[3] = elInfoOld->oppVertex[ochild];
if (fill_opp_coords.isAnySet()) {
for (int j = 0; j < dimOfWorld; j++) {
if (fill_opp_coords.isAnySet())
for (int j = 0; j < dimOfWorld; j++)
oppCoord_[3][j] = elInfoOld->oppCoord_[ochild][j];
}
}
}
}
if (fillFlag__local.isSet(Mesh::FILL_BOUND)) {
for (int i = 0; i < 3; i++) {
for (int i = 0; i < 3; i++)
boundary_[10 + i] = elInfoOld->getBoundary(10 + cv[i]);
}
boundary_[13] = elInfoOld->getBoundary(4);
......@@ -574,9 +560,8 @@ namespace AMDiS {
int geoFace = mesh_->getGeo(FACE);
ce = const_cast<int*>(Tetrahedron::childEdge[el_type_local][ichild]);
for (int iedge = 0; iedge < 4; iedge++) {
boundary_[geoFace + iedge] = elInfoOld->getBoundary(geoFace + ce[iedge]);
}
for (int iedge = 0; iedge < 4; iedge++)
boundary_[geoFace + iedge] = elInfoOld->getBoundary(geoFace + ce[iedge]);
for (int iedge = 4; iedge < 6; iedge++) {
int i = 5 - cv[iedge - 3]; /* old vertex opposite new edge */
boundary_[geoFace + iedge] = elInfoOld->getBoundary(i);
......@@ -592,19 +577,17 @@ namespace AMDiS {
projection_[2] = elInfoOld->getProjection(cv[2]);
projection_[3] = elInfoOld->getProjection(ochild);
for (int iedge = 0; iedge < 4; iedge++) {
for (int iedge = 0; iedge < 4; iedge++)
projection_[geoFace + iedge] = elInfoOld->getProjection(geoFace + ce[iedge]);
}
for (int iedge = 4; iedge < 6; iedge++) {
for (int iedge = 4; iedge < 6; iedge++)
projection_[geoFace + iedge] = elInfoOld->getProjection(5 - cv[iedge - 3]);
}
}
}
if (fillFlag_.isSet(Mesh::FILL_ORIENTATION)) {
orientation =
( dynamic_cast<ElInfo3d*>(const_cast<ElInfo*>(elInfoOld)))->orientation
orientation =
(dynamic_cast<ElInfo3d*>(const_cast<ElInfo*>(elInfoOld)))->orientation
* Tetrahedron::childOrientation[el_type_local][ichild];
}
}
......
AMDiS/src/Global.h
View file @ 08f87122
......@@ -368,7 +368,7 @@ namespace AMDiS {
static const Element *getReferenceElement(int dim)
{
FUNCNAME("Global::getReferenceElement()");
TEST_EXIT((dim > 0) && (dim < 4))("invalid dim: %d\n", dim);
TEST_EXIT(dim > 0 && dim < 4)("invalid dim: %d\n", dim);
return referenceElement[dim];
}
......@@ -388,9 +388,9 @@ namespace AMDiS {
*/
static inline int getGeo(GeoIndex p, int dim)
{
TEST_EXIT_DBG((p >= MINPART) && (p <= MAXPART))
TEST_EXIT_DBG(p >= MINPART && p <= MAXPART)
("Calling for invalid geometry value %d\n",p);
TEST_EXIT_DBG((dim >= 0) && (dim < 4))
TEST_EXIT_DBG(dim >= 0 && dim < 4)
("invalid dim: %d\n", dim);
TEST_EXIT_DBG((dim != 0) || (p == PARTS || p == VERTEX || p == EDGE || p == FACE))
("dim = 0\n");
......
AMDiS/src/Mesh.cc
View file @ 08f87122
......@@ -275,7 +275,6 @@ namespace AMDiS {
typedef std::map<const DegreeOfFreedom*, std::set<MacroElement*> > DofElMap;
typedef std::map<const DegreeOfFreedom*, GeoIndex> DofPosMap;
TEST_EXIT(dim == 2)("Not yet implemented for dim != 2!\n");
TEST_EXIT(admin.size() == 1)("Not yet implemented for multiple admins!\n");
TEST_EXIT(admin[0])("There is something wrong!\n");
......
AMDiS/src/ParallelDomainBase.cc
View file @ 08f87122
......@@ -90,7 +90,7 @@ namespace AMDiS {
#if (DEBUG != 0)
ElementIdxToDofs elMap;