<< KnotTheory`

BeginPackage["KnotTheory`"];

MultivariableAlexander2[PD[Loop[_]]]:={1}&
MultivariableAlexander2[K_]/;Head[K]=!=PD:=MultivariableAlexander2[PD[K]]


MultivariableAlexander2[pd_PD]:=(MultivariableAlexander2[pd]=(l=Length[pd];mat=Table[0, {2*l}, {2*l}];skel=Skeleton[pd];pd1=List@@pd;G=Table[0,{2* l},{l}];pd1//.X[a_, b_, c_, d_]:> If[d==b+1||b-d>1,
{mat[[c,a]] =-t[b]; mat[[c,b]]=t[a]-1;  mat[[c,c]]=1},
{mat[[c,a]]=-1;mat[[c,b]]=1-t[a] ; mat[[c,c]]=t[b]}];c=Times@@pd/.{X[i_,j_,k_,l_]/;(l-j==1||j-l>1):>path[k]path[i] path[j,l],X[i_,j_,k_,l_]/;(j-l==1||l-j>1):>path[k]path[i] path[l,j],P[i_,j_]:>path[i,j]}//.{path[a__,i_]path[i_,b__]:>path[a,i,b],path[a__,i_]path[b__,i_]:>Join[path[a,i],Reverse[path[b]]],path[i_,a__]path[i_,b__]:>Join[Reverse[path[b]],path[i,a]],path[a__,i_]path[i_]:>path[a,i],path[i_,a__]path[i_]:>path[a,i],path[i_]path[i_]:>path[i]};For[i=1,i≤2*l,i++, G=ReplacePart[G,1,{i,First[First[Position[c,i]]]}]];mat=mat/. t[a_]:> t[Position[skel,a][[1,1]]];M=Factor[Simplify[Det[Delete[Transpose[Delete[Transpose[G].mat.G,Position[c,pd1[[1,1]]][[1,1]]]], Position[c,pd1[[1,1]]][[1,1]]]]/( t[ Position[skel,pd1[[1,1]]][[1,1]]]-1)]];emb=Table[Null,{Length[pd]}];done=Table[Null, {2*Length[pd]}];emb[[1]]=0;pd2=pd;rot=Table[0, {Length[skel]}];place[i_, a_] := Module[
{ni, na, arc, dir, oparc},
arc=pd2[[i,a]];
{{ni, na}}=Complement[Position[pd2,arc], {{i,a}}];
If[emb[[ni]]===Null,
emb[[ni]]=3-a+emb[[i]];
pd2[[ni]]=RotateLeft[pd1[[ni]], na-1];
place[ni, #]& /@ {2,3,4},
(* Else *)  oparc=RotateLeft[pd2[[i]], 2][[a]];
If[done[[arc]]===Null,
done[[arc]]=1;
dir=If[arc-oparc==1 || arc-oparc<-1, 1, -1];
rot[[Position[skel, arc][[1,1]]]] += dir*(emb[[ni]]-emb[[i]]+a-na-2)]]];place[1,#]& /@ {1,2,3,4};k=-rot/4;For[j=1,j≤l,j++, k=ReplacePart[k,-1+k[[Position[skel,pd[[j,1]]][[1,1]]]],Position[skel,pd[[j,1]]][[1,1]]]];For[i=1,i≤Length[k],i++, M=t[i]^{(1/2)*k[[i]]}*M];If[pd[[1,4]]==pd[[1,2]]+1||pd[[1,2]]-pd[[1,4]]>1,M=M*t[Position[skel,pd[[1,1]]][[1,1]]]*t[Position[skel,pd[[1,2]]][[1,1]]],
M=M*t[Position[skel,pd[[1,1]]][[1,1]]]];
Evaluate[M/.t->#]&))

MV=MultivariableAlexander; MV2=MultivariableAlexander2

MV2[Link["L5a1"]][t]

MV[Link["L5a1"]][t]

test[L_] :=Head[Expand[Simplify[MV2[L][t][[1]]/MV[L][t]]]]=!=Plus

test1[L_] := (
mv=MV[L][t]; mv2=First[MV2[L][t]];
Or @@ Map[
(mv1=mv /. t[i_] :> t[#[[i]]]; Head[Expand[Simplify[mv2/mv1]]]=!=Plus)&,
Permutations[Range[Length[Skeleton[L]]]]
]
)

test1[Link["L5a1"]]

test1 /@ AllLinks[9]

Timing[MV[PD[#]][t]& /@ AllLinks[10]]

Timing[MV2[PD[#]][t]& /@ AllLinks[10]]

MV3[L_] :=Expand[PowerExpand[MV2[L][t][[1]] * Sqrt[Product[t[i], {i, 1, Length[Skeleton[L]]}]]]]

MV2[AllLinks[7][[5]]][t]

Sqrt[t[1]]/Sqrt[t[1]]

Flip[X[i_,j_,k_,l_]] := If[l==j+1 || j-l>1, X[j,k,l,i], X[l,i,j,k]]

VCube[pd_, l_List] := Module[
{f},
Expand[pd*Times @@ ((1-f[#])& /@ l)] //. pd1_PD*f[i_] :> MapAt[Flip, pd1, i]
]

pd=PD[Knot[3,1]]

Series[VCube[PD[Knot[8,17]], {1,3, 6}] /. pd_PD :> Jones[pd][E^x], {x, 0, 5}]

Series[VCube[PD[#], {1,2, 7}] /. pd_PD :> MV2[pd][t] /. t[i_] -> E^(h x[i]), {h, 0, 3}]& /@AllLinks[8]

Series[VCube[PD[#], {1,2, 7, 9}] /. pd_PD :> MV2[pd][t] /. t[i_] -> E^(h x[i]), {h, 0, 4}]& /@AllLinks[9]

Series[VCube[PD[Link["L11n300"]], {1,2, 7}] /. pd_PD :> MV[pd][t] /. t[i_] -> E^(h x[i]), {h, 0, 3}]