lambda-lanczos 2.0.0
Loading...
Searching...
No Matches
lambda_lanczos::tridiagonal_impl Namespace Reference

Functions

template<typename T >
size_t num_of_eigs_smaller_than (T c, const std::vector< T > &alpha, const std::vector< T > &beta)
 Finds the number of eigenvalues of given tridiagonal matrix smaller than c. More...
 
template<typename T >
tridiagonal_eigen_limit (const std::vector< T > &alpha, const std::vector< T > &beta)
 Computes the upper bound of the absolute value of eigenvalues by Gerschgorin theorem. More...
 
template<typename T >
find_mth_eigenvalue (const std::vector< T > &alpha, const std::vector< T > &beta, const size_t m)
 Finds the mth smaller eigenvalue of given tridiagonal matrix. More...
 
template<typename T >
std::vector< T > compute_tridiagonal_eigenvector_from_eigenvalue (const std::vector< T > &alpha, const std::vector< T > &beta, const size_t index, const T ev)
 Computes an eigenvector corresponding to given eigenvalue for given tri-diagonal matrix. More...
 
template<typename T >
void tridiagonal_eigenpairs_bisection (const std::vector< T > &alpha, const std::vector< T > &beta, std::vector< T > &eigenvalues, std::vector< std::vector< T > > &eigenvectors)
 Computes all eigenpairs (eigenvalues and eigenvectors) for given tri-diagonal matrix. More...
 
template<typename T >
std::pair< T, T > calc_givens_cs (T a, T b)
 Calculates the cosine and sine of Givens rotation that eliminates a specific element (see detail). More...
 
template<typename T >
void isqr_step (std::vector< T > &alpha, std::vector< T > &beta, std::vector< std::vector< T > > &q, size_t offset, size_t nsub, bool rotate_matrix)
 Performs an implicit shift QR step A = Z^T A Z on given sub tridiagonal matrix A. More...
 
template<typename T >
void find_subspace (const std::vector< T > &alpha, std::vector< T > &beta, size_t &submatrix_first, size_t &submatrix_last)
 Find sub-tridiagonal matrix that remains non-diagonal. More...
 
template<typename T >
size_t tridiagonal_eigenpairs (const std::vector< T > &alpha, const std::vector< T > &beta, std::vector< T > &eigenvalues, std::vector< std::vector< T > > &eigenvectors, bool compute_eigenvector=true)
 Computes all eigenpairs (eigenvalues and eigenvectors) for given tri-diagonal matrix using the Implicitly Shifted QR algorithm. More...
 
template<typename T >
size_t tridiagonal_eigenvalues (const std::vector< T > &alpha, const std::vector< T > &beta, std::vector< T > &eigenvalues)
 Computes all eigenvalues for given tri-diagonal matrix using the Implicitly Shifted QR algorithm. More...
 

Function Documentation

◆ calc_givens_cs()

template<typename T >
std::pair< T, T > lambda_lanczos::tridiagonal_impl::calc_givens_cs ( a,
b 
)
inline

Calculates the cosine and sine of Givens rotation that eliminates a specific element (see detail).

This function calculates the cosine and sine that eliminate the element b, i.e. that satisfy

  ( c  s)(a) = (x)
  (-s  c)(b) = (0),
* 

where x >= 0.

Parameters
[in]aInput element to remain non-zero.
[in]bInput element to be eliminated.
Returns
Pair (c, s).

◆ compute_tridiagonal_eigenvector_from_eigenvalue()

template<typename T >
std::vector< T > lambda_lanczos::tridiagonal_impl::compute_tridiagonal_eigenvector_from_eigenvalue ( const std::vector< T > &  alpha,
const std::vector< T > &  beta,
const size_t  index,
const T  ev 
)
inline

Computes an eigenvector corresponding to given eigenvalue for given tri-diagonal matrix.

◆ find_mth_eigenvalue()

template<typename T >
T lambda_lanczos::tridiagonal_impl::find_mth_eigenvalue ( const std::vector< T > &  alpha,
const std::vector< T > &  beta,
const size_t  m 
)
inline

Finds the mth smaller eigenvalue of given tridiagonal matrix.

◆ find_subspace()

template<typename T >
void lambda_lanczos::tridiagonal_impl::find_subspace ( const std::vector< T > &  alpha,
std::vector< T > &  beta,
size_t &  submatrix_first,
size_t &  submatrix_last 
)
inline

Find sub-tridiagonal matrix that remains non-diagonal.

