symbolic expression of eigenvalues for this 3x3 real symmetric matrix

158 views
Skip to first unread message

Junwei Huang

unread,
Jan 21, 2015, 4:29:58 PM1/21/15
to sy...@googlegroups.com
Hello all
Is it too greedy to ask computer to give symbolic expression for this matrix? I can't get the result either from mathematicaTM or sage.
a11,a12,a13,a14,a15,a16,a22,a23,a24,a25,a26,a33,a34,a35,a36,a44,a45,a46,a55,a56,a66 = var('a11,a12,a13,a14,a15,a16,a22,a23,a24,a25,a26,a33,a34,a35,a36,a44,a45,a46,a55,a56,a66',real=True)
p
,q,r=var('p,q,r',real=True)
g11
=a11*p*p+a66*q*q+a55*r*r+2.0*a16*p*q+2.0*a15*p*r+2.0*a56*q*r;
g22
=a66*p*p+a22*q*q+a44*r*r+2.0*a26*p*q+2.0*a46*p*r+2.0*a24*q*r;
g33
=a55*p*p+a44*q*q+a33*r*r+2.0*a45*p*q+2.0*a35*p*r+2.0*a34*q*r;
g12
=a16*p*p+a26*q*q+a45*r*r+(a12+a66)*p*q+(a14+a56)*p*r+(a46+a25)*q*r;
g13
=a15*p*p+a46*q*q+a35*r*r+(a14+a56)*p*q+(a13+a55)*p*r+(a36+a45)*q*r;
g23
=a56*p*p+a24*q*q+a34*r*r+(a46+a25)*p*q+(a36+a45)*p*r+(a23+a44)*q*r;

m
= Matrix([[g11,g12,g13],[g12,g22,g23],[g13,g23,g33]])
m0
= m.subs([(a14,0),(a15,0),(a16,0),(a24,0),(a25,0),(a26,0),(a34,0),(a35,0),(a36,0),(a45,0),(a46,0),(a56,0),(a22,a11),(a23,a13),(a55,a44)])

eigvals0
=m0.eigenvals()
eigvals
=m.eigenvals()

For matrix m0, sage works but sympy fails:
TypeError: cannot determine truth value of
-a11**2*a33*p**2*q**2*r**2 - a11**2*a44*p**4*q**2 - a11**2*a44*p**2*q**4 + 2*a11*a13**2*p**2*q**2*r**2 + 4*a11*a13*a44*p**2*q**2*r**2 - a11*a33*a44*p**2*r**4 - a11*a33*a44*q**2*r**4 - a11*a33*a66*p**4*r**2 - a11*a33*a66*q**4*r**2 - a11*a44**2*p**4*r**2 - a11*a44**2*q**4*r**2 - a11*a44*a66*p**6 - a11*a44*a66*p**4*q**2 - a11*a44*a66*p**2*q**4 - a11*a44*a66*q**6 + a12**2*a33*p**2*q**2*r**2 + a12**2*a44*p**4*q**2 + a12**2*a44*p**2*q**4 - 2*a12*a13**2*p**2*q**2*r**2 - 4*a12*a13*a44*p**2*q**2*r**2 + 2*a12*a33*a66*p**2*q**2*r**2 - 2*a12*a44**2*p**2*q**2*r**2 + 2*a12*a44*a66*p**4*q**2 + 2*a12*a44*a66*p**2*q**4 + a13**2*a44*p**2*r**4 + a13**2*a44*q**2*r**4 + a13**2*a66*p**4*r**2 - 2*a13**2*a66*p**2*q**2*r**2 + a13**2*a66*q**4*r**2 + 2*a13*a44**2*p**2*r**4 + 2*a13*a44**2*q**2*r**4 + 2*a13*a44*a66*p**4*r**2 - 4*a13*a44*a66*p**2*q**2*r**2 + 2*a13*a44*a66*q**4*r**2 - a33*a44**2*r**6 - a33*a44*a66*p**2*r**4 - a33*a44*a66*q**2*r**4 - 4*a44**2*a66*p**2*q**2*r**2 + 2*(-a11*p**2 - a11*q**2 - a33*r**2 - a44*p**2 - a44*q**2 - 2*a44*r**2 - a66*p**2 - a66*q**2)**3/27 - (-a11*p**2 - a11*q**2 - a33*r**2 - a44*p**2 - a44*q**2 - 2*a44*r**2 - a66*p**2 - a66*q**2)*(a11**2*p**2*q**2 + a11*a33*p**2*r**2 + a11*a33*q**2*r**2 + a11*a44*p**4 + 2*a11*a44*p**2*q**2 + a11*a44*p**2*r**2 + a11*a44*q**4 + a11*a44*q**2*r**2 + a11*a66*p**4 + a11*a66*q**4 - a12**2*p**2*q**2 - 2*a12*a66*p**2*q**2 - a13**2*p**2*r**2 - a13**2*q**2*r**2 - 2*a13*a44*p**2*r**2 - 2*a13*a44*q**2*r**2 + 2*a33*a44*r**4 + a33*a66*p**2*r**2 + a33*a66*q**2*r**2 + a44**2*p**2*r**2 + a44**2*q**2*r**2 + a44**2*r**4 + a44*a66*p**4 + 2*a44*a66*p**2*q**2 + a44*a66*p**2*r**2 + a44*a66*q**4 + a44*a66*q**2*r**2)/3 < 0

