# Computer algebra system for calculation of characteristic polynomial of sparse matrix

I have a $n \times n$ matrix, for which i need to calculate the characteristic polynomial. The matrix is over $GF(2)$, and $n \approx 10^4$. However the matrix is very sparse, with around $n$ non zero entries.

Is there a computer algebra system with sparse matrix structures capable of handling the calculation of the characteristic polynomial? What would be the most promising system to try?

-
What do you need the characteristic polynomial for, if I may ask? Maybe you can avoid computing it explicitly. – Federico Poloni Oct 23 '11 at 19:41
@Frederico Polini: The coefficients of the characteristic polynomial will give the next term of a recursion given the n previous terms. – George Arg Oct 23 '11 at 20:45
I don't think sage or gp/pari can do $10^4$ with less than 8G RAM. sage does $10^3$ in about 9 seconds. Probably you don't want to use the fact that in $GF(2)$ $x^n=x$? – joro Oct 24 '11 at 6:05
@unknown: In other words, you're iterating matrix multiplications? I think that this linear recurrence boils down to it... – Federico Poloni Oct 26 '11 at 8:48

I have been experimenting with Fermat to see how this goes. For a 5000x5000 test case I made up, it computed the minimal polynomial (Wiedemann algorithm) in about 30 seconds, then using that, the characteristic poly in about 22 minutes. It used about 45 meg of RAM. This is on a three year old IMac.

I am the author of Fermat. You can find my email online, and the Fermat web site.

To the Maple person, James:

I tried your code in Maple and the kernel crashed, memory allocation failed. (My IMac has 4 gig RAM.)

Does your example really have only 1 entry per row? The test case I quote above with Fermat has about two per row. With one per row, I get the answer in 0.92 seconds.

-
@rhlewis As shown, I used a randomly generated matrix with density about 1/N, so there were, in fact, some zero rows, as well as some rows with 1 or 2 ones, and one row with 3 ones (based on a handful of random samples I looked at). If I change the density to 2/N, the runtime goes up to about 22 seconds. I just now tried a random permutation matrix as well, and got a runtime of about 5 minutes. That seems to make sense, as there are no zero rows to cut down the work. I used a 64-bit Linux machine with 8 Gb of memory. – James Oct 25 '11 at 5:02

In Maple 15 I'm getting times of about 12 seconds on decent, fairly recent hardware. I used

with( LinearAlgebra ):
N := 10^4:
A := RandomMatrix( N, N, 'generator' = 0 .. 1, 'density' = evalf( 1 / N ) ):
time( Modular:-CharacteristicPolynomial( 2, A, t ) );


for testing. The fact that it is sparse is important. Runtimes for dense examples were much, much longer.

EDIT: It seems that if you use a sparse matrix data structure, and smaller (one byte) integers, this can be done much more efficiently. If you have access to Maple 14 or 15, try this (and note that, here, $N = 10^5$):

with( LinearAlgebra ):
N := 10^5:
A := RandomMatrix( N, N, 'generator' = 0 .. 1, 'density' = evalf( 1 / N ), 'storage' = 'sparse', 'datatype' = 'integer'[1] ):
time( CharacteristicPolynomial( A, t ) mod 2 );


On my machine, I get the answer in less then 0.5 seconds. (It took longer - about 0.75 seconds - to construct the matrix!) Memory used was about 45.3 Mb, but going back to size $10^4\times 10^4$ reduces the memory to 2.7 Mb. I did examples with $N = 10^6$ in about 5 seconds and 220 Mb of memory.

Anyway, it's clear from the various answers that there are computer algebra systems out there that can handle the computations you are interested in (at least Maple and Fermat, and likely others), so you should be in a position to choose whichever system is most convenient for you.

-

It seems that Fermat (see http://home.bway.net/lewis/) specializes in these sorts of computations.

-

I have routinely used sage (google "sage math") to compute the characteristic polynomial of adjacency matrices of graphs on around 800 vertices. I believe sage calls linbox for the actual computation. Note that you can read of the characteristic polynomial from the rational normal form and this might be cheaper to compute.

-

Update.

The versions of Fermat available for download were not able to use the best determinant algorithms for matrices with > 5000 rows. I've changed that now. Google "Lewis Fermat"

I experimented some more with examples like this. I created a random 10000 row square matrix with density 1.86. Using the best algorithm over GF(2), which turns out to be just Gaussian elimination for this example, the determinant (characteristic poly) is computed in 18 seconds using 18 megabytes, on a three-year old iMac. For larger densities, at some point the algorithm I mentioned above will be better.

Correction: Apologies for a mistake in the above: Gaussian elimination is effective only for densities < 1.01. However, the other method works well as far as I tried it, up to density 2.0, taking about 100 minutes and only about 68 meg of RAM.

-