(*
Mathematica program by Lenny Ng
to accompany "Braid loops with infinite monodromy on the Legendrian contact DGA" by Roger Casals and Lenny Ng, arXiv:2101.02318.
This version: 18 January 2021.

Evaluate this cell to load the program. Computations in the paper are given after the general code.

IMPORTANT: Right now this is set up for braids with 4 strands. For a different number of strands, change the "n=4;" line.

Braids are denoted as lists as follows:
{i,j} = crossing a[i] between strands j and j+1;
{{t[1],t[2],t[3],t[4]}}=a vertical line of basepoints t[1],t[2],t[3],t[4].

Primary functions:
Differential[braid] gives the DGA differential for the braid.
AugmentationFromFilling[braid,crossings] gives the DGA augmentation coming from resolving the list of crossings in order. 

For this code, it is important to "over-resolve" so that the end result of the resolutions is an n-component unlink (i.e., the resolved braid is a full twist) rather than an unknot. So: resolve to get an unknot, then add n-1 more resolutions to get the unlink back again. E.g. for affine D4, instead of resolving {9,10,11,12,13}, which would give an unknot, add three more resolutions and use {9,10,11,12,13,14,15,16}, which gives the 4-component unlink.
*)
n=4;

P[var_,l_]:=Table[If[i==l,If[j==l+1,1,0],If[i==l+1,If[j==l,1,If[j==l+1,var,0]],If[j==i,1,0]]],{i,n},{j,n}];
Pr[var_,l_]:=Table[If[i==j,If[i==l,-var^(-1),If[i==l+1,var,1]],0],{i,n},{j,n}];
diagonal[vect_]:=Table[If[i==j,vect[[i]],0],{i,n},{j,n}];
(*
general path matrix for a single crossing/resolved crossing/marked point column
*)
Pgeneral[vect_]:=If[vect[[3]]==2,diagonal[vect[[1]]],
If[vect[[3]]==1,P[a[vect[[1]]],vect[[2]]],Pr[s[vect[[1]]],vect[[2]]]]];

listify[list_]:=Table[If[Length[list[[i]]]==2,Flatten[{list[[i]],1}],{list[[i,1]],0,2}],{i,Length[list]}];
allvars[list_]:=Flatten[{Table[If[list[[i,2]]≠0,{a[list[[i,1]]]},{}],{i,Length[list]}],
Table[t[i],{i,n}]}];

diffl[list_]:=
Module[{prodmat=IdentityMatrix[n]},
For[q=1,q≤Length[list],q++,
prodmat=prodmat.Pgeneral[list[[q]]]];
IdentityMatrix[n]+prodmat];
difflmod[list_]:=
Module[{prodmat=IdentityMatrix[n]},
For[q=1,q≤Length[list],q++,
prodmat=prodmat.
If[list[[q,3]]==1,P[0,list[[q,2]]],Pgeneral[list[[q]]]]];
IdentityMatrix[n]+prodmat];
Differential[list_]:=
diffl[listify[list]];

rightcrossingmat[vect_]:=
Table[If[i==vect[[2]]+1,If[j==vect[[2]],s[vect[[1]]]^(-1),If[j==vect[[2]]+1,1,0]],If[i==j,1,0]],{i,n},{j,n}];
leftcrossingmat[vect_]:=
Table[If[i==vect[[2]],If[j==vect[[2]],1,If[j==vect[[2]]+1,s[vect[[1]]]^(-1),0]],If[i==j,1,0]],{i,n},{j,n}];
(*
the position where "k" appears in a list
*)
whereinlist[list_,k_]:=Flatten[Table[If[
list[[i,3]]==2,{},
If[list[[i,1]]==k,{i},{}]],{i,Length[list]}]][[1]];
replaceentryby[list_,new_,i_]:=Table[If[j==i,new,list[[j]]],{j,Length[list]}];