For matrix m, I guess it just takes a bit longer to fail.

I just want the same eigenvalue associated with the same eigenvector or eigensystem as aij, p,q,r vary.  If I use numerical solution I lose the track of the eigensystem. Any alternatives, or suggestions? Thanks very much.

Chris Smith

unread,
Jan 21, 2015, 5:09:19 PM1/21/15
to sy...@googlegroups.com
Not only make your variables real, give them a sign (positive, for example: var('x', positive=True) is real and positive). This might "still" not be enough to allow the sign to be determined, however...I kind of doubt that it will, but that is a step in the right direction.

Junwei Huang

unread,
Jan 21, 2015, 5:23:58 PM1/21/15
to sy...@googlegroups.com
Thanks Chris. At most I can add aij to be real positive, the computation may still take unreasonable of time. I guess this task is unrealistic and numerical solution could be the only way.

Aaron Meurer

unread,
Jan 21, 2015, 6:27:20 PM1/21/15
to sy...@googlegroups.com
On Wed, Jan 21, 2015 at 3:29 PM, Junwei Huang <jwhua...@gmail.com> wrote:
Hello all
Is it too greedy to ask computer to give symbolic expression for this matrix? I can't get the result either from mathematicaTM or sage.
a11,a12,a13,a14,a15,a16,a22,a23,a24,a25,a26,a33,a34,a35,a36,a44,a45,a46,a55,a56,a66 = var('a11,a12,a13,a14,a15,a16,a22,a23,a24,a25,a26,a33,a34,a35,a36,a44,a45,a46,a55,a56,a66',real=True)
p
,q,r=var('p,q,r',real=True)
g11
=a11*p*p+a66*q*q+a55*r*r+2.0*a16*p*q+2.0*a15*p*r+2.0*a56*q*r;
g22
=a66*p*p+a22*q*q+a44*r*r+2.0*a26*p*q+2.0*a46*p*r+2.0*a24*q*r;
g33
=a55*p*p+a44*q*q+a33*r*r+2.0*a45*p*q+2.0*a35*p*r+2.0*a34*q*r;
g12
=a16*p*p+a26*q*q+a45*r*r+(a12+a66)*p*q+(a14+a56)*p*r+(a46+a25)*q*r;
g13
=a15*p*p+a46*q*q+a35*r*r+(a14+a56)*p*q+(a13+a55)*p*r+(a36+a45)*q*r;
g23
=a56*p*p+a24*q*q+a34*r*r+(a46+a25)*p*q+(a36+a45)*p*r+(a23+a44)*q*r;

m
= Matrix([[g11,g12,g13],[g12,g22,g23],[g13,g23,g33]])
m0
= m.subs([(a14,0),(a15,0),(a16,0),(a24,0),(a25,0),(a26,0),(a34,0),(a35,0),(a36,0),(a45,0),(a46,0),(a56,0),(a22,a11),(a23,a13),(a55,a44)])

eigvals0
=m0.eigenvals()
eigvals
=m.eigenvals()

For matrix m0, sage works but sympy fails:
TypeError: cannot determine truth value of

This error typically indicates a bug in SymPy.
 
