Re: Re: Re: Eigenvalue Problem
- To: mathgroup at smc.vnet.net
- Subject: [mg16993] Re: [mg16965] Re: [mg16908] Re: [mg16894] Eigenvalue Problem
- From: Andrzej Kozlowski <andrzej at tuins.ac.jp>
- Date: Sat, 10 Apr 1999 02:13:31 -0400
- Sender: owner-wri-mathgroup at wolfram.com
I must seem mad but I can't resist sending a third reply to this
posting. I got involved in trying to actually prove using Mathematica
that the eigenvalues of m which Mathematica gives are real, without using
the information that the matrix we started from is symmetric. It turns
out that one can do this. Mathematica does not seem to be able to show
this automatically or at least I have not been able to do so. Still,
using one standard package, an auxiliary function, and a little
elementary algebra I have actually been able to prove this. I think the
method may be interesting to people who might want to use Mathematica for
this sort of algebraic computations, and also may perhaps show some ways
in which Mathematica might be improved, so I thought it may be worthwhile
to send it. If somebody can improve this argument I would be very
interested to see it.
(I have removed all the output, to make this message shorter)
Starting from the beginning:
let
In[1]:=
m = {{10 a, 0, b, 0, 0, 0},
{0, -2 a, 0, c, 0, 0},
{b, 0, -8 a, 0,c, 0},
{0, c, 0, -8 a, 0,b},
{0, 0, c, 0, -2 a, 0},
{0, 0, 0,b, 0, 10 a}};
In[2]:=
l=Eigenvalues[m];
We shall only prove that the third eigenvalue is real. All the other
cases can be dealt with in the same way.
In[3]:=
v=l[[3]]
Out[3]=...
We see that the eigenvalue is actually the sum of two expressions, g1 and g2
In[4]:=
g1=v[[1]];g2=v[[2]];
We would like to show that g1 and g2 are always conjugate, when a,b,c are
all real. That will imply that g1+g2 must also be real. Unfortunately
using ComplexExpand (even with various choices for TargetFunctions) does
not work. Mathematica is not able to discover this fact by itself, so we
need to study g1 and g2 directly.
In[5]:=
g1
Out[5]=...
In[6]:=
g2
Out[6]=...
Comparing the two expressions we see that it would be a good idea to
rationalize the denominator in g1, moving all the radicals to the
numerator, so that the expressions would look more similar. Rationalizing
the denominator (getting rid of the radicals in the denominator) is a
useful operation: should there not be a built in function to do this? It
is certainly true that one can do always do this but the algorithm seems
rather complicated(?) In any case, I have written a limited version that
will do what is needed just in this case only. In fact it does not
actually quite rationalize the denominator, just gets rid of the square
root, but I gave it the name which reflects what I want to do:
rationalizeDenominator[p_*Power[a_+(b_:1)*Sqrt[v_],m_?Negative]]:=
Expand[p*(a-b*Sqrt[v])^(-m)]/Expand[(a-b*Sqrt[v])*(a+b*Sqrt[v])]^(-m)//
Simplify;
rationalizeDenominator[Power[a_+(b_:1)*Sqrt[v_],m_?Negative]]:=
Expand[(a-b*Sqrt[v])^(-m)]/Expand[(a-b*Sqrt[v])*(a+b*Sqrt[v])]^(-m)//
Simplify;
(The second definition is not needed, I just added it for completeness.
rationalizeDenominator is supposed to do this sort of thing:
rationalizeDenominator[3/(2 - 3*Sqrt[x - 1])^3]
3 (-46 + 9 Sqrt[-1 + x] + 27 (2 + Sqrt[-1 + x]) x)
-(--------------------------------------------------)
3
(-13 + 9 x)
)
Now:
In[8]:=
rationalizeDenominator[g1]
Out[8]=
Since all the terms in the denominator are positive we can safely apply
PowerExpand, which will now really rationalize the denominator:
In[9]:=
g3=PowerExpand[%]
Out[9]=
We are nearly home. This looks almost like the conjugate of g2. The terms
under the square root sign in the numerator of g3 and in g2 are actually
the same (=h), though for some reason Mathematica arranged them slightly
differently. We can actually check that:
In[11]:=
Expand[h = 2916*a^2*(80*a^2 + b^2 - 5*c^2)^2 + 4*(-252*a^2 - 3*b^2 -
3*c^2)^3] ==
Expand[2916*a^2*(80*a^2 + b^2 - 5*c^2)^2 - 108*(84*a^2 + b^2 + c^2)^3]
Out[11]=
True
The last thing to note is that h is in fact negative, for all real values
of a,b,c. Mathematica can actually show this, if we load the package:
<<Algebra`AlgebraicInequalities`
Now we can check:
In[13]:=
SemialgebraicComponents[h>0,{a,b,c}]
Out[13]=
{}
This actually means that h is always positive, and hence Sqrt[h] pure
imaginary. Wouldn't it be useful if ComplexExpand could call on the
AlgebraicInequalities package to identify those expressions inside
radicals which are always positive or negative?
The numerators of g2 and g3 contain expressions of the form
(x+Sqrt[h])^(1/3) and (x-Sqrt[h])^(1/3). Since Sqrt[h] is pure imaginary
and x=54*a*(80*a^2 + b^2 - 5*c^2) is real we see at once that these
expressions are always conjugate. Comparing the remaining factors of g2
and g3 (=g1) we see at once that they are indeed always conjugate, and
thus g1+g2 is always real.
>>
>>>Hello everybody
>>>
>>>I am trying to solve the eigenvalue problem for the following matrix:
>>>
>>>m = {{10 A, 0, B, 0, 0, 0},
>>> {0, -2 A, 0, C, 0, 0},
>>> {B, 0, -8 A, 0, C, 0},
>>> {0, C, 0, -8 A, 0, B},
>>> {0, 0, C, 0, -2 A, 0},
>>> {0, 0, 0, B, 0, 10 A}}
>>>
>>>which is symmetric. Now mathematica returns some complex eigenvalues
>>>which is not
>>>possible for a real, symmetric matrix. Can anybody help me ? Maybe the
>>>error occurs because
>>>mathematica means that the coefficients are complex but how can I make
>>>them real ?
>>>
>>>Thank's in advance for any help.
>>>
>>>
>>> Peter Huesser
>>>
>>
>>
>>Andrzej Kozlowski
>>Toyama International University
>>JAPAN
>>http://sigma.tuins.ac.jp/
>>http://eri2.tuins.ac.jp/
>>
>