0
$\begingroup$

I swear I recall hearing about a method for computing eigenvalues that takes a vector of known largest eigenvalues (n-1 eigenvalues) and computes the next largest eigenvalue (eigenvalue n) of a matrix more rapidly than re-computing the first n from scratch. I'm computing eigenpairs for a very large sparse Hermetian matrix and it would be really nice if I decide I need, say, 205 eigs rather than 200 to be able to decrease the amount of time computing those last 5 eigenvalues takes.

Anybody know what I'm talking about? Got a citation for me or a link to a C++ or Matlab library? Am I crazy and there's no such algorithm?

Thanks folks.

$\endgroup$
  • 1
    $\begingroup$ Do you have the 200 eigenvectors as well or just the 200 eigenvalues? $\endgroup$ – Carl Christian Feb 25 '16 at 15:50
  • $\begingroup$ @CarlChristian Eigenvectors are available as well, yes. $\endgroup$ – sintax Feb 25 '16 at 15:54
1
$\begingroup$

Since you have the dominant eigenspace you can deflate the problem and find eigenpairs which are less significant using, say, subspace iteration with or without Rayleigh Ritz acceleration.

Deflation is easily understood in the simplest case where $A \in \mathbb{R}^{n \times n}$ be symmetric positive definite and has distinct eigenvalues \begin{equation} \lambda_1 > \lambda_2 > \dotsc > \lambda_n > 0 \end{equation} Let \begin{equation} V_k = \begin{pmatrix} v_1 & v_2 & \dotsc & v_k \end{pmatrix} \end{equation} denote an $n$ by $k$ orthonormal matrix such that \begin{equation} \text{Ran}(V_k) = \text{span}\{v_1,v_2,\dotsc,v_k\} \end{equation} is the $k$th domininant eigenspace. Then the deflated matrix $A_k$ given by \begin{equation} A_k = A - V_k (V_k^T A V_k) V_k^T \end{equation} is symmetric positive semi-definite and has eigenvalues \begin{equation} \lambda_{k+1} > \lambda_{k+2} > \dotsc > \lambda_{n} > 0 \end{equation} and also \begin{equation} \lambda = 0 \end{equation} which has multiplicity $k$.

If you replace the action of $A$ with the action of $A_k$ in the subspace iteration, then you can extract information about the next set of eigen-pairs.

Certainly, it is more costly to compute the action of $A$ in terms of flops, but you will be doing dense arithmetic instead of sparse arithmetic so you flop rate is likely to be much higher.

Subspace iteration is a glorified version of the power method. Instead of a single vector you will use $m$ vectors. The normalization process is replaced by a tall $QR$ factorization. Apart for the cost of computing the action of the matrix, you have to pay $O(nm^2)$ operation per $QR$ factorization. Rayleigh Ritz acceleration is worth your while as it reduces the overall number of iterations, but costs $O(m^3)$ operations per iteration. It is certainly to your advantage to work with the deflated matrix and use $m=200$, instead of the full matrix and $m=400$.

You will find the algorithms describe in Gene Golub and Charles van Loan's book, which is freely available through MIT at this address

http://web.mit.edu/ehliu/Public/sclark/Golub%20G.H.,%20Van%20Loan%20C.F.-%20Matrix%20Computations.pdf

The dense kernels which you will need are part of LAPACK and BLAS, but the sparse matrix vector multiply is usually left to the user. Moreover, when you are doing subspace iteration using $m$ vectors, then the cost of Rayleigh Ritz acceleration is $O(m^3)$ (dense eigenvalue problem)

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.