-a11**2*a33*p**2*q**2*r**2 - a11**2*a44*p**4*q**2 - a11**2*a44*p**2*q**4 + 2*a11*a13**2*p**2*q**2*r**2 + 4*a11*a13*a44*p**2*q**2*r**2 - a11*a33*a44*p**2*r**4 - a11*a33*a44*q**2*r**4 - a11*a33*a66*p**4*r**2 - a11*a33*a66*q**4*r**2 - a11*a44**2*p**4*r**2 - a11*a44**2*q**4*r**2 - a11*a44*a66*p**6 - a11*a44*a66*p**4*q**2 - a11*a44*a66*p**2*q**4 - a11*a44*a66*q**6 + a12**2*a33*p**2*q**2*r**2 + a12**2*a44*p**4*q**2 + a12**2*a44*p**2*q**4 - 2*a12*a13**2*p**2*q**2*r**2 - 4*a12*a13*a44*p**2*q**2*r**2 + 2*a12*a33*a66*p**2*q**2*r**2 - 2*a12*a44**2*p**2*q**2*r**2 + 2*a12*a44*a66*p**4*q**2 + 2*a12*a44*a66*p**2*q**4 + a13**2*a44*p**2*r**4 + a13**2*a44*q**2*r**4 + a13**2*a66*p**4*r**2 - 2*a13**2*a66*p**2*q**2*r**2 + a13**2*a66*q**4*r**2 + 2*a13*a44**2*p**2*r**4 + 2*a13*a44**2*q**2*r**4 + 2*a13*a44*a66*p**4*r**2 - 4*a13*a44*a66*p**2*q**2*r**2 + 2*a13*a44*a66*q**4*r**2 - a33*a44**2*r**6 - a33*a44*a66*p**2*r**4 - a33*a44*a66*q**2*r**4 - 4*a44**2*a66*p**2*q**2*r**2 + 2*(-a11*p**2 - a11*q**2 - a33*r**2 - a44*p**2 - a44*q**2 - 2*a44*r**2 - a66*p**2 - a66*q**2)**3/27 - (-a11*p**2 - a11*q**2 - a33*r**2 - a44*p**2 - a44*q**2 - 2*a44*r**2 - a66*p**2 - a66*q**2)*(a11**2*p**2*q**2 + a11*a33*p**2*r**2 + a11*a33*q**2*r**2 + a11*a44*p**4 + 2*a11*a44*p**2*q**2 + a11*a44*p**2*r**2 + a11*a44*q**4 + a11*a44*q**2*r**2 + a11*a66*p**4 + a11*a66*q**4 - a12**2*p**2*q**2 - 2*a12*a66*p**2*q**2 - a13**2*p**2*r**2 - a13**2*q**2*r**2 - 2*a13*a44*p**2*r**2 - 2*a13*a44*q**2*r**2 + 2*a33*a44*r**4 + a33*a66*p**2*r**2 + a33*a66*q**2*r**2 + a44**2*p**2*r**2 + a44**2*q**2*r**2 + a44**2*r**4 + a44*a66*p**4 + 2*a44*a66*p**2*q**2 + a44*a66*p**2*r**2 + a44*a66*q**4 + a44*a66*q**2*r**2)/3 < 0

For matrix m, I guess it just takes a bit longer to fail.

I just want the same eigenvalue associated with the same eigenvector or eigensystem as aij, p,q,r vary.  If I use numerical solution I lose the track of the eigensystem. Any alternatives, or suggestions? Thanks very much.

The general eigenvalues of a 3x3 matrix involve the roots of a cubic polynomial, which involve a somewhat nasty equation. It will also be difficult for SymPy to manipulate these eigenvalues into finding the eigenvectors.

It might help to make g11, g12, ... into Symbols, find the eigenvectors of the matrix of symbols, and replace them with their full values as late as possible.

Aaron Meurer
 

--
You received this message because you are subscribed to the Google Groups "sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
To post to this group, send email to sy...@googlegroups.com.
Visit this group at http://groups.google.com/group/sympy.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/aae216a8-482e-46b6-9fac-41de10b1cd38%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Junwei Huang

unread,
Jan 22, 2015, 2:15:29 PM1/22/15
to sy...@googlegroups.com
Thanks Aaron. It is worth of trying before turning to full numerical solutions.
Reply all
Reply to author
Forward
0 new messages