resolvenew[list0_,listofres_]:=
Module[{list=list0,imageofvars=allvars[list0],shift,qq,rmat,lmat,replfnright,replfnleft},
For[p=1,p≤Length[listofres],p++,
shift=whereinlist[list0,listofres[[p]]];
imageofvars=imageofvars/.{a[list[[shift,1]]]->s[list[[shift,1]]]};
(*
push rightcrossingmat all the way around to the right
*)
rmat=rightcrossingmat[list0[[shift]]];
For[q=1,q≤Length[list0]-1,q++,
qq=If[q≤Length[list0]-shift,q+shift,q+shift-Length[list0]];
If[list[[qq,3]]==1,
replfnright[q]={a[list[[qq,1]]]->a[list[[qq,1]]]-rmat[[list[[qq,2]]+1,list[[qq,2]]]]},
replfnright[q]={}];
If[list[[qq,3]]==1,
rmat=Inverse[P[a[list[[qq,1]]]+rmat[[list[[qq,2]]+1,list[[qq,2]]]],list[[qq,2]]]].rmat.P[a[list[[qq,1]]],list[[qq,2]]],
rmat=Inverse[Pgeneral[list[[qq]]]].rmat.Pgeneral[list[[qq]]]]];
If[rmat≠IdentityMatrix[n],Print["Error: resolution here possibly disallowed."],];
For[q=1,q≤Length[list0]-1,q++,
imageofvars=imageofvars/.replfnright[Length[list0]-q]];
(*
push leftcrossingmat all the way around to the left
*)
lmat=leftcrossingmat[list0[[shift]]];
For[q=1,q≤Length[list0]-1,q++,
qq=If[q≤shift-1,shift-q,shift-q+Length[list0]];
If[list[[qq,3]]==1,
replfnleft[q]={a[list[[qq,1]]]->a[list[[qq,1]]]-lmat[[list[[qq,2]],list[[qq,2]]+1]]},
replfnleft[q]={}];
If[list[[qq,3]]==1,
lmat=P[a[list[[qq,1]]],list[[qq,2]]].lmat.Inverse[P[a[list[[qq,1]]]+lmat[[list[[qq,2]],list[[qq,2]]+1]],list[[qq,2]]]],
lmat=Pgeneral[list[[qq]]].lmat.Inverse[Pgeneral[list[[qq]]]]]];
If[lmat≠IdentityMatrix[n],Print["Error: resolution here possibly disallowed."],];
For[q=1,q≤Length[list0]-1,q++,
imageofvars=imageofvars/.replfnleft[Length[list0]-q]];
list=replaceentryby[list,{list[[shift,1]],list[[shift,2]],0},shift]];
{list,imageofvars}];

fillin[list_]:=
Flatten[{Table[If[list[[i,3]]==1,{a[list[[i,1]]]->0},{}],{i,Length[list]}],
Solve[LogicalExpand[difflmod[list]==Table[0,{i,n},{j,n}]],Table[t[i],{i,n}]]}];

augfromfillauxnew[list_,listofres_]:=
resolvenew[list,listofres][[2]]/.fillin[resolvenew[list,listofres][[1]]];
assignaux[list1_,list2_]:=Table[list1[[i]]->list2[[i]],{i,Length[list1]}];
augfromfillnew[list_,listofres_]:=assignaux[allvars[list],augfromfillauxnew[list,listofres]];

listify[list_]:=Table[If[Length[list[[i]]]==2,Flatten[{list[[i]],1}],{list[[i,1]],0,2}],{i,Length[list]}];
Differential[list_]:=
diffl[listify[list]];
AugmentationFromFilling[list_,listofres_]:=
augfromfillnew[listify[list],listofres];

(*
some auxiliary functions
*)
(*
2x2 matrix for a purple box with a given list of crossings
*)
PurpleBox[list_]:=
Module[{prodmat=IdentityMatrix[2]},
For[q=1,q≤Length[list],q++,
prodmat = prodmat.{{0,1},{1,a[list[[q]]]}}];
prodmat.{{t[1],0},{0,t[2]}}];
(*
these are the filling augmentation composed with small powers of the monodromy map, applied to the pair of crossings "pair" numbered from right to left. IterateUp is if the pair comes from a strand crossing the double strand from bottom left to top right; IterateDown is if the pair comes from a strand crossing the double strand from top left to bottom right. If the answers get progressively more complicated then that's a good sign that the monodromy has infinite order.
*)
IterateUp[pair_,box_,list_,listofres_,numberofiterates_]:=
Simplify[Table[MatrixPower[PurpleBox[box],p,{a[pair[[2]]],a[pair[[1]]]}],{p,0,numberofiterates}]/.AugmentationFromFilling[list,listofres]];
IterateDown[pair_,box_,list_,listofres_,numberofiterates_]:=
Simplify[Table[MatrixPower[Inverse[Transpose[PurpleBox[box]]],p,{a[pair[[1]]],a[pair[[2]]]}],{p,0,numberofiterates}]/.AugmentationFromFilling[list,listofres]];

