linalg.eigh orders eigenvalues/eigenvectors differently than linalg.eig
Hello all, It seems that the 'eigh' routine from numpy.linalg does not follow the same convention as numpy.linalg.eig in terms of the order of the returned eigenvalues. (And thus eigenvectors as well...) Specifically, eig returns eigenvalues in order from largest to smallest, while eigh returns them from smallest to largest. Example:
a = numpy.array([[21, 28, 35],[28, 38, 48],[35, 48, 61]]) numpy.linalg.eigh(a) (array([ -1.02825542e-14, 7.04131679e-01, 1.19295868e+02]), array([[ 0.40824829, -0.81314396, -0.41488581], [-0.81649658, -0.12200588, -0.56431188], [ 0.40824829, 0.56913221, -0.71373795]]))
numpy.linalg.eig(a) (array([ 1.19295868e+02, 7.04131679e-01, 4.62814557e-15]), array([[-0.41488581, -0.81314396, 0.40824829], [-0.56431188, -0.12200588, -0.81649658], [-0.71373795, 0.56913221, 0.40824829]]))
Is this a bug? If it is, though, fixing it now might break code that depends on this 'wrong' order. (This is also present in scipy.linalg.) If not a bug, or not-fixable-now, then at least some documentation as to the convention regarding ordering of eigenvalues is probably worthwhile... Any thoughts? Zach
Zachary Pincus schrieb:
Hello all,
It seems that the 'eigh' routine from numpy.linalg does not follow the same convention as numpy.linalg.eig in terms of the order of the returned eigenvalues. (And thus eigenvectors as well...)
I was told on this list that the ordering should not be relied upon, and that it might change in the future. So it seems that user code should explicitly re-order the eigenvalues (and corresponding eigenvectors, probably using argsort and fancy indexing -- if those are the right terms). ...
scipy.linalg.) If not a bug, or not-fixable-now, then at least some documentation as to the convention regarding ordering of eigenvalues is probably worthwhile...
true cheers, sven
On Monday 19 February 2007 12:06, Sven Schreiber wrote:
Zachary Pincus schrieb:
Hello all,
It seems that the 'eigh' routine from numpy.linalg does not follow the same convention as numpy.linalg.eig in terms of the order of the returned eigenvalues. (And thus eigenvectors as well...)
I was told on this list that the ordering should not be relied upon, and that it might change in the future. So it seems that user code should explicitly re-order the eigenvalues (and corresponding eigenvectors, probably using argsort and fancy indexing -- if those are the right terms).
Indeed. eig and eigh are wrappers for lapack functions, so the result is whatever those give back. Do not rely on a particular order of eigenvalues, sort yourself. Short example for convenience: #--------- eigvals, eigvecs = eig(some_matrix) ind = argsort(eigvals) eigvals = eigvals[ind] eigvecs = eigvecs[:, ind] # second axis !! # etc. #------------ Johannes
Johannes Loehnert wrote:
On Monday 19 February 2007 12:06, Sven Schreiber wrote:
Zachary Pincus schrieb:
Hello all,
It seems that the 'eigh' routine from numpy.linalg does not follow the same convention as numpy.linalg.eig in terms of the order of the returned eigenvalues. (And thus eigenvectors as well...)
I was told on this list that the ordering should not be relied upon, and that it might change in the future. So it seems that user code should explicitly re-order the eigenvalues (and corresponding eigenvectors, probably using argsort and fancy indexing -- if those are the right terms).
Indeed. eig and eigh are wrappers for lapack functions, so the result is whatever those give back. Do not rely on a particular order of eigenvalues, sort yourself.
Short example for convenience: #--------- eigvals, eigvecs = eig(some_matrix) ind = argsort(eigvals) eigvals = eigvals[ind] eigvecs = eigvecs[:, ind] # second axis !! # etc. #------------
Johannes _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion
eigh is based on *heev (* placeholder for s,d,c,z) http://www.netlib.org/lapack/complex16/zheev.f * W (output) DOUBLE PRECISION array, dimension (N) * If INFO = 0, the eigenvalues in ascending order. Apparently the eigenvalues are in ascending order for eigh. Hence, It would be nice if scipy could offer more wrappers for matrices with certain properties, e.g. http://www.netlib.org/lapack/double/dsygv.f http://www.netlib.org/lapack/double/dsyevr.f http://www.netlib.org/lapack/double/dsygvx.f http://www.netlib.org/lapack/double/dstevr.f Any comments ? Nils Also LAPACK3.1 has many improved eigensolvers. Please try suzuki.py with LAPACK3.0
w_e array([-200. +0.j , -100. +0.j , 700.00307816+0.00103379j, 700.00307816-0.00103379j, 699.9993562 +0.00318262j, 699.9993562 -0.00318262j, 699.99756564+0.00214883j, 699.99756564-0.00214883j])
versus LAPACK3.1.
w_e array([-200. +0.j , -100. +0.j , 699.99891191+0.00188463j, 699.99891191-0.00188463j, 700.00217619+0.j , 700.00294523+0.j , 699.99852738+0.0025507j , 699.99852738-0.0025507j ])
participants (4)
-
Johannes Loehnert -
Nils Wagner -
Sven Schreiber -
Zachary Pincus