1
$\begingroup$

I am considering implementing the accelerated inverse power (AIP) method with Rayleigh quotient to speed up eigendecomposition of a real square symmetric matrix. Halton (1996) gives an example algorithm which is listed below

Step 1. $\quad \mathbf{M}= \mathbf{H} - \mu \mathbf{I} $

Step 2. $\quad \text{Solve } \mathbf{M}\mathbf{y}= \mathbf{z}$

Step 3. $\quad \mathbf{z}= (\mathbf{y}^\top \mathbf{y})^{-1/2} \mathbf{y} $

Step 4. $\quad \rho = z_{k_{h}} / y_{k_{h}} $

Step 5. $\quad \mu = \mu + \rho $

Repeat Steps 1-5 until the delta between the updated eigenvector $\mathbf{z}$ and previous $\mathbf{z}$ are below some value.

Question is, why does he have subscripts for $k$ and $h$ on the ratio elements in the Rayleigh quotient? Which vector elements of $\mathbf{z}$ and $\mathbf{y}$ is he referring to?

I have also seen a matrix form of the Rayleigh quotient used in AIP which looked similar to a conjugate gradient step involving a numerator and denominator with vectors sandwiched around a matrix. What would that be, as an example?

Last, I also need all eigenvalues and their eigenvectors, and AIP is only geared toward finding the closest eigenvalue to $\mu$. Since my matrix $\mathbf{H}$ is a $p \times p$ correlation matrix, I know the sum of eigenvalues is equal to the $p$-dimensions of $\mathbf{H}$. Thus, I believe that if I first set $\mu=p$, then I will obtain the major eigenvalue/eigenvector. After $\lambda_1$ and $\mathbf{e}_1$ are found, dimension both $\mathbf{z}$ and $\mathbf{y}$ to be $p \times 2$, and still set $\mu=p$, but don't update $\lambda_1$ inside the $2 \times 2$ $\boldsymbol{\Lambda}$ or $\mathbf{e}_1$ inside the $p \times 2$ $\mathbf{E}$. The results for $\lambda_2$ and $\mathbf{e}_2$ will be inside $\boldsymbol{\Lambda}$ and $\mathbf{E}$.

Has anyone developed steps for finding all eigenvalue/eigenvector pairs iteratively using AIP?

Improve this question
$\endgroup$
6
  • 1
    $\begingroup$ do you want to implement the method to use in some project or simply to understand the mathematics better? $\endgroup$
    – Thijs Steel
    Sep 14, 2021 at 7:33
  • $\begingroup$ Yes, implement. However, it is looking like any variant of AIP won't be able to outperform QL and QR iteration for determining all pairs at once, which I already use (QL). But they're both $6n^3 + O(n^2)$. Dhillon's MMRR is $O(n^2)$, but I am not sure if it's implemented in LAPACK yet. $\endgroup$
    – user39981
    Sep 15, 2021 at 16:06
  • 1
    $\begingroup$ In that case, just use LAPACK. In fact, the QR algorithm can be seen as a more complicated variant of the AIP method. And the reduction step that precedes MMRR is still $O(n^3)$ $\endgroup$
    – Thijs Steel
    Sep 15, 2021 at 19:19
  • $\begingroup$ Thanks! What's the specific LAPACK package for QR algorithm for a real square symmetric matrix [which is always the correlation matrix having no pathologies, and is always positive-semidefinite (could have a few zero eigenvalues)]? I ask because there's dozens of similar LAPACK algorithms - I know MATLAB studies a matrix, then decides which LAPACK algorithm to use. $\endgroup$
    – user39981
    Sep 16, 2021 at 3:51
  • 1
    $\begingroup$ DSYEV is a good option, you will need to change the routine depending the working precision $\endgroup$
    – Thijs Steel
    Sep 16, 2021 at 7:30

1 Answer 1

Reset to default
0
$\begingroup$

Per comments from @Thijs Steel, this has been answered via suggestion to use the LAPACK DSYEV algorithm.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.