Re: Determinant and Characteristic Polynomial not working properly
- To: mathgroup at smc.vnet.net
- Subject: [mg86421] Re: Determinant and Characteristic Polynomial not working properly
- From: Sebastian Meznaric <meznaric at gmail.com>
- Date: Tue, 11 Mar 2008 02:56:43 -0500 (EST)
- References: <200803091000.FAA12297@smc.vnet.net> <fr2miu$o6s$1@smc.vnet.net>
On 10 Mar, 07:05, Andrzej Kozlowski <a... at mimuw.edu.pl> wrote:
> On 9 Mar 2008, at 11:00, Sebastian Meznaric wrote:
>
>
>
> > Recently I reported a problem with the symbolic evaluation of
> > eigenvalues of a Hermitian matrix. This time I also have a Hermitian
> > matrix, for which Eigenvalues works as expected but the
> > CharacteristicPolynomial and Det do not give the right results. Here
> > goes
>
> > The matrix is
>
> > mat = {{-1, 0, 1/Sqrt[3], 0, 0, 2 Sqrt[2/3], 0, 0, 0, 0, 0, 0, 0,
> > 0}, {0, -1, 0, 1/Sqrt[3], 0, 0, 2 Sqrt[2/3], 0, 0, 0, 0, 0, 0,
> > 0}, {1/Sqrt[3], 0, -5/3, 4/Sqrt[3], 2 Sqrt[2/3], (2 Sqrt[2])/3, 0,
> > 0, 0, 0, 0, 0, 0, 0}, {0, 1/Sqrt[3], 4/Sqrt[3], 1, -(2 Sqrt[2])/3,
> > 0, (2 Sqrt[2])/3, 0, 0, 0, 0, 0, 0, 0}, {0, 0,
> > 2 Sqrt[2/3], -(2 Sqrt[2])/3, 2/3, 0, 0, (5 Sqrt[2])/3, Sqrt[10/3],
> > 0, 0, 0, 0, 0}, {2 Sqrt[2/3], 0, (2 Sqrt[2])/3, 0, 0, -41/6, 4/Sqrt[
> > 3], 2 Sqrt[2/3], 0, Sqrt[15]/2, 0, 0, 0, 0}, {0, 2 Sqrt[2/3], 0, (
> > 2 Sqrt[2])/3, 0, 4/Sqrt[3], -25/6, -(2 Sqrt[2])/3, 0, 0, Sqrt[15]/2,
> > 0, 0, 0}, {0, 0, 0, 0, (5 Sqrt[2])/3,
> > 2 Sqrt[2/3], -(2 Sqrt[2])/3, -1, -Sqrt[5/3]/2, 0, 0, Sqrt[15]/2, 0,
> > 0}, {0, 0, 0, 0, Sqrt[10/3], 0, 0, -Sqrt[5/3]/2, 2, 0, 0, 0, Sqrt[
> > 15]/2, 0}, {0, 0, 0, 0, 0, Sqrt[15]/2, 0, 0, 0, -5/2, 4/Sqrt[3],
> > 2 Sqrt[2/3], 0, 0}, {0, 0, 0, 0, 0, 0, Sqrt[15]/2, 0, 0, 4/Sqrt[3],
> > 1/6, -(2 Sqrt[2])/3, 0, 0}, {0, 0, 0, 0, 0, 0, 0, Sqrt[15]/2, 0,
> > 2 Sqrt[2/3], -(2 Sqrt[2])/3, -2/3, Sqrt[3/5]/2, 3 Sqrt[2/5]}, {0, 0,
> > 0, 0, 0, 0, 0, 0, Sqrt[15]/2, 0, 0, Sqrt[3/5]/2, 9/5, (3 Sqrt[6])/
> > 5}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3 Sqrt[2/5], (3 Sqrt[6])/5, 6/
> > 5}}
>
> > It is Hermitian so there should be no imaginary eigenvalues. Now
> > Factor[CharacteristicPolynomial[ham52, x]]
>
> > 1/69984000(-3439853567779 - 24078974976000 x - 56470929408000 x^2 -
> > 39558316032000 x^3 + 27303837696000 x^4 + 35330162688000 x^5 -
> > 5536014336000 x^6 - 8814624768000 x^7 + 425502720000 x^8 +
> > 954021888000 x^9 + 10077696000 x^10 - 47029248000 x^11 -
> > 2519424000 x^12 + 839808000 x^13 + 69984000 x^14)
>
> > Now we do
> > N[Solve[CharacteristicPolynomial[ham52, x] == 0, x]]
> > (note we are doing Solve first and then N) and we get
> > {{x -> -9.81623}, {x -> -4.34022}, {x -> -2.72532}, {x -> -0.62248}, \
> > {x -> -0.365633}, {x -> 2.}, {x -> 2.96246}, {x ->
> > 4.90719}, {x -> -4.34015 -
> > 0.0000367756 \[ImaginaryI]}, {x -> -4.34015 +
> > 0.0000367756 \[ImaginaryI]}, {x -> -0.622083 -
> > 0.000229308 \[ImaginaryI]}, {x -> -0.622083 +
> > 0.000229308 \[ImaginaryI]}, {x ->
> > 2.96235- 0.0000632158 \[ImaginaryI]}, {x ->
> > 2.96235+ 0.0000632158 \[ImaginaryI]}}
>
> > N[Eigenvalues[mat]]
>
> > {-9.81623, 4.90719, -4.34017, -4.34017, -4.34017, 2.96239, 2.96239, \
> > 2.96239, -2.72532, 2., -0.622216, -0.622216, -0.622216, -0.365633}
>
> > and Eigenvalues[N[mat]]
> > {-9.81623, 4.90719, -4.34017, -4.34017, -4.34017, 2.96239, 2.96239, \
> > 2.96239, -2.72532, 2., -0.622216, -0.622216, -0.622216, -0.365633}
>
> Yes, indeed, there must be a bug in in the way Det of symbolic
> matrices is comuted. One can get the correct charctristic polynomial
> by using significance arithmetic:
>
> f = Factor[Rationalize[CharacteristicPolynomial[N[mat, 10], x], 0]]
> (x - 2)*(x^3 + 2*x^2 - 12*x - 8)^3*(x^4 + 8*x^3 - 32*x^2 - 144*x - 48)
>
> x /. NSolve[f == 0, x]
> {-9.81623254755391, -4.340172973252067,
> -4.340172973252067, -4.340172973252067,
> -2.725323288781692, -0.6222156349319639,
> -0.6222156349319639, -0.6222156349319639,
> -0.3656331685835141, 2., 2.9623886081840314,
> 2.9623886081840314, 2.9623886081840314,
> 4.907189004919114}
>
> Also, as you correclty note, Det gives the wrong answer:
>
> Det[mat]
> -(332576637342150268561733/188956800)
>
> One way to get the right answer is to use the LUDecomposition:
>
> {lu, p, c} = LUDecomposition[mat];
> u = lu SparseArray[{i_, j_} /; j >= i -> 1, Dimensions[mat]];
> Det[u]
> -49152
>
> which agrees with
>
> FullSimplify[Times @@ Eigenvalues[mat]]
> -49152
>
> Andrzej Kozlowski
One thing one has to be careful about when using LUDecomposition is
that the resulting L and U matrices are such that L.U is a permutation
of the rows of the initial matrix. The second element of the result is
the permutation. If this permutation happens to be odd then you need
to multiply your obtained determinant with -1 (in your case it happens
to be even). The Signature function in Mathematica returns the
necessary result. Also you don't need the full triangular matrix but
just the diagonal elements so to save some computation you can also
multiply the matrix with IdentityMatrix. So in general I would use the
following
{lu, p, c} = LUDecomposition[mat];
u = lu * IdentityMatrix[Length[mat]];
Det[u]*Signature[p]
This works great and actually gives an improvement in the evaluation
time of the determinant for all the matrices I've been dealing with.
- References:
- Determinant and Characteristic Polynomial not working properly
- From: Sebastian Meznaric <meznaric@gmail.com>
- Determinant and Characteristic Polynomial not working properly