Code snippets

The following code is for informational use only, with no warranties, express or implied. (Any corrections or improvements to these snippets are welcome!)

1. Fox derivatives and the Alexander ideals

Here is some Mathematica code for calculating Fox derivatives, Wirtinger presentations, and Alexander ideals. There are a few parts: w represents words in the free group algebra, FD calculates Fox derivatives of such elements, wirt calculates Wirtinger presentations of oriented PDs (planar diagrams), and alex calculates elementary ideals of a group presentation, assuming there is a homomorphism to an infinite cyclic group such that all the generators map to $t$ (!) (see Abelianizations and Alexander ideals of groups for how to do multivariate Alexander ideals).

ClearAll[w];
(* w is multiplication in the free group algebra, as a multilinear map *)
w[x___, w[y___], z___] := w[x, y, z];
w[x___, a_ + b_, y___] := w[x, a, y] + w[x, b, y];
w[x___, 0, y___] := 0;
w[x___, (coeff : Except[_w]) a_w, y___] := coeff w[x, a, y];
w /: x_w^-1 := #^-1 & /@ Reverse[x];
w /: x_w^n_Integer := w @@ Table[If[n >= 0, x, x^-1], Abs[n]];
w[x___, g_Symbol^n_Integer, z___] /; n != -1 := w[x, w[g]^n, z];
w[x___, a_, a_^-1, y___] := w[x, y];
w[x___, a_^-1, a_, y___] := w[x, y];

ClearAll[FD];
FD::usage = 
  "FD[word_w,g_Symbol] is the Fox derivative of w with respect to g.";
FD[w[g_Symbol], h_] := If[SameQ[g, h], w[], 0];
FD[w[], h_] := 0;
FD[w[g_Symbol^-1], h_] := -w[g^-1, FD[w[g], h]];
FD[w[u_, v___], h_] := FD[w[u], h] + w[u, FD[w[v], h]];

For example, FD[w[x^3, y^-2], x] is w[] + w[x] + w[x, x] and FD[w[x^3, y^-2], y] is -w[x, x, x, 1/y] - w[x, x, x, 1/y, 1/y]. (Mathematica considers y^-1 and 1/y to be exactly the same.)

ClearAll[wirt, PD, Xr, Xl];
PD::usage = "PD[Xr...,Xl...] is an oriented planar diagram.";
Xr::usage = 
  "Xr[a_,b_,c_,d_] d-to-b is the overstrand, a-to-c is the understrand.";
Xl::usage = 
  "Xl[b_,c_,d_,a_] a-to-c is the overstrand, d-to-b is the understrand.";
wirt[p_PD] := Module[{g, rel},
   g[i_] := Symbol["g" <> ToString[i]];
   rel[Xr[a_, b_, c_, d_]] := {
     w[g[d], g[b]^-1], 
     w[g[c], g[b], g[a]^-1, g[d]^-1]};
   rel[Xl[b_, c_, d_, a_]] := {
     w[g[a], g[c]^-1], 
     w[g[c], g[b], g[a]^-1, g[d]^-1]};
   pres[g /@ Union[Flatten[List @@@ List @@ p]],
    Flatten[rel /@ List @@ p]]
   ];
(* a right-handed trefoil knot: *)
k31 = PD[Xr[1, 2, 5, 4], Xr[2, 3, 6, 5], Xr[3, 1, 4, 6]];
(* a figure-eight knot: *)
k41 = PD[Xl[2, 7, 1, 6], Xl[6, 3, 5, 2], Xr[7, 5, 8, 4], 
   Xr[3, 1, 4, 8]];

Figure 1. X is a crossing for unoriented planar diagrams, and Xr and Xl are crossings for oriented planar diagrams.

Figure 2. The edge numbering used in the planar diagram k41.

This definition of the Wirtinger presentation is not particularly efficient with the number of generators, but the approach makes it so we do not need to figure out the arcs.

ClearAll[alex];
alex::usage = 
  "alex[presentation,t] assumes all generators have image t in the abelianization";
