Macaulay 2 code for computing Hilbert-Mumford index of Hilbert point of a variety

 

This algorithm was worked out in my paper

 

[HHL] Stability computation via Gröbner basis (with Brendan Hassett and Yongnam Lee)
Journal of Korean Mathematical Society (to appear)

 

(1) The following computes the filtered Hilbert function (Proposition 2 of [HHL]) of a homogeneous ideal I with respect to a given nonnegative weight vector w.

 

 mumfordIndex = (I,w) -> (

  S = ring I;

  r = dim Proj(S/I);

  regI = regularity resolution I;

  MUm = (I,w,m) -> (

    S = ring I;

    N = numgens S;

    K = coefficientRing S;

    Sw = K[gens S, Weights => w, MonomialOrder => GLex];

    W = map(Sw, S, vars Sw);

    I = W(I);

    P = hilbertPolynomial I;

    inI = ideal leadTerm I;

    Sbar = Sw/inI;

    F = map(Sbar, Sw, vars Sbar);

    Bm = basis(m, Sw);

    -- Computes a basis of the degree m piece of Sw.

    Bmbar = basis(m, Sbar);

    -- Computes a basis of the degree m piece of Sbar.

    Bm = flatten entries Bm;

    PSm = #Bm;

    monomialWeight = (f) ->

     (expf = flatten exponents f;

      sum(expf, w, times));

    e = apply(0..(PSm-1),i->(if F(Bm_i)===F(0) then 0 else 1));

    TOTALWT = sum for i from 0 to PSm-1 list

    product{monomialWeight(Bm_i), e_i};

    -- Computes the total weight.

    mu = sum{product{N,-TOTALWT}, product{m, P(m), sum w}}

    -- Computes the Hilbert-Mumford index.

    );

  b = transpose matrix(QQ,table(1,r+2,(i,j)->MUm(I,w,regI + j)));

  A = matrix(QQ,table(r+2,r+2, (i,j) -> (regI + i)^j));

  c = A^(-1)*b;

  QQ[m];

  fHilbFun = sum(r+2, i->c_(i,0)*m^i);

  -- Computes the filtered Hilbert polynomial.

  if regI > 2 then

   val = apply(i = 2..(regI-1), i->MUm(I,w,i))

    else val = ();

  print(regI, val, fHilbFun)

  )