Issue in Det[...] when computing determinants of polynomials?

When I form a matrix of at least size $12 \times 12$ of second order polynomials in $w$ and I calculate the determinant of that, I get something that is a rational expression in $w$.

Since the computation of a determinant only involves addition, subtraction and multiplication, I'd expect to get a polynomial, and specifically of order $2N$ for a $N \times N$ matrix. This indeed holds as long as $N < 12$.

MWE:

mat = Table[Table[RandomReal[] + RandomReal[] w + RandomReal[] w^2,
{i, 1, 12}], {j, 1, 12}]
det = Det[mat]
Denominator[det]


The final line produces something that is very much not $1$. In addition, the order of the nominator is not 24, as I'd expect.

Am I missing something here, or does Mathematica have a serious bug?

This is similar to this question, but that one didn't receive an answer, and I have a much simpler MWE, which might help in finding the cause of the problem.

• There are different methods to compute determinants. The one following directly from the definition has factorial complexity, which is monstrous. No reasonable software would use it for a 12×12 matrix. You can calculate determinant in something like n³ using methods based on Gaussian elimination — but that involves division as well. The result is the same, of course, but may be obtained in a different equivalent form. Mathematica just didn't bother dividing the two polynomials for you because they are so bulky. Commented Mar 10, 2016 at 10:34
• Definitely NAB. Changed subject header accordingly. Will comment further later. Commented Mar 10, 2016 at 16:21
• For symbolic matrices, at dimension 12, Det switches from a memoization version of cofactor expansion to one-step row reduction. If the input is exact then denominators can be (and I believe are) removed. The approximate coefficient case is another story entirely; "exact" polynomial division will fail due to round-off error. Commented Mar 10, 2016 at 16:31
• @Daniel Your comment above seems to be a clear answer to this question. Would you care to post it as such? I find the comment very insightful, and I would like it preserved. It would seems likely that this question may end up never receiving a proper answer otherwise. Commented Mar 10, 2016 at 17:39
• Okay, I posted it verbatim as an answer. Short on time at the moment but if I get a chance I'll try to say a bit more later. Commented Mar 10, 2016 at 17:41

For symbolic matrices, at dimension 12, Det switches from a memoization version of cofactor expansion to one-step row reduction. If the input is exact then denominators can be (and I believe are) removed. The approximate coefficient case is another story entirely; "exact" polynomial division will fail due to round-off error.

Here is a 12x12 example with a few variables, mostly on the diagonal. It will help to illustrate the problems of using approximate coefficients in this situation (I will not go into much detail as to why this problem exists however)

n = 6;
mat1 = RandomInteger[{-10, 10}, {2 n, 2 n}];
mat2 = DiagonalMatrix[
Riffle[ConstantArray[x, n], ConstantArray[y, n]]];
mat = mat1 + mat2;
mat[[1, -1]] = z;
mat[[-1, 1]] = w;


We can compute the exact determinant with little fuss.

Timing[det = Det[mat]]

