4
$\begingroup$

I'm trying to solve a generalized eigenvalue problem

$$Ax = \lambda Bx, \quad A = A^\top > 0,\; B = B^\top > 0$$

with $\lambda \approx \sigma$ using Rayleigh Quotient Iteration (RQI) (RQI is applied to $B^{-1/2}AB^{-1/2}$, but the $B^{-1/2}$ cancels from the final formulas): $$ (A - \mu_k B) \tilde x_{k+1} = B x_k\\ x_{k+1} = \frac{\tilde x_{k+1}}{||\tilde x_{k+1}||}\\ \mu_k = \frac{x_m^\top A x_m}{x_m^\top B x_m}\\ \mu_0 = \sigma, \quad x_0 = \text{random vector} $$ I've decided to use some iterative method to solve $$ (A - \mu_k B) \tilde x_{k+1} = B x_k $$ for $\tilde x_{k+1}$, since I do not want to store $A$ and $B$ explicitly. I only wish to code computing $Ax$ and $Bx$ products.

The matrix is symmetric, but not positive definite, otherwise I would have used conjugate gradients method. Which method would you suggest? Are there any pitfalls when using Krylov subspace methods with RQI?

Improve this question
$\endgroup$
2
  • 1
    $\begingroup$ @nicoguaro $A > 0$ is a standard notation to denote a positive definite matrix en.wikipedia.org/wiki/… . $\endgroup$
    – uranix
    Aug 18, 2015 at 6:07
  • 1
    $\begingroup$ $A \geq 0$ is psd, $A > 0$ is pd. Yes, both matrices are positive definite $\endgroup$
    – uranix
    Aug 18, 2015 at 6:12

1 Answer 1

Reset to default
4
+50
$\begingroup$

If your matrices are large, why not use a library like ARPACK? The shift-and-invert mode of ARPACK will help you calculate the eigenvalues close to $\sigma$.

There are interfaces to ARPACK for most high-level programming languages used for numerical computations (Fortran, C, Python, MATLAB, etc.).

For example, a quick tutorial for Python can be found here (reasonably small, dense matrices), or here (large, sparse matrices, available only in operator form).

Later edit: I read your post more carefully and noticed you are looking for a way to solve the shift-invert linear system... so why not try MINRES or SYMMLQ? Both algorithms are based on a Lanczos (3-term Arnoldi) iteration specifically designed for symmetric (but not necessarily positive definite) matrices. SYMMLQ should perform better if $A$ is non-singular.

Improve this answer
$\endgroup$
3
  • $\begingroup$ Since I'm trying to solve near-to-singular system $(A - \mu B)$ when $\mu$ is close to an eigenvalue, I'd better stick with MINRES instead of SYMMLQ. But I have troubles with finding some canonical pseudocode for MINRES, $\endgroup$
    – uranix
    Aug 20, 2015 at 20:58
  • $\begingroup$ Look here. It is a rather concise MATLAB implementation of MINRES with lots of helpful comments (and a "tutorial" PDF). $\endgroup$
    – GoHokies
    Aug 21, 2015 at 7:31
  • $\begingroup$ That implementation is mostly unreadable, already tried. I found one good here $\endgroup$
    – uranix
    Aug 21, 2015 at 15:45

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.

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