15
$\begingroup$

I have two graphs with nearly n~100000 nodes each. In both graphs, each node is connected to exactly 3 other nodes so the adjacency matrix is symmetric and very sparse.

The hard part is I need all eigenvalues of the adjacency matrix but not eigenvectors. To be accurate, this is going to be once in my lifetime (as far as I can see, at least!) so I want to get all eigenvalues and don't mind waiting a few days to get them.

I tried scipy wrappers around ARPACK, but it takes way too long. I found multiple libraries but they work best for getting a subset of largest/smallest eigenvalues. Is there any library which works for symmetric sparse matrices with possibly parallel implementation to get all eigenvalues?

Improve this question
$\endgroup$
4
  • 6
    $\begingroup$ Out of curiosity, why exactly do you need all eigenvalues? Most problems of this size are approximations of even larger (or even infinite dimensional) problems, and so the eigenvalues of the small problems only approximate that of the problem one really want to solve. Most often, the quality of approximation is only good for the smallest or largest eigenvalues, and all others are only poorly approximated and consequently not of much practical interest. $\endgroup$
    – Wolfgang Bangerth
    Sep 19, 2016 at 3:55
  • $\begingroup$ related question and answer $\endgroup$
    – GoHokies
    Sep 19, 2016 at 17:36
  • $\begingroup$ @WolfgangBangerth: (Forgive me if these are obvious to you) The problem is coming from physics of materials. It's related to the tight binding approximation of materials to get band structure, vibrational and electrical properties. To get these, I need the full spectrum of eigenvalues. BTW, this is nothing new and it goes back to 70's and 80's but since my system is amorphous, I needed to have a very large system to get good results. Although most people care about crystals only which is extremely easy compared to my case. $\endgroup$
    – Mahdi
    Sep 19, 2016 at 19:01
  • 4
    $\begingroup$ @Mahdi: Well, what I meant is that the physical properties are determined by the spectrum of some partial differential operator. I suspect (but of course don't know, because you don't describe where the problem comes from) that the large matrix eigenvalue problem you have is just an approximation of the PDE problem. Consequently, your eigenvalues will also only be approximations. $\endgroup$
    – Wolfgang Bangerth
    Sep 19, 2016 at 21:31

2 Answers 2

Reset to default
11
$\begingroup$

You can use the shift-invert spectral transform [1] and compute the spectrum band by band.

The technique is also explained in my article [2]. Besides the implementation in [1], an implementation is available in C++ in my Graphite software [3] (update Jan 17: now everything is ported to geogram/graphite version 3), that I used to compute the eigenfunctions of the Laplace operator for meshes with up to 1 million vertices (a problem that is similar to yours).

How it works:

The idea is that if $A$ is invertible, then if $(V,\lambda)$ is an eigenpair of $A$, $(V, 1/\lambda)$ is an eigenpair of $A^{-1}$. The iterative method in ARPACK is very efficient for computing the large eigenvalues (high frequencies), but much less efficient for the small ones (small frequencies). Thus, when one needs to compute small frequencies, it is a good idea to replace $A$ with $A^{-1}$. Now, since ARPACK just needs to compute matrix-vector products, it is not necessary to really invert $A$: one can instead factor it (using for instance a sparse LU or sparse $LL^t$ factorization), then solve $Ax=b$ whenever ARPACK asks for a matrix-vector product. This is the "invert" transform. Now when the number of eigenvalues becomes large, ARPACK becomes to be slow, but there is another trick/transform that can be used, and one computes the eigenvalues of $A - \sigma Id$ where $\sigma$ is a "shift" that determines which part of the spectrum is explored (this is the "shift" transform). Combining both transforms, one computes a certain number of eigenvalues of $(A-\sigma Id)^{-1}$, and then explores the whole spectrum band-by-band, by increasing $\sigma$. The details are in [1],[2].

[1] http://www.mcs.anl.gov/uploads/cels/papers/P1263.pdf

[2] http://alice.loria.fr/index.php/publications.html?redirect=0&Paper=ManifoldHarmonics@2008

[3] http://alice.loria.fr/software/graphite/doc/html/

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thanks Bruno! These look very promising, I'll look into them! $\endgroup$
    – Mahdi
    Sep 17, 2016 at 6:49
1
$\begingroup$

Another option would be using Jacobi rotations. Since your matrix is already almost diagonal, it shouldn't take much time to converge. Generally it converges in linear rate, but after enough iterations the convergence rate becomes quadratic.

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.

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