Fork me on GitHub
 All Classes Files Functions Variables Groups Pages
Functions
EigenSolvers

Iterative methods for computing eigenpairs of hermitian and non-hermitian linear systems. More...

Functions

template<typename Matrix , typename Array2d >
void cusp::eigen::arnoldi (const Matrix &A, Array2d &H, size_t k=10)
 Approximate spectral radius of A using Arnoldi. More...
 
template<typename Matrix , typename Array1d1 , typename Array2d1 , typename Array1d2 , typename Array2d2 , typename LanczosOptions >
void cusp::eigen::lanczos (const Matrix &A, Array1d1 &eigVals, Array2d1 &eigVecs, LanczosOptions &options)
 Lanczos method. More...
 
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::eigen::lobpcg (LinearOperator &A, Vector &S, Vector &X, Monitor &monitor, Preconditioner &M, bool largest=true)
 LOBPCG method. More...
 
template<typename MatrixType >
double cusp::eigen::disks_spectral_radius (const MatrixType &A)
 Uses Gershgorin disks to approximate spectral radius. More...
 
template<typename MatrixType >
double cusp::eigen::estimate_rho_Dinv_A (const MatrixType &A)
 Approximate spectral radius of (D^-1)A. More...
 
template<typename MatrixType >
double cusp::eigen::estimate_spectral_radius (const MatrixType &A, size_t k=20)
 Approximate spectral radius of A using Lanczos. More...
 
template<typename MatrixType >
double cusp::eigen::ritz_spectral_radius (const MatrixType &A, size_t k=10, bool symmetric=false)
 Approximate spectral radius of A using Lanczos or Arnoldi. More...
 

Function Documentation

template<typename Matrix , typename Array2d >
void cusp::eigen::arnoldi ( const Matrix &  A,
Array2d &  H,
size_t  k = 10 
)

Approximate spectral radius of A using Arnoldi.

Template Parameters
Matrixtype of a sparse or dense matrix
Array2dtype of dense matrix of partials
Parameters
Amatrix of the linear system
Hdense matrix of ritz values
kmaximum number of outer Arnoldi iterations
Returns
spectral radius approximation
Overview
Approximates the spectral radius A using a specified number of Arnoldi iterations.
Example
The following code snippet demonstrates how to use arnoldi to compute the spectral radius of a 16x16 Laplacian matrix.
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// compute the largest eigenpair of A using 20 Arnoldi iterations
float rho = cusp::eigen::arnoldi(A, 20);
std::cout << "Spectral radius of A : " << rho << std::endl;
return 0;
}
template<typename MatrixType >
double cusp::eigen::disks_spectral_radius ( const MatrixType &  A)

Uses Gershgorin disks to approximate spectral radius.

Template Parameters
MatrixTypetype of a sparse or dense matrix
Parameters
Amatrix of the linear system
Returns
spectral radius approximation
Overview
Approximates the spectral radius of a matrix using Gershgorin disks.
See Also
https://en.wikipedia.org/wiki/Gershgorin_circle_theorem
Example
The following code snippet demonstrates how to use disks_spectral_radius to compute the spectral radius of a 16x16 Laplacian matrix.
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// compute the largest eigenpair of A
std::cout << "Spectral radius of A : " << rho << std::endl;
return 0;
}
template<typename MatrixType >
double cusp::eigen::estimate_rho_Dinv_A ( const MatrixType &  A)

Approximate spectral radius of (D^-1)A.

Template Parameters
MatrixTypetype of a sparse or dense matrix
Parameters
Amatrix of the linear system
Returns
spectral radius approximation
Overview
Approximates the spectral radius (D^-1)A, where D is a diagonal matrix containing the diagonal entries of A. The spectral radius of (D^-1)A is computed using either Lanczos or Arnoldi.
Example
The following code snippet demonstrates how to use estimate_rho_Dinv_A to compute the spectral radius of a 16x16 Laplacian matrix.
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// compute the largest eigenpair of A
std::cout << "Spectral radius of (D^-1)A : " << rho << std::endl;
return 0;
}
template<typename MatrixType >
double cusp::eigen::estimate_spectral_radius ( const MatrixType &  A,
size_t  k = 20 
)

Approximate spectral radius of A using Lanczos.

