slanczos
Compute truncated eigen decomposition of a symmetric matrix
Description
Uses Lanczos iteration to find the truncated eigen-decomposition of a symmetric matrix.
Usage
slanczos(A,k=10,kl=-1,tol=.Machine$double.eps^.5,nt=1)
Arguments
A | A symmetric matrix. |
k | Must be non-negative. If |
kl | If |
tol | tolerance to use for convergence testing of eigenvalues. Error in eigenvalues will be less than the magnitude of the dominant eigenvalue multiplied by |
nt | number of threads to use for leading order iterative multiplication of A by vector. May show no speed improvement on two processor machine. |
Details
If kl
is non-negative, returns the highest k
and lowest kl
eigenvalues, with their corresponding eigenvectors. If kl
is negative, returns the largest magnitude k
eigenvalues, with corresponding eigenvectors.
The routine implements Lanczos iteration with full re-orthogonalization as described in Demmel (1997). Lanczos iteraction iteratively constructs a tridiagonal matrix, the eigenvalues of which converge to the eigenvalues of A
, as the iteration proceeds (most extreme first). Eigenvectors can also be computed. For small k
and kl
the approach is faster than computing the full symmetric eigendecompostion. The tridiagonal eigenproblems are handled using LAPACK.
The implementation is not optimal: in particular the inner triadiagonal problems could be handled more efficiently, and there would be some savings to be made by not always returning eigenvectors.
Value
A list with elements values
(array of eigenvalues); vectors
(matrix with eigenvectors in its columns); iter
(number of iterations required).
Author(s)
Simon N. Wood [email protected]
References
Demmel, J. (1997) Applied Numerical Linear Algebra. SIAM
See Also
Examples
require(mgcv) ## create some x's and knots... set.seed(1); n <- 700;A <- matrix(runif(n*n),n,n);A <- A+t(A) ## compare timings of slanczos and eigen system.time(er <- slanczos(A,10)) system.time(um <- eigen(A,symmetric=TRUE)) ## confirm values are the same... ind <- c(1:6,(n-3):n) range(er$values-um$values[ind]);range(abs(er$vectors)-abs(um$vectors[,ind]))
Copyright (©) 1999–2012 R Foundation for Statistical Computing.
Licensed under the GNU General Public License.