1
$\begingroup$

I am trying to write a QR algorithm in Python for eigenvectors and eigenvalues finding for large symmetric matrices,

My initial thought was to use Householder transformation with a Wilkinson shift on a Hessenberg matrix when I ran it on an 800×800 matrix is did not converge within reasonable time. I tested the algorithm on smaller matrices and it did converge and the results were correct

What else can I apply to reduce the convergence time and how can I use the fact that the given matrix is symmetric to my aid?

$\endgroup$
1
  • 1
    $\begingroup$ Do you know that there is a bug? It may simply be that your algorithm has a run time that grows with a certain power of the size of the matrix, and that you ran out of patience when you say "did not converge within a reasonable time". $\endgroup$ Apr 30, 2023 at 23:31

1 Answer 1

5
$\begingroup$

Before you drill into the problem for it's own sake, you might just reach for existing LAPACK algorithms for this exact problem that are very robust, and should be accessible through numpy. But if it's just a learning exercise, carry on.

It's been some time since I looked at this, but if memory serves the symmetric eigenproblem reduces to tridiagonal matrices, not Hessenberg ones. (To be more specific, the projection in question is both symmetric and Hessenberg, which implies it is tridiagonal). The frontend routine in LAPACK to orthogonally reduce a symmetric/input matrix A to this tridiagonal one T is dsytrd/ssytrd. There are a number of backend routines from there, to iteratively reduce T to diagonal/eigenvalue form, but you should probably imitate dstev/sstev .. or perhaps more specifically dsteqr/ssteqr as it sounds most similar to what you've attempted so far (Wilkinson-shifted QR iterations to seek eigenvalues, simultaneously capturing eigenvectors by accumulating the Givens rotations).

$\endgroup$

Your Answer

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

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