From 39d603fd9e9300c54985f3f9e3f814e4d6738337 Mon Sep 17 00:00:00 2001 From: young Date: Wed, 9 Jan 2013 20:52:56 +0000 Subject: [PATCH] Distinguish between Generalized and Jacobi Davidson solvers. Keep backward comparibility for deal.II 7.3 git-svn-id: https://svn.dealii.org/trunk@28000 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/slepc_solver.h | 77 ++++++++++++++++++++-- deal.II/source/lac/slepc_solver.cc | 37 +++++++++-- 2 files changed, 102 insertions(+), 12 deletions(-) diff --git a/deal.II/include/deal.II/lac/slepc_solver.h b/deal.II/include/deal.II/lac/slepc_solver.h index 6aa1c52cb9..ed85982912 100644 --- a/deal.II/include/deal.II/lac/slepc_solver.h +++ b/deal.II/include/deal.II/lac/slepc_solver.h @@ -619,13 +619,74 @@ 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. + * An implementation of the solver interface using the SLEPc + * Generalized Davidson solver. Usage: All problem types. + * + * @deprecated The use of this class has been superceded by the + * class SolverGeneralizedDavidson. * * @ingroup SLEPcWrappers - * @author Toby D. Young 2010 + * @author Toby D. Young 2011 + */ +#define SolverDavidson SolverGeneralizedDavidson; + + /** + * An implementation of the solver interface using the SLEPc + * Generalized Davidson solver. Usage: All problem types. + * + * @ingroup SLEPcWrappers + * @author Toby D. Young 2013 + */ + class SolverGeneralizedDavidson : 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. + */ + SolverGeneralizedDavidson (SolverControl &cn, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, + 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; + + }; + + /** + * An implementation of the solver interface using the SLEPc + * Jacobi-Davidson solver. Usage: All problem types. + * + * @ingroup SLEPcWrappers + * @author Toby D. Young 2013 */ - class SolverDavidson : public SolverBase + class SolverJacobiDavidson : public SolverBase { public: /** @@ -645,9 +706,9 @@ namespace SLEPcWrappers * the PETScWrappers, but you can * change that. */ - SolverDavidson (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, - const AdditionalData &data = AdditionalData()); + SolverJacobiDavidson (SolverControl &cn, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, + const AdditionalData &data = AdditionalData()); protected: @@ -668,6 +729,8 @@ namespace SLEPcWrappers }; #endif + + // --------------------------- inline and template functions ----------- /** * This is declared here to make it possible to take a std::vector diff --git a/deal.II/source/lac/slepc_solver.cc b/deal.II/source/lac/slepc_solver.cc index a0a5c6e865..85345bd9a1 100644 --- a/deal.II/source/lac/slepc_solver.cc +++ b/deal.II/source/lac/slepc_solver.cc @@ -377,18 +377,18 @@ namespace SLEPcWrappers } #if DEAL_II_PETSC_VERSION_GTE(3,1,0) - /* ---------------------- Davidson ----------------------- */ - SolverDavidson::SolverDavidson (SolverControl &cn, - const MPI_Comm &mpi_communicator, - const AdditionalData &data) + /* ---------------- Generalized Davidson ----------------- */ + SolverGeneralizedDavidson::SolverGeneralizedDavidson (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 + SolverGeneralizedDavidson::set_solver_type (EPS &eps) const { int ierr; ierr = EPSSetType (eps, const_cast(EPSGD)); @@ -403,6 +403,33 @@ namespace SLEPcWrappers this->solver_control.max_steps()); AssertThrow (ierr == 0, ExcSLEPcError(ierr)); } + + /* ------------------ Jacobi Davidson -------------------- */ + + SolverJacobiDavidson::SolverJacobiDavidson (SolverControl &cn, + const MPI_Comm &mpi_communicator, + const AdditionalData &data) + : + SolverBase (cn, mpi_communicator), + additional_data (data) + {} + + void + SolverJacobiDavidson::set_solver_type (EPS &eps) const + { + int ierr; + ierr = EPSSetType (eps, const_cast(EPSJD)); + 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 } -- 2.39.5