Bugfix in petsc solver global matrix with periodic bc.

AMDiS/src/DOFMatrix.cc

@@ -238,10 +238,10 @@ namespace AMDiS {
 	if (condition->applyBoundaryCondition()) {
 
 #ifdef HAVE_PARALLEL_DOMAIN_AMDIS
-	   	  if (dofMap->isRankDof(rowIndices[i])) {
+	  if (dofMap->isRankDof(rowIndices[i])) {
  	    applyDBCs.insert(row);
 	    //	    dirichletDofs.push_back(row);
-	    	  }
+	  }
 #else
  	  applyDBCs.insert(row);
 #endif

AMDiS/src/DirichletBC.cc

@@ -76,9 +76,9 @@ namespace AMDiS {
 	if (dofVec)
 	  value = (*dofVec)[dofIndices[i]];
 
-// #ifdef HAVE_PARALLEL_DOMAIN_AMDIS
-// 	vector->setDirichletDofValue(dofIndices[i], value);
-// #else
+	//#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+	//	vector->setDirichletDofValue(dofIndices[i], value);
+	//#else
 	(*vector)[dofIndices[i]] = value;
 	//#endif
       }

AMDiS/src/Element.cc

@@ -92,7 +92,7 @@ namespace AMDiS {
 	for (int i = 0; i < mesh->getGeo(EDGE); i++)
 	  dof[mesh->getNode(EDGE) + i] = mesh->getDof(EDGE);
 
-      if (mesh->getNumberOfDofs(FACE)) 
+      if (mesh->getDim() == 3 && mesh->getNumberOfDofs(FACE)) 
 	for (int i = 0; i < mesh->getGeo(FACE); i++)
 	  dof[mesh->getNode(FACE) + i] = mesh->getDof(FACE);
 

AMDiS/src/Mesh.h

@@ -160,6 +160,7 @@ namespace AMDiS {
     /// Returns nDof[i] of the mesh
     inline const int getNumberOfDofs(int i) const 
     { 
+      TEST_EXIT_DBG(i <= dim)("Wrong index: %d %d\n", i, dim);
       return nDof[i]; 
     }
 
@@ -181,12 +182,10 @@ namespace AMDiS {
       return node[i]; 
     }
 
-    /** \brief
-     * Allocates the number of DOFs needed at position and registers the DOFs
-     * at the DOFAdmins. The number of needed DOFs is the sum over the needed
-     * DOFs of all DOFAdmin objects belonging to this mesh. 
-     * The return value is a pointer to the first allocated DOF. 
-     */
+    /// Allocates the number of DOFs needed at position and registers the DOFs
+    /// at the DOFAdmins. The number of needed DOFs is the sum over the needed
+    /// DOFs of all DOFAdmin objects belonging to this mesh. 
+    /// The return value is a pointer to the first allocated DOF. 
     DegreeOfFreedom *getDof(GeoIndex position);
 
     /// Returns *(\ref admin[i]) of the mesh
@@ -195,26 +194,20 @@ namespace AMDiS {
       return *(admin[i]);
     }
 
-    /** \brief
-     * Creates a DOFAdmin with name lname. nDof specifies how many DOFs 
-     * are needed at the different positions (see \ref DOFAdmin::nrDOF).
-     * A pointer to the created DOFAdmin is returned.
-     */
+    /// Creates a DOFAdmin with name lname. nDof specifies how many DOFs 
+    /// are needed at the different positions (see \ref DOFAdmin::nrDOF).
+    /// A pointer to the created DOFAdmin is returned.
     const DOFAdmin* createDOFAdmin(string lname, DimVec<int> nDof);
 
-    /** \brief
-     * Returns the size of \ref admin which is the number of the DOFAdmins
-     * belonging to this mesh
-     */
+    /// Returns the size of \ref admin which is the number of the DOFAdmins
+    /// belonging to this mesh
     const int getNumberOfDOFAdmin() const 
     {
       return admin.size();
     }
 
-    /** \brief
-     * Returns the size of \ref macroElements which is the number of
-     * of macro elements of this mesh
-     */
+    /// Returns the size of \ref macroElements which is the number of
+    /// of macro elements of this mesh
     const int getNumberOfMacros() const 
     {
       return macroElements.size();

AMDiS/src/parallel/ParallelDofMapping.cc

@@ -340,7 +340,6 @@ namespace AMDiS {
     // DOFs in the first FE spaces is the smalled global index of a DOF that is
     // owned by the rank.
     int offset = rStartDofs;
-    
 
     // === Create the matrix indices for all component FE spaces. ===
 
@@ -358,7 +357,7 @@ namespace AMDiS {
 	    dofToMatIndex.add(i, it->first, globalMatIndex);
 	}
       }
-      
+
       // Increase the offset for the next FE space by the number of DOFs owned 
       // by the rank in the current FE space.
       offset += data[feSpaces[i]].nRankDofs;
@@ -399,7 +398,7 @@ namespace AMDiS {
 
 	stdMpi.recv(rank);
       }
-      
+
       stdMpi.startCommunication();
       
       for (DofComm::Iterator it(dofComm->getRecvDofs(), 0, feSpaces[i]); 
@@ -461,9 +460,13 @@ namespace AMDiS {
 	for (; !it.endDofIter(); it.nextDof())
 	  if (dofMap.count(it.getDofIndex())) {
 	    if (globalIndex) {
+	      TEST_EXIT_DBG(counter + 2 <= stdMpi.getRecvData(it.getRank()).size())
+		("Should not happen!\n");
+
 	      dofToMatIndex.add(i, 
-				stdMpi.getRecvData(it.getRank())[counter++],
-				stdMpi.getRecvData(it.getRank())[counter++]);
+				stdMpi.getRecvData(it.getRank())[counter],
+				stdMpi.getRecvData(it.getRank())[counter + 1]);
+	      counter += 2;
 	    } else {
 	      dofToMatIndex.add(i, it.getDofIndex(),
 				stdMpi.getRecvData(it.getRank())[counter++]);

AMDiS/src/parallel/PetscSolverGlobalMatrix.cc

@@ -989,12 +989,13 @@ namespace AMDiS {
 	    int localPetscRowIdx = petscRowIdx - rankStartIndex;
 	    
 	    TEST_EXIT_DBG(localPetscRowIdx >= 0 && localPetscRowIdx < nRankRows)
-	      ("Should not happen! \n Debug info: localRowIdx = %d   globalRowIndx = %d   petscRowIdx = %d   localPetscRowIdx = %d   rStart = %d   nCompontens = %d   nRankRows = %d\n",
+	      ("Should not happen! \n Debug info: DOF = %d   globalRowIndx = %d   petscRowIdx = %d   localPetscRowIdx = %d   rStart = %d   compontens = %d from %d   nRankRows = %d\n",
 	       *cursor,
 	       (*interiorMap)[feSpaces[i]][*cursor].global,
 	       petscRowIdx, 
 	       localPetscRowIdx, 
 	       rankStartIndex,
+	       i,
 	       nComponents, 
 	       nRankRows);