MathGroup Archive 2009

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Eigensystem sometimes returns non-orthonormal

  • To: mathgroup at
  • Subject: [mg96874] Re: [mg96805] Eigensystem sometimes returns non-orthonormal
  • From: Yen Lee Loh <yloh at>
  • Date: Thu, 26 Feb 2009 07:59:03 -0500 (EST)
  • References: <> <>

Dear Daniel,

Thanks for your suggestion.  Unfortunately it doesn't fix the problem:

    h = [a certain 3x3 hermitian matrix]
    {u, d} = Eigensystem[h];   d=DiagonalMatrix[d]

The eigenvalue equation is satisfied, as per the specifications:
h.Transpose[u] - Transpose[u].d == 0.
The eigenvectors are normalized, but not orthonormal:
u.ConjugateTranspose[u] != IdentityMatrix[3].
If we compute the q part of the QRDecomposition,

    q = First@QRDecomposition[u]  ,

then q.ConjugateTranspose[q] == IdentityMatrix[3], but h.Transpose[q] -
Transpose[q].d ==0,
so the eigenvalue equation is no longer satisfied.   q =Transpose@First
@QRDecomposition[Transpose@u] doesn't work either.

I haven't analysed things in detail, but since the task at hand is to
"orthogonalize only the sets of eigenvectors u corresponding to degenerate
eigenvalues", it may not be safe to use an approach
(Orthogonalize/QRDecomposition) that doesn't "know" about the pattern of
eigenvalues.  Of course one can write a program beginning with
and then calling Orthogonalize within each set, but this is kinda messy.
I'm wondering if it is actually easier to do the calculation in Fortran/C++
and use MathLink, somehow.  Does Wolfram have any plans to add support for
orthonormal eigenvectors?

Out of interest, does anyone know WHY Eigensystem can returns
non-orthonormal eigenvectors?  I was under the impression that Eigensystem[]
first checks to see if the matrix is Hermitian, and if it is, it uses
LAPACK-based algorithms for Hermitian matrices.  I thought LAPACK's
algorithms involved tridiagonalization followed by some form of QR
factorization/divide-and-conquer Givens rotations, in which the eigenvector
matrix started off as the identity matrix and accumulated successive unitary
transformations, thus remaining unitary at every stage;  I'm curious how
Eigensystem[] even manages to break orthonormality at all.

p.s.  DrMajorBob: Thanks for your suggestion, but my eigenvectors are
already normalized --- the problem is that they're not mutually orthogonal,
which seems to be trickier to fix.
p.p.s.  Sorry this post is rather badly written.  Maybe when the final
solution appears I will summarize things properly for the benefit of future

On Wed, Feb 25, 2009 at 10:02 AM, Daniel Lichtblau <danl at> wrote:

> Yen Lee Loh wrote:
>> Dear Mathematica users,
>> In Mathematica 7.0.0 (for Linux), calling Eigenvectors[H] or
>> Eigensystem[H]
>> for a numerical Hermitian matrix H sometimes returns eigenvectors that are
>> not orthonormal.
>> This happens when some eigenvalues are degenerate.  (I can supply example
>> code that illustrates the problem, if necessary.)
>> This is not really a bug -- the documentation for Eigensystem[] doesn't
>> make
>> any guarantees of orthonormality -- but nevertheless it is an annoying
>> part
>> of Mathematica's design.
>> This issue has been raised 11 years ago (
>> ), but
>> that post is corrupted!  (Surely Wolfram isn't resorting to censorship?)
> A non-corrupted version is located at the URL below.
> Note that the behavior in that report was of a more serious nature than
> that which you describe (eigenvectors sometimes had complex values).
>  I
>> was hoping that in Mathematica 7 I would be able to write something like
>>    Eigensystem[H, Method->"LAPACK-ZHEEVR"]
>> or
>>    Eigensystem[H, OrthonormalizeEigenvectors->True]
>> but no such options seem to exist.  One workaround is to apply
>> Orthogonalize[] to the matrix of eigenvectors, but the documentation for
>> Orthogonalize[] doesn't guarantee that the orthonormalization will only
>> occur within the "degenerate subspace".  So one has to resort to
>> complicated
>> fixes (e.g., ).
>> Does anyone have a simpler solution?  (For example, is there an easy way
>> to
>> call LAPACK'S ZHEEVR routine, which guarantees orthornormal eigenvectors,
>> from Mathematica?)
>> Thanks a lot.
>> Yen Lee Loh
> If you apply QRDecomposition to the eigenvectors, then I believe the Q part
> will provide what you are looking for.
> Daniel Lichtblau
> Wolfram Research

Yen Lee Loh
Postdoctoral Associate, The Ohio State University
Home: 544 Stinchcomb Dr Apt 10, Columbus OH 43202-1728, USA
Office: 2043 Physics Research Building, 191 W Woodruff Ave, Columbus OH
43210-1117, USA
Office phone: +1 614 247 4772
Mobile phone: +1 765 532 9457
Email: yloh at

  • Prev by Date: Re: Graphics export for high quality documents
  • Next by Date: Re: Re: listing user defined, independent variables
  • Previous by thread: Re: Eigensystem sometimes returns non-orthonormal
  • Next by thread: Re: Eigensystem[hermitianMatrix]sometimes returns