Why eigs( 'lm') is much faster than eigs('sm')

Question

I use eigs to calculate the eigen vectors of sparse square matrices which are large (tens of thousands). What I want is the smallest set of eigen vectors. But

eigs(A, 10, 'sm')      % Note: A is the matrix

runs very slow.

However, using eigs(A, 10, 'lm') gives me the answer relatively faster. And as I tried, replacing 10 with A_width in eigs(A, 10, 'lm') so that this includes all the eigen vectors, doesn't solve this problem, 'cause this make it the as slow as using 'sm'.

So, I want to know why calculating the smallest vectors(using 'sm') is much slower than calculating the largest?

BTW, if you have any idea about how to use eigs with 'sm' as fast as with 'lm', please tell me that.

Answer 1

Since eigs is actually an m-file function, we can profile it. I have run a couple of basic tests, and it depends very much on the nature of the data in the matrix. If we run the profiler separately on the following two lines of code:

eigs(eye(1000), 10, 'lm'), and
eigs(eye(1000), 10, 'sm'),

then in the first instance it calls arpackc (the main function that does the work - according to the comments in eigs it's probably from here) a total of 22 times. In the second instance it is called 103 times.

On the other hand, trying it with

eigs(rand(1000), 10, 'lm'), and
eigs(rand(1000), 10, 'sm'),

I get results where the 'lm' option consistently calls arpackc many more times than the sm option.

I'm afraid I don't know the details of the algorithm, and so can't explain it in any deeper mathematical sense, but the page that I linked suggests ARPACK is best for matrices with some structure. Since matrices generated by rand have little structure, it is probably safe to assume the latter behaviour I described is not what you'd expect under normal operating conditions.

In short: it simply takes the algorithm more iterations to converge when you ask it for the smallest eigenvalues of a structured matrix. This being an iterative process, however, it very much depends on the actual data you give it.

Edit: There is a wealth of information and references about this method here, and the key to understanding exactly why this happens is surely contained somewhere therein.

Answer 2

The algorithm used in pretty much any standard eigs function is (some variation of) the Lanczos algorithm. It is iterative and the first iterations give you the largest eigenvalues. This explains pretty much every observation you make:

1. Largest eigenvalues take the least amount of iterations,
2. Smallest eigenvalues take the maximum amount of iterations,
3. All eigenvalues also take the maximum amount of iterations.

There are tricks to "fool" eigs into calculating the smallest eigenvalues by actually making them the largest eigenvalues of another problem. This is usually accomplished by a shift parameter. Skimming over the Matlab documentation for eigs, I see that they have a sigma parameter, which might help you. Note the same documentation recommends proper eig if the matrix fits into memory, as eigs has its numerical quirks.

Answer 3

The reason is actually much more simple and due to the basics of solving large sparse eigenvalue problems. These are all based on solving:

(1) A x = lam x

Most solution methods use some power law (e.g. a Krylov subspace spanned in both the Lanczos and Arnoldi methods)

The thing is that the a power series converge to the largest eigenvalue of (1). Therefore we have that the largest eigenvalues are found by the subspace spanned by: K^k = {A*r0,....,A^k*r0}, which requires only matrix vector multiplications (cheap).

To find the smallest, we have to reformulate (1) as follows:

(2) 1/lam x = A^(-1) x or A^(-1) x = invlam x  

Now solving for the largest eigenvalue of (2) is equivalent to finding the smallest eigenvalue of (1). In this case the subspace is spanned by K^k = {A^(-1)*r0,....,A^(-k)*r0}, which requires solving several linear system (expensive!).