bugGNU Octave - Bugs: bug #38577, Eig returns non-unitary...

 
 

bug #38577: Eig returns non-unitary transformation matrix

Submitted by:  None
Submitted on:  Sat 23 Mar 2013 04:44:18 PM UTC  
 
Category: NoneSeverity: 3 - Normal
Priority: 5 - NormalItem Group: Incorrect Result
Status: InvalidAssigned to: None
Originator Name: PeterOriginator Email: -unavailable-
Open/Closed: ClosedRelease: 3.6.2
Operating System: GNU/Linux

Add a New Comment (Rich MarkupRich Markup):
   

You are not logged in

Please log in, so followups can be emailed to you.

 

(Jump to the original submission Jump to the original submission)

Thu 28 Mar 2013 08:49:33 AM UTC, comment #10:

Hi,

indeed, I have LaPack 3.4.1 installed.

Using LD_LIBRARY_PATH to make octave using
ACML-5.3 I get correct results.

Amazing.

Thanks for the clarification,
Peter

Anonymous
Thu 28 Mar 2013 03:42:45 AM UTC, comment #9:

I also see this bug with LAPACK 3.4.1 but not with LAPACK 3.4.2 from Debian experimental. Looks resolved.

Mike Miller <mtmiller>
Project Member
Wed 27 Mar 2013 01:52:18 PM UTC, comment #8:

I managed to get access to my development machine, and I confirm that Octave is calling LAPACK's dsyev correctly. The problem is likely therefore to be in LAPACK:

http://www.netlib.org/lapack/bug_list.html

This seems like bug 0098, or at least related to it.

I am running LAPACK 3.4.1, where this bug is present. What LAPACK version are you using? I am tentatively closing this report. If you still believe that the bug is in Octave and not in LAPACK, we may reopen it.

Jordi Gutiérrez Hermoso <jordigh>
Project Administrator
Wed 27 Mar 2013 01:19:10 PM UTC, comment #7:

> Hermitian matrices can never be spectrally defective.


But they can be very close to defective. A defective matrix is somewhat numerically invisible, since any perturbation destroys its defect. This is why the Jordan form is somewhat useless, since it's numerically unstable.

> Eig returning a zero eigenvector for a hermitian matrix is a sign of a bug.


Not necessarily, a defective matrix can be close to being Hermitian, e.g. [1 eps; 0 1].

The eigenvalue routine that Octave uses isn't implemented in Octave. It's just calling LAPACK's dgeev family of functions:

http://hg.savannah.gnu.org/hgweb/octave/file/a2f65b8f1955/liboctave/numeric/EIG.cc#l31

I'm not sure how your matrix is getting close to being defective, and it seems to depend on the algorithm since the Schur factorisation doesn't exhibit this problem. I have had reports that different Octave builds on Windows don't have this problem.

I don't have a development machine right now to test on, but if you're able to build Octave and step through it with a debugger, you should check if indeed the LAPACK dseev symmetric eigenvalue algorithm is being called for your matrix. If it isn't, then there's a possible bug in Octave. If it is, then the problem is in LAPACK.

Jordi Gutiérrez Hermoso <jordigh>
Project Administrator
Wed 27 Mar 2013 09:08:59 AM UTC, comment #6:

Hermitian matrices can never be spectrally defective. They are guaranteed to have an orthogonal basis of eigenvectors. Eig returning a zero eigenvector for a hermitian matrix is a sign of a bug.

For hermitian matrices eig and schur decompositions coincide. The schur decomposition
S = U' * A * U
implies that the triangular matrix S is hermitian if A is. Being both triangular and hermitian forces it to be diagonal.

Anonymous
Mon 25 Mar 2013 04:06:27 PM UTC, comment #5:

I can confirm this. My first instinct is that the matrix is numerically defective, so that v not being of full rank would make sense:

That's of rank 199, not 200. It is akin to the following situation:

However, this doesn't seem to explain it. The original 200x200 matrix is not apparently numerically defective:

Compare with

which is 1, showing the spectral defect.

At any rate, eig returning a zero eigenvector is a sign of a defective matrix. I would need to study deeper the Schur decomposition method to understand this... Curiously enough, as a workaround, you can use schur() instead of eig, since in this case it does produce the eigenvectors (the upper triangular S matrix ends up being not just upper triangular, but diagonal).

