]> https://gitweb.dealii.org/ - dealii.git/commitdiff
A simple slepc 3.1 compatibility patch to the SLPEcWrappers: Also added working David...
authorToby D. Young <tyoung@ippt.pan.pl>
Wed, 4 Aug 2010 14:06:34 +0000 (14:06 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Wed, 4 Aug 2010 14:06:34 +0000 (14:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@21609 0785d39b-7218-0410-832d-ea1e28bc413d

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

index bdc08e0408510c31e1d8d72777c84906cf4ca28b..3accb6229a71cede3a0214f6d2575da618694533 100644 (file)
@@ -206,14 +206,14 @@ namespace SLEPcWrappers
                                     */
       void
        set_initial_vector 
-       (const PETScWrappers::VectorBase &set_initial_vector);
+       (const PETScWrappers::VectorBase &this_initial_vector);
 
                                    /**
                                    * Set the spectral transformation
                                    * to be used.
                                     */
       void
-       set_transformation (SLEPcWrappers::TransformationBase &set_transformation);
+       set_transformation (SLEPcWrappers::TransformationBase &this_transformation);
 
                                    /**
                                    * Indicate which part of the
index f48a0b30cf28ad840413f8de8b64cc7e2ca581a7..4d13adb12477a17b45b206da3bc9d8b460c26ecb 100644 (file)
@@ -71,15 +71,15 @@ namespace SLEPcWrappers
   }
 
   void
-  SolverBase::set_initial_vector (const PETScWrappers::VectorBase &set_initial_vector)
+  SolverBase::set_initial_vector (const PETScWrappers::VectorBase &this_initial_vector)
   {
-    initial_vector = (&set_initial_vector);
+    initial_vector = (&this_initial_vector);
   }
 
   void
-  SolverBase::set_transformation (SLEPcWrappers::TransformationBase &set_transformation)
+  SolverBase::set_transformation (SLEPcWrappers::TransformationBase &this_transformation)
   {
-    transformation = &set_transformation;
+    transformation = &this_transformation;
   }
 
   void
@@ -175,8 +175,8 @@ namespace SLEPcWrappers
     AssertThrow (solver_data.get() != 0, ExcSLEPcWrappersUsageError());
 
                                     // get converged eigenpair
-    int ierr = EPSGetEigenpair(solver_data->eps, index,
-                              &kr, PETSC_NULL, vr, PETSC_NULL);
+    int ierr = EPSGetEigenpair (solver_data->eps, index,
+                               &kr, PETSC_NULL, vr, PETSC_NULL);
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 
@@ -193,8 +193,9 @@ namespace SLEPcWrappers
   EPS *
   SolverBase::get_EPS ()
   {
-    if( solver_data.get() == 0 )
+    if (solver_data.get() == 0)
       return NULL;
+
     return &solver_data->eps;
   }
 
@@ -321,8 +322,8 @@ namespace SLEPcWrappers
                                     // number of iteration steps to
                                     // the SLEPc convergence
                                     // criterion.
-    ierr = EPSSetTolerances(eps, this->solver_control.tolerance(),
-                           this->solver_control.max_steps());
+    ierr = EPSSetTolerances (eps, this->solver_control.tolerance(),
+                            this->solver_control.max_steps());
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 
@@ -348,8 +349,8 @@ namespace SLEPcWrappers
                                     // number of iteration steps to
                                     // the SLEPc convergence
                                     // criterion.
-    ierr = EPSSetTolerances(eps, this->solver_control.tolerance(),
-                           this->solver_control.max_steps());
+    ierr = EPSSetTolerances (eps, this->solver_control.tolerance(),
+                            this->solver_control.max_steps());
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 
@@ -368,7 +369,7 @@ namespace SLEPcWrappers
   SolverDavidson::set_solver_type (EPS &eps) const
   {
     int ierr;
-    ierr = EPSSetType (eps, const_cast<char *>(EPSDAVIDSON));
+    ierr = EPSSetType (eps, const_cast<char *>(EPSGD));
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
 
                                     // hand over the absolute
@@ -376,8 +377,8 @@ namespace SLEPcWrappers
                                     // number of iteration steps to
                                     // the SLEPc convergence
                                     // criterion.
-    ierr = EPSSetTolerances(eps, this->solver_control.tolerance(),
-                           this->solver_control.max_steps());
+    ierr = EPSSetTolerances (eps, this->solver_control.tolerance(),
+                            this->solver_control.max_steps());
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 #endif
index 7545234a8a9e1ca1aa5cdf4d3c91bcd9063025cf..beb45d6bfd65a35e63e18c0368559fad5ae4ab19 100644 (file)
@@ -32,16 +32,13 @@ DEAL_II_NAMESPACE_OPEN
 namespace SLEPcWrappers
 {
   TransformationBase::TransformationData::~TransformationData ()
-  {
-  }
+  {}
 
   TransformationBase::TransformationBase ()
-  {
-  }
+  {}
 
   TransformationBase::~TransformationBase ()
-  {
-  }
+  {}
 
   void TransformationBase::set_context (EPS &eps)
   {
@@ -56,6 +53,7 @@ namespace SLEPcWrappers
   }
 
   /* ------------------- TransformationShift --------------------- */
+
   TransformationShift::AdditionalData::
   AdditionalData (const double shift_parameter)
                   :
@@ -79,6 +77,7 @@ namespace SLEPcWrappers
   }
 
   /* ---------------- TransformationShiftInvert ------------------ */
+
   TransformationShiftInvert::AdditionalData::
   AdditionalData (const double shift_parameter)
                   :
@@ -94,7 +93,11 @@ namespace SLEPcWrappers
   TransformationShiftInvert::set_transformation_type (ST &st) const
   {
     int ierr;
+#if DEAL_II_PETSC_VERSION_LT(3,1,0)
     ierr = STSetType (st, const_cast<char *>(STSINV));
+#else
+    ierr = STSetType (st, const_cast<char *>(STSINVERT));
+#endif
     AssertThrow (ierr == 0, SolverBase::ExcSLEPcError(ierr));
 
     ierr = STSetShift (st, additional_data.shift_parameter);
@@ -102,6 +105,7 @@ namespace SLEPcWrappers
   }
 
   /* --------------- TransformationSpectrumFolding ----------------- */
+
   TransformationSpectrumFolding::AdditionalData::
   AdditionalData (const double shift_parameter)
                   :
@@ -126,6 +130,7 @@ namespace SLEPcWrappers
   }
 
   /* ------------------- TransformationCayley --------------------- */
+
   TransformationCayley::AdditionalData::
   AdditionalData (const double shift_parameter,
                  const double antishift_parameter)

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.