]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update for SLEPcWrappers to slepc-dev (to be slepc-3.1)
authorToby D. Young <tyoung@ippt.pan.pl>
Wed, 21 Jul 2010 15:07:15 +0000 (15:07 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Wed, 21 Jul 2010 15:07:15 +0000 (15:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@21552 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/slepc_solver.h
deal.II/lac/include/lac/slepc_spectral_transformation.h
deal.II/lac/source/slepc_solver.cc

index bdffad50d3263cb1d541ef95ffb9b0442e7985ae..bdc08e0408510c31e1d8d72777c84906cf4ca28b 100644 (file)
@@ -546,7 +546,7 @@ namespace SLEPcWrappers
     };
 
 /**
- * An implementation of the solver interface using the SLEPc Lanczos
+ * An implementation of the solver interface using the SLEPc Power
  * solver. Usage: Largest values of spectrum only, all problem types,
  * complex.
  *
@@ -595,6 +595,56 @@ namespace SLEPcWrappers
 
     };
 
+#if DEAL_II_PETSC_VERSION_GTE(3,1,0)
+/**
+ * An implementation of the solver interface using the SLEPc Davidson
+ * solver. Usage (incomplete/untested): All problem types.
+ *
+ * @ingroup SLEPcWrappers
+ * @author Toby D. Young 2010
+ */
+  class SolverDavidson : public SolverBase
+    {
+    public:
+                                   /**
+                                    * Standardized data struct to pipe
+                                    * additional data to the solver,
+                                    * should it be needed.
+                                    */
+      struct AdditionalData
+      {};
+
+                                   /**
+                                    * SLEPc solvers will want to have
+                                    * an MPI communicator context over
+                                    * which computations are
+                                    * parallelized. By default, this
+                                    * carries the same behaviour has
+                                    * the PETScWrappers, but you can
+                                    * change that.
+                                    */
+      SolverDavidson (SolverControl        &cn,
+                     const MPI_Comm       &mpi_communicator = PETSC_COMM_WORLD,
+                     const AdditionalData &data = AdditionalData());
+
+    protected:
+
+                                   /**
+                                    * Store a copy of the flags for
+                                    * this particular solver.
+                                    */
+      const AdditionalData additional_data;
+
+                                   /**
+                                    * Function that takes a Eigenvalue
+                                    * Problem Solver context object,
+                                    * and sets the type of solver that
+                                    * is appropriate for this class.
+                                    */
+      virtual void set_solver_type (EPS &eps) const;
+
+    };
+#endif
 
   // --------------------------- inline and template functions -----------
   /**
@@ -638,7 +688,7 @@ namespace SLEPcWrappers
     {
       unsigned int n_converged;
 
-      set_matrices (A,B);
+      set_matrices (A, B);
       solve (n_eigenvectors, &n_converged);
 
       if (n_converged >= n_eigenvectors)
index 736d024e15ab1fa17fe5d6cb03052cc13a2324a0..82ebf671623f5c9e91307b58d509ffa62f5ac652 100644 (file)
@@ -45,9 +45,7 @@ namespace SLEPcWrappers
     public:
 
                                    /**
-                                    * Constructor. Takes the MPI
-                                    * communicator over which parallel
-                                    * computations are to happen.
+                                    * Constructor. 
                                     */
       TransformationBase ();
 
index 930e874065e295a8d999b2364c006ed8cc3a6e23..10a0b223d0e3fc816e46631ecde2df9ccc01200e 100644 (file)
@@ -246,6 +246,7 @@ namespace SLEPcWrappers
   }
 
   /* ---------------------- SolverKrylovSchur ------------------------ */
+
   SolverKrylovSchur::SolverKrylovSchur (SolverControl        &cn,
                                        const MPI_Comm       &mpi_communicator,
                                        const AdditionalData &data)
@@ -272,6 +273,7 @@ namespace SLEPcWrappers
   }
 
   /* ---------------------- SolverArnoldi ------------------------ */
+
   SolverArnoldi::SolverArnoldi (SolverControl        &cn,
                                const MPI_Comm       &mpi_communicator,
                                const AdditionalData &data)
@@ -351,6 +353,35 @@ namespace SLEPcWrappers
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 
+  /* ---------------------- Davidson ----------------------- */
+
+#if DEAL_II_PETSC_VERSION_LT(3,1,0)
+  SolverDavidson::SolverDavidson (SolverControl        &cn,
+                                 const MPI_Comm       &mpi_communicator,
+                                 const AdditionalData &data)
+  :
+    SolverBase (cn, mpi_communicator),
+    additional_data (data)
+  {}
+
+  void
+  SolverDavidson::set_solver_type (EPS &eps) const
+  {
+    int ierr;
+    ierr = EPSSetType (eps, const_cast<char *>(EPSDAVIDSON));
+    AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+
+                                    // hand over the absolute
+                                    // tolerance and the maximum
+                                    // number of iteration steps to
+                                    // the SLEPc convergence
+                                    // criterion.
+    ierr = EPSSetTolerances(eps, this->solver_control.tolerance(),
+                           this->solver_control.max_steps());
+    AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+  }
+#endif
+
 }
 
 DEAL_II_NAMESPACE_CLOSE

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.