Mathematica 9 is now available
Services & Resources / Wolfram Forums / MathGroup Archive
-----

MathGroup Archive 2008

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

Search the Archive

Re: Mathematica 6 obtains imaginary eigenvalues for a Hermitian

  • To: mathgroup at smc.vnet.net
  • Subject: [mg86174] Re: [mg86010] Mathematica 6 obtains imaginary eigenvalues for a Hermitian
  • From: Daniel Lichtblau <danl at wolfram.com>
  • Date: Tue, 4 Mar 2008 02:10:30 -0500 (EST)
  • References: <200803010947.EAA23583@smc.vnet.net>

Sebastian Meznaric wrote:
> I have a 14x14 Hermitian matrix, posted at the bottom of this message.
> The eigenvalues that Mathematica obtains using the
> N[Eigenvalues[matrix]] include non-real numbers:
> {-9.41358 + 0.88758 \[ImaginaryI], -9.41358 -
>   0.88758 \[ImaginaryI], -7.37965 + 2.32729 \[ImaginaryI], -7.37965 -
>   2.32729 \[ImaginaryI], -4.46655 + 2.59738 \[ImaginaryI], -4.46655 -
>   2.59738 \[ImaginaryI], 4.36971, 3.21081, -2.32456 +
>   2.10914 \[ImaginaryI], -2.32456 - 2.10914 \[ImaginaryI],
>  2.04366+ 0.552265 \[ImaginaryI],
>  2.04366- 0.552265 \[ImaginaryI], -0.249588 +
>   1.29034 \[ImaginaryI], -0.249588 - 1.29034 \[ImaginaryI]}.
> However, if you do Eigenvalues[N[matrix]] it obtains different results
> {-9.09122, -7.41855, -7.41855, -7.2915, 4.33734, -4., -4., 3.2915, \
> -3.24612, -2.38787, -2.38787, 1.80642, 1.80642, 0}.
> 
> These results agree with Solve[CharacteristicPolynomial[matrix,x],x].
> Therefore I assume that the latter are correct. Has anyone seen this?
> I am using 6.0.0.
> 
> 
> Here is the matrix:
> {{-6, 0, -Sqrt[3], 0, 0, Sqrt[3], 0, 0, 0, 0, 0, 0, 0, 0}, {0, -6,
>   0, -Sqrt[3], 0, 0, Sqrt[3], 0, 0, 0, 0, 0, 0, 0}, {-Sqrt[3], 0, -4,
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), 2 Sqrt[2/3], 0, 0, Sqrt[3], 0,
>   0, 0, 0, 0, 0}, {0, -Sqrt[3],
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), -4/3, -(2 Sqrt[2])/3, 0, 0, 0,
>   Sqrt[3], 0, 0, 0, 0, 0}, {0, 0, 2 Sqrt[2/3], -(2 Sqrt[2])/3, 7/3, 0,
>    0, 0, 0, Sqrt[3], 0, 0, 0, 0}, {Sqrt[3], 0, 0, 0, 0, -4, 0,
>   2 (-1/(4 Sqrt[3]) + Sqrt[3]/4), 0, 0, 2 Sqrt[2/3], 0, 0, 0}, {0,
>   Sqrt[3], 0, 0, 0, 0, -4, 0, 2 (-1/(4 Sqrt[3]) + Sqrt[3]/4), 0, 0,
>   2 Sqrt[2/3], 0, 0}, {0, 0, Sqrt[3], 0, 0,
>   2 (-1/(4 Sqrt[3]) + Sqrt[3]/4), 0, -14/3,
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), 2 Sqrt[2/3], (2 Sqrt[2])/3, 0,
>   0, 0}, {0, 0, 0, Sqrt[3], 0, 0, 2 (-1/(4 Sqrt[3]) + Sqrt[3]/4),
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), -2, -(2 Sqrt[2])/3, 0, (
>   2 Sqrt[2])/3, 0, 0}, {0, 0, 0, 0, Sqrt[3], 0, 0,
>   2 Sqrt[2/3], -(2 Sqrt[2])/3, -7/3, 0, 0,
>   2 (1/(3 Sqrt[2]) + (2 Sqrt[2])/3), Sqrt[10/3]}, {0, 0, 0, 0, 0,
>   2 Sqrt[2/3], 0, (2 Sqrt[2])/3, 0, 0, -16/3,
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), 2 Sqrt[2/3], 0}, {0, 0, 0, 0, 0,
>    0, 2 Sqrt[2/3], 0, (2 Sqrt[2])/3, 0,
>   2 (-1/(4 Sqrt[3]) + (3 Sqrt[3])/4), -8/3, -(2 Sqrt[2])/3, 0}, {0, 0,
>    0, 0, 0, 0, 0, 0, 0, 2 (1/(3 Sqrt[2]) + (2 Sqrt[2])/3),
>   2 Sqrt[2/3], -(2 Sqrt[2])/3, 1/2,
>   2 (-Sqrt[5/3]/16 - Sqrt[15]/16)}, {0, 0, 0, 0, 0, 0, 0, 0, 0, Sqrt[
>   10/3], 0, 0, 2 (-Sqrt[5/3]/16 - Sqrt[15]/16), 7/2}}

