]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Prepare for slepc-3.1 release by making necessary trivial changes to slepcwrappers
authorToby D. Young <tyoung@ippt.pan.pl>
Fri, 9 Jul 2010 19:53:28 +0000 (19:53 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Fri, 9 Jul 2010 19:53:28 +0000 (19:53 +0000)
git-svn-id: https://svn.dealii.org/trunk@21467 0785d39b-7218-0410-832d-ea1e28bc413d

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

index afd0a37427fbd229d76f99a7fe00511d353230f0..a5ae859aeb06dbe6fee217996c8dfef9ac52d2f2 100644 (file)
@@ -44,14 +44,10 @@ DEAL_II_NAMESPACE_OPEN
  * where $A$ is a system matrix, $M$ 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 and therefore, most of the
- * methods offered by the SLEPc library are projection methods or
- * other methods with similar properties. SLEPc implements these basic
- * methods as well as more sophisticated algorithms. On the other
- * hand, SLEPc is a general library in the sense that it covers
- * standard and generalized eigenvalue problems, and wrappers are
- * provided to interface to SLEPc solvers that handle both of these
- * problem sets.
+ * which the associated matrices are sparse. Most of the methods
+ * offered by the SLEPc library are projection methods or other
+ * methods with similar properties; and wrappers are provided to
+ * interface to SLEPc solvers that handle both of these problem sets.
  *
  * SLEPcWrappers can be implemented in application codes in the
  * following way:
@@ -59,14 +55,12 @@ DEAL_II_NAMESPACE_OPEN
   SolverControl solver_control (1000, 1e-9);
   SolverArnoldi system (solver_control,
                         mpi_communicator);
-  system.solve (A, B, eigenvalues, eigenvectors,
-                size_of_spectrum);
+  system.solve (A, M, lambda, x, size_of_spectrum);
  @endverbatim
-
- * for the generalized eigenvalue problem $Ax=B\lambda x$, where the
+ * for the generalized eigenvalue problem $Ax=M\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 step-36 for a hands-on example.
+ * also <code>step-36</code> for a hands-on example.
  *
  * An alternative implementation to the one above is to use the API
  * internals directly within the application code. In this way the
@@ -91,11 +85,13 @@ DEAL_II_NAMESPACE_OPEN
  * "PETScWrappers", on which they depend.
  *
  * @ingroup SLEPcWrappers
- * @author Toby D. Young 2008, 2009, 2010; and Rickard Armiento 2008
  *
- * Various tweaks to the SLEPcWrappers have been contributed by Eloy
- * Romeoro and Jose Roman.
+ * @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.
  */
+
 namespace SLEPcWrappers
 {
 
@@ -210,14 +206,15 @@ namespace SLEPcWrappers
                                     * solver.
                                     */
       void
-       set_initial_vector (const PETScWrappers::VectorBase &initial_vec);
+       set_initial_vector 
+       (const PETScWrappers::VectorBase &set_initial_vector);
 
                                    /**
                                    * Set the spectral transformation
                                    * to be used.
                                     */
       void
-       set_transformation (SLEPcWrappers::TransformationBase &trans);
+       set_transformation (SLEPcWrappers::TransformationBase &set_transformation);
 
                                    /**
                                    * Indicate which part of the
@@ -332,9 +329,9 @@ namespace SLEPcWrappers
 
       const PETScWrappers::MatrixBase   *opA;
       const PETScWrappers::MatrixBase   *opB;
-      const PETScWrappers::VectorBase   *ini_vec;
+      const PETScWrappers::VectorBase   *initial_vector;
 
-      SLEPcWrappers::TransformationBase *transform;
+      SLEPcWrappers::TransformationBase *transformation;
 
     private:
 
@@ -594,8 +591,7 @@ namespace SLEPcWrappers
       vr.resize (n_converged, vr.front());
       kr.resize (n_converged);
 
-      for (unsigned int index=0; index < n_converged;
-          ++index)
+      for (unsigned int index=0; index < n_converged; ++index)
        get_eigenpair (index, kr[index], vr[index]);
     }
 }
index 49e60c9f56bfbc426cf1bf21b6a652dbc20e13d1..d8da4f5b9bbc3cad356b7e9ec74613cbcd0dc4e9 100644 (file)
@@ -46,10 +46,9 @@ namespace SLEPcWrappers
     mpi_communicator (mpi_communicator),
     set_which (EPS_LARGEST_MAGNITUDE),
     opA (NULL), opB (NULL),
-    ini_vec (NULL),
-    transform (NULL)
-  {
-  }
+    initial_vector (NULL),
+    transformation (NULL)
+  {}
 
   SolverBase::~SolverBase ()
   {}
@@ -70,15 +69,15 @@ namespace SLEPcWrappers
   }
 
   void
-  SolverBase::set_initial_vector (const PETScWrappers::VectorBase &initial_vec)
+  SolverBase::set_initial_vector (const PETScWrappers::VectorBase &set_initial_vector)
   {
-    ini_vec = &initial_vec;
+    initial_vector = (&set_initial_vector);
   }
 
   void
-  SolverBase::set_transformation (SLEPcWrappers::TransformationBase &trans)
+  SolverBase::set_transformation (SLEPcWrappers::TransformationBase &set_transformation)
   {
-    transform = &trans;
+    transformation = &set_transformation;
   }
 
   void
@@ -107,14 +106,21 @@ namespace SLEPcWrappers
       ierr = EPSSetOperators (solver_data->eps, *opA, PETSC_NULL);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
-    if (ini_vec && ini_vec->size() != 0)
+    if (initial_vector && initial_vector->size() != 0)
       {
-       ierr = EPSSetInitialVector(solver_data->eps, *ini_vec);
+
+#if DEAL_II_PETSC_VERSION_LT(3,1,0)
+       ierr = EPSSetInitialVector (solver_data->eps, *initial_vector);
+#else
+       Vec this_vector = *initial_vector;
+       ierr = EPSSetInitialSpace (solver_data->eps, 1, &this_vector);
+#endif
+
        AssertThrow (ierr == 0, ExcSLEPcError(ierr));
       }
 
-    if (transform)
-      transform->set_context(solver_data->eps);
+    if (transformation)
+      transformation->set_context(solver_data->eps);
 
                                     // set runtime options.
     set_solver_type (solver_data->eps);
@@ -128,6 +134,9 @@ namespace SLEPcWrappers
                             PETSC_DECIDE, PETSC_DECIDE);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
+                                    // set the solve options to the
+                                    // eigenvalue problem solver
+                                    // context
     ierr = EPSSetFromOptions (solver_data->eps);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 

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.