18 #ifndef __deal2__slepc_solver_h
19 #define __deal2__slepc_solver_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_SLEPC
25 # include <deal.II/base/std_cxx1x/shared_ptr.h>
26 # include <deal.II/lac/exceptions.h>
27 # include <deal.II/lac/solver_control.h>
28 # include <deal.II/lac/petsc_matrix_base.h>
29 # include <deal.II/lac/slepc_spectral_transformation.h>
31 # include <petscconf.h>
32 # include <petscksp.h>
33 # include <slepceps.h>
35 DEAL_II_NAMESPACE_OPEN
108 namespace SLEPcWrappers
154 template <
typename OutputVector>
157 std::vector<PetscScalar> &eigenvalues,
158 std::vector<OutputVector> &eigenvectors,
159 const unsigned int n_eigenpairs = 1);
166 template <
typename OutputVector>
170 std::vector<PetscScalar> &eigenvalues,
171 std::vector<OutputVector> &eigenvectors,
172 const unsigned int n_eigenpairs = 1);
179 template <
typename OutputVector>
183 std::vector<double> &real_eigenvalues,
184 std::vector<double> &imag_eigenvalues,
185 std::vector<OutputVector> &real_eigenvectors,
186 std::vector<OutputVector> &imag_eigenvectors,
187 const unsigned int n_eigenpairs = 1);
247 <<
" An error with error number " << arg1
248 <<
" occurred while calling a SLEPc function");
255 <<
" The number of converged eigenvectors is " << arg1
256 <<
" but " << arg2 <<
" were requested. ");
302 solve (
const unsigned int n_eigenpairs,
303 unsigned int *n_converged);
312 PetscScalar &eigenvalues,
322 double &real_eigenvalues,
323 double &imag_eigenvalues,
392 PetscScalar real_eigenvalue,
393 PetscScalar imag_eigenvalue,
394 PetscReal residual_norm,
395 PetscReal *estimated_error,
748 template <
typename OutputVector>
751 std::vector<PetscScalar> &eigenvalues,
752 std::vector<OutputVector> &eigenvectors,
753 const unsigned int n_eigenpairs)
756 AssertThrow ((n_eigenpairs > 0) && (n_eigenpairs <= A.
m ()),
757 ExcSLEPcWrappersUsageError());
763 unsigned int n_converged = 0;
764 solve (n_eigenpairs, &n_converged);
766 if (n_converged > n_eigenpairs)
767 n_converged = n_eigenpairs;
769 ExcSLEPcEigenvectorConvergenceMismatchError(n_converged, n_eigenpairs));
771 AssertThrow (eigenvectors.size() != 0, ExcSLEPcWrappersUsageError());
772 eigenvectors.resize (n_converged, eigenvectors.front());
773 eigenvalues.resize (n_converged);
775 for (
unsigned int index=0; index<n_converged; ++index)
776 get_eigenpair (index, eigenvalues[index], eigenvectors[index]);
779 template <
typename OutputVector>
783 std::vector<PetscScalar> &eigenvalues,
784 std::vector<OutputVector> &eigenvectors,
785 const unsigned int n_eigenpairs)
792 AssertThrow ((n_eigenpairs>0) && (n_eigenpairs<=A.
m ()),
793 ExcSLEPcWrappersUsageError());
799 unsigned int n_converged = 0;
800 solve (n_eigenpairs, &n_converged);
802 if (n_converged>=n_eigenpairs)
803 n_converged = n_eigenpairs;
806 ExcSLEPcEigenvectorConvergenceMismatchError(n_converged, n_eigenpairs));
807 AssertThrow (eigenvectors.size() != 0, ExcSLEPcWrappersUsageError());
809 eigenvectors.resize (n_converged, eigenvectors.front());
810 eigenvalues.resize (n_converged);
812 for (
unsigned int index=0; index<n_converged; ++index)
813 get_eigenpair (index, eigenvalues[index], eigenvectors[index]);
816 template <
typename OutputVector>
820 std::vector<double> &real_eigenvalues,
821 std::vector<double> &imag_eigenvalues,
822 std::vector<OutputVector> &real_eigenvectors,
823 std::vector<OutputVector> &imag_eigenvectors,
824 const unsigned int n_eigenpairs)
831 AssertThrow (real_eigenvalues.size() == imag_eigenvalues.size(),
833 AssertThrow (real_eigenvectors.size() == imag_eigenvectors.size (),
837 AssertThrow ((n_eigenpairs>0) && (n_eigenpairs<=A.
m ()),
838 ExcSLEPcWrappersUsageError());
844 unsigned int n_converged = 0;
845 solve (n_eigenpairs, &n_converged);
847 if (n_converged>=n_eigenpairs)
848 n_converged = n_eigenpairs;
851 ExcSLEPcEigenvectorConvergenceMismatchError(n_converged, n_eigenpairs));
852 AssertThrow ((real_eigenvectors.size()!=0) && (imag_eigenvectors.size()!=0),
853 ExcSLEPcWrappersUsageError());
855 real_eigenvectors.resize (n_converged, real_eigenvectors.front());
856 imag_eigenvectors.resize (n_converged, imag_eigenvectors.front());
857 real_eigenvalues.resize (n_converged);
858 imag_eigenvalues.resize (n_converged);
860 for (
unsigned int index=0; index<n_converged; ++index)
862 real_eigenvalues[index], imag_eigenvalues[index],
863 real_eigenvectors[index], imag_eigenvectors[index]);
868 DEAL_II_NAMESPACE_CLOSE
870 #endif // DEAL_II_WITH_SLEPC
virtual void set_solver_type(EPS &eps) const
EPSProblemType set_problem
PetscScalar target_eigenvalue
bool delayed_reorthogonalization
const AdditionalData additional_data
SolverJacobiDavidson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverPower(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
void set_initial_vector(const PETScWrappers::VectorBase &this_initial_vector)
static int convergence_test(EPS eps, PetscScalar real_eigenvalue, PetscScalar imag_eigenvalue, PetscReal residual_norm, PetscReal *estimated_error, void *solver_control)
const PETScWrappers::MatrixBase * opB
DeclException0(ExcSLEPcWrappersUsageError)
#define AssertThrow(cond, exc)
const PETScWrappers::VectorBase * initial_vector
void solve(const PETScWrappers::MatrixBase &A, std::vector< PetscScalar > &eigenvalues, std::vector< OutputVector > &eigenvectors, const unsigned int n_eigenpairs=1)
const MPI_Comm mpi_communicator
const AdditionalData additional_data
SolverKrylovSchur(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator)
AdditionalData(const bool delayed_reorthogonalization=false)
void set_problem_type(EPSProblemType set_problem)
DeclException2(ExcSLEPcEigenvectorConvergenceMismatchError, int, int,<< " The number of converged eigenvectors is "<< arg1<< " but "<< arg2<< " were requested. ")
virtual void set_solver_type(EPS &eps) const =0
void set_target_eigenvalue(const PetscScalar &this_target)
void set_matrices(const PETScWrappers::MatrixBase &A)
void set_which_eigenpairs(EPSWhich set_which)
SolverGeneralizedDavidson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
virtual void set_solver_type(EPS &eps) const
const AdditionalData additional_data
SolverControl & control() const
virtual void set_solver_type(EPS &eps) const
const AdditionalData additional_data
const PETScWrappers::MatrixBase * opA
void get_eigenpair(const unsigned int index, PetscScalar &eigenvalues, PETScWrappers::VectorBase &eigenvectors)
SolverLanczos(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(EPS &eps) const
virtual void set_solver_type(EPS &eps) const
SolverLAPACK(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
virtual void set_solver_type(EPS &eps) const
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
const AdditionalData additional_data
std_cxx1x::shared_ptr< SolverData > solver_data
void set_transformation(SLEPcWrappers::TransformationBase &this_transformation)
virtual void set_solver_type(EPS &eps) const
SolverControl & solver_control
DeclException1(ExcSLEPcError, int,<< " An error with error number "<< arg1<< " occurred while calling a SLEPc function")
EPSConvergedReason reason
SolverArnoldi(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SLEPcWrappers::TransformationBase * transformation
void get_solver_state(const SolverControl::State state)