MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
mfem::SlepcEigenSolver Class Reference

#include <slepc.hpp>

Public Types

enum  Which {
  LARGEST_MAGNITUDE , SMALLEST_MAGNITUDE , LARGEST_REAL , SMALLEST_REAL ,
  LARGEST_IMAGINARY , SMALLEST_IMAGINARY , TARGET_MAGNITUDE , TARGET_REAL
}
 Target spectrum for the eigensolver. More...
 
enum  SpectralTransformation { SHIFT , SHIFT_INVERT }
 Spectral transformations that can be used by the solver in order to accelerate the convergence to the target eignevalues. More...
 

Public Member Functions

 SlepcEigenSolver (MPI_Comm comm, const std::string &prefix=std::string())
 Constructors.
 
virtual ~SlepcEigenSolver ()
 
void SetTol (real_t tol)
 Set solver convergence tolerance relative to the magnitude of the eigenvalue.
 
void SetMaxIter (int max_iter)
 Set maximum number of iterations allowed in the call to SlepcEigenSolver::Solve.
 
void SetNumModes (int num_eigs)
 Set the number of eigenmodes to compute.
 
void SetOperator (const PetscParMatrix &op)
 Set operator for standard eigenvalue problem.
 
void SetOperators (const PetscParMatrix &op, const PetscParMatrix &opB)
 Set operators for generalized eigenvalue problem.
 
void Customize (bool customize=true) const
 Customize object with options set.
 
void Solve ()
 Solve the eigenvalue problem for the specified number of eigenvalues.
 
int GetNumConverged ()
 Get the number of converged eigenvalues after the call to SlepcEigenSolver::Solve.
 
void GetEigenvalue (unsigned int i, real_t &lr) const
 Get the ith eigenvalue after the system has been solved.
 
void GetEigenvalue (unsigned int i, real_t &lr, real_t &lc) const
 Get the ith eigenvalue after the system has been solved.
 
void GetEigenvector (unsigned int i, Vector &vr) const
 Get the ith eigenvector after the system has been solved.
 
void GetEigenvector (unsigned int i, Vector &vr, Vector &vc) const
 Get the ith eigenvector after the system has been solved.
 
void SetWhichEigenpairs (Which which)
 Set the which eigenvalues the solver will target and the order they will be indexed in.
 
void SetTarget (real_t target)
 Set the target value for the eigenpairs you want when using SlepcEigenSolver::TARGET_MAGNITUDE or SlepcEigenSolver::TARGET_REAL in the SlepcEigenSolver::SetWhichEigenpairs method.
 
void SetSpectralTransformation (SpectralTransformation transformation)
 Set the spectral transformation strategy for acceletating convergenvce. Both SlepcEigenSolver::SHIFT and SlepcEigenSolver::SHIFT_INVERT are available.
 
 operator slepc::EPS () const
 Conversion function to SLEPc's EPS type.
 
 operator PetscObject () const
 Conversion function to PetscObject.
 

Detailed Description

Definition at line 36 of file slepc.hpp.

Member Enumeration Documentation

◆ SpectralTransformation

Spectral transformations that can be used by the solver in order to accelerate the convergence to the target eignevalues.

Enumerator
SHIFT 

Utilize the shift of origin strategy.

SHIFT_INVERT 

Utilize the shift and invert strategy.

Definition at line 152 of file slepc.hpp.

◆ Which

Target spectrum for the eigensolver.

This will define the order in which the eigenvalues/eigenvectors are indexed after the call to SlepcEigenSolver::Solve.

Note
Target imaginary is not supported without complex support in SLEPc, and intervals are not implemented.
Enumerator
LARGEST_MAGNITUDE 

The eigenvalues with the largest complex magnitude (default)

SMALLEST_MAGNITUDE 

The eigenvalues with the smallest complex magnitude.

LARGEST_REAL 

The eigenvalues with the largest real component.

SMALLEST_REAL 

The eigenvalues with the smallest real component.

LARGEST_IMAGINARY 

The eigenvalues with the largest imaginary component.

SMALLEST_IMAGINARY 

The eigenvalues with the smallest imaginary component.

TARGET_MAGNITUDE 

The eigenvalues with complex magnitude closest to the target value.

TARGET_REAL 

The eigenvalues with the real component closest to the target value.

Definition at line 129 of file slepc.hpp.

Constructor & Destructor Documentation

◆ SlepcEigenSolver()

mfem::SlepcEigenSolver::SlepcEigenSolver ( MPI_Comm comm,
const std::string & prefix = std::string() )

Constructors.

Definition at line 60 of file slepc.cpp.

◆ ~SlepcEigenSolver()

mfem::SlepcEigenSolver::~SlepcEigenSolver ( )
virtual

Definition at line 70 of file slepc.cpp.

Member Function Documentation

◆ Customize()

void mfem::SlepcEigenSolver::Customize ( bool customize = true) const

Customize object with options set.

Definition at line 137 of file slepc.cpp.

◆ GetEigenvalue() [1/2]

void mfem::SlepcEigenSolver::GetEigenvalue ( unsigned int i,
real_t & lr ) const

Get the ith eigenvalue after the system has been solved.

