Re: Re: Eigensystem sometimes returns
- To: mathgroup at smc.vnet.net
- Subject: [mg96905] Re: [mg96874] Re: [mg96805] Eigensystem sometimes returns
- From: DrMajorBob <btreat1 at austin.rr.com>
- Date: Fri, 27 Feb 2009 06:10:29 -0500 (EST)
- References: <200902250904.EAA15732@smc.vnet.net>
- Reply-to: drmajorbob at bigfoot.com
"Out of interest, does anyone know WHY Eigensystem can returns
non-orthonormal eigenvectors?"
The eigenvectors of a matrix are generally NOT orthogonal (since that's
what you meant, rather than orthonormal).
The eigenvectors ARE usually independent (if the matrix is full rank), so
you can transform all of space until they are orthogonal... but that has
little to do with the original matrix and its eigensystem.
Bobby
On Thu, 26 Feb 2009 06:59:03 -0600, Yen Lee Loh <yloh at mps.ohio-state.edu>
wrote:
> 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
> degenEigvalIndices=Split[Sort[d]]
> 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
> readers.
>
>
>
>
>
>
> On Wed, Feb 25, 2009 at 10:02 AM, Daniel Lichtblau <danl at wolfram.com>
> 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 (
>>> http://forums.wolfram.com/mathgroup/archive/1998/Mar/msg00418.html ),
>>> but
>>> that post is corrupted! (Surely Wolfram isn't resorting to
>>> censorship?)
>>>
>>
>> A non-corrupted version is located at the URL below.
>>
>> http://forums.wolfram.com/mathgroup/archive/1998/Mar/msg00425.html
>>
>> 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., http://arxiv.org/pdf/hep-ph/9607313 ).
>>> 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
>>
>>
>>
>
>
--
DrMajorBob at bigfoot.com
- References:
- Eigensystem[hermitianMatrix] sometimes returns non-orthonormal
- From: Yen Lee Loh <yloh@mps.ohio-state.edu>
- Eigensystem[hermitianMatrix] sometimes returns non-orthonormal