From 5c29fa0c1699ffeda1fa9b1b613cdbec8854d95e Mon Sep 17 00:00:00 2001 From: "Toby D. Young" Date: Sat, 20 Apr 2013 12:18:30 +0000 Subject: [PATCH] Interface to SLEPc's interface to a LAPACK direct solver. git-svn-id: https://svn.dealii.org/trunk@29349 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/slepc_solver.h | 44 ++++++++++++++++++++++ deal.II/source/lac/slepc_solver.cc | 35 +++++++++++++++++ 2 files changed, 79 insertions(+) diff --git a/deal.II/include/deal.II/lac/slepc_solver.h b/deal.II/include/deal.II/lac/slepc_solver.h index 3d1eb2ff10..7985623114 100644 --- a/deal.II/include/deal.II/lac/slepc_solver.h +++ b/deal.II/include/deal.II/lac/slepc_solver.h @@ -655,6 +655,50 @@ namespace SLEPcWrappers virtual void set_solver_type (EPS &eps) const; }; + + /** + * An implementation of the solver interface using the SLEPc LAPACK + * direct solver. + * + * @ingroup SLEPcWrappers + * + * @author Toby D. Young 2013 + */ + class SolverLAPACK : 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. + */ + SolverLAPACK (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; + }; + // --------------------------- 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 a7bdf48670..c76b6c4331 100644 --- a/deal.II/source/lac/slepc_solver.cc +++ b/deal.II/source/lac/slepc_solver.cc @@ -460,6 +460,41 @@ namespace SLEPcWrappers #endif } + + /* ---------------------- LAPACK ------------------------- */ + SolverLAPACK::SolverLAPACK (SolverControl &cn, + const MPI_Comm &mpi_communicator, + const AdditionalData &data) + : + SolverBase (cn, mpi_communicator), + additional_data (data) + {} + + void + SolverLAPACK::set_solver_type (EPS &eps) const + { + // 'Tis overwhelmingly likely that PETSc/SLEPc *always* has + // BLAS/LAPACK, but let's be defensive. +#if PETSC_HAVE_BLASLAPACK + int ierr; + ierr = EPSSetType (eps, const_cast(EPSLAPACK)); + 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)); +#else + // Supress compiler warnings about unused paameters. + (void) eps; + + Assert ((false), + ExcMessage ("Your PETSc/SLEPc installation was not configured with BLAS/LAPACK " + "but this is needed to use the LAPACK solver.")); +#endif + } + } DEAL_II_NAMESPACE_CLOSE -- 2.39.5