
Classes | |
| struct | SolverData |
Public Member Functions | |
| SolverBase (SolverControl &cn, const MPI_Comm &mpi_communicator) | |
| virtual | ~SolverBase () |
| template<typename OutputVector > | |
| void | solve (const PETScWrappers::MatrixBase &A, std::vector< double > &kr, std::vector< OutputVector > &vr, const unsigned int n_eigenvectors) |
| template<typename OutputVector > | |
| void | solve (const PETScWrappers::MatrixBase &A, const PETScWrappers::MatrixBase &B, std::vector< double > &kr, std::vector< OutputVector > &vr, const unsigned int n_eigenvectors) |
| void | set_matrices (const PETScWrappers::MatrixBase &A) |
| void | set_matrices (const PETScWrappers::MatrixBase &A, const PETScWrappers::MatrixBase &B) |
| void | set_initial_vector (const PETScWrappers::VectorBase &this_initial_vector) |
| void | set_transformation (SLEPcWrappers::TransformationBase &this_transformation) |
| void | set_which_eigenpairs (EPSWhich set_which) |
| void | solve (const unsigned int n_eigenvectors, unsigned int *n_converged) |
| void | get_eigenpair (const unsigned int index, double &kr, PETScWrappers::VectorBase &vr) |
| void | reset () |
| EPS * | get_EPS () |
| SolverControl & | control () const |
Static Public Member Functions | |
| ::ExceptionBase & | ExcSLEPcWrappersUsageError () |
| ::ExceptionBase & | ExcSLEPcError (int arg1) throw (errortext << " An error with error number " << arg1 << " occured while calling a SLEPc function" ) |
| ::ExceptionBase & | ExcSLEPcEigenvectorConvergenceMismatchError (int arg1, int arg2) throw (errortext << " The number of converged eigenvectors is " << arg1 << " but " << arg2 << " were requested. " ) |
Protected Member Functions | |
| virtual void | set_solver_type (EPS &eps) const =0 |
Protected Attributes | |
| SolverControl & | solver_control |
| const MPI_Comm | mpi_communicator |
| EPSWhich | set_which |
| const PETScWrappers::MatrixBase * | opA |
| const PETScWrappers::MatrixBase * | opB |
| const PETScWrappers::VectorBase * | initial_vector |
| SLEPcWrappers::TransformationBase * | transformation |
Static Private Member Functions | |
| static int | convergence_test (EPS eps, const int iteration, const PetscReal residual_norm, EPSConvergedReason *reason, void *solver_control) |
Private Attributes | |
| std_cxx1x::shared_ptr< SolverData > | solver_data |
Base class for solver classes using the SLEPc solvers. Since solvers in SLEPc are selected based on flags passed to a generic solver object, basically all the actual solver calls happen in this class, and derived classes simply set the right flags to select one solver or another, or to set certain parameters for individual solvers.
Definition at line 110 of file slepc_solver.h.
| SLEPcWrappers::SolverBase::SolverBase | ( | SolverControl & | cn, |
| const MPI_Comm & | mpi_communicator | ||
| ) |
Constructor. Takes the MPI communicator over which parallel computations are to happen.
| virtual SLEPcWrappers::SolverBase::~SolverBase | ( | ) | [virtual] |
Destructor.
| void SLEPcWrappers::SolverBase::solve | ( | const PETScWrappers::MatrixBase & | A, |
| std::vector< double > & | kr, | ||
| std::vector< OutputVector > & | vr, | ||
| const unsigned int | n_eigenvectors = 1 |
||
| ) |
Composite method that solves the eigensystem
. The eigenvector sent in has to have at least one element that we can use as a template when resizing, since we do not know the parameters of the specific vector class used (i.e. local_dofs for MPI vectors). However, while copying eigenvectors, at least twice the memory size of vr is being used (and can be more). To avoid doing this, the fairly standard calling sequence excecuted here is used: Initialise; Set up matrices for solving; Actually solve the system; Gather the solution(s); and reset.
Note that the number of converged eigenvectors can be larger than the number of eigenvectors requested; this is due to a round off error (success) of the eigenproblem solver context. If this is found to be the case we simply do not bother with more eigenpairs than requested, but handle that it may be more than specified by ignoring any extras. By default one eigenvector/eigenvalue pair is computed.
This is declared here to make it possible to take a std::vector of different PETScWrappers vector types
Definition at line 678 of file slepc_solver.h.
References AssertThrow, ExcSLEPcEigenvectorConvergenceMismatchError(), ExcSLEPcWrappersUsageError(), get_eigenpair(), and set_matrices().
Referenced by solve().
| void SLEPcWrappers::SolverBase::solve | ( | const PETScWrappers::MatrixBase & | A, |
| const PETScWrappers::MatrixBase & | B, | ||
| std::vector< double > & | kr, | ||
| std::vector< OutputVector > & | vr, | ||
| const unsigned int | n_eigenvectors = 1 |
||
| ) |
Same as above, but here a composite method for solving the system
.
Definition at line 703 of file slepc_solver.h.
References AssertThrow, ExcSLEPcEigenvectorConvergenceMismatchError(), ExcSLEPcWrappersUsageError(), get_eigenpair(), set_matrices(), and solve().
| void SLEPcWrappers::SolverBase::set_matrices | ( | const PETScWrappers::MatrixBase & | A ) |
Initialize solver for the linear system
. (Note: this is required before calling solve ())
Referenced by solve().
| void SLEPcWrappers::SolverBase::set_matrices | ( | const PETScWrappers::MatrixBase & | A, |
| const PETScWrappers::MatrixBase & | B | ||
| ) |
Same as above, but here initialize solver for the linear system
.
| void SLEPcWrappers::SolverBase::set_initial_vector | ( | const PETScWrappers::VectorBase & | this_initial_vector ) |
Set the initial vector for the solver.
| void SLEPcWrappers::SolverBase::set_transformation | ( | SLEPcWrappers::TransformationBase & | this_transformation ) |
Set the spectral transformation to be used.
| void SLEPcWrappers::SolverBase::set_which_eigenpairs | ( | EPSWhich | set_which ) |
Indicate which part of the spectrum is to be computed. By default largest magnitude eigenvalues are computed. For other allowed values see the SLEPc documentation.
| void SLEPcWrappers::SolverBase::solve | ( | const unsigned int | n_eigenvectors, |
| unsigned int * | n_converged | ||
| ) |
Solve the linear system for n_eigenvectors eigenstates. Parameter n_converged contains the actual number of eigenstates that have . converged; this can be both fewer or more than n_eigenvectors, depending on the SLEPc eigensolver used.
| void SLEPcWrappers::SolverBase::get_eigenpair | ( | const unsigned int | index, |
| double & | kr, | ||
| PETScWrappers::VectorBase & | vr | ||
| ) |
Access the solutions for a solved eigenvector problem, pair index solutions,
.
Referenced by solve().
| void SLEPcWrappers::SolverBase::reset | ( | ) |
Reset the solver, and return memory for eigenvectors
| EPS* SLEPcWrappers::SolverBase::get_EPS | ( | ) |
Retrieve the SLEPc solver object that is internally used.
| SolverControl& SLEPcWrappers::SolverBase::control | ( | ) | const |
Access to the object that controls convergence.
| ::ExceptionBase& SLEPcWrappers::SolverBase::ExcSLEPcWrappersUsageError | ( | ) | [static] |
Exceptions.
Referenced by solve().
| ::ExceptionBase& SLEPcWrappers::SolverBase::ExcSLEPcError | ( | int | arg1 ) | throw (errortext << " An error with error number " << arg1 << " occured while calling a SLEPc function" ) [static] |
| ::ExceptionBase& SLEPcWrappers::SolverBase::ExcSLEPcEigenvectorConvergenceMismatchError | ( | int | arg1, |
| int | arg2 | ||
| ) | throw (errortext << " The number of converged eigenvectors is " << arg1 << " but " << arg2 << " were requested. " ) [static] |
Referenced by solve().
| virtual void SLEPcWrappers::SolverBase::set_solver_type | ( | EPS & | eps ) | const [protected, pure virtual] |
Function that takes an Eigenvalue Problem Solver context object, and sets the type of solver that is requested by the derived class.
Implemented in SLEPcWrappers::SolverKrylovSchur, SLEPcWrappers::SolverArnoldi, SLEPcWrappers::SolverLanczos, and SLEPcWrappers::SolverPower.
| static int SLEPcWrappers::SolverBase::convergence_test | ( | EPS | eps, |
| const int | iteration, | ||
| const PetscReal | residual_norm, | ||
| EPSConvergedReason * | reason, | ||
| void * | solver_control | ||
| ) | [static, private] |
A function that is used in SLEPc as a callback to check on convergence. It takes the information provided from SLEPc and checks it against deal.II's own SolverControl objects to see if convergence has been reached.
SolverControl& SLEPcWrappers::SolverBase::solver_control [protected] |
Reference to the object that controls convergence of the iterative solver.
Definition at line 303 of file slepc_solver.h.
const MPI_Comm SLEPcWrappers::SolverBase::mpi_communicator [protected] |
Copy of the MPI communicator object to be used for the solver.
Definition at line 310 of file slepc_solver.h.
EPSWhich SLEPcWrappers::SolverBase::set_which [protected] |
Which portion of the spectrum to solve from.
Definition at line 325 of file slepc_solver.h.
const PETScWrappers::MatrixBase* SLEPcWrappers::SolverBase::opA [protected] |
The matrix
of the generalized eigenvalue problem
, or the standard eigenvalue problem
.
Definition at line 334 of file slepc_solver.h.
const PETScWrappers::MatrixBase* SLEPcWrappers::SolverBase::opB [protected] |
The matrix
of the generalized eigenvalue problem
.
Definition at line 341 of file slepc_solver.h.
const PETScWrappers::VectorBase* SLEPcWrappers::SolverBase::initial_vector [protected] |
Definition at line 343 of file slepc_solver.h.
Pointer to an an object that describes transformations that can be applied to the eigenvalue problem.
Definition at line 351 of file slepc_solver.h.
std_cxx1x::shared_ptr<SolverData> SLEPcWrappers::SolverBase::solver_data [private] |
Definition at line 404 of file slepc_solver.h.
documentation generated on Fri Feb 3 2012 06:04:16 by
doxygen
1.7.2