ClearAll[pmax,Qmax,h1,GVTab,PTTabnew,DTTabnew,PTTabGv,h1fun,h1int,h1ser,JTab];
SetDirectory[NotebookDirectory[]];
(* choose your favorite CICY *)
<<X6.m
{CY,chi,kappa,c2}
(* constructing the table of coefficients *)
JTab=CoefficientList[h1int,q];

(* Castelnuovo bound on the D0-charge for PT invariants *)
MMin[Q_]:=Ceiling[-((Q^2)/(2kappa))-(Q/2)];
(* Castelnovo bound on the genus of GV invariants *)
gMax[Q_]:=Floor[((Q^2)/(2 kappa))+(Q/2)+1];

(* spectral flow transformation bringing D2 brane charge to the interval [0,κ/2] *)
SFmu[κ_,μ_]:=Module[{x= Round[μ/κ]},If[κ x > μ,κ x-μ ,μ-κ x]];
SFM[κ_,μ_,M_]:=Module[{x= Round[μ/κ]},If[κ x > μ,M+μ x+(κ/2)x(1-x)-(κ x-μ),M+μ x+(κ/2)x(1-x)]];
(* PT from J invariants using Soheyla's formula for r=1 *)
(* Conditions for application of Soheyla's Theorem *)
fcon[x_]:=If[x<1, x+(1/2),If[x<(15/8),Sqrt[2x+(1/4)],If[x<(9/4),(2/3)x+(3/4),If[x<3,(1/3)x+(3/2),(1/2)x+1
]]]];
Con[beta_,m_]:=((2 beta)/3)fcon[(beta/kappa)]+m;
(* limits of summation *)
ULbeta[beta_,m_]:=beta+((3 kappa m)/(2 beta))+(kappa/2);
ULm[beta_,m_,betap_]:=((((beta-betap))^2)/(2 kappa))+(1/2)(beta+betap)+m;
LLm[betap_]:=-((((betap))^2)/(2 kappa))-(1/2)betap;
ChiJ[beta_,m_]:=- beta-m+(kappa/6)+(c2/12);
(* answer from Soheyla's Thm 2 *)
SJ[mu_,M_]:=Sum[Sum[(((-1))^(ChiJ[mu,M]-mup+mp+1))(ChiJ[mu,M]-mup+mp)PTfromJ[mup,mp] J[SFmu[kappa,mu-mup],SFM[kappa,mu-mup,M-mp+mup]],{mp,Ceiling[LLm[mup]],Floor[ULm[mu,M,mup]]}],{mup,0,IntegerPart[ULbeta[mu,M]]}];

(* rule for substitution of trank 0 DT invariants *)
repJ={J[mu_,n_]:>If[n<MMin[mu],0,If[mu+1<=Length[JTab]&& n+1<=Length[JTab[[mu+1]]]+MMin[mu],JTab[[mu+1,n-MMin[mu]+1]],J[mu,n]]]};
(* evaluation of a PT invariant *)
PTfromJ[mu_,n_]:=If[n<MMin[mu],0,If[mu+1<=Length[PTTabnew]&& n-MMin[mu]+2<=Length[PTTabnew[[mu+1]]],PTTabnew[[mu+1,n-MMin[mu]+2]],If[Con[mu,n]<0 , SJ[mu,n]/.repJ,PT[mu,n]]]];

QMax=12;
MMax=12;
(* qq[x_]:=-((((x/2)+Sqrt[((x^2)/4)+1]))^2);
ZGV[QMax_,x_,y_]:=Sum[((MoebiusMu[k])/k)Log[1+Sum[Sum[PTfromJ[Q,M] ((qq[x])^(k M)) (y^(k Q)),{M,MMin[Q],MMax}],{Q,1,QMax}]],{k,1,QMax}]; *)
ZGV[QMax_,x_,y_]:=Sum[((MoebiusMu[k])/k)Log[1+Sum[Sum[PTfromJ[Q,M] (((-1))^M)(((-x))^(k M)) (y^(k Q)),{M,MMin[Q],MMax}],{Q,1,QMax}]],{k,1,QMax}]
(* y-expand the generating function of PT invariants *)
ZGVexpy=Simplify[Series[ ZGV[QMax,x,y],{y,0,QMax}]];

GVTabnew=Array[0 &, {gMax[QMax]+1,QMax}];
Do[ZGVexp=ZGVexpy[[3,Q]];
(* ZGVexp=CoefficientList[Series[(x^(gMax[Q]-1))ZGVexpy[[3,Q]],{x,0,gMax[Qmax]}],x]; *)
Do[GVTabnew[[g+1,Q]]=Limit[(x^(g-1)) ZGVexp,x->0];
ZGVexp=Series[ZGVexp-GVTabnew[[g+1,Q]] (x^(1-g))(((1+x))^(2(g-1))),{x,0,MMax}],{g,gMax[Q],0,-1}];
GVTabnew[[gMax[QMax]+1,Q]]=ZGVexp,{Q,1,QMax}];
MatrixForm[GVTabnew]

With[{mu=15,M=-29},
{MMin[mu]-M,PTfromJ[mu,M],Con[mu,M]}
]