Pensieve header:  A fresh implementation of baby DoPeGDO. Continues pensieve://2020-09/, pensieve://2020-03/Testing123.nb, and pensieve://People/VanDerVeen/TimidHeisenbergRGeneralForm@.nb.

𝔼[ω,Q,P_ϵSeries] represents ω (e^(Q+P)), where ω is a scalar, Q is an ϵ-free quadratic, and P=(∑_(k=0)^$k)P[[k]](ϵ^k) is a perturbation (it is ill-advised to include ω in P because then it will have log terms).

Scheme: (Subscript[𝔼, _])[_]//(Subscript[𝔼, _])[_] calls FZip or Zip, which are functionally the same. Zip works by handling the quadratic part and calling PZip for the perturbation-only part. PZip works by iteratively solving the synthesis equation. FZip works by encapsulating coefficients, calling Zip, and back-substituting.

Initialization, minor utilities, and “Define” Code

SetDirectory["C:\\drorbn\\AcademicPensieve\\Projects\\BabyDoPeGDO"];
Once[<<KnotTheory`];
Once[Get@"../Profile/Profile.m"];

$k=1;

CCF[ℰ_] := ExpandDenominator@ExpandNumerator@Together[ℰ];
CF[ℰ_List]:=CF/@ℰ;
CF[ℰ_ϵSeries]:=CF/@ℰ;
CF[ℰ_] := (Subscript[PP, CF])@Module[
{vs=Cases[ℰ,(Subscript[(y|x|η|ξ), _]),∞]⋃{y|x|η|ξ}},
Total[CoefficientRules[Expand[ℰ], vs] /. (ps_->c_) :> CCF[c](Times@@(vs^ps))]
];
(*CF[ℰ_] := (Subscript[PP, CF])@CCF[ℰ];*)
CF[ℰ_𝔼]:=CF/@ℰ; 
CF[(Subscript[𝔼, sp___])[ℰs___]]:=CF/@(Subscript[𝔼, sp])[ℰs];

ϵSeries/: S1_ϵSeries≡S2_ϵSeries:=Length[S1]==Length[S2]∧Inner[CF[#1]==CF[#2]&,S1,S2,And];
ϵSeries[0]:=ϵSeries@@Table[0,$k+1];
ϵSeries/: S1_ϵSeries+S2_ϵSeries := ϵSeries@@Table[S1[[k]]+S2[[k]],{k,Min[Length@S1,Length@S2]}];
ϵSeries/: S1_ϵSeries*S2_ϵSeries := ϵSeries@@Table[Sum[S1[[j+1]]*S2[[k-j+1]], {j,0,k}],{k,0,Min[Length@S1,Length@S2]-1}];
ϵSeries /:c_*S_ϵSeries:=(c#)& /@ S;
ϵSeries/:(∂_(vs___))S_ϵSeries:=(s|->(∂_vs)s)/@S;

Define[lhs=rhs, ...] defines the lhs to be rhs, except that rhs is computed only once for each value of $k. Fancy Mathematica not for the faint of heart. Most readers should ignore.

SetAttributes[Define, HoldAll];
Define[def_,defs__] := (Define[def]; Define[defs];);
Define[(Subscript[op_, is__])=ℰ_] := Module[{SD,ii,jj,kk,isp,nis,nisp,sis},Block[{i,j,k},
ReleaseHold[Hold[
SD[(Subscript[op, nisp,$k_Integer]),Block[{i,j,k},(Subscript[op, isp,$k])=ℰ;(Subscript[op, nis,$k])]];
SD[(Subscript[op, isp]),(Subscript[op, {is},$k])]; SD[(Subscript[op, sis__]),(Subscript[op, {sis}])];
]/. {SD->SetDelayed,
isp->{is}/.{i->i_,j->j_,k->k_},
nis->{is}/.{i->ii,j->jj,k->kk},
nisp->{is}/.{i->ii_,j->jj_,k->kk_}
}]  ]]

The Basic Tensors

Define[(Subscript[m, i,j->k])=(Subscript[𝔼, {i,j}->{k}])[1,-(Subscript[ξ, i])(Subscript[η, j])+((Subscript[η, i])+(Subscript[η, j]))(Subscript[y, k])+((Subscript[ξ, i])+(Subscript[ξ, j]))(Subscript[x, k]),ϵSeries[0]]]

Define[
(Subscript[R, i,j])=(Subscript[𝔼, {}->{i,j}])[1,(-1+T) (Subscript[x, j]) ((Subscript[y, i])-(Subscript[y, j])),ϵSeries[0,-(1/2) (1-T) (Subscript[x, j]^2) (Subscript[y, i]^2)+(Subscript[x, i]) (Subscript[x, j]) (Subscript[y, i]) (Subscript[y, j])+(1/2) (1-3 T) (Subscript[x, j]^2) (Subscript[y, i]) (Subscript[y, j])]],
(Subscript[(OverBar[R]), i,j])=(Subscript[𝔼, {}->{i,j}])[1,(-1+(1/T)) (Subscript[x, j]) ((Subscript[y, i])-(Subscript[y, j])),ϵSeries[0,-(((-1+T) (Subscript[x, i]) (Subscript[x, j]) (Subscript[y, i]^2))/(T^2))-(((1-T) (Subscript[x, j]^2) (Subscript[y, i]^2))/(2 (T^3)))-(((Subscript[x, i]) (Subscript[x, j]) (Subscript[y, i]) (Subscript[y, j]))/(T^2))-(((-1-T) (Subscript[x, j]^2) (Subscript[y, i]) (Subscript[y, j]))/(2 (T^3)))]],
(Subscript[CC, i])=(Subscript[𝔼, {}->{i}])[(Sqrt[T]),0,ϵSeries[0,-(((Subscript[x, i]) (Subscript[y, i]))/T)]],
(Subscript[(OverBar[CC]), i])=(Subscript[𝔼, {}->{i}])[(1/(Sqrt[T])),0,ϵSeries[0,(((Subscript[x, i]) (Subscript[y, i]))/T)]]
]

Define[(Subscript[Kink, i])=(Subscript[CC, 3])(Subscript[R, 1,2])//(Subscript[m, 2,3->2])//(Subscript[m, 2,1->i]), (Subscript[(OverBar[Kink]), i])=(Subscript[CC, 3])(Subscript[(OverBar[R]), 1,2])//(Subscript[m, 1,3->1])//(Subscript[m, 1,2->i])]

The Main Program

Variables and their duals:

{(SuperStar[y]),(SuperStar[x]),(SuperStar[η]),(SuperStar[ξ])}={η,ξ,y,x}; 
(SuperStar[(vs_List)]):=(v|->(SuperStar[v]))/@vs;
(SuperStar[((Subscript[u_, i_]))]):=(Subscript[((SuperStar[u])), i]);

𝔼 operations:

𝔼/:𝔼[ω1_,Q1_,P1_]≡𝔼[ω2_,Q2_,P2_]:=CF[ω1==ω2]∧CF[Q1==Q2]∧(P1≡P2);
𝔼/:𝔼[ω1_,Q1_,P1_]𝔼[ω2_,Q2_,P2_]:=𝔼[ω1 ω2,Q1+Q2,P1+P2];
(Subscript[𝔼, d1_->r1_])[ℰ1s___]≡(Subscript[𝔼, d2_->r2_])[ℰ2s___] ^:= (d1==d2)∧(r1==r2)∧(𝔼[ℰ1s]≡𝔼[ℰ2s]);
(Subscript[𝔼, d1_->r1_])[ℰ1s___](Subscript[𝔼, d2_->r2_])[ℰ2s___] ^:= (Subscript[𝔼, (d1⋃d2)->(r1⋃r2)])@@(𝔼[ℰ1s]𝔼[ℰ2s]);
(Subscript[(Subscript[𝔼, dr_])[ℰs___], $k_]):=(Subscript[𝔼, dr])@@(Subscript[𝔼[ℰs], $k]);

(Subscript[𝔼, d1_->r1_])[ℰ1s___]//(Subscript[𝔼, d2_->r2_])[ℰ2s___] :=Module[{is=r1⋂d2,lvs},
lvs=Flatten@Table[{(Subscript[x, $@i]),(Subscript[y, $@i])},{i,is}];
(Subscript[𝔼, (d1⋃Complement[d2,is])->(r2⋃Complement[r1,is])])@@((Subscript[Zip, lvs⋃(SuperStar[lvs])])[(SuperStar[lvs]).lvs,Times[
𝔼[ℰ1s]/.Table[(Subscript[(v:x|y), i])->(Subscript[v, $@i]),{i,is}],
𝔼[ℰ2s]/.Table[(Subscript[(v:ξ|η), i])->(Subscript[v, $@i]),{i,is}]
]])
]

(Subscript[Zip, vs_])[ℱ_,ℰ_]:=⟨ℱ,ℰ⟩//(Subscript[Zip1, vs])//(Subscript[Zip2, vs])//(Subscript[Zip3, vs])

Getting rid of the quadratic.

(Subscript[Zip1, vs_])@⟨ℱ_,𝔼[ω_,Q_,P_]⟩ := (Subscript[PP, Zip1])@Module[{ℐ,F,G,u,v},
ℐ=IdentityMatrix@Length@vs;
F=Table[(∂_(u,v))ℱ,{u,(SuperStar[vs])},{v,(SuperStar[vs])}];
G=Table[(∂_(u,v))Q,{u,vs},{v,vs}];
CF/@⟨(SuperStar[vs]).F.Inverse[ℐ-G.F].(SuperStar[vs])/2, 𝔼[PowerExpand@Factor[ω  ((Det[ℐ- G.F])^(-1/2))],Q-vs.G.vs/2,P]⟩
]

Getting rid of linear terms.

(Subscript[Zip2, vs_])@⟨ℱ_,𝔼[ω_,Q_,P_]⟩ := (Subscript[PP, Zip2])@Module[{F,Y,u,v},
F=Table[(∂_(u,v))ℱ,{u,(SuperStar[vs])},{v,(SuperStar[vs])}];
Y=Table[(∂_v)Q,{v,vs}];
CF/@⟨ℱ, 𝔼[ω,Q-Y.vs+Y.F.Y/2, P/.Thread[vs->vs+F.Y]]⟩
]

Dealing with Feynman diagrams.

Note that the power m of λ is at most k-1+((2k+2)/2)=2k. We write (Subscript[Z, λ])=∑Z[m](λ^m).

(Subscript[Zip3, vs_])@⟨ℱ_,𝔼[ω_,Q_,P_]⟩ := (Subscript[PP, Zip3])@Module[{u,v,m,j},
Z[0]=P;
For[m=0,m<2$k,++m,
Z[m+1]=CF[(1/(2(m+1)))Sum[(∂_((SuperStar[u]),(SuperStar[v])))ℱ((∂_(u,v))Z[m]+Sum[((∂_u)Z[j])((∂_v)Z[m-j]),{j,0,m}]),{u,vs},{v,vs}]]
];
𝔼[ω,Q,CF[Sum[Z[m],{m,0,2$k}]/.Table[v->0,{v,vs}]]]
]

Some Knot Theory

RVK[pd_PD]:=Module[{n,xs,x,rots,front={0},k},
n=Length@pd;rots=Table[0,{2n}];
xs=Cases[pd,x_X:>(Piecewise[{{Xp[x[[4]],x[[1]]], PositiveQ@x}, {Xm[x[[2]],x[[1]]], True}}])];
For[k=0,k<2n,++k,If[k==0∨FreeQ[front,-k],
front=Flatten[front/.k->(xs/.{
Xp[k+1,l_]|Xm[l_,k+1]:>{l,k+1,1-l},
Xp[l_,k+1]|Xm[k+1,l_]:>(++rots[[l]];{1-l,k+1,l})
})],
Cases[front,k|-k]/.{k,-k}:>--rots[[k+1]];
]];
RVK[xs,rots] ];
RVK[K_]:=RVK[PD[K]];

rot[i_,0]:=(Subscript[𝔼, {}->{i}])[1,0,ϵSeries@0];
rot[i_,n_]:=Module[{j},
rot[i,n]=If[n>0,rot[i,n-1](Subscript[CC, j]),rot[i,n+1](Subscript[(OverBar[CC]), j])]//(Subscript[m, i,j->i])];

Z[K_] := Z[RVK@K];
Z[rvk_RVK] := (*Z[rvk] =*)Module[{todo,n,rots,ζ,done,st,cx, ζ1, i,j,k,k1,k2,k3},
{todo,rots}=List@@rvk;
AppendTo[rots,0];
n=Length[todo];
ζ=(Subscript[𝔼, {}->{0}])[1,0,ϵSeries@0];
done={0};
st=Range[0,2n+1];
While[{}=!=($M=todo),
{cx}=MaximalBy[todo, Length[done⋂{#[[1]], #[[2]], #[[1]]-1,#[[2]]-1}]&,1];
{i,j}=List@@cx;
ζ1=Switch[Head[cx],
 Xp, ((Subscript[R, i,j])(Subscript[(OverBar[Kink]), k]))//(Subscript[m, j,k->j]),
 Xm, ((Subscript[(OverBar[R]), i,j])(Subscript[Kink, k]))//(Subscript[m, j,k->j])
];
ζ1=(rot[k,rots[[i]]]ζ1)//(Subscript[m, k,i->i]); rots[[i]]=0;
ζ1=(ζ1 rot[k, rots[[i+1]]])//(Subscript[m, i,k->i]); rots[[i+1]]=0;
ζ1=(rot[k,rots[[j]]]ζ1)//(Subscript[m, k,j->j]); rots[[j]]=0;
ζ1=(ζ1 rot[k, rots[[j+1]]])//(Subscript[m, j,k->j]); rots[[j+1]]=0;
ζ*=ζ1;
If[MemberQ[done, i], ζ=ζ//(Subscript[m, i,i+1->i]); st=st/.st[[i+2]]->st[[i+1]]];
If[MemberQ[done, i-1],  ζ=ζ//(Subscript[m, st[[i]],i->st[[i]]]); st=st/.st[[i+1]]->st[[i]]]; 
If[MemberQ[done, j],  ζ=ζ//(Subscript[m, j,j+1->j]); st=st/.st[[j+2]]->st[[j+1]]]; 
If[MemberQ[done, j-1],  ζ=ζ//(Subscript[m, st[[j]],j->st[[j]]]); st=st/.st[[j+1]]->st[[j]]];
done=done⋃{i-1,i,j-1,j};
todo=DeleteCases[todo, cx]
];
CF/@ (ζ (*/. {(Subscript[x, 0])->x,(Subscript[y, 0])->y,(Subscript[a, 0])->a}*))
]

BeginProfile[];
PopupWindow[Button["Show Profile Monitor"],Dynamic[PrintProfile[], UpdateInterval->3, TrackedSymbols->{}]]

NewBit[K_]:=Module[{Alex=Alexander[K][T]},
(T^3)((Alex^2)/(T-1))Z[K][[3,2]]//Factor]

NewBit/@AllKnots[{3,5}]

(*Two knots with equal Alexander, new bit does not agree*)
Alexander[Knot[6,1]]==Alexander[Knot[9,46]]
Timing[NewBit[Knot[6,1]]==NewBit[Knot[9,46]]]

PrintProfile[]