(*
Affine D4 = β22
*)
(*
The first output is the augmentation computed in the paper; the second output is the verification that this actually is an augmentation (it's the differential on degree 1 generators, composed with the augmentation); the second output is the matrix M1 from the paper, which has all positive coefficients
*)
AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,14,15,16}]
Differential[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}}]/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,14,15,16}]//Simplify
(*
PurpleBox[{19,20}] is the path matrix M in the paper; N1 is the matrix N in the paper (or more precisely N is the image of N1 under the augmentation)
*)
N1={{a[11],1},{a[9],0}};
Inverse[N1].PurpleBox[{19,20}].N1/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,14,15,16}]/.
{s[11]->-s[11],s[12]->-s[12],s[15]->-s[15],s[16]->-s[16]}//Expand

(*
In each of the following, the first output is the images of several relevant variables under the filling augmentation; the second output is the matrix M1, which has all positive coefficients.
*)
(*
β12
*)
{t[1],t[2],a[9],a[11],a[19],a[20]}/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,15,16}]
N1={{a[11],1},{a[9],0}};
Inverse[N1].PurpleBox[{19,20}].N1/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{19,1},{20,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,15,16}]/.{s[11]->-s[11],s[12]->-s[12],s[15]->-s[15],s[16]->-s[16]}//Expand

(*
β21
*)
{t[1],t[2],a[9],a[11],a[19]}/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,14,16}]
N1={{a[11],1},{a[9],0}};
Inverse[N1].PurpleBox[{19}].N1/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{18,3},{19,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,14,16}]/.
{s[11]->-s[11],s[12]->-s[12],s[16]->-s[16]}//Expand

(*
β11
*)
{t[1],t[2],a[9],a[11],a[19]}/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{19,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,16}]
N1={{a[11],1},{a[9],0}};
Inverse[N1].PurpleBox[{19}].N1/.AugmentationFromFilling[{{16,2},{15,1},{14,3},{13,2},{12,2},{11,1},{10,3},{9,2},{8,2},{7,1},{6,3},{5,2},{4,2},{3,1},{2,3},{1,2},{17,3},{19,1},{{t[1],t[2],t[3],t[4]}}},{9,10,11,12,13,16}]/.
{s[11]->-s[11],s[12]->-s[12],s[16]->-s[16]}//Expand

(*
Λ1: note n=3
*)
{t[1],t[2],a[9],a[10],a[13]}/.AugmentationFromFilling[{{12,2},{11,1},{10,1},{9,2},{8,2},{7,1},{6,1},{5,2},{4,2},{3,1},{2,1},{1,2},{13,1},{{t[1],t[2],t[3]}}},
{10,9,8,7,6,5,12}]
N1={{a[10],1},{a[9],0}};
Inverse[N1].PurpleBox[{13}].N1/.AugmentationFromFilling[{{12,2},{11,1},{10,1},{9,2},{8,2},{7,1},{6,1},{5,2},{4,2},{3,1},{2,1},{1,2},{13,1},{{t[1],t[2],t[3]}}},
{10,9,8,7,6,5,12}]/.{s[5]->-s[5],s[7]->-s[7],s[8]->-s[8],s[10]->-s[10]}//Expand

(*
Λ2: note n=3
*)
{t[1],t[2],a[9],a[10],a[13],a[14]}/.AugmentationFromFilling[{{12,2},{11,1},{10,1},{9,2},{8,2},{7,1},{6,1},{5,2},{4,2},{3,1},{2,1},{1,2},{13,1},{14,1},{{t[1],t[2],t[3]}}},
{10,9,8,7,6,5,12,11}]
N1={{a[10],1},{a[9],0}};
Inverse[N1].PurpleBox[{13,14}].N1/.AugmentationFromFilling[{{12,2},{11,1},{10,1},{9,2},{8,2},{7,1},{6,1},{5,2},{4,2},{3,1},{2,1},{1,2},{13,1},{14,1},{{t[1],t[2],t[3]}}},
{10,9,8,7,6,5,12,11}]/.{s[5]->-s[5],s[7]->-s[7],s[8]->-s[8],s[10]->-s[10]}//Expand