[SciPy-dev] Comments on API for Matlab's eigs equivalent (computing a few eigenvalues only)

David Cournapeau david at silveregg.co.jp
Thu Feb 4 01:49:42 EST 2010


Charles R Harris wrote:
> 
> 
> On Wed, Feb 3, 2010 at 11:27 PM, David Cournapeau <david at silveregg.co.jp 
> <mailto:david at silveregg.co.jp>> wrote:
> 
>     Hi,
> 
>     I have played a bit with adding support for computing a few eigenvalues
>     of full symmetric matrices. I would like some comments on the
>     current API:
> 
>     import numpy as np
>     from scipy.linalg import eigs
>     x = np.random.randn(10, 10)
>     x = np.dot(x.T, x)
>     # Retrieve the 3 biggest eigenvalues
>     eigs(x, 3)[0]
>     # Retrieve the 3 smallest eigenvalues
>     eigs(x, -3)[0]
>     # Retrieve the 3 biggest eigenvalues
>     eigs(x, [0, 3], mode="index")[0]
>     # Retrieve the 2nd and 3rd biggest
>     eigs(x, [1, 3], mode="index")[0]
>     # Retrieve all the eigenvalues in the range [1.5, 3.5[
>     eigs(x, [1.5, 3.5], mode="range")[0]
> 
>     One thing which does not feel right is that that the range in the
>     "index" mode is exactly inverted compared to the output (i.e. if you ask
> 
> 
> Why not use separate keywords for index and range.

Then I am not sure how to handle the "give me the k biggest/smallest" 
case, which is the most common I think (that's the only one I care 
personally :) ).

Maybe this just warrants several functions ? But then there is the issue 
of naming (supporting non-symmetric matrices would be nice, but requires 
a totally different implementation, as LAPACK does not support it AFAIK 
- the easiest way would be to use ARPACK ATM).

> Changing the meaning 
> of an argument using another keyword is just weird. 

Agreed - this is just the best I could came up with one function to 
support all cases while keeping the common case simple. But I feel a bit 
like I outsmarted myself here.

> Why the index on the end? Is the return a list?

Yes, it also returns eigenvectors

> 
>     for the range [0, 3], you get the last three items from what you would
>     get if you asked for the full range [0, 10]), but this is because I kept
>     compatibility with Octave (always showing from biggest to smallest). It
> 
> 
> Why follow Octave, isn't Octave like Matlab? Matlab functions are a mess.

Actually, there is another reason I forgot to mention: that's how eigen 
in scipy.sparse does as well (in decreasing order). I think that the 
ARPACK wrappers are not that well thought out, though.

cheers,

David



More information about the SciPy-Dev mailing list