Parameters
[in]iThe index for the eigenvalue you want ordered by SlepcEigenSolver::SetWhichEigenpairs
[out]lrThe real component of the eigenvalue
Note
the index i must be between 0 and SlepcEigenSolver::GetNumConverged - 1

Definition at line 147 of file slepc.cpp.

◆ GetEigenvalue() [2/2]

void mfem::SlepcEigenSolver::GetEigenvalue ( unsigned int i,
real_t & lr,
real_t & lc ) const

Get the ith eigenvalue after the system has been solved.

Parameters
[in]iThe index for the eigenvalue you want ordered by SlepcEigenSolver::SetWhichEigenpairs
[out]lrThe real component of the eigenvalue
[out]lcThe imaginary component of the eigenvalue
Note
the index i must be between 0 and SlepcEigenSolver::GetNumConverged - 1

Definition at line 152 of file slepc.cpp.

◆ GetEigenvector() [1/2]

void mfem::SlepcEigenSolver::GetEigenvector ( unsigned int i,
Vector & vr ) const

Get the ith eigenvector after the system has been solved.

Parameters
[in]iThe index for the eigenvector you want ordered by SlepcEigenSolver::SetWhichEigenpairs
[out]vrThe real components of the eigenvector
Note
the index i must be between 0 and SlepcEigenSolver::GetNumConverged - 1

Definition at line 158 of file slepc.cpp.

◆ GetEigenvector() [2/2]

void mfem::SlepcEigenSolver::GetEigenvector ( unsigned int i,
Vector & vr,
Vector & vc ) const

Get the ith eigenvector after the system has been solved.

Parameters
[in]iThe index for the eigenvector you want ordered by SlepcEigenSolver::SetWhichEigenpairs
[out]vrThe real components of the eigenvector
[out]vcThe imaginary components of the eigenvector
Note
the index i must be between 0 and SlepcEigenSolver::GetNumConverged - 1

Definition at line 171 of file slepc.cpp.

◆ GetNumConverged()

int mfem::SlepcEigenSolver::GetNumConverged ( )

Get the number of converged eigenvalues after the call to SlepcEigenSolver::Solve.

Definition at line 188 of file slepc.cpp.

◆ operator PetscObject()

mfem::SlepcEigenSolver::operator PetscObject ( ) const
inline

Conversion function to PetscObject.

Definition at line 185 of file slepc.hpp.

◆ operator slepc::EPS()

mfem::SlepcEigenSolver::operator slepc::EPS ( ) const
inline

Conversion function to SLEPc's EPS type.

Definition at line 182 of file slepc.hpp.

◆ SetMaxIter()

void mfem::SlepcEigenSolver::SetMaxIter ( int max_iter)

Set maximum number of iterations allowed in the call to SlepcEigenSolver::Solve.

Definition at line 116 of file slepc.cpp.

◆ SetNumModes()

void mfem::SlepcEigenSolver::SetNumModes ( int num_eigs)

Set the number of eigenmodes to compute.

Definition at line 124 of file slepc.cpp.

◆ SetOperator()

void mfem::SlepcEigenSolver::SetOperator ( const PetscParMatrix & op)

Set operator for standard eigenvalue problem.

Definition at line 81 of file slepc.cpp.

◆ SetOperators()

void mfem::SlepcEigenSolver::SetOperators ( const PetscParMatrix & op,
const PetscParMatrix & opB )

Set operators for generalized eigenvalue problem.

Definition at line 93 of file slepc.cpp.

◆ SetSpectralTransformation()

void mfem::SlepcEigenSolver::SetSpectralTransformation ( SlepcEigenSolver::SpectralTransformation transformation)

Set the spectral transformation strategy for acceletating convergenvce. Both SlepcEigenSolver::SHIFT and SlepcEigenSolver::SHIFT_INVERT are available.

Definition at line 234 of file slepc.cpp.

◆ SetTarget()

void mfem::SlepcEigenSolver::SetTarget ( real_t target)

Set the target value for the eigenpairs you want when using SlepcEigenSolver::TARGET_MAGNITUDE or SlepcEigenSolver::TARGET_REAL in the SlepcEigenSolver::SetWhichEigenpairs method.

Definition at line 229 of file slepc.cpp.

◆ SetTol()

void mfem::SlepcEigenSolver::SetTol ( real_t tol)

Set solver convergence tolerance relative to the magnitude of the eigenvalue.

Note
Default value is 1e-8

Definition at line 106 of file slepc.cpp.

◆ SetWhichEigenpairs()

void mfem::SlepcEigenSolver::SetWhichEigenpairs ( SlepcEigenSolver::Which which)

Set the which eigenvalues the solver will target and the order they will be indexed in.

For SlepcEigenSolver::TARGET_MAGNITUDE or SlepcEigenSolver::TARGET_REAL you will also need to set the target value with SlepcEigenSolver::SetTarget.

Definition at line 195 of file slepc.cpp.

◆ Solve()

void mfem::SlepcEigenSolver::Solve ( )

Solve the eigenvalue problem for the specified number of eigenvalues.

Definition at line 130 of file slepc.cpp.


The documentation for this class was generated from the following files: