Speeding up computation of symbolic determinant in SymPy

Question

I have a 4x4 matrix A with rather long but simple symbolic expressions in each of its entries. About 30 different symbols are involved. By "simple" I mean that these symbols are combined using only addition/subtraction, multiplication/division, and integer powers. By "long" I mean that if I print out the matrix, it covers three or four screens worth.

I need the determinant of this matrix. Or, to be more specific, I know that the determinant is a fourth-order polynomial in one particular symbol, and I need the coefficients of this polynomial. A.det() does not terminate after hours and hours of running, so I need a different approach. Any ideas? So far I've tried to throw various simplify functions at each element of A without any success.

Is there some strategy I can employ to let SymPy be aware of the simple structure of my expressions, or that I know that the result is a polynomial in one of the symbols?

Answer 1

Maybe it would work to create the general expression for a 4x4 determinant

In [30]: A = Matrix(4, 4, symbols('A:4:4'))

In [31]: A
Out[31]:
⎡A₀₀  A₀₁  A₀₂  A₀₃⎤
⎢                  ⎥
⎢A₁₀  A₁₁  A₁₂  A₁₃⎥
⎢                  ⎥
⎢A₂₀  A₂₁  A₂₂  A₂₃⎥
⎢                  ⎥
⎣A₃₀  A₃₁  A₃₂  A₃₃⎦

In [32]: A.det()
Out[32]:
A₀₀⋅A₁₁⋅A₂₂⋅A₃₃ - A₀₀⋅A₁₁⋅A₂₃⋅A₃₂ - A₀₀⋅A₁₂⋅A₂₁⋅A₃₃ + A₀₀⋅A₁₂⋅A₂₃⋅A₃₁ + A₀₀⋅A₁₃⋅A₂₁⋅A₃₂ - A₀₀⋅A₁₃⋅A₂₂⋅A₃₁ - A₀₁⋅A₁₀⋅A₂₂⋅A₃₃ + A₀₁⋅A₁₀⋅A₂₃⋅A₃₂ + A₀₁⋅A₁₂⋅A₂₀⋅
A₃₃ - A₀₁⋅A₁₂⋅A₂₃⋅A₃₀ - A₀₁⋅A₁₃⋅A₂₀⋅A₃₂ + A₀₁⋅A₁₃⋅A₂₂⋅A₃₀ + A₀₂⋅A₁₀⋅A₂₁⋅A₃₃ - A₀₂⋅A₁₀⋅A₂₃⋅A₃₁ - A₀₂⋅A₁₁⋅A₂₀⋅A₃₃ + A₀₂⋅A₁₁⋅A₂₃⋅A₃₀ + A₀₂⋅A₁₃⋅A₂₀⋅A₃₁ - A₀₂⋅A₁
₃⋅A₂₁⋅A₃₀ - A₀₃⋅A₁₀⋅A₂₁⋅A₃₂ + A₀₃⋅A₁₀⋅A₂₂⋅A₃₁ + A₀₃⋅A₁₁⋅A₂₀⋅A₃₂ - A₀₃⋅A₁₁⋅A₂₂⋅A₃₀ - A₀₃⋅A₁₂⋅A₂₀⋅A₃₁ + A₀₃⋅A₁₂⋅A₂₁⋅A₃₀

and then substitute in the entries with something like

A.det().subs(zip(list(A), list(your_matrix)))

SymPy being slow to generate a 4x4 determinant is a bug, though. You should report it at https://github.com/sympy/sympy/issues/new.

EDIT (this wouldn't fit in a comment)

It looks like Matrix.det is calling a simplification function. For matrices 3x3 and smaller, the determinant formula is written out explicitly, but for larger matrices, it is computed using the Bareis algorithm. You can see where the simplification function (cancel) is called here, which is necesssary as part of the computation, but end up doing a lot of work because it tries to simplify your very large expressions. It would probably be smarter to only do the simplifications that are needed to cancel terms of the determinant itself. I opened an issue for this.

Another possibility to speed this up, which I'm not sure will work or not, would be to select a different determinant algorithm. The options are Matrix.det(method=alg) where alg is one of "bareis" (the default), "berkowitz", or "det_LU".