Documentation of parallel AMDiS version.

AMDiS/src/Debug.cc

@@ -476,12 +476,20 @@ namespace AMDiS {
 	    for (int i = 0; i < (elInfo->getLevel() - elLevel); i++)
 	      oss << "   ";
 	    oss << "|--" << elInfo->getElement()->getIndex();
-	    MSG("%s\n", oss.str().c_str());
+	    if (oss.str().length() >= 255) {
+	      WARNING("String is to long! Element index is %d.\n", elInfo->getElement()->getIndex());
+	    } else {
+	      MSG("%s\n", oss.str().c_str());
+	    }
 	  }
 	}
 
 	elInfo = stack.traverseNext(elInfo);
       }
+
+      if (!printInfo) {
+	MSG("Could not find element with index %d\n", elIndex);
+      }
     }
 
 

AMDiS/src/parallel/MeshDistributor.cc

@@ -593,12 +593,23 @@ namespace AMDiS {
 	if (elCode.getCode() != recvCodes[i].getCode()) {
 	  TEST_EXIT_DBG(refineManager)("Refinement manager is not set correctly!\n");
 
-	  bool b = fitElementToMeshCode(recvCodes[i], 
-					boundIt->rankObj.el,
-					boundIt->rankObj.subObj,
-					boundIt->rankObj.ithObj, 
-					boundIt->rankObj.elType,
-					boundIt->rankObj.reverseMode);  
+	  MSG("CHECK EL %d %d %d %d WITH El %d %d %d %d in RANK %d\n",
+	      boundIt->rankObj.elIndex,
+	      boundIt->rankObj.subObj,
+	      boundIt->rankObj.ithObj, 
+	      boundIt->rankObj.elType,	      
+	      boundIt->neighObj.elIndex,
+	      boundIt->neighObj.subObj,
+	      boundIt->neighObj.ithObj, 
+	      boundIt->neighObj.elType,	      
+	      it->first);
+
+	  bool b = startFitElementToMeshCode(recvCodes[i], 
+					     boundIt->rankObj.el,
+					     boundIt->rankObj.subObj,
+					     boundIt->rankObj.ithObj, 
+					     boundIt->rankObj.elType,
+					     boundIt->rankObj.reverseMode);
 
 	  if (b)
 	    meshFitTogether = false;
@@ -612,35 +623,45 @@ namespace AMDiS {
   }
 
 
-  bool MeshDistributor::fitElementToMeshCode(MeshStructure &code, 
-					     Element *el, 
-					     GeoIndex subObj,
-					     int ithObj, 
-					     int elType,
-					     bool reverseMode)
+  bool MeshDistributor::startFitElementToMeshCode(MeshStructure &code, 
+						  Element *el, 
+						  GeoIndex subObj,
+						  int ithObj, 
+						  int elType,
+						  bool reverseMode)
   {
-    FUNCNAME("MeshDistributor::fitElementToMeshCode()");
+    FUNCNAME("MeshDistributor::startFitElementToMeshCode()");
 
     TEST_EXIT_DBG(el)("No element given!\n");
 
+    // If the code is empty, the element does not matter and the function can
+    // return without chaning the mesh.
     if (code.empty())
       return false;
 
-    int s1 = el->getSubObjOfChild(0, subObj, ithObj, elType);
-    int s2 = el->getSubObjOfChild(1, subObj, ithObj, elType);
+    // s0 and s1 are the number of the edge/face in both child of the element,
+    // which contain the edge/face the function has to traverse through. If the
+    // edge/face is not contained in one of the children, s0 or s1 is -1.
+    int s0 = el->getSubObjOfChild(0, subObj, ithObj, elType);
+    int s1 = el->getSubObjOfChild(1, subObj, ithObj, elType);
 
-    TEST_EXIT_DBG(s1 != -1 || s2 != -1)("This should not happen!\n");
+    TEST_EXIT_DBG(s0 != -1 || s1 != -1)("This should not happen!\n");
 
     bool meshChanged = false;
     Flag traverseFlag = 
       Mesh::CALL_EVERY_EL_PREORDER | Mesh::FILL_NEIGH | Mesh::FILL_BOUND;
 
-    if (reverseMode) {
-      std::swap(s1, s2);
-      traverseFlag |= Mesh::CALL_REVERSE_MODE;
-    }
+    // Test for reverse mode, in which the left and right children of elements
+    // are flipped.
+    if (reverseMode)
+      traverseFlag |= Mesh::CALL_REVERSE_MODE;    
+
 
-    if (s1 != -1 && s2 != -1) {
+    // === If the edge/face is contained in both children. ===
+
+    if (s0 != -1 && s1 != -1) {
+      // Create traverse stack and traverse within the mesh until the element,
+      // which should be fitted to the mesh structure code, is reached.
       TraverseStack stack;
       ElInfo *elInfo = stack.traverseFirst(el->getMesh(), -1, traverseFlag);
       while (elInfo && elInfo->getElement() != el)
@@ -648,78 +669,88 @@ namespace AMDiS {
 
       TEST_EXIT_DBG(elInfo->getElement() == el)("This should not happen!\n");
 
-      meshChanged = 
-	fitElementToMeshCode2(code, stack, subObj, ithObj, elType, reverseMode);
+      meshChanged = fitElementToMeshCode(code, stack, subObj, ithObj, reverseMode);
       return meshChanged;
     }
 
+
+    // === The edge/face is contained in only on of the both children. ===
+
     if (el->isLeaf()) {
+
+      // If element is leaf and code contains only one leaf element, we are finished.
       if (code.getNumElements() == 1 && code.isLeafElement())
 	return false;     
 
-      meshChanged = true;
-
+      // Create traverse stack and traverse the mesh to the element.
       TraverseStack stack;
       ElInfo *elInfo = stack.traverseFirst(el->getMesh(), -1, traverseFlag);
-
-      while (elInfo && elInfo->getElement() != el) {
-	elInfo = stack.traverseNext(elInfo);
-      }
+      while (elInfo && elInfo->getElement() != el)
+	elInfo = stack.traverseNext(elInfo);      
 
       TEST_EXIT_DBG(elInfo)("This should not happen!\n");
 
+      // Code is not leaf, therefore refine the element.
       el->setMark(1);
       refineManager->setMesh(el->getMesh());
       refineManager->setStack(&stack);
       refineManager->refineFunction(elInfo);
+
+      MSG("A-NEW ELs %d %d\n", el->getChild(0)->getIndex(), el->getChild(1)->getIndex());
+
+      meshChanged = true;
     }
 
     Element *child0 = el->getFirstChild();
     Element *child1 = el->getSecondChild();
-    if (reverseMode)
-      std::swap(child0, child1);
+    if (reverseMode) {
+      std::swap(s0, s1);
+      std::swap(child0, child1);    
+    }
 
-    if (s1 != -1) {
-      TraverseStack stack;
-      ElInfo *elInfo = stack.traverseFirst(el->getMesh(), -1, traverseFlag);
+    // === We know that the edge/face is contained in only one of the children. ===
+    // === Therefore, traverse the mesh to this children and fit this element   ===
+    // === To the mesh structure code.                                          ===
 
-      while (elInfo && elInfo->getElement() != child0) {
-	elInfo = stack.traverseNext(elInfo);
-      }
+    TraverseStack stack;
+    ElInfo *elInfo = stack.traverseFirst(el->getMesh(), -1, traverseFlag);
 
-      meshChanged |= 
-	fitElementToMeshCode2(code, stack, subObj, s1, el->getChildType(elType), reverseMode);
-    } else {
-      TraverseStack stack;
-      ElInfo *elInfo = stack.traverseFirst(el->getMesh(), -1, traverseFlag);
+    if (s0 != -1) {
+      while (elInfo && elInfo->getElement() != child0)
+	elInfo = stack.traverseNext(elInfo);     
 
-      while (elInfo && elInfo->getElement() != child1) {
-	elInfo = stack.traverseNext(elInfo);
-      }
+      meshChanged |= fitElementToMeshCode(code, stack, subObj, s0, reverseMode);
+    } else {
+      while (elInfo && elInfo->getElement() != child1) 
+	elInfo = stack.traverseNext(elInfo);      
 
-      meshChanged |= 
-	fitElementToMeshCode2(code, stack, subObj, s2, el->getChildType(elType), reverseMode);
+      meshChanged |= fitElementToMeshCode(code, stack, subObj, s1, reverseMode);
     }
 
+
     return meshChanged;
   }
 
 
-  bool MeshDistributor::fitElementToMeshCode2(MeshStructure &code, 
-					      TraverseStack &stack,
-					      GeoIndex subObj,
-					      int ithObj, 
-					      int elType,
-					      bool reverseMode)
+  bool MeshDistributor::fitElementToMeshCode(MeshStructure &code, 
+					     TraverseStack &stack,
+					     GeoIndex subObj,
+					     int ithObj, 
+					     bool reverseMode)
   {
-    FUNCNAME("MeshDistributor::fitElementToMeshCode2()");
+    FUNCNAME("MeshDistributor::fitElementToMeshCode()");
+
+
+    // === Test if there are more elements in stack to check with the code. ===
 
     ElInfo *elInfo = stack.getElInfo();
-    bool value = false;
     if (!elInfo)
-      return value;
+      return false;
 
-    Element *el = elInfo->getElement();
+
+    // === Test if code contains a leaf element. If this is the case, the ===
+    // === current element is finished. Traverse the mesh to the next     ===
+    // === coarser element.                                               ===
 
     if (code.isLeafElement()) {
       int level = elInfo->getLevel();
@@ -728,46 +759,80 @@ namespace AMDiS {
 	elInfo = stack.traverseNext(elInfo);
       } while (elInfo && elInfo->getLevel() > level);
 
-      return value;
+      return false;
     }
 
-    if (!elInfo)
-      return value;
+
+    bool meshChanged = false;
+    Element *el = elInfo->getElement();
+
+
+    // === If element is leaf (and the code is not), refine the element. ===
 
     if (el->isLeaf()) {
+      TEST_EXIT_DBG(elInfo->getLevel() < 255)("This should not happen!\n");
+
       el->setMark(1);
       refineManager->setMesh(el->getMesh());
       refineManager->setStack(&stack);
       refineManager->refineFunction(elInfo);
-      value = true;
+      meshChanged = true;
     }
 
-    int s1 = el->getSubObjOfChild(0, subObj, ithObj, elType);
-    int s2 = el->getSubObjOfChild(1, subObj, ithObj, elType);
+
+    // === Continue fitting the mesh structure code to the children of the ===
+    // === current element.                                                ===
+
+    int s0 = el->getSubObjOfChild(0, subObj, ithObj, elInfo->getType());
+    int s1 = el->getSubObjOfChild(1, subObj, ithObj, elInfo->getType());
     Element *child0 = el->getFirstChild();
     Element *child1 = el->getSecondChild();
     if (reverseMode) {
-      std::swap(s1, s2);
+      std::swap(s0, s1);
       std::swap(child0, child1);
     }
 
-    if (s1 != -1) {
+    
+    // === Traverse left child. ===
+
+    if (s0 != -1) {
+      // The edge/face is contained in the left child, therefore fit this
+      // child to the mesh structure code.
+
       stack.traverseNext(elInfo);
       code.nextElement();
-      value |= fitElementToMeshCode2(code, stack, subObj, s1, el->getChildType(elType), reverseMode);
+      meshChanged |= fitElementToMeshCode(code, stack, subObj, s0, reverseMode);
       elInfo = stack.getElInfo();
     } else {
+      // The edge/face is not contained in the left child. Hence we need
+      // to traverse through all subelements of the left child until we get
+      // the second child of the current element.
+
       do {
 	elInfo = stack.traverseNext(elInfo);
-      } while (elInfo && elInfo->getElement() != child1);  
+      } while (elInfo && elInfo->getElement() != child1); 
+
+      TEST_EXIT_DBG(elInfo != NULL)("This should not happen!\n");
     }  
 
-    TEST_EXIT_DBG(elInfo->getElement() == child1)("This should not happen!\n");
+    TEST_EXIT_DBG(elInfo->getElement() == child1)
+      ("Got wrong child with idx = %d! Searched for child idx = %d\n",
+       elInfo->getElement()->getIndex(), child1->getIndex());
+
+
+    // === Traverse right child. ===
+
+    if (s1 != -1) {
+      // The edge/face is contained in the right child, therefore fit this
+      // child to the mesh structure code.
 
-    if (s2 != -1) {
       code.nextElement();
-      value |= fitElementToMeshCode2(code, stack, subObj, s2, el->getChildType(elType), reverseMode);
+      meshChanged |= fitElementToMeshCode(code, stack, subObj, s1, reverseMode);
     } else {
+      // The edge/face is not contained in the right child. Hence we need
+      // to traverse through all subelements of the right child until we have
+      // finished traversing the current element with all its subelements.
+
       int level = elInfo->getLevel();
 
       do {
@@ -775,7 +840,8 @@ namespace AMDiS {
       } while (elInfo && elInfo->getLevel() > level);
     }
 
-    return value;
+
+    return meshChanged;
   }
 
   
@@ -1003,18 +1069,24 @@ namespace AMDiS {
   {
     FUNCNAME("MeshDistributor::createBoundaryDataStructure()");
 
+    // Data structure to store all sub-objects of all elements of the mesh.
     ElementObjects elObjects(partitionVec);
-
+    // Maps to each element index a pointer to the corresponding element.
     std::map<int, Element*> elIndexMap;
+    // Maps to each element index the type of this element.
     std::map<int, int> elIndexTypeMap;
 
+    // The following three data structures store periodic DOFs, edges and faces,
+    // i.e., 
     std::map<std::pair<DegreeOfFreedom, DegreeOfFreedom>, BoundaryType> periodicDofs;
-    std::map<DegreeOfFreedom, std::set<BoundaryType> > periodicDofAssoc;
-
     std::map<std::pair<DofEdge, DofEdge>, BoundaryType> periodicEdges;
     std::map<std::pair<DofFace, DofFace>, BoundaryType> periodicFaces;
 
-    // === PHASE 1 ===
+    // Stores to each DOF all its periodic associations.
+    std::map<DegreeOfFreedom, std::set<BoundaryType> > periodicDofAssoc;
+
+
+    // === Phase 1, fills the data structures defined above. ===
 
     TraverseStack stack;
     ElInfo *elInfo = 
@@ -1024,7 +1096,10 @@ namespace AMDiS {
       Element *el = elInfo->getElement();
       elIndexMap[el->getIndex()] = el;
       elIndexTypeMap[el->getIndex()] = elInfo->getType();
+
       
+      // === Add all sub object of the element to the variable elObjects. ===
+
       for (int i = 0; i < el->getGeo(VERTEX); i++)
 	elObjects.addVertex(el->getDOF(i, 0), el->getIndex(), i);
 
@@ -1035,6 +1110,8 @@ namespace AMDiS {
 	elObjects.addFace(el->getFace(i), el->getIndex(), i);
 
 
+      // === === 
+
       switch (mesh->getDim()) {
       case 2:
 	for (int i = 0; i < el->getGeo(EDGE); i++) {

AMDiS/src/parallel/MeshDistributor.h

@@ -326,19 +326,49 @@ namespace AMDiS {
     // Removes all periodic boundaries from a given boundary map.
     void removePeriodicBoundaryConditions(BoundaryIndexMap& boundaryMap);
 
+    /** \brief
+     * Starts the procedure to fit a given edge/face of an element with a mesh
+     * structure code. This functions prepares some data structures and call
+     * then \ref fitElementToMeshCode, that mainly refines the element such that
+     * it fits to the mesh structure code.
+     *
+     * \param[in] code         The mesh structure code to which the edge/face of
+     *                         an element must be fitted.
+     * \param[in] el           Pointer to the element.
+     * \param[in] subObj       Defines whether an edge or a face must be fitted.
+     * \param[in] ithObj       Defines which edge/face must be fitted.
+     * \param[in] elType       Element type of the element (only important in 3D).
+     * \param[in] reverseMode  Defines, whether the mesh structure code is given
+     *                         in reverse mode, i.e., left and right children where
+     *                         changed when the code was created.
+     */
+    bool startFitElementToMeshCode(MeshStructure &code, 
+				   Element *el, 
+				   GeoIndex subObj,
+				   int ithObj, 
+				   int elType,
+				   bool reverseMode);
+    
+    /** \brief
+     * Recursively fits a given mesh structure code to an edge/face of an element.
+     * This function is always initialy called from \ref startFitElementToMeshCode.
+     *
+     * \param[in] code         The mesh structure code which is used to fit an
+     *                         edge/face of an element.
+     * \param[in] stack        A traverse stack object. The upper most element in this
+     *                         stack must be used for fitting the mesh structure code
+     *                         at the current position.
+     * \param[in] subObj       Defines whether an edge or a face must be fitted.
+     * \param[in] ithObj       Defines which edge/face must be fitted.
+     * \param[in] reverseMode  Defines, whether the mesh structure code is given
+     *                         in reverse mode, i.e., left and right children where
+     *                         changed when the code was created.
+     */
     bool fitElementToMeshCode(MeshStructure &code, 
-			      Element *el, 
+			      TraverseStack &stack,
 			      GeoIndex subObj,
-			      int ithObj, 
-			      int elType,
+			      int ithObj,
 			      bool reverseMode);
-    
-    bool fitElementToMeshCode2(MeshStructure &code, 
-			       TraverseStack &stack,
-			       GeoIndex subObj,
-			       int ithObj, 
-			       int elType,
-			       bool reverseMode);
 
     /// Writes a vector of dof pointers to an output stream.
     void serialize(std::ostream &out, DofContainer &data);