From 27b9eae66d80748d7d3086a58b7e09ff24e97b16 Mon Sep 17 00:00:00 2001 From: "Toby D. Young" Date: Fri, 9 Jul 2010 19:53:28 +0000 Subject: [PATCH] Prepare for slepc-3.1 release by making necessary trivial changes to slepcwrappers git-svn-id: https://svn.dealii.org/trunk@21467 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/slepc_solver.h | 40 ++++++++++++-------------- deal.II/lac/source/slepc_solver.cc | 33 +++++++++++++-------- 2 files changed, 39 insertions(+), 34 deletions(-) diff --git a/deal.II/lac/include/lac/slepc_solver.h b/deal.II/lac/include/lac/slepc_solver.h index afd0a37427..a5ae859aeb 100644 --- a/deal.II/lac/include/lac/slepc_solver.h +++ b/deal.II/lac/include/lac/slepc_solver.h @@ -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 const unsigned int size_of_spectrum tells * SLEPc the number of eigenvector/eigenvalue pairs to solve for: See - * also step-36 for a hands-on example. + * also step-36 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]); } } diff --git a/deal.II/lac/source/slepc_solver.cc b/deal.II/lac/source/slepc_solver.cc index 49e60c9f56..d8da4f5b9b 100644 --- a/deal.II/lac/source/slepc_solver.cc +++ b/deal.II/lac/source/slepc_solver.cc @@ -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)); -- 2.39.5