--- edit ---

More recent versions of the work containing the code mentioned below are found at these links: SymbolicFAQ and PCwGB

--- end edit ---

This will not scale to dimension 100 but will be an improvement on what you now have. It is cribbed from the section "Linear Algebra over Galois Fields here as well as the section "Groebner bases over modules and related computations" in this notebook.

```
deg = 12;
flen = 3;
j = 0;
While[flen > 2 && j++ < 100,
defpoly =
x^deg + 1 + RandomInteger[{0, 1}, deg - 1].x^Range[1, deg - 1];
flen = Length[FactorList[defpoly, Modulus -> 2]]
]
dim = 20;
mat = Table[
RandomInteger[{0, 1}, deg].x^Range[0, deg - 1], {dim}, {dim}];
newvars = Array[z, 2*dim];
augmat = Transpose[
ArrayFlatten[{mat, IdentityMatrix[Length[mat]]}, 1]];
auxpolys = Union[Flatten[Outer[Times, newvars, newvars]]];
linpolys = Join[augmat.newvars, {defpoly}, auxpolys];
allvars = Append[newvars, x];
```

Here is the bulk of the computation.

```
Timing[gb = GroebnerBasis[linpolys, allvars, Modulus -> 2];]
(* {168.240000, Null} *)
```

Now extract the inverse matrix.

```
modulegb = Complement[gb, Join[auxpolys, {defpoly}]];
redmat = Reverse[Sort[Outer[D, modulegb, newvars]]];
invmat = Transpose[redmat[[All, dim + 1 ;; 2*dim]]];
```

We'll check the result.

```
PolynomialMod[invmat.mat, {defpoly, 2}] === IdentityMatrix[dim]
(* True *)
```

A serious bottleneck is the need, using this approach, to have $O(\tt{dim}^2)$ auxiliary polynomials.

One can improve on this by writing a direct row reduction using a list representation form the field elements. This would require working with undocumented functionality in the `Algebra``

context, such as `Algebra`PolynomialTimesModList`

. This TMJ article may give slightly more information should you opt to go that route.

--- edit ---

Here is an approach that will sometimes work. Treat the elements as algebraics over Q instead of GF[2]. Now you can use some built in functionality. If you are lucky you get substantially the same row reduction and can recover the mod-2 result. I'll show what I mean with the example above.

We first set up a matrix of AlgebraicNumber opjects.

```
alg = algnum[root[defpoly], CoefficientList[#1, x]] &;
mat2 = Map[alg, mat, {2}] /. x -> # /. root[a_] :> Root[a &, 1] /.
algnum -> AlgebraicNumber;
augmat = Transpose[
ArrayFlatten[{mat2, IdentityMatrix[Length[mat2]]}, 1]];
```

Now row reduce it and hope we do not get any denominators divisible by 2.

```
Timing[
rred = RowReduce[augmat, Method -> "OneStepRowReduction"];]
(* {36.53499999999994, Null} *)
```

Extract the inverse part.

```
invmat = Transpose[
rred[[All, dim + 1 ;; 2*dim]] /.
AlgebraicNumber[aa_, bb_] :>
PolynomialMod[bb, 2].x^Range[0, deg - 1]];
```

Check:

```
PolynomialMod[invmat.mat, {defpoly, 2}] === IdentityMatrix[dim]
(* True *)
```

--- end edit ---

1Thank you for great references! It seems that I has to reimplement finite fields stuff in Mathematica to improve performance. This is bad news for me :( – Piotr Semenov – 2012-11-25T11:28:31.640

I think the key is

`Algebra`PolynomialExtendedGCDModList`

because it allows you to find field inverses. You invoke it using polynomials represented as their coefficient lists (constant term on the left and ascending from there). If defpoly is the list for your field defining polynomial, then to invert the field element`r1`

you would do`r1inv=Algebra`PolynomialExtendedGCDModList[r1,defpoly,p][[2,1]]`

. In this notation`p`

is the prime charatceristic, so it would be 2 for the cases of apparent interest. – Daniel Lichtblau – 2012-11-26T19:04:29.940@Daniel Lichtblau Thank you for the point! I see, the easiest ways are to implement finite field extensions myself or do the patch for

`FiniteFields`

. If I get some progress, I will share it in the community. – Piotr Semenov – 2012-11-27T10:11:09.803Yes, please do share with community. Here and/or library.wolfram.com would be good places. Is your work only with linear algebra on field elements, or does it involve e.g. polynomials over extension fields? if the latter, there may be some existing functionality of interest I could locate. – Daniel Lichtblau – 2012-11-27T15:35:40.620

@DanielLichtblau My work includes both, the linear algebra and polynomials over extension fields. For example, I have to compute the determinant of big matrices over field extensions and to search roots of polynomials by their factorization. Sure, it will be great! Thanks. – Piotr Semenov – 2012-12-03T12:32:14.653

Okay. Please pose a separate question on working with polynomials over finite fields. That way whatever I (or others) can find won't get buried in a response to a somewhat different query. – Daniel Lichtblau – 2012-12-03T15:52:54.520

@DanielLichtblau Thank you! I have opened a new question (please, see http://mathematica.stackexchange.com/questions/16159/how-to-do-the-polynomial-stuff-over-finite-fields-extensions-fast) especially for you answer :)

– Piotr Semenov – 2012-12-12T09:22:48.113