Thu Mar 28 08:49:33 2013, comment #10:
Hi,
indeed, I have LaPack 3.4.1 installed.
Using LD_LIBRARY_PATH to make octave using
ACML5.3 I get correct results.
Amazing.
Thanks for the clarification,
Peter

Thu Mar 28 03:42:45 2013, comment #9:
I also see this bug with LAPACK 3.4.1 but not with LAPACK 3.4.2 from Debian experimental. Looks resolved.

Wed Mar 27 13:52:18 2013, 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.

Wed Mar 27 13:19:10 2013, 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.

Wed Mar 27 09:08:59 2013, 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.

Mon Mar 25 16:06:27 2013, 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?

Mon Mar 25 09:24:46 2013, comment #4:
Hi,
< http://en.wikipedia.org/wiki/Eigenvalues_and_eigenvectors > :
An eigenvector of a square matrix A is a nonzero _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

Mon Mar 25 07:26:44 2013, 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 bugtracker documentation to understand the markup.

Sun Mar 24 18:54:52 2013, 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?

Sun Mar 24 15:38:41 2013, 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')(end1).
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?

Sat Mar 23 16:44:18 2013, original submission:
Dear All,
the following code
M=200;
A=zeros(M,M);
Index1= ( M  mod(M,2) ) / 2;
A(M,MIndex1) = 1.0;
A(MIndex1,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
