Uses Lanczos iteration to find the truncated eigen-decomposition of a symmetric matrix.
slanczos(A,k=10,kl=-1,tol=.Machine$double.eps^.5,nt=1)A list with elements values (array of eigenvalues); vectors (matrix with eigenvectors in its columns);
iter (number of iterations required).
A symmetric matrix.
Must be non-negative. If kl is negative, then the k largest magnitude eigenvalues
are found, together with the corresponding eigenvectors. If kl is non-negative then the k
highest eigenvalues are found together with their eigenvectors and the kl lowest eigenvalues with
eigenvectors are also returned.
If kl is non-negative then the kl lowest eigenvalues are returned together with their
corresponding eigenvectors (in addition to the k highest eignevalues + vectors).
negative kl signals that the k largest magnitude eigenvalues should be returned, with eigenvectors.
tolerance to use for convergence testing of eigenvalues. Error in eigenvalues will be less
than the magnitude of the dominant eigenvalue multiplied by tol (or the machine precision!).
number of threads to use for leading order iterative multiplication of A by vector. May show no speed improvement on two processor machine.
Simon N. Wood simon.wood@r-project.org
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.
Demmel, J. (1997) Applied Numerical Linear Algebra. SIAM
cyclic.p.spline
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]))
Run the code above in your browser using DataLab