I have heard reports that this problem doesn't exist in some Windows builds. Which version of LAPACK and which BLAS are you using?

Jordi Gutiérrez Hermoso <jordigh>
Project Administrator
Mon 25 Mar 2013 09:24:46 AM UTC, comment #4:

Hi,

< http://en.wikipedia.org/wiki/Eigenvalues_and_eigenvectors > :
An eigenvector of a square matrix A is a non-zero _vecto_ v that, when multiplied by A, yields the original vector multiplied by a single number .

<http://mathworld.wolfram.com/Eigenvector.html > :
Eigenvectors may not be equal to the zero vector.

Peter

Anonymous
Mon 25 Mar 2013 07:26:44 AM UTC, comment #3:

Hi,

> The zero vector for x solves this equation and is orthogonal
> to the other eigenvectors which are all zeros except for a single
> '1' on the diagonal. Octave appears to place this as the very
> last eigenvector >every time. You can check this with
> [...]
> So is this issue merely definitional?


Actually it is not merely definitional.
A * 0 = 0 is true for any matrix. The main point is that the
eigen vectors of a hermitian matrix form a complete basis.
With one zero vector this is not true, one vector is missing.

Actually, I didn't construct this example on purpose, but found it
while debugging a script, where I perform base transformations into the basis of eigen vectors of a hermitian matrix. At least
in quantum mechanics this is essential.

Peter

P.S. I'm sorry, but I still have to read the bug-tracker documentation to understand the mark-up.

Anonymous
Sun 24 Mar 2013 06:54:52 PM UTC, comment #2:

Apparently the markup language ate some of my last post. Here it is again:

The eigenvalues in E are correct ([1, -1, lots of zeros]). This can be checked with sum (E(:)) which is always 0.

For the eigenvectors with zero eigenvalues one needs to solve

The zero vector for x solves this equation and is orthogonal to the other eigenvectors which are all zeros except for a single '1' on the diagonal. Octave appears to place this as the very last eigenvector every time. You can check this with

So is this issue merely definitional? Octave believes the zero vector is a valid eigenvector because it meets the criteria above, but ordinary usage does not include the zero vector as an eigenvector?

Rik <rik5>
Project Administrator
Sun 24 Mar 2013 03:38:41 PM UTC, comment #1:

The eigenvalues in E are correct ([1, -1, lots of zeros]). This can be checked with sum (E(:)) which is always 0.

For the eigenvectors with zero eigenvalues one needs to solve Ax = 0x. The zero vector for x solves this equation and is orthogonal to the other eigenvectors which are all zeros except for a single '1' on the diagonal. Octave appears to place this as the very last eigenvector every time. You can check this with diag (UU')(end-1).

So in this issue merely definitional? Octave believes the zero vector is a valid eigenvector because it meets the criteria above, but ordinary usage does not include the zero vector as an eigenvector?

Rik <rik5>
Project Administrator
Sat 23 Mar 2013 04:44:18 PM UTC, original submission:

Dear All,

the following code

M=200;
A=zeros(M,M);

Index1= ( M - mod(M,2) ) / 2;

A(M,M-Index1) = -1.0;
A(M-Index1,M) = -1.0;

[ U, E ] = eig( A );

z = max ( abs( diag( U * U' ) - ones (M,1)))

produces a '1' as result, instead of a numerical zero.
Tested on octave 3.2.4 and 3.6.2.
It looks like one eigen vector is just a zero vector.
On 3.6.2 it work for M<130, and fails for M >=130.
At least I checked until M=1000.

Peter

Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach File(s):
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -unavailable- added by mtmiller (Posted a comment)
  • -unavailable- added by jordigh (Posted a comment)
  • -unavailable- added by rik5 (Posted a comment)
  •  

    Do you think this task is very important?
    If so, you can click here to add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

    Please enter the title of George Orwell's famous dystopian book (it's a date):

     

     

    Follow 3 latest changes.

    Date Changed By Updated Field Previous Value => Replaced By
    Wed 27 Mar 2013 01:52:18 PM UTCjordighStatusConfirmed=>Invalid
      Open/ClosedOpen=>Closed
    Mon 25 Mar 2013 04:06:27 PM UTCjordighStatusNone=>Confirmed

    Back to the top


    Powered by Savane 3.1-cleanup