]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add solver to SLEPcWrappers and trivially update documentation
authorToby D. Young <tyoung@ippt.pan.pl>
Wed, 21 Jul 2010 10:21:35 +0000 (10:21 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Wed, 21 Jul 2010 10:21:35 +0000 (10:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@21550 0785d39b-7218-0410-832d-ea1e28bc413d

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

index a5ae859aeb06dbe6fee217996c8dfef9ac52d2f2..bdffad50d3263cb1d541ef95ffb9b0442e7985ae 100644 (file)
@@ -35,13 +35,13 @@ DEAL_II_NAMESPACE_OPEN
  * Base class for solver classes using the SLEPc solvers which are
  * selected based on flags passed to the eigenvalue problem solver
  * context. Derived classes set the right flags to set the right
- * solver. On the other hand, note that: the
- * <code>AdditionalData</code> structure is a dummy structure that
- * currently exists for backward/forward compatibility only.
+ * solver. Note that: the <code>AdditionalData</code> structure is a
+ * dummy structure that currently exists for backward/forward
+ * compatibility only.
  *
  * The SLEPc solvers are intended to be used for solving the
- * generalized eigenspectrum problem $(A-\lambda M)x=0$, for $x\neq0$;
- * where $A$ is a system matrix, $M$ is a mass matrix, and $\lambda,
+ * generalized eigenspectrum problem $(A-\lambda B)x=0$, for $x\neq0$;
+ * where $A$ is a system matrix, $B$ is a mass matrix, and $\lambda,
  * x$ are a set of eigenvalues and eigenvectors respectively. The
  * emphasis is on methods and techniques appropriate for problems in
  * which the associated matrices are sparse. Most of the methods
@@ -53,11 +53,10 @@ DEAL_II_NAMESPACE_OPEN
  * following way:
  @verbatim
   SolverControl solver_control (1000, 1e-9);
-  SolverArnoldi system (solver_control,
-                        mpi_communicator);
-  system.solve (A, M, lambda, x, size_of_spectrum);
+  SolverArnoldi system (solver_control, mpi_communicator);
+  system.solve (A, B, lambda, x, size_of_spectrum);
  @endverbatim
- * for the generalized eigenvalue problem $Ax=M\lambda x$, where the
+ * for the generalized eigenvalue problem $Ax=B\lambda x$, where the
  * variable <code>const unsigned int size_of_spectrum</code> tells
  * SLEPc the number of eigenvector/eigenvalue pairs to solve for: See
  * also <code>step-36</code> for a hands-on example.
@@ -89,7 +88,7 @@ DEAL_II_NAMESPACE_OPEN
  * @author Toby D. Young 2008, 2009, 2010; and Rickard Armiento 2008.
  *
  * @note Various tweaks and enhancments contributed by Eloy Romero and
- * Jose E.  Roman 2009, 2010.
+ * Jose E. Roman 2009, 2010.
  */
 
 namespace SLEPcWrappers
@@ -320,10 +319,15 @@ namespace SLEPcWrappers
                                     */
       virtual void set_solver_type (EPS &eps) const = 0;
 
-                                   /**
+                                   /*
                                     * Attributes that store the
                                    * relevant information for the
-                                   * eigenproblem solver context.
+                                   * eigenproblem solver
+                                   * context. These are: which
+                                   * portion of the spectrum to solve
+                                   * from; and the matrices and
+                                   * vectors that discribe the
+                                   * problem.
                                     */
       EPSWhich                           set_which;
 
@@ -331,6 +335,12 @@ namespace SLEPcWrappers
       const PETScWrappers::MatrixBase   *opB;
       const PETScWrappers::VectorBase   *initial_vector;
 
+                                   /**
+                                    * Pointer to an an object that
+                                    * describes transformations that
+                                    * can be applied to the eigenvalue
+                                    * problem.
+                                    */
       SLEPcWrappers::TransformationBase *transformation;
 
     private:
@@ -535,6 +545,56 @@ namespace SLEPcWrappers
 
     };
 
+/**
+ * An implementation of the solver interface using the SLEPc Lanczos
+ * solver. Usage: Largest values of spectrum only, all problem types,
+ * complex.
+ *
+ * @ingroup SLEPcWrappers
+ * @author Toby D. Young 2010
+ */
+  class SolverPower : 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.
+                                    */
+      SolverPower (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;
+
+    };
+
 
   // --------------------------- inline and template functions -----------
   /**
index d8da4f5b9bbc3cad356b7e9ec74613cbcd0dc4e9..930e874065e295a8d999b2364c006ed8cc3a6e23 100644 (file)
@@ -56,6 +56,7 @@ namespace SLEPcWrappers
   void
   SolverBase::set_matrices (const PETScWrappers::MatrixBase &A)
   {
+                                   // standard eigenspectrum problem
     opA = &A;
     opB = NULL;
   }
@@ -64,6 +65,7 @@ namespace SLEPcWrappers
   SolverBase::set_matrices (const PETScWrappers::MatrixBase &A,
                            const PETScWrappers::MatrixBase &B)
   {
+                                   // generalized eigenspectrum problem
     opA = &A;
     opB = &B;
   }
@@ -95,10 +97,12 @@ namespace SLEPcWrappers
     solver_data.reset (new SolverData());
 
                                     // create eigensolver context and
-                                    // set operators.
+                                    // set operators
     ierr = EPSCreate (mpi_communicator, &solver_data->eps);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
+                                    // set eigenspectrum problem type
+                                    // (general/standard)
     AssertThrow (opA, ExcSLEPcWrappersUsageError());
     if (opB)
       ierr = EPSSetOperators (solver_data->eps, *opA, *opB);
@@ -106,6 +110,7 @@ namespace SLEPcWrappers
       ierr = EPSSetOperators (solver_data->eps, *opA, PETSC_NULL);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
+                                    // set the initial vector(s) if any
     if (initial_vector && initial_vector->size() != 0)
       {
 
@@ -119,12 +124,15 @@ namespace SLEPcWrappers
        AssertThrow (ierr == 0, ExcSLEPcError(ierr));
       }
 
+                                    // set transformation type if any
     if (transformation)
       transformation->set_context(solver_data->eps);
 
-                                    // set runtime options.
+                                    // set runtime options
     set_solver_type (solver_data->eps);
 
+                                    // set which portion of the
+                                    // eigenspectrum to solve for
     ierr = EPSSetWhichEigenpairs (solver_data->eps, set_which);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
@@ -316,6 +324,33 @@ namespace SLEPcWrappers
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 
+  /* ----------------------- Power ------------------------- */
+
+  SolverPower::SolverPower (SolverControl        &cn,
+                           const MPI_Comm       &mpi_communicator,
+                           const AdditionalData &data)
+    :
+    SolverBase (cn, mpi_communicator),
+    additional_data (data)
+  {}
+
+  void
+  SolverPower::set_solver_type (EPS &eps) const
+  {
+    int ierr;
+    ierr = EPSSetType (eps, const_cast<char *>(EPSPOWER));
+    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));
+  }
+
 }
 
 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.