From bef826f882d16167438177eca3fa3c5d23535fe4 Mon Sep 17 00:00:00 2001 From: "Toby D. Young" Date: Wed, 21 Jul 2010 15:07:15 +0000 Subject: [PATCH] Update for SLEPcWrappers to slepc-dev (to be slepc-3.1) git-svn-id: https://svn.dealii.org/trunk@21552 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/slepc_solver.h | 54 ++++++++++++++++++- .../lac/slepc_spectral_transformation.h | 4 +- deal.II/lac/source/slepc_solver.cc | 31 +++++++++++ 3 files changed, 84 insertions(+), 5 deletions(-) diff --git a/deal.II/lac/include/lac/slepc_solver.h b/deal.II/lac/include/lac/slepc_solver.h index bdffad50d3..bdc08e0408 100644 --- a/deal.II/lac/include/lac/slepc_solver.h +++ b/deal.II/lac/include/lac/slepc_solver.h @@ -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) diff --git a/deal.II/lac/include/lac/slepc_spectral_transformation.h b/deal.II/lac/include/lac/slepc_spectral_transformation.h index 736d024e15..82ebf67162 100644 --- a/deal.II/lac/include/lac/slepc_spectral_transformation.h +++ b/deal.II/lac/include/lac/slepc_spectral_transformation.h @@ -45,9 +45,7 @@ namespace SLEPcWrappers public: /** - * Constructor. Takes the MPI - * communicator over which parallel - * computations are to happen. + * Constructor. */ TransformationBase (); diff --git a/deal.II/lac/source/slepc_solver.cc b/deal.II/lac/source/slepc_solver.cc index 930e874065..10a0b223d0 100644 --- a/deal.II/lac/source/slepc_solver.cc +++ b/deal.II/lac/source/slepc_solver.cc @@ -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(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 -- 2.39.5