Template Parameters
MatrixTypetype of a sparse or dense matrix
Parameters
Amatrix of the linear system
Returns
spectral radius approximation
Overview
Approximates the spectral radius A using a specified number of Lanczos iterations.
Example
The following code snippet demonstrates how to use estimate_spectral_radius to compute the spectral radius of a 16x16 Laplacian matrix.
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// compute the largest eigenpair of A using 20 Lanczos iterations
std::cout << "Spectral radius of A : " << rho << std::endl;
return 0;
}
template<typename Matrix , typename Array1d1 , typename Array2d1 , typename Array1d2 , typename Array2d2 , typename LanczosOptions >
void cusp::eigen::lanczos ( const Matrix &  A,
Array1d1 &  eigVals,
Array2d1 &  eigVecs,
LanczosOptions &  options 
)

Lanczos method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
Seigenvalues
Xeigenvectors
monitormonitors iteration and determines stopping conditions
Mpreconditioner for A
largestIf true compute the eigenpair corresponding to the largest eigenvalue otherwise compute the smallest.
Overview
Computes the extreme eigenpairs of hermitian linear systems A x = s x with preconditioner M.
Note
A and M must be symmetric.
Example
The following code snippet demonstrates how to use lobpcg to solve a 10x10 Poisson problem.
#include <cusp/csr_matrix.h>
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
cusp::gallery::poisson5pt(A, 1024, 1024);
// allocate storage and initialize eigenpairs
// initialize Lanzcos option
cusp::eigen::lanczos_options<float> options;
s
options.tol = 1e-6;
options.maxIter = 100;
options.verbose = true;
options.computeEigVecs = false;
options.reorth = cusp::eigen::Full;
// compute the largest eigenpair of A
cusp::eigen::lanczos(A, S, V, options);
std::cout << "Largest eigenvalue : " << S[4] << std::endl;
return 0;
}
template<class LinearOperator , class Vector , class Monitor , class Preconditioner >
void cusp::eigen::lobpcg ( LinearOperator &  A,
Vector &  S,
Vector &  X,
Monitor &  monitor,
Preconditioner &  M,
bool  largest = true 
)

LOBPCG method.

Template Parameters
LinearOperatoris a matrix or subclass of linear_operator
Vectorvector
Monitoris a monitor
Preconditioneris a matrix or subclass of linear_operator
Parameters
Amatrix of the linear system
Seigenvalues
Xeigenvectors
monitormonitors iteration and determines stopping conditions
Mpreconditioner for A
largestIf true compute the eigenpair corresponding to the largest eigenvalue otherwise compute the smallest.
Overview
Computes the extreme eigenpairs of hermitian linear systems A x = s x using LOBPCG.
Note
A and M must be symmetric.
See Also
https://en.wikipedia.org/wiki/LOBPCG
Example
The following code snippet demonstrates how to use lobpcg to solve a 10x10 Poisson problem.
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// allocate storage and initialize eigenpairs
cusp::random_array<double> randx(A.num_rows);
// set stopping criteria:
// iteration_limit = 100
// relative_tolerance = 1e-6
// absolute_tolerance = 0
// verbose = true
cusp::monitor<double> monitor(X, 100, 1e-6, 0, true);
// set preconditioner (identity)
// Compute the largest eigenpair of A
cusp::eigen::lobpcg(A, S, X, monitor, M, true);
std::cout << "Largest eigenvalue : " << S[0] << std::endl;
return 0;
}
See Also
monitor
template<typename MatrixType >
double cusp::eigen::ritz_spectral_radius ( const MatrixType &  A,
size_t  k = 10,
bool  symmetric = false 
)

Approximate spectral radius of A using Lanczos or Arnoldi.

Template Parameters
MatrixTypetype of a sparse or dense matrix
Parameters
Amatrix of the linear system
Returns
spectral radius approximation
Overview
Approximates the spectral radius A using a specified number of Lanczos or Arnoldi iterations.
Example
The following code snippet demonstrates how to use ritz_spectral_radius to compute the spectral radius of a 16x16 Laplacian matrix.
#include <cusp/monitor.h>
int main(void)
{
// create an empty sparse matrix structure (CSR format)
// initialize matrix
// compute the largest eigenpair of A using 20 Arnoldi iterations
float rho = cusp::eigen::ritz_spectral_radius(A, 20, false);
std::cout << "Spectral radius of A : " << rho << std::endl;
return 0;
}