Evaluating large determinants with multivariate polynomial entries

Question

I have some large (n~100) square matrices with entries two variable polynomials of bounded degree (roughly <20, but many entries are smaller) and integer coefficients, and I'd like to be able to compute their determinants (exactly, as two variable polynomials with integer coefficients).

What method should I use? Is there an efficient software implementation already?

In smaller instances of the problem, I've just been using Mathematica's builtin Det function and hoping for the best, but I suspect it's just doing row reduction with two variable rational functions, which quickly becomes inefficient.

I've also tried (thanks Emily!) substituting integer values, calculating the determinant using row reduction with (much faster) arbitrary precision rational arithmetic, and then Newton interpolation. Sadly my implementation seems even slower than Mathematica.

At this point I'm happy with 'good guesses', where if I check a few new values and consistently get the right answers I'm satisfied --- if proving the answer is really correct is much slower that's okay.

Finally, in case this suggests a good trick, there appears to some structure to the answers for my problem instances: the determinant factors, with a large fraction of the total degree coming from high powers of a few low degree polynomials, and then a handful of large factors without multiplicity. I can even predict some of the small factors and (approximately) their multiplicity.

Answer 1

Division-free algorithms for the determinant are the Samuelson-Berkowitz algorithm (see german wikipedia) and the Leverrier-Faddeev algorithm (see german wikipedia). The latter requires division by integers.

There is also a variant of the Gauß algorithm, the Gauß-Bareiss algorithm, that does book-keeping on numerators, i.e., keeping it in a factored form and reducing the fractions if possible. Wikipedia only gives integer matrices as example, but the method applies to integral (euclidean?) rings in general, so the polynomial ring over the integers is covered.

Answer 2

Computing the determinant of matrices of size 100 is an expensive problem to begin with, even with just numbers. I would say it's still within the feasible range, but it's definitely at the limit of what can conveniently be done in a stable way.

It's even more expensive if you want to use symbolic computations for polynomials. The question, of course, is what happens even if you managed to do it: you'd get a polynomial of degree ~2,000 that you will find is impossible to evaluate stably unless it has a very particular structure.

My recommendation is not to wonder how it can be done, but whether you really need to do it and what you want to do with the result. As I say, whatever you get will be rather unreliable and you should try to come up with a method for whatever you want to do that does not involve computing the determinant of these matrices at all.

Answer 3

I've actually come up with a nice trick!

Compute the determinants obtained by evaluating one of the variables (call it $a$) at a small number of large primes $\{p_i\}$. These are rational functions in the other variable (call it $b$). Assume for simplicity they are actually irreducible polynomials $f_i(b)$ (but there's an easy variation of this scheme for the general case).

Now, the constant terms satisfy $f_i(0) = Det(p_i, 0) \pmod{p_i}$, and using the Chinese remainder theorem we can compute $Det(a, 0)$ as a polynomial in $a$, up to an ambiguity in each coefficient of $\prod_i p_i$. Ignore that ambiguity for now, subtract off $Det(a,0)$, divide by $b$, and continue until you're done.

At this point, we have a guess for the determinant as a two-variable polynomial. Use the Schwartz-Zippel lemma, which only requires evaluating the determinant at rational values of the two variables, to efficiently obtain arbitrarily good bounds on the probability you've got the wrong answer.