On Fri, Jan 15, 2010 at 12:07 PM, <josef.p...@gmail.com> wrote: > On Fri, Jan 15, 2010 at 11:32 AM, Sebastian Walter > <sebastian.wal...@gmail.com> wrote: >> numpy.linalg.eig guarantees to return right eigenvectors. >> evec is not necessarily an orthonormal matrix when there are >> eigenvalues with multiplicity >1. >> For symmetrical matrices you'll have mutually orthogonal eigenspaces >> but each eigenspace might be spanned by >> vectors that are not orthogonal to each other. >> >> Your omega has eigenvalue 1 with multiplicity 3. > > Yes, I thought about the multiplicity. However, even for random > symmetric matrices, I don't get the result > I change the example matrix to > omega0 = np.random.randn(20,8) > omega = np.dot(omega0.T, omega0) > print np.max(np.abs(omega == omega.T)) > > I have been playing with left and right eigenvectors, but I cannot > figure out how I could compose my original matrix with them either. > > I checked with wikipedia, to make sure I remember my (basic) linear algebra > http://en.wikipedia.org/wiki/Eigendecomposition_(matrix)#Symmetric_matrices > > The left and right eigenvectors are almost orthogonal > ev, evecl, evecr = sp.linalg.eig(omega, left=1, right=1) >>>> np.abs(np.dot(evecl.T, evecl) - np.eye(8))>1e-10 >>>> np.abs(np.dot(evecr.T, evecr) - np.eye(8))>1e-10 > > shows three non-orthogonal pairs
This doesn't seem to be correct. I think, I had an old omega with multiplicity of eigenvalues in the interpreter. Writing it as a clean script, I get orthogonal left and right eigenvectors. Thanks for the reply, Josef > >>>> ev > array([ 6.27688862, 8.45055356, 15.03789945, 19.55477818, > 20.33315408, 24.58589363, 28.71796764, 42.88603728]) > > > I always thought eigenvectors are always orthogonal, at least in the > case without multiple roots > > I had assumed that eig will treat symmetric matrices in the same way as eigh. > Since I'm mostly or always working with symmetric matrices, I will > stick to eigh which does what I expect. > > Still, I'm currently not able to reproduce any of the composition > result on the wikipedia page with linalg.eig which is puzzling. > > Josef > >> >> >> >> >> On Fri, Jan 15, 2010 at 4:31 PM, <josef.p...@gmail.com> wrote: >>> I had a problem because linal.eig doesn't rebuild the original matrix, >>> linalg.eigh does, see script below >>> >>> Whats the trick with linalg.eig to get the original (or the inverse) >>> back ? None of my variations on the formulas worked. >>> >>> Thanks, >>> Josef >>> >>> >>> import numpy as np >>> import scipy as sp >>> import scipy.linalg >>> >>> omega = np.array([[ 6., 2., 2., 0., 0., 3., 0., 0.], >>> [ 2., 6., 2., 3., 0., 0., 3., 0.], >>> [ 2., 2., 6., 0., 3., 0., 0., 3.], >>> [ 0., 3., 0., 6., 2., 0., 3., 0.], >>> [ 0., 0., 3., 2., 6., 0., 0., 3.], >>> [ 3., 0., 0., 0., 0., 6., 2., 2.], >>> [ 0., 3., 0., 3., 0., 2., 6., 2.], >>> [ 0., 0., 3., 0., 3., 2., 2., 6.]]) >>> >>> for fun in [np.linalg.eig, np.linalg.eigh, sp.linalg.eig, sp.linalg.eigh]: >>> print fun.__module__, fun >>> ev, evec = fun(omega) >>> omegainv = np.dot(evec, (1/ev * evec).T) >>> omegainv2 = np.linalg.inv(omega) >>> omegacomp = np.dot(evec, (ev * evec).T) >>> print 'composition', >>> print np.max(np.abs(omegacomp - omega)) >>> print 'inverse', >>> print np.max(np.abs(omegainv - omegainv2)) >>> >>> this prints: >>> >>> numpy.linalg.linalg <function eig at 0x017EDDF0> >>> composition 0.405241032278 >>> inverse 0.405241032278 >>> >>> numpy.linalg.linalg <function eigh at 0x017EDE30> >>> composition 3.5527136788e-015 >>> inverse 7.21644966006e-016 >>> >>> scipy.linalg.decomp <function eig at 0x01DB14F0> >>> composition 0.238386662463 >>> inverse 0.238386662463 >>> >>> scipy.linalg.decomp <function eigh at 0x01DB1530> >>> composition 3.99680288865e-015 >>> inverse 4.99600361081e-016 >>> _______________________________________________ >>> NumPy-Discussion mailing list >>> NumPy-Discussion@scipy.org >>> http://mail.scipy.org/mailman/listinfo/numpy-discussion >>> >> _______________________________________________ >> NumPy-Discussion mailing list >> NumPy-Discussion@scipy.org >> http://mail.scipy.org/mailman/listinfo/numpy-discussion >> > _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion