Petsc-solvers added to Creator-Map, PetscSolver is now based on OEMSolver

AMDiS/src/CreatorMap.cc

@@ -29,6 +29,16 @@
 #include "time/RosenbrockMethod.h"
 #include "nonlin/NonLinSolver.h"
 
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+  #include "parallel/PetscSolverFeti.h"
+  #include "parallel/PetscSolverGlobalBlockMatrix.h"
+  #include "parallel/PetscSolverGlobalMatrix.h"
+  #include "parallel/PetscSolverNavierStokes.h"
+  #include "parallel/PetscSolverSchur.h"
+  #include "parallel/BddcMlSolver.h"
+#endif
+
+
 namespace AMDiS {
 
   template<>
@@ -73,6 +83,17 @@ namespace AMDiS {
     creator = new UmfPackSolver::Creator;
     addCreator("umfpack", creator);
 #endif
+    
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+    addCreator("petsc", new PetscSolverGlobalMatrix::Creator);
+    addCreator("petsc-schur", new PetscSolverSchur::Creator);
+    addCreator("petsc-block", new PetscSolverGlobalBlockMatrix::Creator);
+    addCreator("petsc-feti", new PetscSolverFeti::Creator);
+    #ifdef HAVE_BDDC_ML
+      addCreator("bddcml", new BddcMlSolver::Creator);
+    #endif
+    addCreator("petsc-stokes", new PetscSolverNavierStokes::Creator);
+#endif
   }
 
 

AMDiS/src/DOFVector.cc

@@ -770,7 +770,7 @@ namespace AMDiS {
     Quadrature* quad =
       Quadrature::provideQuadrature(vecs[0]->getFeSpace()->getMesh()->getDim(), deg);
     FastQuadrature *fastQuad =
-      FastQuadrature::provideFastQuadrature(vecs[0]->getFeSpace()->getBasisFcts(), *quad, INIT_PHI);
+      FastQuadrature::provideFastQuadrature(vecs[0]->getFeSpace()->getBasisFcts(), *quad,  INIT_PHI | INIT_GRD_PHI);
 
     std::vector<mtl::dense_vector<double> > qp(vecs.size());
     std::vector<mtl::dense_vector<WorldVector<double> > > qpGrd(vecs.size());
@@ -782,7 +782,7 @@ namespace AMDiS {
       qpGrd[i].change_dim(fastQuad->getNumPoints());
 
     double value = 0.0;
-    Flag traverseFlag = Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET;
+    Flag traverseFlag = Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET | Mesh::FILL_GRD_LAMBDA;
     TraverseStack stack;
     ElInfo *elInfo = stack.traverseFirst(vecs[0]->getFeSpace()->getMesh(), -1, traverseFlag);
     while (elInfo) {

AMDiS/src/ElInfo2d.cc

@@ -270,7 +270,7 @@ namespace AMDiS {
 		  neighbourCoord[i][0] = macroNeighbour->coord[0];
 		  neighbourCoord[i][1] = macroNeighbour->coord[2];		 
 		} else {
-		  ERROR_EXIT("Should not happen!\n");
+		  ERROR_EXIT("Should not happen! Non-conforming AMDiS-mesh? Periodic mesh with corrected index-circles?\n");
 		}
 
 		// I've deleted here some code, be I think that this case is not

AMDiS/src/FixVec.h

@@ -638,6 +638,11 @@ namespace AMDiS {
     typedef T1 type;
   };
   
+  template<> struct ProductType<double, int> { typedef double type; };
+  template<> struct ProductType<int, double> { typedef double type; };
+  template<> struct ProductType<float, int> { typedef float type; };
+  template<> struct ProductType<int, float> { typedef float type; };
+  
   // scalar times vector
   template<typename T1>
   struct ProductType<T1, WorldVector<T1> >

AMDiS/src/ProblemStat.cc

@@ -507,7 +507,6 @@ namespace AMDiS {
   {
     FUNCNAME("ProblemStat::createSolver()");
 
-#ifndef HAVE_PARALLEL_DOMAIN_AMDIS
     // === create solver ===
     string solverType("0");
     string initFileStr = name + "->solver";
@@ -518,7 +517,6 @@ namespace AMDiS {
     solverCreator->setName(initFileStr);
     solver = solverCreator->create();
     solver->initParameters();
-#endif
   }
 
 

AMDiS/src/parallel/BddcMlSolver.h

@@ -35,6 +35,19 @@ namespace AMDiS {
   class BddcMlSolver : public PetscSolver
   {
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new BddcMlSolver(this->name); 
+      }
+    };
+    
     BddcMlSolver(string name)
       : PetscSolver(name),
 	rhsVec(NULL),

AMDiS/src/parallel/ParallelCoarseSpaceSolver.cc

@@ -18,8 +18,10 @@ namespace AMDiS {
 
   using namespace std;
 
-  ParallelCoarseSpaceSolver::ParallelCoarseSpaceSolver()
-    : interiorMap(NULL),
+  ParallelCoarseSpaceSolver::ParallelCoarseSpaceSolver(string name)
+    : OEMSolver(name), 
+      initFileStr(name),
+      interiorMap(NULL),
       lastMeshNnz(-1),
       alwaysCreateNnzStructure(false),
       meshDistributor(NULL),

AMDiS/src/parallel/ParallelCoarseSpaceSolver.h

@@ -31,9 +31,10 @@
 #include "parallel/ParallelDofMapping.h"
 #include "parallel/MeshDistributor.h"
 #include "parallel/MatrixNnzStructure.h"
+#include "OEMSolver.h"
 
 namespace AMDiS {
-
+  
   /**
    * This class implements a block structured PETSc matrix/vec which seperates
    * the discretization of the interior of subdomains and the discretization 
@@ -48,10 +49,10 @@ namespace AMDiS {
    * - This class also manages the creation of the corresponding non zero 
    *   structure of the matrices.
    */
-  class ParallelCoarseSpaceSolver {
+  class ParallelCoarseSpaceSolver : public OEMSolver {
   public:
     /// Constructor
-    ParallelCoarseSpaceSolver();
+    ParallelCoarseSpaceSolver(string name);
 
     /// Set parallel DOF mapping for the interior DOFs.
     void setDofMapping(ParallelDofMapping *interiorDofs)
@@ -305,6 +306,9 @@ namespace AMDiS {
     bool alwaysCreateNnzStructure;
 
   protected:
+    /// Prefix string for parameters in init file.
+    string initFileStr;
+    
     /// Pointer to a mesh distributor object.
     MeshDistributor *meshDistributor;
 
@@ -340,6 +344,7 @@ namespace AMDiS {
     /// different coarse spaces for different components.
     map<int, ParallelDofMapping*> coarseSpaceMap;   
   };
+  
 }
 
 #endif

AMDiS/src/parallel/PetscProblemStat.cc

@@ -32,33 +32,6 @@ namespace AMDiS {
   {
     FUNCNAME("PetscProblemStat::PetscProblemStat()");
 
-    string tmp("");
-    string initFileStr = nameStr + "->solver";
-    Parameters::get(initFileStr, tmp);
-   
-    if (tmp == "petsc-schur") {
-      petscSolver = new PetscSolverSchur(initFileStr);
-    } else if (tmp == "petsc-feti") {
-      petscSolver = new PetscSolverFeti(initFileStr);
-    } else if (tmp == "petsc-block") {
-      petscSolver = new PetscSolverGlobalBlockMatrix(initFileStr);
-    } else if (tmp == "petsc") {
-      petscSolver = new PetscSolverGlobalMatrix(initFileStr);
-    } else if (tmp == "bddcml") {
-#ifdef HAVE_BDDC_ML
-      petscSolver = new BddcMlSolver(initFileStr);
-#else
-      ERROR_EXIT("AMDiS was compiled without BDDC-ML support!\n");
-#endif
-    } else if (tmp == "petsc-navierstokes") {
-      petscSolver = new PetscSolverNavierStokes(initFileStr);
-    } else {
-      ERROR_EXIT("No parallel solver %s available!\n", tmp.c_str());
-    }
-
-    tmp = "";
-    Parameters::get(nameStr + "->solver->petsc prefix", tmp);
-    petscSolver->setKspPrefix(tmp);    
   }
 
 
@@ -76,6 +49,10 @@ namespace AMDiS {
   {
     ParallelProblemStatBase::initialize(initFlag, adoptProblem, adoptFlag);
     
+    string tmp("");
+    Parameters::get(getName() + "->solver->petsc prefix", tmp);
+    petscSolver->setKspPrefix(tmp);    
+    
     meshDistributor->setBoundaryDofRequirement(petscSolver->getBoundaryDofRequirement());
 
     petscSolver->setMeshDistributor(meshDistributor, 
@@ -83,6 +60,16 @@ namespace AMDiS {
 				    PETSC_COMM_SELF);   
     petscSolver->init(this->getComponentSpaces(), this->getFeSpaces());
   }
+  
+    
+  void PetscProblemStat::createSolver()
+  {
+    FUNCNAME("createSolver()");
+    ParallelProblemStatBase::createSolver();
+    
+    petscSolver = dynamic_cast<PetscSolver*>(solver);
+    TEST_EXIT(petscSolver)("wrong petsc-solver type\n");
+  }
 
 
   void PetscProblemStat::solve(AdaptInfo *adaptInfo,

AMDiS/src/parallel/PetscProblemStat.h

@@ -52,6 +52,8 @@ namespace AMDiS {
     void initialize(Flag initFlag,
 		    ProblemStatSeq *adoptProblem = NULL,
 		    Flag adoptFlag = INIT_NOTHING);
+    
+    virtual void createSolver();
 
     void solve(AdaptInfo *adaptInfo,
 	       bool createMatrixData = true,

AMDiS/src/parallel/PetscSolver.cc

@@ -21,8 +21,7 @@ namespace AMDiS {
   using namespace std;
 
   PetscSolver::PetscSolver(string name)
-    : ParallelCoarseSpaceSolver(),
-      initFileStr(name),
+    : ParallelCoarseSpaceSolver(name),
       dofMap(FESPACE_WISE, true),
       dofMapSd(FESPACE_WISE, true),
       kspPrefix(""),

AMDiS/src/parallel/PetscSolver.h

@@ -49,7 +49,7 @@ namespace AMDiS {
    */
   class PetscSolver : public ParallelCoarseSpaceSolver
   {
-  public:
+  public:    
     PetscSolver(string name);
 
     virtual ~PetscSolver() {}
@@ -173,9 +173,6 @@ namespace AMDiS {
     /// Run test, if matrix is symmetric.
     bool testMatrixSymmetric(Mat mat, bool advancedTest = false);
   protected:
-    /// Prefix string for parameters in init file.
-    string initFileStr;
-
     /// FE spaces of all components for the stationary problem the specific
     /// solver object is registered to.
     vector<const FiniteElemSpace*> componentSpaces;
@@ -215,7 +212,6 @@ namespace AMDiS {
     vector<int> constNullspaceComponent;
   };
 
-
 } // namespace AMDiS
 
 #endif

AMDiS/src/parallel/PetscSolverFeti.h

@@ -41,6 +41,19 @@ namespace AMDiS {
   class PetscSolverFeti : public PetscSolver
   {
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new PetscSolverFeti(this->name); 
+      }
+    };
+    
     /// Constructor of FETI-DP solver class.
     PetscSolverFeti(string name);
 

AMDiS/src/parallel/PetscSolverGlobalBlockMatrix.cc

@@ -92,6 +92,8 @@ namespace AMDiS {
 
     matAssembly();
 
+    
+    /// initPreconditioner(...)
 
     // === Init PETSc solver and preconditioner objects. ===
 

AMDiS/src/parallel/PetscSolverGlobalBlockMatrix.h

@@ -34,6 +34,19 @@ namespace AMDiS {
   class PetscSolverGlobalBlockMatrix : public PetscSolver
   {
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new PetscSolverGlobalBlockMatrix(this->name); 
+      }
+    };
+    
     PetscSolverGlobalBlockMatrix(string name)
       : PetscSolver(name),
 	nComponents(0),

AMDiS/src/parallel/PetscSolverGlobalMatrix.h

@@ -39,6 +39,19 @@ namespace AMDiS {
   class PetscSolverGlobalMatrix : public PetscSolver
   {
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new PetscSolverGlobalMatrix(this->name); 
+      }
+    };
+    
     PetscSolverGlobalMatrix(string name);
 
     void fillPetscMatrix(Matrix<DOFMatrix*> *mat);

AMDiS/src/parallel/PetscSolverNavierStokes.h

@@ -39,7 +39,7 @@ namespace AMDiS {
 
   class PetscSolverNavierStokes : public PetscSolverGlobalMatrix
   {
-  private:
+  private:    
     class IdFct : public AbstractFunction<double, double>
     {
     public:
@@ -108,6 +108,18 @@ namespace AMDiS {
 
 
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new PetscSolverNavierStokes(this->name); 
+      }
+    };
     PetscSolverNavierStokes(string name);
 
     void setStokesData(double *nuPtr, double *invTauPtr, SystemVector *vec)

AMDiS/src/parallel/PetscSolverSchur.h

@@ -33,6 +33,19 @@ namespace AMDiS {
   class PetscSolverSchur : public PetscSolver
   {
   public:
+    /// Creator class
+    class Creator : public OEMSolverCreator
+    {
+    public:
+      virtual ~Creator() {}
+
+      /// Returns a new PetscSolver object.
+      OEMSolver* create() 
+      { 
+	return new PetscSolverSchur(this->name); 
+      }
+    };
+    
     PetscSolverSchur(string name) 
       : PetscSolver(name)
     {}