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

josef.pktd@gmai... josef.pktd@gmai...
Fri Feb 5 05:59:05 CST 2010


On Thu, Feb 4, 2010 at 9:24 PM,  <josef.pktd@gmail.com> wrote:
> On Thu, Feb 4, 2010 at 8:27 PM, David Cournapeau <cournape@gmail.com> wrote:
>> On Fri, Feb 5, 2010 at 9:43 AM, Warren Weckesser
>> <warren.weckesser@enthought.com> wrote:
>>
>>> David C,
>>>
>>> You said this was for symmetric matrices, but do you envision later
>>> allowing nonsymmetric matrices?
>>
>> Yes. I have only implemented symmetric/Hermitian because that's the
>> only solver that handles getting only a few eigenvalues in LAPACK.
>> Matlab own eigs function uses ARPACK for both symmetric/unsymmetric
>> cases, which is not good according to one of the Lapack developer
>> (http://mail.scipy.org/pipermail/scipy-dev/2007-March/006778.html).
>> But we could use ARPACK for the general case in eigs.
>>
>>>
>>> If not, then perhaps the name of the function should be 'eigsh',
>>> following the precedent set by numpy.linalg and scipy.linalg.
>>
>> I wonder whether those eigh functions are a good idea: I fear that
>> most people will always use eig - maybe one could use the underlying
>> eig*h solver in eig* if the matrix is detected as being symmetric ? I
>> am not really knowledgeable about those issues, though. For example, I
>> don't know whether the symmetric aspect should be checked exactly, or
>> if it is better to use a symmetric solver  even if there are very
>> small errors in say A'-A.
>>
>>>
>>> In particular, the choice of ordering by magnitude or by real part is
>>> convenient.
>>
>> It seems that one would need to implement the non-symmetric
>> capabilities to sort this out. I fear that those options are solver
>> specific, though - maybe the solution is to have two levels of API,
>> one low-level and one high level.
>>
>> cheers,
>>
>> David
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev@scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
>
> the current version of scipy.linalg.eigh sorts from smallest to largest
>
>>>> scipy.linalg.eigh(x)[0]
> array([  4.04457427e-03,   1.84073286e-01,   6.74875960e-01,
>         3.23328824e+00,   4.00741304e+00,   6.98333680e+00,
>         1.45314842e+01,   1.49260377e+01,   2.47166702e+01,
>         2.98955755e+01])
>>>> scipy.linalg.eigh(x, eigvals=(0,3))[0]
> array([ 0.00404457,  0.18407329,  0.67487596,  3.23328824])
>>>> scipy.linalg.eigh(x, eigvals=(len(x)-3, len(x)-1))[0]
> array([ 14.92603769,  24.71667021,  29.89557547])
>
>
>>>> x = np.random.randn(10, 10) + 1j*np.random.randn(10, 10)
>>>> x = np.dot(x.T, x)
>>>> scipy.linalg.eigh(x, eigvals=(0,3))[0]
> array([-36.82039662, -18.20362967, -12.21716065,  -9.79752274])
>>>> scipy.linalg.eigh(x, eigvals=(len(x)-3, len(x)-1))[0]
> array([ 14.44815917,  28.8394026 ,  42.45471058])

How does the new implementation relate to the existing implementation
of selecting just a few eigenvalues in a range that is possible with
the current scipy.linalg.eigh ?

Tiziano added this with the integration of symeig into scipy.linalg.

Josef


> Josef
>


More information about the SciPy-Dev mailing list