diff --git a/AMDiS/src/ElInfo3d.cc b/AMDiS/src/ElInfo3d.cc index 83a05c56bc70767474794b0046865d2d19124a15..62deae090b8dffd2bbf9d12ddeb77a88da3c2ddb 100644 --- a/AMDiS/src/ElInfo3d.cc +++ b/AMDiS/src/ElInfo3d.cc @@ -501,18 +501,18 @@ namespace AMDiS { const int (*cvg)[4] = NULL; /* cvg = child_vertex[el_type] */ int *ce; /* ce = child_edge[el_type][ichild] */ Element *nb, *nbk; - Element *el_old = elInfoOld->element; + Element *elOld = elInfoOld->element; Flag fillFlag_local = elInfoOld->fillFlag; DegreeOfFreedom *dof; int ov = -1; Mesh *mesh = elInfoOld->getMesh(); - TEST_EXIT_DBG(el_old->getChild(0))("missing child?\n"); + TEST_EXIT_DBG(elOld->getChild(0))("missing child?\n"); - element = const_cast<Element*>(el_old->getChild(ichild)); + element = const_cast<Element*>(elOld->getChild(ichild)); macroElement = elInfoOld->macroElement; fillFlag = fillFlag_local; - parent = el_old; + parent = elOld; level = elInfoOld->level + 1; iChild = ichild; int el_type_local = 0; @@ -539,8 +539,8 @@ namespace AMDiS { for (int j = 0; j < dimOfWorld; j++) coord[i][j] = elInfoOld->coord[cv[i]][j]; } - if (el_old->getNewCoord()) { - coord[3] = *(el_old->getNewCoord()); + if (elOld->getNewCoord()) { + coord[3] = *(elOld->getNewCoord()); } else { for (int j = 0; j < dimOfWorld; j++) coord[3][j] = (elInfoOld->coord[0][j] + elInfoOld->coord[1][j]) / 2; @@ -556,11 +556,10 @@ namespace AMDiS { Flag fill_opp_coords; fill_opp_coords.setFlags(fillFlag_local & Mesh::FILL_OPP_COORDS); - /*----- nb[0] is other child --------------------------------------------*/ + // === nb[0] is other child === - /* if (nb = el_old->child[ochild]) old version */ - if (el_old->getChild(0) && - (nb = const_cast<Element*>(el_old->getChild(ochild)))) { + if (elOld->getChild(0) && + (nb = const_cast<Element*>(elOld->getChild(ochild)))) { if (nb->getChild(0)) { /* go down one level for direct neighbour */ if (fill_opp_coords.isAnySet()) { @@ -588,16 +587,43 @@ namespace AMDiS { } - /*----- nb[1],nb[2] are childs of old neighbours nb_old[cv[i]] ----------*/ + /*----- nb[1], nb[2] are childs of old neighbours nb_old[cv[i]] ----------*/ for (int i = 1; i < 3; i++) { if ((nb = const_cast<Element*>((*neigh_old)[cv[i]]))) { TEST_EXIT_DBG(nb->getChild(0))("nonconforming triangulation\n"); + + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("NEIGH %d OF EL is %d\n", i, nb->getIndex()); + MSG("EL HAS DOFS: %d %d %d %d\n", + elOld->getDof(0, 0), + elOld->getDof(1, 0), + elOld->getDof(2, 0), + elOld->getDof(3, 0)); + MSG("NEIGH-EL HAS DOFS: %d %d %d %d\n", + nb->getDof(0, 0), + nb->getDof(1, 0), + nb->getDof(2, 0), + nb->getDof(3, 0)); + } + int k; for (k = 0; k < 2; k++) { /* look at both childs of old neighbour */ nbk = const_cast<Element*>(nb->getChild(k)); - if (nbk->getDof(0) == el_old->getDof(ichild)) { + + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("TEST k = %d, nbk = %d, el = %d\n", k, + nbk->getDof(0,0), elOld->getDof(ichild,0)); + MSG("TEST PTR nbk = %p, el = %p\n", + nbk->getDof(0), elOld->getDof(ichild)); + } + + if (nbk->getDof(0) == elOld->getDof(ichild)) { + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("DRIN!\n"); + } + /* opp. vertex */ dof = const_cast<int*>(nb->getDof(elInfoOld->oppVertex[cv[i]])); @@ -615,11 +641,17 @@ namespace AMDiS { } (*neigh_local)[i] = nbk->getChild(0); oppVertex[i] = 3; + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("BREAK 1\n"); + } break; } } else { if (dof != nbk->getDof(2)) { ov = -1; + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("BREAK 2\n"); + } break; } ov = 2; @@ -631,18 +663,30 @@ namespace AMDiS { (*neigh_local)[i] = nbk; oppVertex[i] = ov; + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("BREAK 3\n"); + } break; } } /* end for k */ - // periodic ? + if (elOld->getIndex() == 15376 && ichild == 0) { + MSG("k = %d ov = %d\n", k, ov); + } + + + // === Test if periodic. === + if (k == 2 || ov == -1) { - for (k = 0; k < 2; k++) { /* look at both childs of old neighbour */ + + // Look at both childs of old neighbour. + for (k = 0; k < 2; k++) { nbk = const_cast<Element*>(nb->getChild(k)); - if (nbk->getDof(0) == el_old->getDof(ichild) || - mesh->associated(nbk->getDof(0, 0), el_old->getDof(ichild, 0))) { - /* opp. vertex */ + if (nbk->getDof(0) == elOld->getDof(ichild) || + mesh->associated(nbk->getDof(0, 0), elOld->getDof(ichild, 0))) { + + // opp. vertex dof = const_cast<int*>(nb->getDof(elInfoOld->oppVertex[cv[i]])); if (dof == nbk->getDof(1) || @@ -680,7 +724,9 @@ namespace AMDiS { } /* end for k */ - TEST_EXIT_DBG(k < 2)("child not found with vertex\n"); + TEST_EXIT_DBG(k < 2) + ("Fill %d child of element %d (0-Level: %d): Child not found with vertex!\n", + ichild, elOld->getIndex(), elInfoOld->getMacroElement()->getIndex()); } } else { (*neigh_local)[i] = NULL; diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc index 0b3f4428819fb364badaeac39d60ad45442c45ec..72a57542fe1110b35316b4fedf1a8adbfc4eccfa 100644 --- a/AMDiS/src/parallel/MeshDistributor.cc +++ b/AMDiS/src/parallel/MeshDistributor.cc @@ -204,11 +204,15 @@ namespace AMDiS { for (int i = 0; i < mesh->getGeo(NEIGH); i++) { if ((*it)->getNeighbour(i) && mesh->isPeriodicAssociation((*it)->getBoundary(i))) { + + int neighIndex = (*it)->getNeighbour(i)->getIndex(); + (*it)->getNeighbour(i)->setNeighbour((*it)->getOppVertex(i), NULL); (*it)->setNeighbour(i, NULL); (*it)->setBoundary(i, 0); macroElementNeighbours[(*it)->getIndex()][i] = -1; + macroElementNeighbours[neighIndex][(*it)->getOppVertex(i)] = -1; } } } @@ -218,11 +222,15 @@ namespace AMDiS { for (int i = 0; i < mesh->getGeo(NEIGH); i++) { if ((*it)->getNeighbour(i) && mesh->isPeriodicAssociation((*it)->getBoundary(i))) { + + int neighIndex = (*it)->getNeighbour(i)->getIndex(); + (*it)->getNeighbour(i)->setNeighbour((*it)->getOppVertex(i), NULL); (*it)->setNeighbour(i, NULL); (*it)->setBoundary(i, 0); macroElementNeighbours[(*it)->getIndex()][i] = -1; + macroElementNeighbours[neighIndex][(*it)->getOppVertex(i)] = -1; } } } @@ -974,7 +982,7 @@ namespace AMDiS { } - // === Run ParMETiS to calculate a new mesh partitioning. === + // === Run mesh partitioner to calculate a new mesh partitioning. === partitioner->setLocalGlobalDofMap(&mapLocalGlobalDofs); bool partitioningSucceed = partitioner->partition(elemWeights, ADAPTIVE_REPART); @@ -983,6 +991,7 @@ namespace AMDiS { return; } + // === Create map that maps macro element indices to pointers to the === // === macro elements. === @@ -1060,6 +1069,7 @@ namespace AMDiS { } } + StdMpi<MeshCodeVec> stdMpi(mpiComm, true); stdMpi.send(sendCodes); for (std::map<int, std::vector<int> >::iterator it = partitioner->getRecvElements().begin(); @@ -1067,6 +1077,7 @@ namespace AMDiS { stdMpi.recv(it->first); stdMpi.startCommunication(); + if (interchangeVectors.size() == 0) WARNING("There are no interchange vectors defined!\n"); diff --git a/AMDiS/src/parallel/MeshManipulation.cc b/AMDiS/src/parallel/MeshManipulation.cc index e17f001194de3a3b01f6d01e5ea24d34f66b4934..bd27f89ac2e1254558384bc151e631b26950be5d 100644 --- a/AMDiS/src/parallel/MeshManipulation.cc +++ b/AMDiS/src/parallel/MeshManipulation.cc @@ -141,7 +141,7 @@ namespace AMDiS { el1->getNonVertexDofs(feSpace, b1, dofs1); #if (DEBUG != 0) - debug::testDofsByCoords(feSpace, dofs0, dofs1); + // debug::testDofsByCoords(feSpace, dofs0, dofs1); #endif for (unsigned int i = 0; i < dofs0.size(); i++) @@ -183,7 +183,7 @@ namespace AMDiS { el1->getNonVertexDofs(feSpace, b1, dofs1); #if (DEBUG != 0) - debug::testDofsByCoords(feSpace, dofs0, dofs1); + // debug::testDofsByCoords(feSpace, dofs0, dofs1); #endif for (unsigned int i = 0; i < dofs0.size(); i++)