From cd07b6394137af02e13c91d2324ce55db23c62af Mon Sep 17 00:00:00 2001 From: young Date: Sat, 10 Sep 2011 18:01:38 +0000 Subject: [PATCH] Additional data for SLEPcWrappers Arnoldi solver git-svn-id: https://svn.dealii.org/trunk@24301 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/slepc_solver.h | 64 +++++++++++++++------- deal.II/source/lac/slepc_solver.cc | 15 +++++ 2 files changed, 58 insertions(+), 21 deletions(-) diff --git a/deal.II/include/deal.II/lac/slepc_solver.h b/deal.II/include/deal.II/lac/slepc_solver.h index e0e7b88d57..cc857d4cec 100644 --- a/deal.II/include/deal.II/lac/slepc_solver.h +++ b/deal.II/include/deal.II/lac/slepc_solver.h @@ -318,21 +318,29 @@ namespace SLEPcWrappers */ virtual void set_solver_type (EPS &eps) const = 0; - /* - * Attributes that store the - * relevant information for the - * eigenproblem solver - * context. These are: which - * portion of the spectrum to solve - * from; and the matrices and - * vectors that discribe the - * problem. + /** + * Which portion of the spectrum to + * solve from. + */ + EPSWhich set_which; + + /** + * The matrix $A$ of the + * generalized eigenvalue problem + * $Ax=B\lambda x$, or the standard + * eigenvalue problem $Ax=\lambda + * x$. */ - EPSWhich set_which; + const PETScWrappers::MatrixBase *opA; - const PETScWrappers::MatrixBase *opA; - const PETScWrappers::MatrixBase *opB; - const PETScWrappers::VectorBase *initial_vector; + /** + * The matrix $B$ of the + * generalized eigenvalue problem + * $Ax=B\lambda x$. + */ + const PETScWrappers::MatrixBase *opB; + + const PETScWrappers::VectorBase *initial_vector; /** * Pointer to an an object that @@ -425,8 +433,8 @@ namespace SLEPcWrappers * the PETScWrappers, but you can * change that. */ - SolverKrylovSchur (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD, + SolverKrylovSchur (SolverControl &cn, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, const AdditionalData &data = AdditionalData()); protected: @@ -451,7 +459,7 @@ namespace SLEPcWrappers * solver. Usage: All spectrum, all problem types, complex. * * @ingroup SLEPcWrappers - * @author Toby D. Young 2008 + * @author Toby D. Young 2008, 2011 */ class SolverArnoldi : public SolverBase { @@ -462,7 +470,21 @@ namespace SLEPcWrappers * should it be needed. */ struct AdditionalData - {}; + { + /** + * Constructor. By default, set the + * option of delayed + * reorthogonalization to false, + * i.e. don't do it. + */ + AdditionalData (const bool delayed_reorthogonalization = false); + + /** + * Flag for delayed + * reorthogonalization. + */ + bool delayed_reorthogonalization; + }; /** * SLEPc solvers will want to have @@ -474,7 +496,7 @@ namespace SLEPcWrappers * change that. */ SolverArnoldi (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, const AdditionalData &data = AdditionalData()); protected: @@ -523,7 +545,7 @@ namespace SLEPcWrappers * change that. */ SolverLanczos (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, const AdditionalData &data = AdditionalData()); protected: @@ -573,7 +595,7 @@ namespace SLEPcWrappers * change that. */ SolverPower (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, const AdditionalData &data = AdditionalData()); protected: @@ -623,7 +645,7 @@ namespace SLEPcWrappers * change that. */ SolverDavidson (SolverControl &cn, - const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD, + const MPI_Comm &mpi_communicator = PETSC_COMM_SELF, const AdditionalData &data = AdditionalData()); protected: diff --git a/deal.II/source/lac/slepc_solver.cc b/deal.II/source/lac/slepc_solver.cc index 7ddd9b8797..97e33ab3c9 100644 --- a/deal.II/source/lac/slepc_solver.cc +++ b/deal.II/source/lac/slepc_solver.cc @@ -282,6 +282,12 @@ namespace SLEPcWrappers /* ---------------------- SolverArnoldi ------------------------ */ + SolverArnoldi::AdditionalData:: + AdditionalData (const bool delayed_reorthogonalization) + : + delayed_reorthogonalization (delayed_reorthogonalization) + {} + SolverArnoldi::SolverArnoldi (SolverControl &cn, const MPI_Comm &mpi_communicator, const AdditionalData &data) @@ -305,6 +311,15 @@ namespace SLEPcWrappers ierr = EPSSetTolerances(eps, this->solver_control.tolerance(), this->solver_control.max_steps()); AssertThrow (ierr == 0, ExcSLEPcError(ierr)); + + // if requested, set delayed + // reorthogonalization in the + // Arnoldi iteration. + if (additional_data.delayed_reorthogonalization) + { + ierr = EPSArnoldiSetDelayed (eps, PETSC_TRUE); + AssertThrow (ierr == 0, ExcSLEPcError(ierr)); + } } /* ---------------------- Lanczos ------------------------ */ -- 2.39.5