Re: eig returns nonorthogonal eigenvectors

[Nicholas Jankowski]

On Mon, Nov 21, 2016 at 5:33 PM, chris2 <address@hidden> wrote:

Hello,
  The following commands show a case where octave returns non-orthogonal
eigenvectors.  Mathematically, the eigenvalue is degenerate, so nothing
prohibits this, even for a symmetric matrix, but it would be very useful if
it did.  Is there any way to force octave to orthogonalize the basis?
--thanks, Chris

octave> [q,r]=qr([1, -1, 1, 0; 2, 3, 0, 1; 4, 3, 0, 0; -1, 2, 0, 0])
octave> a=q*diag([5,5,-7,9])*q'
octave> [va,da]=eig(a)
octave> va'*va

not a subject I'm intimately familiar with, but following a similar discussion regarding matlab's eig:
http://stackoverflow.com/questions/33258139/matlab-not-returning-orthonormal-matrix-of-eigenvectors

it seems your non-orthogonal results are the result of numerical precision errors.  From that discussion they mention that eig will return orthogonal eigenvectors if the input is exactly symmetric.  However, from your input after creating a:

>> a
a =
   -3.7972   -1.4947    1.7509   -4.7829
   -1.4947    7.0956   -1.1977   -2.0942
    1.7509   -1.1977    5.6030    1.7675
   -4.7829   -2.0942    1.7675    3.0986
>> a-a'
ans =
   1.0e-15 *
         0    0.2220         0         0
   -0.2220         0         0         0
         0         0         0   -0.2220
         0         0    0.2220         0

(not exactly symmetric)

>> [va,da]=eig(a)
va =
    0.8562   -0.5166   -0.0000    0.1073
    0.1455    0.2411   -0.7664   -0.6148
   -0.1704   -0.2824    0.4877   -0.7320
    0.4655    0.7715    0.4180   -0.2732
da =
   -7.0000         0         0         0
         0    5.0000         0         0
         0         0    9.0000         0
         0         0         0    5.0000
>> va'*va
ans =
    1.0000    0.0000    0.0000   -0.0000
    0.0000    1.0000   -0.0000   -0.2077
    0.0000   -0.0000    1.0000   -0.0000
   -0.0000   -0.2077   -0.0000    1.0000

Now, if you force a to be symmetric:

>> a = (a+a')/2
a =
   -3.7972   -1.4947    1.7509   -4.7829
   -1.4947    7.0956   -1.1977   -2.0942
    1.7509   -1.1977    5.6030    1.7675
   -4.7829   -2.0942    1.7675    3.0986
>> a-a'
ans =
     0     0     0     0
     0     0     0     0
     0     0     0     0
     0     0     0     0
>> [va,da]=eig(a)
va =
    0.8562    0.5068   -0.1003    0.0000
    0.1455   -0.1244    0.6132    0.7664
   -0.1704    0.4340    0.7381   -0.4877
    0.4655   -0.7344    0.2631   -0.4180
da =
   -7.0000         0         0         0
         0    5.0000         0         0
         0         0    5.0000         0
         0         0         0    9.0000
>> va'*va
ans =
    1.0000    0.0000   -0.0000   -0.0000
    0.0000    1.0000   -0.0000    0.0000
   -0.0000   -0.0000    1.0000    0.0000
   -0.0000    0.0000    0.0000    1.0000

>> va'*va-eye(4)
ans =
   1.0e-15 *
    0.4441    0.0555   -0.0278   -0.0555
    0.0555    0.2220   -0.1110    0.0555
   -0.0278   -0.1110    0.4441    0.1388
   -0.0555    0.0555    0.1388    0.6661

So it's orthogonal to within 1e-15.  there was some discussion at that link about how such small errors in a producing large differences in the result.