The algorithm described in Section 4.3 of “The Hilbert series of symplectic quotients by the 2-torus” by Hans-Christian Herbig, Daniel Herden, and Christopher Seaton to compute the first Laurent coefficient of the Hilbert series of the real regular functions on the symplectic quotient associated to a representation of the 2-torus

herbighc@gmail.com
daniel herden@baylor.edu
seatonc@rhodes.edu

Definition of the algorithm function gam0

The function gam0 takes a single 2xn matrix (the weight matrix) as an argument, checks that it is effective and in standard form, and if it is, computes the first Laurent coefficient (Subscript[gamma, 0]) of the Hilbert series.

gam0[mat_]:=Module[{a=mat,nn,i,j,d,ds,gamPre,sub,num,den,subs,vars,sings,red,term1,term2,sCount,t,count,test},
n=Length[a[[1]]];
nn=Table[count,{count,1,n}];
(* the (Subscript[d, i,j]) are determinants of the matrix with the ith column perturbed by (Subscript[t, i]) *)
(Subscript[d, i_,j_]):=a[[1,i]]a[[2,j]](Subscript[t, i])(Subscript[t, j])-a[[2,i]]a[[1,j]](Subscript[t, i]) (Subscript[t, j]);
(* the (Subscript[ds, i,j]) are the actual values of the determinants, unperturbed *)
(Subscript[ds, i_,j_]):=a[[1,i]]a[[2,j]]-a[[2,i]]a[[1,j]];

(* check to see if the representation is effective *)

(* make a sequence of the d's that isn't a list for the GCD function *)
dSeq=Flatten[Table[(Subscript[ds, i,j]),{i,nn},{j,nn}]]/.{stuff___}:>(stuff);
gcd=GCD[dSeq];
rank=MatrixRank[a];
If[rank<2||gcd>1,
Print["The matrix is not effective; the rank is ",rank," and the gcd of the 2x2 minors is ",gcd],
(* the remainder of the module is for the effective cases *)

(* check to see if the representation is in standard form *)
count=1;
test=True;
While[count<=n,
If[a[[1,count]]≤0, test=False];
count++];
If[test==False,
Print["The matrix is not in standard form. Every entry of the first row needs to be positive. You can fix this by adding multiples of the second row to the first and/or multiplying columns by -1."],
(* the remainder of the module is for the cases in standard form *)


(* compute the formula for gamma0 *)
gamPre=Together[Sum[Sum[
(* there is a term for each i,j such that (Subscript[d, i,j])>0 *)
If[(Subscript[ds, i,j])>0,
(* check to see if the denominator of the ij term is 0.
 if it is, then use the perturbed values.
if it is not, use the actual values. *)
If[Product[((Subscript[ds, i,j])+(Subscript[ds, i,k])+(Subscript[ds, j,k]))((Subscript[ds, i,j])-(Subscript[ds, i,k])-(Subscript[ds, j,k])),{k,Complement[nn,{i,j}]}]==0,
((Subscript[d, i,j]) )^(2n-5)/Product[
((Subscript[d, i,j])+(Subscript[d, i,k])+(Subscript[d, j,k]))((Subscript[d, i,j])-(Subscript[d, i,k])-(Subscript[d, j,k])),
{k,Complement[nn,{i,j}]}],
((Subscript[ds, i,j]) )^(2n-5)/Product[
((Subscript[ds, i,j])+(Subscript[ds, i,k])+(Subscript[ds, j,k]))((Subscript[ds, i,j])-(Subscript[ds, i,k])-(Subscript[ds, j,k])),
{k,Complement[nn,{i,j}]}]]
,0],
{j,1,n}],
{i,1,n}]];

(* print the matrix *)
Print[MatrixForm[a]];
(* determine the numerator and denominator *)
num=Numerator[gamPre];
den=Factor[Denominator[gamPre]];
(* define the substitutions to set each t = 1 *)
subs=Table[(Subscript[t, ii])->1,{ii,nn}];
(* check if the denominator after setting each t to 1 is 0.
if it is not, set each t = 1 and compute gamma0. *)
If[(den//.subs)!=0,
Print["gamma0 = "];
Print[gamPre//.subs],
(* if the denominator is zero, display that we used the second process. *)
Print["(degeneracies were removed using the second process)"];
(* define sings to be the product of all factors of the denominator that are 0. *)
sings=Product[
If[(den[[i]]//.subs)==0, den[[i]], 1],
{i,1,Length[den]}];
vars=Table[(Subscript[t, ii]),{ii,1,n}];
(* divide sings into the numerator. *)
red=PolynomialReduce[num,{ sings},vars];
(* term 1 is the quotient, which we can now cancel sings from and substitute *)
term1=Cancel[red[[1]][[1]]sings/den]//.subs;
(* term2 is the remainder, which still has a removable singularity, but is much simpler. take the limit to evaluate. *)
term2=Limit[red[[2]]/den,subs];
Print["gamma0 = "];
Print[term1+term2]
]]]
];

Examples

gam0[{{2,2,2},{1,1,1}}]//Timing
NotebookSave[]

gam0[{{-5,1,1,1,1,1,2,1,1},{1,1,1,1,1,1,1,1,1}}]//Timing
NotebookSave[]

gam0[{{1,1,2},{1,0,1}}]//Timing
(* comparing the Laurent expansion of the Hilbert series, computed using other methods *)
Series[((1 - T)^6(1 - T^2)^2)*(1+T^3)/((1-T^3)*(1-T^2)^3*(1-T)^6),{T,1,1}]