Revision history [back]
 | 1 | initial version |
You can change the base ring of the matrix to the ring of rational functions (as you did), factor the characteristic polynomial, and set the factors equal to zero symbolically again (and solve).
R = PolynomialRing(QQ, names='a,b,c,d')
M_poly = M.change_ring(R.fraction_field())
S.<x> = R[]
f = M_rat.characteristic_polynomial(x)
f_factors = f.factor()
linear_factors = [(g,e) for (g,e) in f_factors if g.degree() == 1]
rational_roots = [(solve(SR(g)==0,SR(x),solution_dict=True)[0][SR(x)],e) for (g,e) in linear_factors]
nonlinear_factors = [(g,e) for (g,e) in f_factors if g.degree() > 1]
nonrational_roots = sum([map(lambda eqn: (eqn.rhs(), e), solve(SR(g)==0,SR(x))) for (g,e) in nonlinear_factors], [])
eigenvalues = rational_roots + nonrational_roots
This runs in 33 seconds on my machine. The result is
[(3*(a^2 - b^2 + 2*b*c - c^2)*d + 6*a + 1, 1),
(-3*(a^2 - b^2 - 2*a*c + c^2)*d + 6*b + 1, 1),
(-3*(a^2 - 2*a*b + b^2 - c^2)*d + 6*c + 1, 1),
((3*a^2 + 2*a*b + 3*b^2 + 2*(a + b)*c + 3*c^2)*d - 2*a - 2*b - 2*c + 1, 2),
((3*a^2 - 3*b^2 - 2*b*c - 3*c^2)*d - 2*a + 1, 3),
(-(3*a^2 - 3*b^2 + 2*a*c + 3*c^2)*d - 2*b + 1, 3),
(-(3*a^2 + 2*a*b + 3*b^2 - 3*c^2)*d - 2*c + 1, 3),
((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c - 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1,
1),
((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c + 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1,
1)]
Maybe an algorithm like this could be included into Sage in the future.
(The hard part is to detect when this can be applied.)
| | 2 | No.2 Revision |
You can change the base ring of the matrix to the ring of rational functions (as you did), factor the characteristic polynomial, and set the factors equal to zero symbolically again (and solve).
R = PolynomialRing(QQ, names='a,b,c,d')
M_poly M_rat = M.change_ring(R.fraction_field())
S.<x> = R[]
f = M_rat.characteristic_polynomial(x)
f_factors = f.factor()
linear_factors = [(g,e) for (g,e) in f_factors if g.degree() == 1]
rational_roots = [(solve(SR(g)==0,SR(x),solution_dict=True)[0][SR(x)],e) for (g,e) in linear_factors]
nonlinear_factors = [(g,e) for (g,e) in f_factors if g.degree() > 1]
nonrational_roots = sum([map(lambda eqn: (eqn.rhs(), e), solve(SR(g)==0,SR(x))) for (g,e) in nonlinear_factors], [])
eigenvalues = rational_roots + nonrational_roots
This runs in 33 seconds on my machine. The result is
[(3*(a^2 - b^2 + 2*b*c - c^2)*d + 6*a + 1, 1),
(-3*(a^2 - b^2 - 2*a*c + c^2)*d + 6*b + 1, 1),
(-3*(a^2 - 2*a*b + b^2 - c^2)*d + 6*c + 1, 1),
((3*a^2 + 2*a*b + 3*b^2 + 2*(a + b)*c + 3*c^2)*d - 2*a - 2*b - 2*c + 1, 2),
((3*a^2 - 3*b^2 - 2*b*c - 3*c^2)*d - 2*a + 1, 3),
(-(3*a^2 - 3*b^2 + 2*a*c + 3*c^2)*d - 2*b + 1, 3),
(-(3*a^2 + 2*a*b + 3*b^2 - 3*c^2)*d - 2*c + 1, 3),
((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c - 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1,
1),
((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c + 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1,
1)]
Maybe an algorithm like this could be included into Sage in the future.
(The hard part is to detect when this can be applied.)
