9

I was wondering if there is a Python package, numpy or otherwise, that has a function that computes the first eigenvalue and eigenvector of a small matrix, say 2x2. I could use the linalg package in numpy as follows.

import numpy as np

def whatever():
    A = np.asmatrix(np.rand(2, 2))
    evals, evecs = np.linalg.eig(A)
    #Assume that the eigenvalues are ordered from large to small and that the
    #eigenvectors are ordered accordingly.
    return evals[0], evecs[:, 0]

But this takes a really long time. I suspect that it's because numpy computes eigenvectors through some sort of iterative process. So I was wondering if there were a much faster algorithm that only returns the first (largest) eigenvalue and eigenvector, since I only need the first.

For 2x2 matrices of course I can write a function myself, that computes the eigenvalue and eigenvector analytically, but then there are problems with floating point computations, for example when I divide a very big number by a very small number, I get infinity or NaN. Does anyone know anything about this? Please help! Thank you in advance!

Ray
  • 7,833
  • 13
  • 57
  • 91
  • Umm, what do you mean by "long time"? `%timeit np.linalg.eig(np.random.rand(2,2))` gives `10000 loops, best of 3: 208 us per loop`. That's reasonably fast, I must say. – Avaris Oct 20 '11 at 17:17
  • 1
    For 2x2 Matrix numpy is slow due the overhead, just use: http://www.math.harvard.edu/archive/21b_fall_04/exhibits/2dmatrices/index.html – tillsten Oct 21 '11 at 14:00
  • @Avaris, don't include a call to a random-number generator in your timing tests! – Ahmed Fasih Oct 01 '14 at 16:44

3 Answers3

7

Use this: http://docs.scipy.org/doc/scipy/reference/sparse.linalg.html

http://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.eigs.html#scipy.sparse.linalg.eigs

Find k eigenvalues and eigenvectors of the square matrix A.
Jack Twain
  • 6,273
  • 15
  • 67
  • 107
  • Selecting / upvoting link-only answers doesn't help in the long term (and apparently not in the short term either, as the [eigs function](http://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.eigs.html#scipy.sparse.linalg.eigs) doesn't seem to be a valid solution). The site encourages true answers, like [this one](https://stackoverflow.com/a/7839553/774575). – mins Jan 14 '21 at 12:14
2

According to the docs:

http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eig.html

and also to my own experience, numpy.linalg.eig(A) does NOT sort the eigenvectors in any particular order, which is what the OP and subsequent seem to be assuming. I suggest something like:

rearrangedEvalsVecs = sorted(zip(evals,evecs.T),\
                                    key=lambda x: x[0].real, reverse=True)
Josh Berryman
  • 21
  • 3
  • [Greg Willden](http://stackoverflow.com/users/3385672/user3385672) [comments](http://stackoverflow.com/review/suggested-edits/4240532): It should be `evecs.T` here to add transpose to evecs in the zip statement. The zip was associating the eigenvalues with the rows of evecs but should have been associated with the columns of evecs. – Rup Mar 05 '14 at 22:12
1

There doesn't appear to be a numpy equivalent of Matlab's eigs(A,B,k) for finding the k largest eigenvectors.

If you're interested, Enthought has compiled a table showing the differences between Matlab and numpy. That should be helpful for answering questions like this one: Link

One other thought, for 2x2 matrices, I don't think eigs(A,B,1) would help anyway. The effort involved in computing the first eigenpair leaving the matrix transformed to where the second emerges directly. There is only benefit for 3x3 and larger.

Glorfindel
  • 21,988
  • 13
  • 81
  • 109
Raymond Hettinger
  • 216,523
  • 63
  • 388
  • 485
  • 1
    There is an eigs in scipy.sparse.linalg.eigen.arpack (might be under different names depending on the version of scipy), but it is much slower than using numpy.linalg.eig for 3x3 matrices which isn't surprising. It is great when you have huge sparse matrices though. – Justin Peel Oct 20 '11 at 18:34
  • I see, thanks Mr. Hettinger. I guess I'll just have to program my own analytic algorithm using some sort of arbitrary-precision floating point package so I don't get those annoying NaN errors. – Ray Oct 22 '11 at 10:41
  • Matlab's `eigs` uses ARPACK, which is exactly what `scipy.sparse.linalg.eigs` uses, even down to the arguments to specify largest/smallest, real/imaginary, etc. (Don't let the "sparse" fool you, it works for dense matrixes too.) Also, `eigs(r, 1)` is slower than `eig(r)` even for 10x10 `r`: the algorithm is just slower. – Ahmed Fasih Oct 01 '14 at 16:42