5
$\begingroup$

I am trying to calculate all of the eigenvectors/eigenvalues of large (40000x40000), sparse, asymmetric matrix. I am using MATLAB and have 3 GB of working memory. The way I am calculating them is by sweeping through (as suggested here -http://www.mathworks.com/support/solutions/en/data/1-A0QRGR/index.html?product=SL&solution=1-A0QRGR). All of my eigenvalues are negative. I start from near the most negative eigenvalue, use eigs, then select the maximum. I repeat eigs on that maximum plus 1E-5 to avoid using eigs with a shift equal to an eigenvalue and then collect the eigenvalues (and the corresponding eigenvectors) that are greater than that maximum plus 1E-14 to avoid duplicate eigenvalues due to numerical error. I repeat this until the maximum of a new set of eigenvalues reaches close to 0. However, I always seem to fall a few (~5) eigenvalues short of the total. I have checked for the possibility of a maximum being a complex number and making sure that the complex conjugate has been tabulated, but that does not do anything. Moreover, some (~2) of the zero eigenvalues have corresponding eigenvectors that are entirely NaN. I am currently at a complete loss of how to proceed, and any suggestions regarding useful texts, etc., would be greatly appreciated.

$\endgroup$
  • 3
    $\begingroup$ I generally consider a problem of your size intractable (unless you have access to code of academic research-grade algorithms and also a supercomputer). I would like to ask why you need to compute ALL of the eigenpairs, and what is the application of this? $\endgroup$ – Victor Liu Aug 5 '13 at 18:54
  • $\begingroup$ You have $1.6\cdot 10^9$ matrix elements, which requires 12.8GB to store just the matrix. Regardless of the algorithm you choose, it will not yield any results in any reasonable amount of time if you can't even store the matrix in memory. I agree with Victor that it is impossible in practical terms to compute all eigenvalues of such a matrix. $\endgroup$ – Wolfgang Bangerth Aug 5 '13 at 20:12
  • $\begingroup$ I am working on a simple Master equation (link) in which the transition matrix contains rate constants do not change in time. I have calculated the rate matrix, which is ~40000x40000. I would like to apply this rate matrix to some initial population vector and look at the time evolution. This can be calculated exactly by calculating the eigenvalues/eigenvectors of the rate matrix. Using the algorithm I described, I can recapitulate the calculation for a smaller rate matrix for which the eigenpairs can be calculated exactly. $\endgroup$ – Kapil Aug 5 '13 at 21:35
  • $\begingroup$ Why should it be intractable? Using the algorithm described above, it takes ~2 hrs to calculate (nearly) all of the eigen-pairs. $\endgroup$ – Kapil Aug 5 '13 at 21:40
  • 2
    $\begingroup$ If the matrix is sparse, why does it require so much storage for the eigenvalues? (The eigenvectors, if not sparse, are another matter.) There's nothing in the problem statement that doesn't permit the input A(i,j)=-\delta_{i,j}+0.001*\delta_{i,j+1}, although I recognize that this is a stupid case. $\endgroup$ – Jeff Aug 6 '13 at 2:16
2
$\begingroup$

Over the course of 40,000 eigenvalues/eigenvectors, I wouldn't be surprised at all if at least a couple were closer to each other than the 1e-14 tolerance you're allowing, and that's why you're missing a few. If you really need all of these eigenvectors, the easiest possible solution would be to just find a computer with enough memory to do this in matlab. As Wolfgang pointed out, this will probably take ~20gb+ of memory; pushes the limit of tractability on a single computer but you can try and find one available that works. Another alternative is to set up a lot of virtual memory (swap space), keeping in mind that this will likely be very slow. So either get a bigger computer, use the eigenpairs you do already have, or just use swap space and wait.

Also, you used the phrase "matrix for which the eigenpairs can be calculated exactly." Keep in mind you are never computing these "exactly". They're being determined by iterative methods, and are subject to all floating point errors the same as any other method.

$\endgroup$
  • $\begingroup$ You're right. Even if I use the method described and instead of 1E-14 use 0 and try to remove duplicate eigenvalues by checking for duplicate eigenvectors, I will have to check the results. Simply using the eig function on a computer with more memory is the easiest way to go. Thanks. $\endgroup$ – Kapil Aug 6 '13 at 2:36

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.