This function does the following two things:

  1. Overwrite small elements with zero,
  2. Find subspace that remains non-diagonal. The dimension of the resulting sub-matrix is submatrix_last - submatrix_first + 1.
Parameters
[in]alphaDiagonal elements of the full tridiagonal matrix.
[in,out]betaSub-diagonal elements of the full tridiagonal matrix.
[out]submatrix_firstIndex to point the first element of the sub-tridiagonal matrix.
[in,out]submatrix_lastIndex to point the last element of the sub-tridiagonal matrix.

◆ isqr_step()

template<typename T >
void lambda_lanczos::tridiagonal_impl::isqr_step ( std::vector< T > &  alpha,
std::vector< T > &  beta,
std::vector< std::vector< T > > &  q,
size_t  offset,
size_t  nsub,
bool  rotate_matrix 
)
inline

Performs an implicit shift QR step A = Z^T A Z on given sub tridiagonal matrix A.

Parameters
[in,out]alphaDiagonal elements of the full tridiagonal matrix.
[in,out]betaSubdiagonal elements of the full tridiagonal matrix.
[in,out]qA "matrix" Q that will be overwritten as Q = QZ if rotate_matrix is true. See note for details.
offsetThe first index of the submatrix.
nsubThe size of the submatrix.
rotate_matrixTrue to apply Given's rotations to the matrix q. If false, the matrix q won't be accessed.
Note
A series of the QR steps produces an eigenvector "matrix" q that stores the k-th eigenvector as q[k][:]. This definition differs from a usual mathematical sense (i.e., Q_{:,k} specifies the k-th eigenvector).

◆ num_of_eigs_smaller_than()

template<typename T >
size_t lambda_lanczos::tridiagonal_impl::num_of_eigs_smaller_than ( c,
const std::vector< T > &  alpha,
const std::vector< T > &  beta 
)
inline

Finds the number of eigenvalues of given tridiagonal matrix smaller than c.

Algorithm from Peter Arbenz et al. / "High Performance Algorithms for Structured Matrix Problems" / Nova Science Publishers, Inc.

◆ tridiagonal_eigen_limit()

template<typename T >
T lambda_lanczos::tridiagonal_impl::tridiagonal_eigen_limit ( const std::vector< T > &  alpha,
const std::vector< T > &  beta 
)
inline

Computes the upper bound of the absolute value of eigenvalues by Gerschgorin theorem.

This routine gives a rough upper bound, but it is sufficient because the bisection routine using the upper bound converges exponentially.

◆ tridiagonal_eigenpairs()

template<typename T >
size_t lambda_lanczos::tridiagonal_impl::tridiagonal_eigenpairs ( const std::vector< T > &  alpha,
const std::vector< T > &  beta,
std::vector< T > &  eigenvalues,
std::vector< std::vector< T > > &  eigenvectors,
bool  compute_eigenvector = true 
)
inline

Computes all eigenpairs (eigenvalues and eigenvectors) for given tri-diagonal matrix using the Implicitly Shifted QR algorithm.

Parameters
[in]alphaDiagonal elements of the full tridiagonal matrix.
[in]betaSub-diagonal elements of the full tridiagonal matrix.
[out]eigenvaluesEigenvalues.
[out]eigenvectorsEigenvectors. The k-th eigenvector will be stored in eigenvectors[k].
[in]compute_eigenvectorTrue to calculate eigenvectors. If false, eigenvectors won't be accessed.
Returns
Count of forced breaks due to unconvergence.

◆ tridiagonal_eigenpairs_bisection()

template<typename T >
void lambda_lanczos::tridiagonal_impl::tridiagonal_eigenpairs_bisection ( const std::vector< T > &  alpha,
const std::vector< T > &  beta,
std::vector< T > &  eigenvalues,
std::vector< std::vector< T > > &  eigenvectors 
)
inline

Computes all eigenpairs (eigenvalues and eigenvectors) for given tri-diagonal matrix.

◆ tridiagonal_eigenvalues()

template<typename T >
size_t lambda_lanczos::tridiagonal_impl::tridiagonal_eigenvalues ( const std::vector< T > &  alpha,
const std::vector< T > &  beta,
std::vector< T > &  eigenvalues 
)
inline

Computes all eigenvalues for given tri-diagonal matrix using the Implicitly Shifted QR algorithm.

Parameters
[in]alphaDiagonal elements of the full tridiagonal matrix.
[in]betaSub-diagonal elements of the full tridiagonal matrix.
[out]eigenvaluesEigenvalues.
Returns
Count of forced breaks due to unconvergence.