(* Out[101]= {0.15625, -34796927854880 - 740164577174 w -
2983120104655 x - 84187702890 w x - 264375018691 x^2 -
11690392910 w x^2 - 14814118187 x^3 - 1053612211 w x^3 +
220708431 x^4 - 37197707 w x^4 + 47954613 x^5 - 257537 w x^5 +
1004649 x^6 - 8706213336912 y - 123404401991 w y -
777789619275 x y - 11892286072 w x y - 71298944066 x^2 y -
1647030923 w x^2 y - 2853939052 x^3 y - 142866207 w x^3 y +
101631272 x^4 y - 5396450 w x^4 y + 12577417 x^5 y -
44451 w x^5 y + 248390 x^6 y - 820562923550 y^2 -
14758208944 w y^2 - 35848751817 x y^2 + 459140990 w x y^2 -
4524395785 x^2 y^2 + 70093624 w x^2 y^2 - 297116521 x^3 y^2 -
4036426 w x^3 y^2 - 4088112 x^4 y^2 - 508478 w x^4 y^2 +
532204 x^5 y^2 - 8357 w x^5 y^2 + 5048 x^6 y^2 - 10914223054 y^3 -
627458161 w y^3 + 1127946475 x y^3 + 26513411 w x y^3 -
464613021 x^2 y^3 + 2721759 w x^2 y^3 - 36290593 x^3 y^3 +
150247 w x^3 y^3 - 544385 x^4 y^3 - 11557 w x^4 y^3 +
51376 x^5 y^3 - 675 w x^5 y^3 + 352 x^6 y^3 + 3207192229 y^4 +
41191335 w y^4 + 90317803 x y^4 - 4328216 w x y^4 -
38243012 x^2 y^4 - 462794 w x^2 y^4 - 1714285 x^3 y^4 +
29272 w x^3 y^4 + 19046 x^4 y^4 + 369 w x^4 y^4 + 5310 x^5 y^4 -
72 w x^5 y^4 + 93 x^6 y^4 + 58114354 y^5 + 929340 w y^5 +
4146890 x y^5 - 33075 w x y^5 + 1129902 x^2 y^5 + 3290 w x^2 y^5 +
42431 x^3 y^5 + 1286 w x^3 y^5 - 5011 x^4 y^5 + 59 w x^4 y^5 -
126 x^5 y^5 + 4 x^6 y^5 + 156402 y^6 + 130756 x y^6 +
29417 x^2 y^6 - 1836 x^3 y^6 - 296 x^4 y^6 + 4 x^5 y^6 + x^6 y^6 +
205578743596 z - 19125819290 w z + 16651877167 x z -
195705766 w x z + 5048041852 x^2 z - 112008142 w x^2 z +
639859284 x^3 z - 2921296 w x^3 z + 29169839 x^4 z +
787279 w x^4 z + 531048 x^5 z + 36942 w x^5 z + 71581353222 y z -
10897406915 w y z + 11524235070 x y z - 945342091 w x y z +
871114255 x^2 y z - 139561682 w x^2 y z + 65672596 x^3 y z -
10059156 w x^3 y z + 3303412 x^4 y z - 122125 w x^4 y z +
119739 x^5 y z + 12648 w x^5 y z + 4220674228 y^2 z -
1193381245 w y^2 z + 1787246988 x y^2 z + 32746354 w x y^2 z +
273674422 x^2 y^2 z + 9454739 w x^2 y^2 z + 26097361 x^3 y^2 z +
462295 w x^3 y^2 z + 1218875 x^4 y^2 z + 23103 w x^4 y^2 z +
14245 x^5 y^2 z + 1580 w x^5 y^2 z - 28146406 y^3 z +
47737560 w y^3 z - 151104780 x y^3 z + 4719985 w x y^3 z -
7418138 x^2 y^3 z - 543584 w x^2 y^3 z + 482439 x^3 y^3 z -
111517 w x^3 y^3 z + 33460 x^4 y^3 z - 7347 w x^4 y^3 z -
424 x^5 y^3 z - 91 w x^5 y^3 z - 36963408 y^4 z + 507638 w y^4 z -
1392072 x y^4 z + 34135 w x y^4 z + 133988 x^2 y^4 z +
9873 w x^2 y^4 z - 853 x^3 y^4 z + 2880 w x^3 y^4 z +
1456 x^4 y^4 z + 91 w x^4 y^4 z + 27 x^5 y^4 z - w x^5 y^4 z +
1037744 y^5 z + 70650 w y^5 z + 72490 x y^5 z - 5416 w x y^5 z -
49162 x^2 y^5 z + 245 w x^2 y^5 z - 4830 x^3 y^5 z +
31 w x^3 y^5 z - 138 x^4 y^5 z - 13 w x^4 y^5 z - w x^5 y^5 z} *)


Things are not so nice when we work with approximate coefficients. The underlying issue is that symbolic polynomial algebra methods no longer really apply, so we only use them in cases where results are "obvious or nearly so" (e.g. dividing out exact powers of a variable).

Timing[detN = Det[N@mat];]

(* Out[95]= {0.65625, Null} *)


In contrast to the exact case, we now have something of a monstrosity in terms of size.

detN // LeafCount

(* Out[97]= 751874 *)


This fellow does not have a nice denominator, even in factored form. Trying to divide it out might make for a good exercise in misery.

Denominator[detN]

(* Out[98]= (0. - 329. w)^8 (40. - 1. w)^9 w^10 (0. + 112000. w^2 -
7440. w^3 + 116. w^4)^7 (0. - 6.98059*10^9 w^4 +
1.01177*10^9 w^5 - 3.74997*10^7 w^6 + 414211. w^7)^6 (0. -
1.63506*10^20 w^8 - 1.43991*10^19 w^9 + 4.01636*10^18 w^10 -
2.60547*10^17 w^11 + 7.68752*10^15 w^12 - 1.09171*10^14 w^13 +
6.07244*10^11 w^14 + 3.83203*10^18 w^8 y - 6.676*10^18 w^9 y +
8.97774*10^17 w^10 y - 4.92499*10^16 w^11 y +
1.34776*10^15 w^12 y - 1.83409*10^13 w^13 y +
9.92674*10^10 w^14 y)^5 *)

• Thanks for posting this! Looking forward to more if you have the time. (+1) Commented Mar 10, 2016 at 17:45
• Is there a deep reason for choosing 12 as the cutoff? Commented Mar 10, 2016 at 17:55
• No deep reason, just a memory heuristic from the olden days. If I were to bump it up I might go to 14 or so. Not that I want to get deep into that code again any time soon. Commented Mar 10, 2016 at 21:25
• Ok, thanks for the explanation. I also found a much better way to solve the original problem where the apparent rational terms were giving me a headache (as they say, it's almost never a good idea to compute the determinant for numerical problems!). I'll be traveling next week without internet access, but I'll accept this answer when I'm back, unless another answer with something essentially more in-depth appears (unlikely, as I take it from the "getting deep into that code" that you work at Wolfram?).
– Timo
Commented Mar 12, 2016 at 8:00