alex[pres[gens_, rels_], t_] := 
  Module[{jacobi, abjacobi, ideals},
   jacobi = Outer[FD[#2, #1] &, gens, rels];
   abjacobi = jacobi /. (# -> t & /@ gens) /. val_w :> Times @@ val;
   ideals = 
    Table[GroebnerBasis[Flatten[Minors[abjacobi, k]], t, CoefficientDomain -> Integers],
      {k, 1, Min[Length@gens, Length@rels]}];
   Replace[
    ideals, {{1} ..., good : (Except[{} | {1}] ..), {} ...} :> 
     Reverse[{good}]]
   ];

This function returns the list of elementary ideals, in ascending order of inclusion. The first ideal is generated by a single element, the Alexander polynomial.

For example, alex[wirt[k31], t] is {{1 - t + t^2}}, and alex[wirt[k41], t] is {{1 - 3 t + t^2}}.

As it is, $8_{18}$ is out of reach of this algorithm due to the sizes of the intermediate matrices.

k8i18 = PD[Xl[8, 13, 7, 12], Xr[13, 3, 14, 2], Xr[1, 7, 2, 6], 
   Xl[4, 9, 3, 8], Xl[12, 1, 11, 16], Xr[9, 15, 10, 14], 
   Xr[5, 11, 6, 10], Xl[16, 5, 15, 4]];

But, with some simple Tietze moves:

ClearAll[simplifyPres];
simplifyPres[pres[gens_, {x___, w[], y___}]] := 
  simplifyPres[pres[gens, {x, y}]];
simplifyPres[pres[gens_, {x___, w[g_, h_^-1] | w[h_^-1, g_], y___}]] := 
  simplifyPres[pres[DeleteCases[gens, h], {x, y} /. h -> g]];
simplifyPres[p_pres] := p;

we can get

In[1]:= alex[simplifyPres[wirt[k8i18]], t]

Out[1]= {{1 - 5 t + 10 t^2 - 13 t^3 + 10 t^4 - 5 t^5 + t^6},
         {1 - t + t^2}}

2. The Kauffman bracket polynomial

This Mathematica code is basically from the Knot Atlas. The idea is that the local relation on diagrams for the Kauffman bracket polynomial can be thought of as an expansion in the Temperley-Lieb planar algebra.

The first definition is for calculations in $\mathrm{TL}^{-(A^2+A^{-2})}$ . Loops evaluate to $-(A^2+A^{-2})$ , and paths contract. P[a,b] is an unoriented path between boundary half-edges a and b.

ClearAll[P, A];
SetAttributes[P, Orderless];
P /: P[a_, b_] P[b_, c_] := P[a, c];
P /: P[a_, b_]^2 := P[a, a];
P /: P[a_, a_] := -(A^2 + A^-2);

We need that rule for squaring because of the way Mathematica can simplify expressions.

Then, like in Fox derivatives and the Alexander ideals, we have a definition for PDs of unoriented planar diagrams.

ClearAll[PD, X];
(*a right-handed trefoil*)
k31 = PD[X[1, 2, 5, 4], X[2, 3, 6, 5], X[3, 1, 4, 6]];

The normalized Kauffman bracket can be given by the following expansion:

ClearAll[kb];
kb[diagram_PD] := FullSimplify[Times @@ (diagram /. {
        X[a_, b_, c_, d_] :> A P[a, b] P[c, d] + A^-1 P[a, d] P[b, c]
        })/-(A^2 + A^-2)];

Since every crossing is given two terms, the running time is exponential in the number of crossings. The Knot Atlas gives a heuristic for multiplication order based on minimizing the size of the frontier as one expands out the Kauffman bracket. Using planar separators, one can even find a

$O(\operatorname{poly}(n)2^{C\sqrt{n}})$ time-and-space algorithm (for some

$C>0$ )^[1]

For oriented planar diagrams, in the sense of Fox derivatives and the Alexander ideals, we can obtain the Jones polynomial since the writhe is easy to calculate.

ClearAll[jones, writhe, mirror];
writhe[diagram_PD] := Plus @@ (diagram /. {_Xr -> 1, _Xl -> -1});
mirror[diagram_PD] := diagram /. {
    Xr[a_, b_, c_, d_] :> Xl[b, c, d, a],
    Xl[b_, c_, d_, a_] :> Xr[a, b, c, d]};
jones[diagram_PD] := 
  Collect[FullSimplify[
    kb[diagram /. Xr | Xl -> X]*(-A^-3)^writhe[diagram] /. 
     A -> t^(-1/4)], t];

So, with the right-handed trefoil:

k31 = PD[Xr[1, 2, 5, 4], Xr[2, 3, 6, 5], Xr[3, 1, 4, 6]];

we have

In[1]:= jones[k31]

Out[1]= t + t^3 - t^4

In[2]:= jones[mirror[k31]]

Out[2]= -(1/t^4) + 1/t^3 + 1/t

3. Abelianizations and Alexander ideals of groups

This is Mathematica code for dealing with group words, abelianizations, and Alexander ideals of groups. The first part is group words, which are w terms of symbols raised to powers.

ClearAll[w];
w[x___, 1, y___] := w[x, y];
w[x___, w[y___], z___] := w[x, y, z];
w /: word_w^-1 := #^-1 & /@ Reverse[word];
w /: word_w^n_ := w @@ Table[Sequence @@ If[n >= 0, word, word^-1], Abs[n]];
w[x___, a_, a_, y___] := w[x, a^2, y];
w[x___, a_, a_^n_, y___] := w[x, a^(1 + n), y];
w[x___, a_^n_, a_, y___] := w[x, a^(n + 1), y];
w[x___, a_^n_, a_^m_, y___] := w[x, a^(n + m), y];

We also have group presentations. For example:

s3 = pres[{x, y}, {w[x, x, x], w[y, y], w[x, y, x, y]}];
trefwirt = pres[{a, b}, {w[a, b, a, w[b, a, b]^-1]}];
treftorus = pres[{x, y}, {w[x^2, y^-3]}];
quaternions = pres[{i, j, k}, {w[i, j, k^-1], w[j, k, i^-1], w[k, i, j^-1]}];

Abelianizations are a matter of computing Smith normal form of the linearized presentation matrix.

ClearAll[ab];
ab[pres[gens_, rels_]] := Module[{smod, vec, linRel, u, r, v, ords},
   (*smod is Mod but allows mod by 0*)
   SetAttributes[smod, {Listable}];
   smod[m_, 0] := m;
   smod[m_, n_] := Mod[m, n];

   (*an indicator vector for a generator*)
   vec[g_Symbol] := Replace[gens, {g -> 1, h_ :> 0}, 1];

   (*linearize a relationship*)
   linRel[word_w] := Plus @@ Replace[List @@ word,
      x_Symbol^(n_: 1) :> n vec[x]
      , 1];
   {u, r, v} = SmithDecomposition[linRel /@ rels];

   (*get the diagonal with the correct length*)
   ords = MapThread[Plus, r];

   With[{nontriv = Position[ords, Except[1], 1, Heads -> False]},
    ords = Extract[ords, nontriv];
    abhom[ords, 
     MapThread[#1 -> smod[Extract[#2, nontriv], ords] &, {gens, v}]]
    ]
   ];

The result is a pair of (1) the orders for the direct summands of the abelianization and (2) a homomorphism from the group to the abelianization. For example,

In[1]:= ab[treftorus]

Out[1]= abhom[{0}, {x -> {3}, y -> {2}}]

means that the abelianization of

$\langle x,y\mid x^2y^{-3}\rangle$ is

$\mathbb{Z}$ (since this is

$\mathbb{Z}/0\mathbb{Z}$ ), and the abelianization homomorphism is given by

$x\mapsto 3$ and

$y\mapsto 2$ .

The Alexander ideals of a group can be thought of as the elementary ideals of the presentation matrix for the first homology of the maximal abelian cover of a $K(G,1)$ , with the homology group as a module over the group ring of the abelianization. This can be calculated using Fox derivatives. We deal with torsion by inserting polynomials such as $t^n-1$ into the ideals (as rideal), and so the resulting ideals should be thought of as being in the free Laurent polynomial ring.

ClearAll[alex];
alex[p : pres[gens_, rels_]] := 
  Module[{ords, hom, vars, ghom, rideal, jacobi, abjacobi, fw, FD, ideals},
   
   {ords, hom} = List @@ ab[p];
   vars = Table[Symbol["t" <> ToString[i]], {i, Length@ords}];
   ghom = Replace[hom, {(g_ -> vec_) :> (g -> Times @@ (vars^vec))}, 1];
   rideal = MapThread[#2^#1 - 1 &, {ords, vars}];
   
   (*fw is multiplication in the free group algebra, as a multilinear map*)
   fw[x___, fw[y___], z___] := fw[x, y, z];
   fw[x___, a_ + b_, y___] := fw[x, a, y] + fw[x, b, y];
   fw[x___, 0, y___] := 0;
   fw[x___, (coeff : Except[_fw]) a_fw, y___] := coeff fw[x, a, y];
   fw /: x_fw^-1 := #^-1 & /@ Reverse[x];
   fw /: x_fw^n_Integer := fw @@ Table[If[n >= 0, x, x^-1], Abs[n]];
   fw[x___, g_Symbol^n_Integer, z___] /; n != -1 := fw[x, fw[g]^n, z];

   (*fw is multiplication in the free group algebra,as a multilinear map*)
   fw[x___, fw[y___], z___] := fw[x, y, z];
   fw[x___, a_ + b_, y___] := fw[x, a, y] + fw[x, b, y];
   fw[x___, 0, y___] := 0;
   fw[x___, (coeff : Except[_fw]) a_fw, y___] := coeff fw[x, a, y];
   fw[x___, g_Symbol^n_Integer, z___] /; n != -1 :=
     fw[x, Sequence @@ Table[If[n >= 0, g, g^-1], Abs[n]], z];
   
   FD::usage = "FD[word_w,g_Symbol] is the Fox derivative of w with respect to g.";
   FD[fw[g_Symbol], h_] := If[SameQ[g, h], fw[], 0];
   FD[fw[], h_] := 0;
   FD[fw[g_Symbol^-1], h_] := -fw[g^-1, FD[fw[g], h]];
   FD[fw[u_, v___], h_] := FD[fw[u], h] + fw[u, FD[fw[v], h]];
   
   jacobi = Outer[FD[#2 /. a_w :> fw @@ a, #1] &, gens, rels];
   abjacobi = jacobi /. fw :> Times /. ghom;
   
   ideals = 
    Table[GroebnerBasis[Join[rideal, Flatten[Minors[abjacobi, k]]], 
      vars,
      CoefficientDomain -> Integers],
     {k, 1, Min[Length@gens, Length@rels]}];
   
   {ghom, 
    Replace[ideals, {{1} ..., good : (Except[{} | {1}] ..), {} ...} :>
       Reverse[{good}]]}
   ];

The result is a map from group generators to polynomials (thought as elements of the group ring for the abelianization) and a list of the Alexander ideals in ascending order of inclusion.

In[1]:= alex[s3]
        alex[trefwirt]
        alex[treftorus]
        alex[quaternions]

Out[1]= {{x -> 1, y -> t1},
         {{1 + t1}, {3, 1 + t1}}}

Out[2]= {{a -> t1, b -> t1},
         {{1 - t1 + t1^2}}}

Out[3]= {{x -> t1^3, y -> t1^2},
         {{1 - t1 + t1^2}}}

Out[4]= {{i -> t1 t2, j -> t1, k -> t2},
         {{-1 + t2^2, 1 + t1 + t2 + t1 t2, -1 + t1^2},
          {2, 1 + t2, 1 + t1}}}

A more exciting example is the Borromean rings.

In[1]:= borrom = simplifyPres[
 wirt[PD[Xr[1, 5, 2, 8], Xr[5, 9, 6, 12], Xr[9, 1, 10, 4], 
   Xl[8, 11, 7, 10], Xl[12, 3, 11, 2], Xl[4, 7, 3, 6]]]]

Out[1]= pres[{g2, g4, g6, g8, g10, g12}, {w[g2, g8, 1/g4, 1/g8], 
  w[g6, g12, 1/g8, 1/g12], w[g10, g4, 1/g12, 1/g4], 
  w[g10, g8, 1/g10, 1/g6], w[g2, g12, 1/g2, 1/g10], 
  w[g6, g4, 1/g6, 1/g2]}]

In[2]:= alex[borrom] // Factor

Out[2]= {{g2 -> t2, g4 -> t2, g6 -> t1, g8 -> t1, g10 -> t3, g12 -> t3},
         {{(-1 + t1) (-1 + t2) (-1 + t3)^2,
           (-1 + t1) (-1 + t2)^2 (-1 + t3),
           (-1 + t1)^2 (-1 + t2) (-1 + t3)},
          {(-1 + t2) (-1 + t3), (-1 + t1) (-1 + t3), (-1 + t1) (-1 + t2)}}}

^[1] Lauren Ellenberg, Gabriella Newman, Stephen Sawin, and Jonathan Shi. “Efficient Computation of the Kauffman Bracket.” Journal of Knot Theory and Its Ramifications. 23 (2013). 10.1142/S0218216514500266.