|
| ArPackPlusPlus_Algorithms (const BCRSMatrix &m, const unsigned int nIterationsMax=100000, const unsigned int verbosity_level=0) |
| Construct from required parameters. More...
|
|
void | computeSymMaxMagnitude (const Real &epsilon, BlockVector &x, Real &lambda) const |
| Assume the matrix to be square, symmetric and perform IRLM to compute an approximation lambda of its dominant (i.e. largest magnitude) eigenvalue and the corresponding approximation x of an associated eigenvector. More...
|
|
void | computeSymMinMagnitude (const Real &epsilon, BlockVector &x, Real &lambda) const |
| Assume the matrix to be square, symmetric and perform IRLM to compute an approximation lambda of its least dominant (i.e. smallest magnitude) eigenvalue and the corresponding approximation x of an associated eigenvector. More...
|
|
void | computeSymCond2 (const Real &epsilon, Real &cond_2) const |
| Assume the matrix to be square, symmetric and perform IRLM to compute an approximation of its spectral condition number which, for symmetric matrices, can be expressed as the ratio of the dominant eigenvalue's magnitude and the least dominant eigenvalue's magnitude. More...
|
|
void | computeNonSymMax (const Real &epsilon, BlockVector &x, Real &sigma) const |
| Assume the matrix to be nonsymmetric and perform IRLM to compute an approximation sigma of its largest singular value and the corresponding approximation x of an associated singular vector. More...
|
|
void | computeNonSymMin (const Real &epsilon, BlockVector &x, Real &sigma) const |
| Assume the matrix to be nonsymmetric and perform IRLM to compute an approximation sigma of its smallest singular value and the corresponding approximation x of an associated singular vector. More...
|
|
void | computeNonSymCond2 (const Real &epsilon, Real &cond_2) const |
| Assume the matrix to be nonsymmetric and perform IRLM to compute an approximation of its spectral condition number which can be expressed as the ratio of the largest singular value and the smallest singular value. More...
|
|
unsigned int | getIterationCount () const |
| Return the number of iterations in last application of an algorithm.
|
|
template<typename
BCRSMatrix, typename
BlockVector>
class Dune::ArPackPlusPlus_Algorithms< BCRSMatrix, BlockVector >
Wrapper to use a range of ARPACK++ eigenvalue solvers.
A class template for performing some eigenvalue algorithms provided by the ARPACK++ library which is based on the implicitly restarted Arnoldi/Lanczos method (IRAM/IRLM), a synthesis of the Arnoldi/Lanczos process with the implicitily shifted QR technique. The method is designed to compute eigenvalue-eigenvector pairs of large scale sparse nonsymmetric/symmetric matrices. This class template uses the algorithms to compute the dominant (i.e. largest magnitude) and least dominant (i.e. smallest magnitude) eigenvalue as well as the spectral condition number of square, symmetric matrices and to compute the largest and smallest singular value as well as the spectral condition number of nonsymmetric matrices.
- Note
- For a recent version of the ARPACK++ library working with recent compiler versions see "http://reuter.mit.edu/software/arpackpatch/" or the git repository "https://github.com/m-reuter/arpackpp.git".
-
Note that the Arnoldi/Lanczos process currently is initialized using a vector which is randomly generated by ARPACK++. This could be changed in a future version, since ARPACK++ supports manual initialization of this vector.
- Todo:
- The current implementation is limited to DUNE-ISTL BCRSMatrix types with blocklevel 2. An extension to blocklevel >= 2 might be provided in a future version.
- Todo:
- Maybe make ARPACK++ parameter ncv available to the user.
- Template Parameters
-
BCRSMatrix | Type of a DUNE-ISTL BCRSMatrix whose eigenvalues respectively singular values shall be considered; is assumed to have blocklevel 2. |
BlockVector | Type of the associated vectors; compatible with the rows of a BCRSMatrix object (if #rows >= #ncols) or its columns (if #rows < #ncols). |
- Author
- Sebastian Westerheide.