I wanted to confirm that this is a bug in the symbolic eigenvalue code. 
It was fixed in our development version back in December; the fix was 
deemed too risky to try to backport to version 6.0.2 (which was very 
near to completion at that time).

The problem is a bit obscure and has to do with some details of row 
reduction when a matrix contains algebraics. It is not obvious that the 
problem lies there, because CharacteristicPolynomial gets the correct 
polynomial. What happens is that Eigenvalues first shuffles rows and 
columns (symmetrically) according to some heuristic, and then attempts 
to compute the characteristic polynomial of this modified matrix. It 
should agree with CharacteristicPolynomial of the original, up to a 
factor of +-1. But it does not, because with the shuffled matrix the row 
reduction used to find that char poly runs afoul of said bug.

Here we have IdentityMatrix[14]*x - shuffledmatrix.

m2 = {{6 + x, Sqrt[3], 0, 0, 0, 0, 0, -Sqrt[3], 0, 0, 0, 0, 0, 0},
   {Sqrt[3], 4 + x, -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), 0, 0,
     0, -Sqrt[3], 0, 0, 0, 0, 0, 0, -2*Sqrt[2/3]},
   {0, -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), 4/3 + x, Sqrt[3], 0,
     -Sqrt[3], 0, 0, 0, 0, 0, 0, 0, (2*Sqrt[2])/3},
   {0, 0, Sqrt[3], 6 + x, -Sqrt[3], 0, 0, 0, 0, 0, 0, 0, 0, 0},
   {0, 0, 0, -Sqrt[3], 4 + x, -2*(-1/(4*Sqrt[3]) + Sqrt[3]/4), 0, 0,
     0, -2*Sqrt[2/3], 0, 0, 0, 0},
   {0, 0, -Sqrt[3], 0, -2*(-1/(4*Sqrt[3]) + Sqrt[3]/4), 2 + x,
     -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), 0, 0, (-2*Sqrt[2])/3, 0,
	(2*Sqrt[2])/3, 0, 0},
   {0, -Sqrt[3], 0, 0, 0, -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), 14/3 + x,
     -2*(-1/(4*Sqrt[3]) + Sqrt[3]/4), (-2*Sqrt[2])/3, 0, 0, 
-2*Sqrt[2/3], 0, 0},
   {-Sqrt[3], 0, 0, 0, 0, 0, -2*(-1/(4*Sqrt[3]) + Sqrt[3]/4), 4 + x,
     -2*Sqrt[2/3], 0, 0, 0, 0, 0},
   {0, 0, 0, 0, 0, 0, (-2*Sqrt[2])/3, -2*Sqrt[2/3], 16/3 + x,
     -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), -2*Sqrt[2/3], 0, 0, 0},
   {0, 0, 0, 0, -2*Sqrt[2/3], (-2*Sqrt[2])/3, 0, 0,
     -2*(-1/(4*Sqrt[3]) + (3*Sqrt[3])/4), 8/3 + x, (2*Sqrt[2])/3, 0, 0, 0},
   {0, 0, 0, 0, 0, 0, 0, 0, -2*Sqrt[2/3], (2*Sqrt[2])/3, -1/2 + x,
     -2*(1/(3*Sqrt[2]) + (2*Sqrt[2])/3), -2*(-Sqrt[5/3]/16 - 
Sqrt[15]/16), 0},
   {0, 0, 0, 0, 0, (2*Sqrt[2])/3, -2*Sqrt[2/3], 0, 0, 0,
     -2*(1/(3*Sqrt[2]) + (2*Sqrt[2])/3), 7/3 + x, -Sqrt[10/3], -Sqrt[3]},
   {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2*(-Sqrt[5/3]/16 - Sqrt[15]/16), 
-Sqrt[10/3],
     -7/2 + x, 0},
   {0, -2*Sqrt[2/3], (2*Sqrt[2])/3, 0, 0, 0, 0, 0, 0, 0,
     0, -Sqrt[3], 0, -7/3 + x}};

This is the computation that gets the characteristic polynomial wrong.

Expand[Det[m2]]

This has been wrong since version 5.1, although I suspect the underlying 
bug predates that and some small change (pivot selection, perhaps) 
brought it out in this example at that time.

For the record, this was one of the longest bug hunts I ever set out on. 
It involved a subtlety in a method known as "one-step" row reduction. 
Subduing the creature required a wooden stake, a silver bullet, and 
several incantations from bug-killing lore that are not repeatable.

Daniel Lichtblau
Wolfram Research


  • Prev by Date: Re: how insert hyperlink to doc paclet in ordinary notebook?
  • Next by Date: Re: Pathological list for ListLinePlot?
  • Previous by thread: Re: Mathematica 6 obtains imaginary eigenvalues for a Hermitian matrix
  • Next by thread: Re: Mathematica 6 obtains imaginary eigenvalues for a Hermitian matrix