Work on preconditioners for schur complement.

AMDiS/src/Global.h

@@ -41,6 +41,7 @@
 
 #include <string>
 #include <vector>
+#include <set>
 #include <fstream>
 #include <math.h>
 #include <iostream>
@@ -80,6 +81,8 @@ namespace AMDiS {
   /// Defines type for a vector of DOF pointers.
   typedef std::vector<const DegreeOfFreedom*> DofContainer;
 
+  typedef std::set<const DegreeOfFreedom*> DofContainerSet;
+
   /// Defines a type for global edge identification via its DOFs.
   typedef std::pair<DegreeOfFreedom, DegreeOfFreedom> DofEdge;
 

AMDiS/src/parallel/MeshDistributor.cc

@@ -1649,20 +1649,21 @@ namespace AMDiS {
     recvDofs.clear();
 
     if (createBoundaryDofFlag.isSet(BOUNDARY_SUBOBJ_SORTED)) {
-      MSG("WITH BOUNDARY SUBOBJ SORTED!\n");
       DofContainer dofs;
 
       for (int geo = FACE; geo >= VERTEX; geo--) {
-	boundaryDofInfo.nGeoDofs[static_cast<GeoIndex>(geo)] = 0;
+	std::set<const DegreeOfFreedom*> &dofSet =
+	  boundaryDofInfo.geoDofs[static_cast<GeoIndex>(geo)];
+	dofSet.clear();
 
 	for (InteriorBoundary::iterator it(myIntBoundary); !it.end(); ++it) {
 	  if (it->rankObj.subObj == geo) {
 	    dofs.clear();
 	    it->rankObj.el->getAllDofs(feSpace, it->rankObj, dofs);
-	   
-	    boundaryDofInfo.nGeoDofs[static_cast<GeoIndex>(geo)] += dofs.size();
-	    DofContainer &tmp = sendDofs[it.getRank()];
+
+	    DofContainer& tmp = sendDofs[it.getRank()];
 	    tmp.insert(tmp.end(), dofs.begin(), dofs.end());
+	    dofSet.insert(dofs.begin(), dofs.end());
 	  }
 	}
       }
@@ -1673,7 +1674,6 @@ namespace AMDiS {
 	    it->rankObj.el->getAllDofs(feSpace, it->rankObj, 
 				       recvDofs[it.getRank()]);
     } else {
-      MSG("GAAAAAAANZ NORMAL!\n");
       for (InteriorBoundary::iterator it(myIntBoundary); !it.end(); ++it)
 	it->rankObj.el->getAllDofs(feSpace, it->rankObj, 
 				   sendDofs[it.getRank()]);

AMDiS/src/parallel/MeshDistributor.h

@@ -48,7 +48,7 @@ namespace AMDiS {
 
   struct BoundaryDofInfo
   {
-    map<GeoIndex, int> nGeoDofs;
+    map<GeoIndex, DofContainerSet> geoDofs;
   };
   
 
@@ -305,6 +305,11 @@ namespace AMDiS {
       createBoundaryDofFlag = flag;
     }
 
+    BoundaryDofInfo& getBoundaryDofInfo()
+    {
+      return boundaryDofInfo;
+    }
+
   protected:
     /** \brief
      * Determines the interior boundaries, i.e. boundaries between ranks, and stores

AMDiS/src/parallel/PetscSolverSchur.cc

@@ -44,15 +44,57 @@ namespace AMDiS {
       }
     }
       
+
+
     nBoundaryDofs = boundaryDofs.size();
     mpi::getDofNumbering(mpiComm, nBoundaryDofs, 
 			 rStartBoundaryDofs, nOverallBoundaryDofs);
 
-    int counter = rStartBoundaryDofs;
+
+    DofContainerSet& edgeDofs = meshDistributor->getBoundaryDofInfo().geoDofs[EDGE];
+    DofContainerSet& vertexDofs = meshDistributor->getBoundaryDofInfo().geoDofs[VERTEX];
+    int nEdgeDofs = edgeDofs.size();
+    int nVertexDofs = vertexDofs.size();
+
+    TEST_EXIT_DBG(nEdgeDofs + nVertexDofs == nBoundaryDofs)
+      ("Should not happen!\n");
+
+    int rStartEdgeDofs, nOverallEdgeDofs;
+    mpi::getDofNumbering(mpiComm, nEdgeDofs, 
+			 rStartEdgeDofs, nOverallEdgeDofs);
+
+    int rStartVertexDofs, nOverallVertexDofs;
+    mpi::getDofNumbering(mpiComm, nVertexDofs, 
+			 rStartVertexDofs, nOverallVertexDofs);
+
+    TEST_EXIT_DBG(nOverallEdgeDofs + nOverallVertexDofs == nOverallBoundaryDofs)
+      ("Should not happen!\n");
+    
+
     mapGlobalBoundaryDof.clear();
+#if 1
+    {
+      int counter = rStartEdgeDofs;
+      for (DofContainerSet::iterator it = edgeDofs.begin();
+	   it != edgeDofs.end(); ++it)
+	mapGlobalBoundaryDof[meshDistributor->mapLocalToGlobal(**it)] = 
+	  counter++;
+    }
+    {
+      int counter = nOverallEdgeDofs + rStartVertexDofs;
+      for (DofContainerSet::iterator it = vertexDofs.begin();
+	   it != vertexDofs.end(); ++it)
+	mapGlobalBoundaryDof[meshDistributor->mapLocalToGlobal(**it)] = 
+	  counter++;
+    }
+#else
+ {
+    int counter = rStartBoundaryDofs;
     for (std::set<DegreeOfFreedom>::iterator it = boundaryDofs.begin();
 	 it != boundaryDofs.end(); ++it)
       mapGlobalBoundaryDof[*it] = counter++;
+ }
+#endif
 
 
 
@@ -78,11 +120,13 @@ namespace AMDiS {
     mpi::getDofNumbering(mpiComm, nInteriorDofs, 
 			 rStartInteriorDofs, nOverallInteriorDofs);
 
-    counter = rStartInteriorDofs;
-    mapGlobalInteriorDof.clear();
-    for (std::set<DegreeOfFreedom>::iterator it = interiorDofs.begin();
-	 it != interiorDofs.end(); ++it)
-      mapGlobalInteriorDof[*it] = counter++;
+    {
+      int counter = rStartInteriorDofs;
+      mapGlobalInteriorDof.clear();
+      for (std::set<DegreeOfFreedom>::iterator it = interiorDofs.begin();
+	   it != interiorDofs.end(); ++it)
+	mapGlobalInteriorDof[*it] = counter++;
+    }
 
 
     TEST_EXIT_DBG(nInteriorDofs > 0)("Should not happen!\n");
@@ -394,7 +438,7 @@ namespace AMDiS {
       if (rankOnly && !meshDistributor->getIsRankDof(dofIt.getDOFIndex()))
 	continue;
 
-      // Calculate global row index of the dof.
+      // Calculate global row index of the DOF.
       DegreeOfFreedom globalRowDof = 
 	meshDistributor->mapLocalToGlobal(dofIt.getDOFIndex());