#include <stdlib.h>
#include "kernel/mod2.h"
#include "misc/mylimits.h"
#include "misc/intvec.h"
#include "kernel/combinatorics/hilb.h"
#include "kernel/combinatorics/stairc.h"
#include "kernel/combinatorics/hutil.h"
#include "polys/monomials/ring.h"
#include "polys/monomials/p_polys.h"
#include "polys/simpleideals.h"
#include "polys/weight.h"
#include "polys/flintconv.h"
#include "polys/flint_mpoly.h"
#include "polys/clapconv.h"
#include "Singular/ipid.h"
#include "kernel/ideals.h"
#include "polys/ext_fields/transext.h"
#include "coeffs/coeffs.h"
#include "kernel/linear_algebra/linearAlgebra.h"
#include "coeffs/numbers.h"
#include <vector>
#include "Singular/ipshell.h"
#include <ctime>
#include <iostream>

Macros
#define	OVERFLOW_MAX LONG_MAX

#define	OVERFLOW_MIN LONG_MIN

#define	omsai 1

Functions
intvec *	hSecondSeries (intvec *hseries1)

void	hDegreeSeries (intvec s1, intvec s2, int co, int mu)

poly	hFirst2Second (poly h, const ring Qt, int &co)

static void	hPrintHilb (poly hseries, const ring Qt, intvec *modul_weight)

static ring	makeQt ()

void	hLookSeries (ideal S, intvec modulweight, ideal Q, intvec wdegree)

static void	idInsertMonomial (ideal I, poly p)

static int	comapreMonoIdBases (ideal J, ideal Ob)

static int	CountOnIdUptoTruncationIndex (ideal I, int tr)

static int	comapreMonoIdBases_IG_Case (ideal J, int JCount, ideal Ob, int ObCount)

static int	positionInOrbit_IG_Case (ideal I, poly w, std::vector< ideal > idorb, std::vector< poly > polist, int trInd, int)

static int	positionInOrbit_FG_Case (ideal I, poly, std::vector< ideal > idorb, std::vector< poly >, int, int)

static int	positionInOrbitTruncationCase (ideal I, poly w, std::vector< ideal > idorb, std::vector< poly > polist, int, int trunDegHs)

static int	monCompare (const void m, const void n)

static void	sortMonoIdeal_pCompare (ideal I)

static ideal	minimalMonomialGenSet (ideal I)

static poly	shiftInMon (poly p, int i, int lV, const ring r)

static poly	deleteInMon (poly w, int i, int lV, const ring r)

static void	TwordMap (poly p, poly w, int lV, int d, ideal Jwi, bool &flag)

static ideal	colonIdeal (ideal S, poly w, int lV, ideal Jwi, int trunDegHs)

void	HilbertSeries_OrbitData (ideal S, int lV, bool IG_CASE, bool mgrad, bool odp, int trunDegHs)

ideal	RightColonOperation (ideal S, poly w, int lV)

static BOOLEAN	p_Div_hi (poly p, const int *exp_q, const ring src)

static int	compare_rp_currRing (const void pp1, const void pp2)

static void	id_DelDiv_hi (ideal id, BOOLEAN *bad, const ring r)

static poly	hilbert_series (ideal A, const ring src, const intvec *wdegree, const ring Qt)

poly	hFirstSeries0p (ideal A, ideal Q, intvec *wdegree, const ring src, const ring Qt)

poly	hFirstSeries0m (ideal A, ideal Q, intvec wdegree, intvec shifts, const ring src, const ring Qt)

intvec *	hFirstSeries0 (ideal A, ideal Q, intvec *wdegree, const ring src, const ring Qt)

static ideal	getModuleComp (ideal A, int c, const ring src)

intvec *	hFirstSeries (ideal A, intvec module_w, ideal Q, intvec wdegree)

static int	hMinModulweight (intvec *modulweight)

static void	hWDegree (intvec *wdegree)

static int64 *	hAddHilb (int Nv, int x, int64 pol, int lp)

static void	hLastHilb (scmon pure, int Nv, varset var, int64 *pol, int lp)

static void	hHilbEst (scfmon stc, int Nstc, varset var, int Nvar)

static void	hHilbStep (scmon pure, scfmon stc, int Nstc, varset var, int Nvar, int64 *pol, int Lpol)

static intvec *	hSeries (ideal S, intvec modulweight, intvec wdegree, ideal Q)

intvec *	hFirstSeries1 (ideal S, intvec modulweight, ideal Q, intvec wdegree)

bigintmat *	hPoly2BIV (poly h, const ring Qt, const coeffs biv_cf)

poly	hBIV2Poly (bigintmat *b, const ring Qt, const coeffs biv_cf)

bigintmat *	hFirstSeries0b (ideal I, ideal Q, intvec wdegree, intvec shifts, const ring src, const coeffs biv_cf)

bigintmat *	hSecondSeries0b (ideal I, ideal Q, intvec wdegree, intvec shifts, const ring src, const coeffs biv_cf)

void	scDegree (ideal S, intvec *modulweight, ideal Q)

Variables
STATIC_VAR int64 **	Qpol

STATIC_VAR int64 *	Q0

STATIC_VAR int64 *	Ql

STATIC_VAR int	hLength

STATIC_VAR ring	hilb_Qt =NULL

Macro Definition Documentation

◆ omsai

#define omsai 1

Definition at line 51 of file hilb.cc.

◆ OVERFLOW_MAX

#define OVERFLOW_MAX LONG_MAX

Definition at line 38 of file hilb.cc.

◆ OVERFLOW_MIN

#define OVERFLOW_MIN LONG_MIN

Definition at line 39 of file hilb.cc.

Function Documentation

◆ colonIdeal()

static ideal colonIdeal	(	ideal	S,
		poly	w,
		int	lV,
		ideal	Jwi,
		int	trunDegHs
	)

static

Definition at line 735 of file hilb.cc.

{
  /*
   * It computes the right colon ideal of a two-sided ideal S
   * w.r.t. word w and save it in a new object Jwi.
   * It keeps S and w unchanged.
   */
 
  if(idIs0(S))
  {
    return(S);
  }
 
  int i, d;
  d = p_Totaldegree(w, currRing);
  if(trunDegHs !=0 && d >= trunDegHs)
  {
   idInsertMonomial(Jwi, p_One(currRing));
   return(Jwi);
  }
  bool flag = FALSE;
  int SCount = IDELEMS(S);
  for(i = 0; i < SCount; i++)
  {
    TwordMap(S->m[i], w, lV, d, Jwi, flag);
    if(flag)
    {
      break;
    }
  }
 
  Jwi = minimalMonomialGenSet(Jwi);
  return(Jwi);
}

◆ comapreMonoIdBases()

static int comapreMonoIdBases	(	ideal	J,
		ideal	Ob
	)

static

Definition at line 279 of file hilb.cc.

{
  /*
   * Monomials of J and Ob are assumed to
   * be already sorted. J and Ob are
   * represented by the minimal generating set.
   */
  int i, s;
  s = 1;
  int JCount = IDELEMS(J);
  int ObCount = IDELEMS(Ob);
 
  if(idIs0(J))
  {
    return(1);
  }
  if(JCount != ObCount)
  {
    return(0);
  }
 
  for(i = 0; i < JCount; i++)
  {
    if(!(p_LmEqual(J->m[i], Ob->m[i], currRing)))
    {
      return(0);
    }
  }
  return(s);
}

◆ comapreMonoIdBases_IG_Case()

static int comapreMonoIdBases_IG_Case	(	ideal	J,
		int	JCount,
		ideal	Ob,
		int	ObCount
	)

static

Definition at line 337 of file hilb.cc.

{
  /*
   * Monomials of J and Ob are assumed to
   * be already sorted in increasing degrees.
   * J and Ob are represented by the minimal
   * generating set.  It checks if J and Ob have
   * same monomials  up to deg <=tr.
   */
 
  int i, s;
  s = 1;
  //when J is null
  //
  if(JCount != ObCount)
  {
    return(0);
  }
 
  if(JCount == 0)
  {
    return(1);
  }
 
  for(i = 0; i< JCount; i++)
  {
    if(!(p_LmEqual(J->m[i], Ob->m[i], currRing)))
    {
      return(0);
    }
  }
 
  return(s);
}

◆ compare_rp_currRing()

static int compare_rp_currRing	(	const void *	pp1,
		const void *	pp2
	)

static

Definition at line 1191 of file hilb.cc.

{
  poly p1=*(poly*)pp1;
  poly p2=*(poly*)pp2;
  for(int i=currRing->N;i>0;i--)
  {
    int e1=p_GetExp(p1,i,currRing);
    int e2=p_GetExp(p2,i,currRing);
    if(e1<e2) return -1;
    if(e1>e2) return 1;
  }
  return 0;
}

◆ CountOnIdUptoTruncationIndex()

static int CountOnIdUptoTruncationIndex	(	ideal	I,
		int	tr
	)

static

Definition at line 310 of file hilb.cc.

{
  /*
   * The ideal I must be sorted in increasing total degree.
   * It counts the number of monomials in I up to
   * degree less than or equal to tr.
   */
 
  //case when I=1;
  if(p_Totaldegree(I->m[0], currRing) == 0)
  {
    return(1);
  }
 
  int count = 0;
  for(int i = 0; i < IDELEMS(I); i++)
  {
    if(p_Totaldegree(I->m[i], currRing) > tr)
    {
      return (count);
    }
    count = count + 1;
  }
 
  return(count);
}

◆ deleteInMon()

static poly deleteInMon	(	poly	w,
		int	i,
		int	lV,
		const ring	r
	)

static

Definition at line 640 of file hilb.cc.

{
  /*
   * deletes the variables up to i^th layer of monomial w
   * w remains unchanged
   * creates new poly and returns it for the colon ideal
   */
 
  poly dw = p_One(currRing);
  int *e = (int *)omAlloc((r->N+1)*sizeof(int));
  int *s=(int *)omAlloc0((r->N+1)*sizeof(int));
  p_GetExpV(w, e, r);
  int j, cnt;
  cnt = i*lV;
  /*
  for(j=1;j<=cnt;j++)
  {
    e[j]=0;
  }*/
  for(j = (cnt+1); j < (r->N+1); j++)
  {
    s[j] = e[j];
  }
 
  p_SetExpV(dw, s, currRing);//new exponents
  omFree(e);
  omFree(s);
 
  p_SetComp(dw, p_GetComp(w, currRing), currRing);
  p_Setm(dw, currRing);
 
  return(dw);
}

◆ getModuleComp()

static ideal getModuleComp	(	ideal	A,
		int	c,
		const ring	src
	)

static

Definition at line 1503 of file hilb.cc.

{
  ideal res=idInit(IDELEMS(A),A->rank);
  for (int i=0;i<IDELEMS(A);i++)
  {
    if ((A->m[i]!=NULL) && (p_GetComp(A->m[i],src)==c))
      res->m[i]=p_Head(A->m[i],src);
  }
  return res;
}

◆ hAddHilb()

static int64 * hAddHilb	(	int	Nv,
		int	x,
		int64 *	pol,
		int *	lp
	)

static

Definition at line 1618 of file hilb.cc.

{
  int  l = *lp, ln, i;
  int64  *pon;
  *lp = ln = l + x;
  pon = Qpol[Nv];
  memcpy(pon, pol, l * sizeof(int64));
  if (l > x)
  {/*pon[i] -= pol[i - x];*/
    for (i = x; i < l; i++)
    {
      #ifndef __SIZEOF_INT128__
      int64 t=pon[i];
      int64 t2=pol[i - x];
      t-=t2;
      if ((t>=OVERFLOW_MIN)&&(t<=OVERFLOW_MAX)) pon[i]=t;
      else if (!errorreported) WerrorS("int overflow in hilb 1");
      #else
      __int128 t=pon[i];
      __int128 t2=pol[i - x];
      t-=t2;
      if ((t>=LONG_MIN)&&(t<=LONG_MAX)) pon[i]=t;
      else if (!errorreported) WerrorS("long int overflow in hilb 1");
      #endif
    }
    for (i = l; i < ln; i++)
    { /*pon[i] = -pol[i - x];*/
      #ifndef __SIZEOF_INT128__
      int64 t= -pol[i - x];
      if ((t>=OVERFLOW_MIN)&&(t<=OVERFLOW_MAX)) pon[i]=t;
      else if (!errorreported) WerrorS("int overflow in hilb 2");
      #else
      __int128 t= -pol[i - x];
      if ((t>=LONG_MIN)&&(t<=LONG_MAX)) pon[i]=t;
      else if (!errorreported) WerrorS("long int overflow in hilb 2");
      #endif
    }
  }
  else
  {
    for (i = l; i < x; i++)
      pon[i] = 0;
    for (i = x; i < ln; i++)
      pon[i] = -pol[i - x];
  }
  return pon;
}

◆ hBIV2Poly()

poly hBIV2Poly	(	bigintmat *	b,
		const ring	Qt,
		const coeffs	biv_cf
	)

Definition at line 2005 of file hilb.cc.

{
  poly p=NULL;
  nMapFunc f=n_SetMap(biv_cf,Qt->cf);
  for(int d=0;d<b->rows()-1;d++)
  {
    poly h=p_New(Qt);
    p_SetExp(h,1,d,Qt);p_Setm(h,Qt);
    pSetCoeff0(h,f(BIMATELEM(*b,1,d+1),biv_cf,Qt->cf));
    p=p_Add_q(p,h,Qt);
  }
  return p;
}

◆ hDegreeSeries()

void hDegreeSeries	(	intvec *	s1,
		intvec *	s2,
		int *	co,
		int *	mu
	)

Definition at line 106 of file hilb.cc.

{
  int i, j, k;
  int m;
  *co = *mu = 0;
  if ((s1 == NULL) || (s2 == NULL))
    return;
  i = s1->length();
  j = s2->length();
  if (j > i)
    return;
  m = 0;
  for(k=j-2; k>=0; k--)
    m += (*s2)[k];
  *mu = m;
  *co = i - j;
}

◆ hFirst2Second()

poly hFirst2Second	(	poly	h,
		const ring	Qt,
		int &	co
	)

Definition at line 124 of file hilb.cc.

{
  poly o_t=p_One(Qt);p_SetExp(o_t,1,1,Qt);p_Setm(o_t,Qt);
  o_t=p_Neg(o_t,Qt);
  o_t=p_Add_q(p_One(Qt),o_t,Qt);
  poly di1=p_Copy(h,Qt);
  co=0;
#if defined(HAVE_FLINT) && (__FLINT_RELEASE >= 20503)
  poly di2;
  fmpq_mpoly_ctx_t ctx;
  convSingRFlintR(ctx,Qt);
  loop
  {
    di2=Flint_Divide_MP(di1,0,o_t,0,ctx,Qt);
    if (di2==NULL) break;
    co++;
    p_Delete(&di1,Qt);
    di1=di2;
  }
#else
  if (di1!=NULL)
  {
    CanonicalForm  Di1=convSingPFactoryP(di1,Qt);
    CanonicalForm  O_t=convSingPFactoryP(o_t,Qt);
    CanonicalForm Di2,dummy;
    loop
    {
      Di2=Di1/O_t;
      dummy=Di2*O_t;
      if (dummy!=Di1) break;
      else Di1=Di2;
      co++;
    }
    p_Delete(&di1,Qt);
    di1=convFactoryPSingP(Di1,Qt);
  }
#endif
  p_Delete(&o_t,Qt);
  return di1;
}

◆ hFirstSeries()

intvec * hFirstSeries	(	ideal	A,
		intvec *	module_w,
		ideal	Q,
		intvec *	wdegree
	)

Definition at line 1514 of file hilb.cc.

{
  intvec* res;
  #if 0
  // find degree bound
  int a,b,prod;
  a=rVar(currRing);
  b=1;
  prod=a;
  while(prod<(1<<15) && (a>1))
  {
    a--;b++;
    prod*=a;
    prod/=b;
  }
  if (a==1) b=(1<<15);
  // check degree bound
  BOOLEAN large_deg=FALSE;
  int max=0;
  for(int i=IDELEMS(A)-1;i>=0;i--)
  {
    if (A->m[i]!=NULL)
    {
      int mm=p_Totaldegree(A->m[i],currRing);
      if (mm>max)
      {
        max=mm;
        if (max>=b)
        {
          large_deg=TRUE;
          break;
        }
      }
    }
  }
  if (!large_deg)
  {
    void (*WerrorS_save)(const char *s) = WerrorS_callback;
    WerrorS_callback=WerrorS_dummy;
    res=hFirstSeries1(A,module_w,Q,wdegree);
    WerrorS_callback=WerrorS_save;
    if (errorreported==0)
    {
      return res;
    }
    else errorreported=0;// retry with other alg.:
  }
  #endif
 
  if (hilb_Qt==NULL) hilb_Qt=makeQt();
  if (!id_IsModule(A,currRing))
    return hFirstSeries0(A,Q,wdegree,currRing,hilb_Qt);
  res=NULL;
  int w_max=0,w_min=0;
  if (module_w!=NULL)
  {
    w_max=module_w->max_in();
    w_min=module_w->min_in();
  }
  for(int c=1;c<=A->rank;c++)
  {
    ideal Ac=getModuleComp(A,c,currRing);
    intvec *res_c=hFirstSeries0(Ac,Q,wdegree,currRing,hilb_Qt);
    id_Delete(&Ac,currRing);
    intvec *tmp=NULL;
    if (res==NULL)
      res=new intvec(res_c->length()+(w_max-w_min));
    if ((module_w==NULL) || ((*module_w)[c-1]==0)) tmp=ivAdd(res,res_c);
    else tmp=ivAddShift(res, res_c,(*module_w)[c-1]-w_min);
    delete res_c;
    if (tmp!=NULL)
    {
      delete res;
      res=tmp;
    }
  }
  (*res)[res->length()-1]=w_min;
  return res;
}

◆ hFirstSeries0()

intvec * hFirstSeries0	(	ideal	A,
		ideal	Q,
		intvec *	wdegree,
		const ring	src,
		const ring	Qt
	)

Definition at line 1478 of file hilb.cc.

{
  poly s=hFirstSeries0p(A,Q,wdegree,src,Qt);
  intvec *ss;
  if (s==NULL)
    ss=new intvec(2);
  else
  {
    ss=new intvec(p_Totaldegree(s,Qt)+2);
    while(s!=NULL)
    {
      int i=p_Totaldegree(s,Qt);
      long l=n_Int(pGetCoeff(s),Qt->cf);
      (*ss)[i]=n_Int(pGetCoeff(s),Qt->cf);
      if((l==0)||(l<=-INT_MAX)||(l>INT_MAX))
      {
        if(!errorreported) Werror("overflow in hilb at t^%d\n",i);
      }
      else (*ss)[i]=(int)l;
      p_LmDelete(&s,Qt);
    }
  }
  return ss;
}

◆ hFirstSeries0b()

bigintmat * hFirstSeries0b	(	ideal	I,
		ideal	Q,
		intvec *	wdegree,
		intvec *	shifts,
		const ring	src,
		const coeffs	biv_cf
	)

Definition at line 2019 of file hilb.cc.

{
  if (hilb_Qt==NULL) hilb_Qt=makeQt();
  poly h;
  int m=0;
  if (id_IsModule(I,src))
  {
    h=hFirstSeries0m(I,Q,wdegree,shifts,src,hilb_Qt);
    if (shifts!=NULL) m=shifts->min_in();
  }
  else
    h=hFirstSeries0p(I,Q,wdegree,src,hilb_Qt);
  bigintmat *biv=hPoly2BIV(h,hilb_Qt,biv_cf);
  if (m!=0)
  {
    n_Delete(&BIMATELEM(*biv,1,biv->cols()),biv_cf);
    BIMATELEM(*biv,1,biv->cols())=n_Init(m,biv_cf);
  }
  p_Delete(&h,hilb_Qt);
  return biv;
}

◆ hFirstSeries0m()

poly hFirstSeries0m	(	ideal	A,
		ideal	Q,
		intvec *	wdegree,
		intvec *	shifts,
		const ring	src,
		const ring	Qt
	)

Definition at line 1427 of file hilb.cc.

{
  int rk=A->rank;
  poly h=NULL;
  for(int i=1;i<=rk;i++)
  {
    ideal AA=id_Head(A,src);
    BOOLEAN have_terms=FALSE;
    for(int ii=0;ii<IDELEMS(AA);ii++)
    {
      if (AA->m[ii]!=NULL)
      {
        if(p_GetComp(AA->m[ii],src)!=i)
          p_Delete(&AA->m[ii],src);
        else
        {
          p_SetComp(AA->m[ii],0,src);
          p_Setm(AA->m[ii],src);
          have_terms=TRUE;
        }
      }
    }
    poly h_i=NULL;
    //int sh=0;
    if (have_terms)
    {
      idSkipZeroes(AA);
      h_i=hFirstSeries0p(AA,Q,wdegree,src,Qt);
    }
    else
    {
      h_i=p_One(Qt);
    }
    id_Delete(&AA,src);
    poly s=p_One(Qt);
    if (shifts!=NULL)
    {
      int m=shifts->min_in();
      int sh=(*shifts)[i-1]-m;
      if (sh!=0)
      {
        p_SetExp(s,1,sh,Qt);
        p_Setm(s,Qt);
      }
    }
    h_i=p_Mult_q(h_i,s,Qt);
    h=p_Add_q(h,h_i,Qt);
  }
  return h;
}

◆ hFirstSeries0p()

poly hFirstSeries0p	(	ideal	A,
		ideal	Q,
		intvec *	wdegree,
		const ring	src,
		const ring	Qt
	)

Definition at line 1384 of file hilb.cc.

{
  A=id_Head(A,src);
  id_Test(A,src);
  ideal AA;
  if (Q!=NULL)
  {
    ideal QQ=id_Head(Q,src);
    AA=id_SimpleAdd(A,QQ,src);
    id_Delete(&QQ,src);
    id_Delete(&A,src);
    idSkipZeroes(AA);
    int c=p_GetComp(AA->m[0],src);
    if (c!=0)
    {
      for(int i=0;i<IDELEMS(AA);i++)
        if (AA->m[i]!=NULL) p_SetComp(AA->m[i],c,src);
    }
  }
  else AA=A;
  id_DelDiv(AA,src);
  IDELEMS(AA)=idSkipZeroes0(AA);
   /* sort */
  if (IDELEMS(AA)>1)
  #ifdef HAVE_QSORT_R
    #if defined(__APPLE__) || defined(__FreeBSD__) || defined(__OpenBSD__) || defined(__NetBSD__) || defined(__CYGWIN__)
      qsort_r(AA->m,IDELEMS(AA),sizeof(poly),src,compare_rp);
    #else
      qsort_r(AA->m,IDELEMS(AA),sizeof(poly),compare_rp,src);
    #endif
  #else
  {
    ring r=currRing;
    currRing=src;
    qsort(AA->m,IDELEMS(AA),sizeof(poly),compare_rp_currRing);
    currRing=r;
  }
  #endif
  poly s=hilbert_series(AA,src,wdegree,Qt);
  id_Delete0(&AA,src);
  return s;
}

◆ hFirstSeries1()

intvec * hFirstSeries1	(	ideal	S,
		intvec *	modulweight,
		ideal	Q,
		intvec *	wdegree
	)

Definition at line 1974 of file hilb.cc.

{
  id_LmTest(S, currRing);
  if (Q!= NULL) id_LmTest(Q, currRing);
 
  intvec *hseries1= hSeries(S, modulweight,wdegree, Q);
  if (errorreported) { delete hseries1; hseries1=NULL; }
  return hseries1;
}

◆ hHilbEst()

static void hHilbEst	(	scfmon	stc,
		int	Nstc,
		varset	var,
		int	Nvar
	)

static

Definition at line 1722 of file hilb.cc.

{
  int  i, j;
  int  x, y, z = 1;
  int64  *p;
  for (i = Nvar; i>0; i--)
  {
    x = 0;
    for (j = 0; j < Nstc; j++)
    {
      y = stc[j][var[i]];
      if (y > x)
        x = y;
    }
    z += x;
    j = i - 1;
    if (z > Ql[j])
    {
      if (z>(MAX_INT_VAL)/2)
      {
       WerrorS("internal arrays too big");
       return;
      }
      p = (int64 *)omAlloc((unsigned long)z * sizeof(int64));
      if (Ql[j]!=0)
      {
        if (j==0)
          memcpy(p, Qpol[j], Ql[j] * sizeof(int64));
        omFreeSize((ADDRESS)Qpol[j], Ql[j] * sizeof(int64));
      }
      if (j==0)
      {
        for (x = Ql[j]; x < z; x++)
          p[x] = 0;
      }
      Ql[j] = z;
      Qpol[j] = p;
    }
  }
}

◆ hHilbStep()

static void hHilbStep	(	scmon	pure,
		scfmon	stc,
		int	Nstc,
		varset	var,
		int	Nvar,
		int64 *	pol,
		int	Lpol
	)

static

Definition at line 1763 of file hilb.cc.

{
  int  iv = Nvar -1, ln, a, a0, a1, b, i;
  int  x, x0;
  scmon pn;
  scfmon sn;
  int64  *pon;
  if (Nstc==0)
  {
    hLastHilb(pure, iv, var, pol, Lpol);
    return;
  }
  x = a = 0;
  pn = hGetpure(pure);
  sn = hGetmem(Nstc, stc, stcmem[iv]);
  hStepS(sn, Nstc, var, Nvar, &a, &x);
  Q0[iv] = Q0[Nvar];
  ln = Lpol;
  pon = pol;
  if (a == Nstc)
  {
    x = pure[var[Nvar]];
    if (x!=0)
      pon = hAddHilb(iv, x, pon, &ln);
    hHilbStep(pn, sn, a, var, iv, pon, ln);
    return;
  }
  else
  {
    pon = hAddHilb(iv, x, pon, &ln);
    hHilbStep(pn, sn, a, var, iv, pon, ln);
  }
  b = a;
  x0 = 0;
  loop
  {
    Q0[iv] += (x - x0);
    a0 = a;
    x0 = x;
    hStepS(sn, Nstc, var, Nvar, &a, &x);
    hElimS(sn, &b, a0, a, var, iv);
    a1 = a;
    hPure(sn, a0, &a1, var, iv, pn, &i);
    hLex2S(sn, b, a0, a1, var, iv, hwork);
    b += (a1 - a0);
    ln = Lpol;
    if (a < Nstc)
    {
      pon = hAddHilb(iv, x - x0, pol, &ln);
      hHilbStep(pn, sn, b, var, iv, pon, ln);
    }
    else
    {
      x = pure[var[Nvar]];
      if (x!=0)
        pon = hAddHilb(iv, x - x0, pol, &ln);
      else
        pon = pol;
      hHilbStep(pn, sn, b, var, iv, pon, ln);
      return;
    }
  }
}

◆ hilbert_series()

static poly hilbert_series	(	ideal	A,
		const ring	src,
		const intvec *	wdegree,
		const ring	Qt
	)

static

Definition at line 1289 of file hilb.cc.

{
  int r=id_Elem(A,src);
  poly h=NULL;
  if (r==0)
    return p_One(Qt);
  if (wdegree!=NULL)
  {
    int* exp=(int*)omAlloc((src->N+1)*sizeof(int));
    for(int i=IDELEMS(A)-1; i>=0;i--)
    {
      if (A->m[i]!=NULL)
      {
        p_GetExpV(A->m[i],exp,src);
        for(int j=src->N;j>0;j--)
        {
          int w=(*wdegree)[j-1];
          if (w<=0)
          {
            WerrorS("weights must be positive");
            return NULL;
          }
          exp[j]*=w; /* (*wdegree)[j-1] */
        }
        p_SetExpV(A->m[i],exp,src);
        #ifdef PDEBUG
        p_Setm(A->m[i],src);
        #endif
      }
    }
    omFreeSize(exp,(src->N+1)*sizeof(int));
  }
  h=p_Init(Qt); pSetCoeff0(h,n_Init(-1,Qt->cf));
  p_SetExp(h,1,p_Totaldegree(A->m[0],src),Qt);
  //p_Setm(h,Qt);
  h=p_Add_q(h,p_One(Qt),Qt); // 1-t
  int *exp_q=(int*)omAlloc((src->N+1)*sizeof(int));
  BOOLEAN *bad=(BOOLEAN*)omAlloc0(r*sizeof(BOOLEAN));
  for (int i=1;i<r;i++)
  {
    ideal J=id_CopyFirstK(A,i,src);
    for(int ii=src->N;ii>0;ii--)
      exp_q[ii]=p_GetExp(A->m[i],ii,src);
    memset(bad,0,i*sizeof(BOOLEAN));
    for(int ii=0;ii<i;ii++)
    {
      bad[ii]=p_Div_hi(J->m[ii],exp_q,src);
    }
    id_DelDiv_hi(J,bad,src);
    // variant A
    // search linear elems:
    int k=0;
    for (int ii=IDELEMS(J)-1;ii>=0;ii--)
    {
      if((J->m[ii]!=NULL) && (bad[ii]) && (p_Totaldegree(J->m[ii],src)==1))
      {
        k++;
        p_LmDelete(&J->m[ii],src);
      }
    }
    IDELEMS(J)=idSkipZeroes0(J);
    poly h_J=hilbert_series(J,src,NULL,Qt);// J_1
    poly tmp;
    if (k>0)
    {
      // hilbert_series of unmodified J:
      tmp=p_Init(Qt); pSetCoeff0(tmp,n_Init(-1,Qt->cf));
      p_SetExp(tmp,1,1,Qt);
      //p_Setm(tmp,Qt);
      tmp=p_Add_q(tmp,p_One(Qt),Qt); // 1-t
      if (k>1)
      {
        tmp=p_Power(tmp,k,Qt); // (1-t)^k
      }
      h_J=p_Mult_q(h_J,tmp,Qt);
    }
    // forget about J:
    id_Delete0(&J,src);
    // t^|A_i|
    tmp=p_Init(Qt); pSetCoeff0(tmp,n_Init(-1,Qt->cf));
    p_SetExp(tmp,1,p_Totaldegree(A->m[i],src),Qt);
    //p_Setm(tmp,Qt);
    tmp=p_Mult_q(tmp,h_J,Qt);
    h=p_Add_q(h,tmp,Qt);
  }
  omFreeSize(bad,r*sizeof(BOOLEAN));
  omFreeSize(exp_q,(src->N+1)*sizeof(int));
  //Print("end hilbert_series, r=%d\n",r);
  return h;
}

◆ HilbertSeries_OrbitData()

void HilbertSeries_OrbitData	(	ideal	S,
		int	lV,
		bool	IG_CASE,
		bool	mgrad,
		bool	odp,
		int	trunDegHs
	)

Definition at line 770 of file hilb.cc.

{
 
        /* new story:
        no lV is needed, i.e.  it is to be determined
        the rest is extracted from the interface input list in extra.cc and makes the input of this proc
        called from extra.cc
        */
 
  /*
   * This is based on iterative right colon operations on a
   * two-sided monomial ideal of the free associative algebra.
   * The algorithm terminates for those monomial ideals
   * whose monomials define  "regular formal languages",
   * that is, all monomials of the input ideal can be obtained
   * from finite languages by applying finite number of
   * rational operations.
   */
 
  int trInd;
  S = minimalMonomialGenSet(S);
  if( !idIs0(S) && p_Totaldegree(S->m[0], currRing)==0)
  {
    PrintS("Hilbert Series:\n 0\n");
    return;
  }
  int (*POS)(ideal, poly, std::vector<ideal>, std::vector<poly>, int, int);
  if(trunDegHs != 0)
  {
    Print("\nTruncation degree = %d\n",trunDegHs);
    POS=&positionInOrbitTruncationCase;
  }
  else
  {
    if(IG_CASE)
    {
      if(idIs0(S))
      {
        WerrorS("wrong input: it is not an infinitely gen. case");
        return;
      }
      trInd = p_Totaldegree(S->m[IDELEMS(S)-1], currRing);
      POS = &positionInOrbit_IG_Case;
    }
    else
    POS = &positionInOrbit_FG_Case;
  }
  std::vector<ideal > idorb;
  std::vector< poly > polist;
 
  ideal orb_init = idInit(1, 1);
  idorb.push_back(orb_init);
 
  polist.push_back( p_One(currRing));
 
  std::vector< std::vector<int> > posMat;
  std::vector<int> posRow(lV,0);
  std::vector<int> C;
 
  int ds, is, ps;
  unsigned long lpcnt = 0;
 
  poly w, wi;
  ideal Jwi;
 
  while(lpcnt < idorb.size())
  {
    w = NULL;
    w = polist[lpcnt];
    if(lpcnt >= 1 && idIs0(idorb[lpcnt]) == FALSE)
    {
      if(p_Totaldegree(idorb[lpcnt]->m[0], currRing) != 0)
      {
        C.push_back(1);
      }
      else
        C.push_back(0);
    }
    else
    {
      C.push_back(1);
    }
 
    ds = p_Totaldegree(w, currRing);
    lpcnt++;
 
    for(is = 1; is <= lV; is++)
    {
      wi = NULL;
      //make new copy 'wi' of word w=polist[lpcnt]
      //and update it (for the colon operation).
      //if corresponding to wi, right colon operation gives
      //a new (right colon) ideal of S,
      //keep 'wi' in the polist else delete it
 
      wi = pCopy(w);
      p_SetExp(wi, (ds*lV)+is, 1, currRing);
      p_Setm(wi, currRing);
      Jwi = NULL;
      //Jwi stores (right) colon ideal of S w.r.t. word
      //wi if colon operation gives a new ideal place it
      //in the vector of ideals  'idorb'
      //otherwise delete it
 
      Jwi = idInit(1,1);
 
      Jwi = colonIdeal(S, wi, lV, Jwi, trunDegHs);
      ps = (*POS)(Jwi, wi, idorb, polist, trInd, trunDegHs);
 
      if(ps == 0)  // finds a new ideal
      {
        posRow[is-1] = idorb.size();
 
        idorb.push_back(Jwi);
        polist.push_back(wi);
      }
      else // ideal is already there  in the set
      {
        posRow[is-1]=ps-1;
        idDelete(&Jwi);
        pDelete(&wi);
      }
    }
    posMat.push_back(posRow);
    posRow.resize(lV,0);
  }
  int lO = C.size();//size of the orbit
  PrintLn();
  Print("maximal length of words = %ld\n", p_Totaldegree(polist[lO-1], currRing));
  Print("\nlength of the Orbit = %d", lO);
  PrintLn();
 
  if(odp)
  {
    Print("words description of the Orbit: \n");
    for(is = 0; is < lO; is++)
    {
      pWrite0(polist[is]);
      PrintS("    ");
    }
    PrintLn();
    PrintS("\nmaximal degree,  #(sum_j R(w,w_j))");
    PrintLn();
    for(is = 0; is < lO; is++)
    {
      if(idIs0(idorb[is]))
      {
        PrintS("NULL\n");
      }
      else
      {
        Print("%ld,  %d \n",p_Totaldegree(idorb[is]->m[IDELEMS(idorb[is])-1], currRing),IDELEMS(idorb[is]));
      }
    }
  }
 
  for(is = idorb.size()-1; is >= 0; is--)
  {
    idDelete(&idorb[is]);
  }
  for(is = polist.size()-1; is >= 0; is--)
  {
    pDelete(&polist[is]);
  }
 
  idorb.resize(0);
  polist.resize(0);
 
  int adjMatrix[lO][lO];
  memset(adjMatrix, 0, lO*lO*sizeof(int));
  int rowCount, colCount;
  int tm = 0;
  if(!mgrad)
  {
    for(rowCount = 0; rowCount < lO; rowCount++)
    {
      for(colCount = 0; colCount < lV; colCount++)
      {
        tm = posMat[rowCount][colCount];
        adjMatrix[rowCount][tm] = adjMatrix[rowCount][tm] + 1;
      }
    }
  }
 
  ring r = currRing;
  int npar;
  char** tt;
  TransExtInfo p;
  if(!mgrad)
  {
    tt=(char**)omAlloc(2*sizeof(char*));
    tt[0] = omStrDup("t");
    npar = 1;
  }
  else
  {
    tt=(char**)omalloc(lV*sizeof(char*));
    for(is = 0; is < lV; is++)
    {
      tt[is] = (char*)omAlloc(12*sizeof(char)); //if required enlarge it later
      snprintf (tt[is],12, "t%d", is+1);
    }
    npar = lV;
  }
 
  p.r = rDefault(0, npar, tt);
  coeffs cf = nInitChar(n_transExt, &p);
  char** xx = (char**)omAlloc(sizeof(char*));
  xx[0] = omStrDup("x");
  ring R = rDefault(cf, 1, xx);
  rChangeCurrRing(R);//rWrite(R);
  /*
   * matrix corresponding to the orbit of the ideal
   */
  matrix mR = mpNew(lO, lO);
  matrix cMat = mpNew(lO,1);
  poly rc;
 
  if(!mgrad)
  {
    for(rowCount = 0; rowCount < lO; rowCount++)
    {
      for(colCount = 0; colCount < lO; colCount++)
      {
        if(adjMatrix[rowCount][colCount] != 0)
        {
          MATELEM(mR, rowCount + 1, colCount + 1) = p_ISet(adjMatrix[rowCount][colCount], R);
          p_SetCoeff(MATELEM(mR, rowCount + 1, colCount + 1), n_Mult(pGetCoeff(mR->m[lO*rowCount+colCount]),n_Param(1, R->cf), R->cf), R);
        }
      }
    }
  }
  else
  {
     for(rowCount = 0; rowCount < lO; rowCount++)
     {
       for(colCount = 0; colCount < lV; colCount++)
       {
          rc=NULL;
          rc=p_One(R);
          p_SetCoeff(rc, n_Mult(pGetCoeff(rc), n_Param(colCount+1, R->cf),R->cf), R);
          MATELEM(mR, rowCount +1, posMat[rowCount][colCount]+1)=p_Add_q(rc,MATELEM(mR, rowCount +1, posMat[rowCount][colCount]+1), R);
       }
     }
  }
 
  for(rowCount = 0; rowCount < lO; rowCount++)
  {
    if(C[rowCount] != 0)
    {
      MATELEM(cMat, rowCount + 1, 1) = p_ISet(C[rowCount], R);
    }
  }
 
  matrix u;
  unitMatrix(lO, u); //unit matrix
  matrix gMat = mp_Sub(u, mR, R);
 
  char* s;
 
  if(odp)
  {
    PrintS("\nlinear system:\n");
    if(!mgrad)
    {
      for(rowCount = 0; rowCount < lO; rowCount++)
      {
        Print("H(%d) = ", rowCount+1);
        for(colCount = 0; colCount < lV; colCount++)
        {
          StringSetS(""); nWrite(n_Param(1, R->cf));
          s = StringEndS(); PrintS(s);
          Print("*"); omFree(s);
          Print("H(%d) + ", posMat[rowCount][colCount] + 1);
        }
        Print(" %d\n", C[rowCount] );
      }
      PrintS("where H(1) represents the series corresp. to input ideal\n");
      PrintS("and i^th summand in the rhs of an eqn. is according\n");
      PrintS("to the right colon map corresp. to the i^th variable\n");
    }
    else
    {
      for(rowCount = 0; rowCount < lO; rowCount++)
      {
        Print("H(%d) = ", rowCount+1);
        for(colCount = 0; colCount < lV; colCount++)
        {
           StringSetS(""); nWrite(n_Param(colCount+1, R->cf));
           s = StringEndS(); PrintS(s);
           Print("*");omFree(s);
           Print("H(%d) + ", posMat[rowCount][colCount] + 1);
        }
        Print(" %d\n", C[rowCount] );
      }
      PrintS("where H(1) represents the series corresp. to input ideal\n");
    }
  }
  PrintLn();
  posMat.resize(0);
  C.resize(0);
  matrix pMat;
  matrix lMat;
  matrix uMat;
  matrix H_serVec = mpNew(lO, 1);
  matrix Hnot;
 
  //std::clock_t start;
  //start = std::clock();
 
  luDecomp(gMat, pMat, lMat, uMat, R);
  luSolveViaLUDecomp(pMat, lMat, uMat, cMat, H_serVec, Hnot);
 
  //to print system solving time
  //if(odp){
  //std::cout<<"solving time of the system = "<<(std::clock()-start)/(double)(CLOCKS_PER_SEC / 1000)<<" ms"<<std::endl;}
 
  mp_Delete(&mR, R);
  mp_Delete(&u, R);
  mp_Delete(&pMat, R);
  mp_Delete(&lMat, R);
  mp_Delete(&uMat, R);
  mp_Delete(&cMat, R);
  mp_Delete(&gMat, R);
  mp_Delete(&Hnot, R);
  //print the Hilbert series and length of the Orbit
  PrintLn();
  Print("Hilbert series:");
  PrintLn();
  pWrite(H_serVec->m[0]);
  if(!mgrad)
  {
    omFree(tt[0]);
  }
  else
  {
    for(is = lV-1; is >= 0; is--)
 
      omFree( tt[is]);
  }
  omFree(tt);
  omFree(xx[0]);
  omFree(xx);
  rChangeCurrRing(r);
  rKill(R);
}

◆ hLastHilb()

static void hLastHilb	(	scmon	pure,
		int	Nv,
		varset	var,
		int64 *	pol,
		int	lp
	)

static

Definition at line 1666 of file hilb.cc.

{
  int  l = lp, x, i, j;
  int64  *pl;
  int64  *p;
  p = pol;
  for (i = Nv; i>0; i--)
  {
    x = pure[var[i + 1]];
    if (x!=0)
      p = hAddHilb(i, x, p, &l);
  }
  pl = *Qpol;
  j = Q0[Nv + 1];
  for (i = 0; i < l; i++)
  { /* pl[i + j] += p[i];*/
    #ifndef __SIZEOF_INT128__
    int64 t=pl[i+j];
    int64 t2=p[i];
    t+=t2;
    if ((t>=OVERFLOW_MIN)&&(t<=OVERFLOW_MAX)) pl[i+j]=t;
    else if (!errorreported) WerrorS("int overflow in hilb 3");
    #else
    __int128 t=pl[i+j];
    __int128 t2=p[i];
    t+=t2;
    if ((t>=LONG_MIN)&&(t<=LONG_MAX)) pl[i+j]=t;
    else if (!errorreported) WerrorS("long int overflow in hilb 3");
    #endif
  }
  x = pure[var[1]];
  if (x!=0)
  {
    j += x;
    for (i = 0; i < l; i++)
    { /* pl[i + j] -= p[i];*/
      #ifndef __SIZEOF_INT128__
      int64 t=pl[i+j];
      int64 t2=p[i];
      t-=t2;
      if ((t>=OVERFLOW_MIN)&&(t<=OVERFLOW_MAX)) pl[i+j]=t;
      else if (!errorreported) WerrorS("int overflow in hilb 4");
      #else
      __int128 t=pl[i+j];
      __int128 t2=p[i];
      t-=t2;
      if ((t>=LONG_MIN)&&(t<=LONG_MAX)) pl[i+j]=t;
      else if (!errorreported) WerrorS("long int overflow in hilb 4");
      #endif
    }
  }
  j += l;
  if (j > hLength)
    hLength = j;
}

◆ hLookSeries()

void hLookSeries	(	ideal	S,
		intvec *	modulweight,
		ideal	Q,
		intvec *	wdegree
	)

Definition at line 220 of file hilb.cc.

{
  id_LmTest(S, currRing);
 
  if (!id_IsModule(S,currRing))
  {
    if (hilb_Qt==NULL) hilb_Qt=makeQt();
    poly hseries=hFirstSeries0p(S,Q,wdegree,currRing,hilb_Qt);
 
    hPrintHilb(hseries,hilb_Qt,wdegree);
    p_Delete(&hseries,hilb_Qt);
  }
  else
  {
    if (hilb_Qt==NULL) hilb_Qt=makeQt();
    poly hseries=hFirstSeries0m(S,Q,wdegree,modulweight,currRing,hilb_Qt);
    if ((modulweight!=NULL)&&(modulweight->compare(0)!=0))
    {
      char *s=modulweight->ivString(1,0,1);
      Print("module weights:%s\n",s);
      omFree(s);
    }
    hPrintHilb(hseries,hilb_Qt,wdegree);
    p_Delete(&hseries,hilb_Qt);
  }
}

◆ hMinModulweight()

static int hMinModulweight ( intvec * modulweight )

static

Definition at line 1595 of file hilb.cc.

{
  if(modulweight==NULL) return 0;
  return modulweight->min_in();
}

◆ hPoly2BIV()

bigintmat * hPoly2BIV	(	poly	h,
		const ring	Qt,
		const coeffs	biv_cf
	)

Definition at line 1984 of file hilb.cc.

{
  int td=0;
  nMapFunc f;
  if (h!=NULL)
  {
    td=p_Totaldegree(h,Qt);
    h=p_Copy(h,Qt);
    f=n_SetMap(Qt->cf,biv_cf);
  }
  bigintmat* biv=new bigintmat(1,td+2,biv_cf);
  while(h!=NULL)
  {
    int d=p_Totaldegree(h,Qt);
    n_Delete(&BIMATELEM(*biv,1,d+1),biv_cf);
    BIMATELEM(*biv,1,d+1)=f(pGetCoeff(h),Qt->cf,biv_cf);
    p_LmDelete(&h,Qt);
  }
  return biv;
}

◆ hPrintHilb()

static void hPrintHilb	(	poly	hseries,
		const ring	Qt,
		intvec *	modul_weight
	)

static

Definition at line 164 of file hilb.cc.

{
  if ((modul_weight!=NULL)&&(modul_weight->compare(0)!=0))
  {
    char *s=modul_weight->ivString(1,0,1);
    Print("module weights:%s\n",s);
    omFree(s);
  }
  PrintS("(");p_Write0(hseries,Qt);Print(") / (1-%s)^%d\n",Qt->names[0],currRing->N);
  int co;
  poly h2=hFirst2Second(hseries,Qt,co);
  int di = (currRing->N)-co;
  if (hseries==NULL) di=0;
  PrintS("("); p_Write0(h2,Qt); Print(") / (1-%s)^%d\n",Qt->names[0],di);
  int mu=0;
  poly p=h2;
  while(p!=NULL)
  {
    mu+=n_Int(pGetCoeff(p),Qt->cf);
    p_LmDelete(&p,Qt);
  }
  if (currRing->OrdSgn == 1)
  {
    if (di>0)
      Print("// dimension (proj.)  = %d\n// degree (proj.)   = %d\n", di-1, mu);
    else
      Print("// dimension (affine) = 0\n// degree (affine)  = %d\n",       mu);
  }
  else
    Print("// dimension (local)   = %d\n// multiplicity = %d\n", di, mu);
}

◆ hSecondSeries()

intvec * hSecondSeries ( intvec * hseries1 )

Definition at line 71 of file hilb.cc.

{
  intvec *work, *hseries2;
  int i, j, k, t, l;
  int s;
  if (hseries1 == NULL)
    return NULL;
  work = new intvec(hseries1);
  k = l = work->length()-1;
  s = 0;
  for (i = k-1; i >= 0; i--)
    s += (*work)[i];
  loop
  {
    if ((s != 0) || (k == 1))
      break;
    s = 0;
    t = (*work)[k-1];
    k--;
    for (i = k-1; i >= 0; i--)
    {
      j = (*work)[i];
      (*work)[i] = -t;
      s += t;
      t += j;
    }
  }
  hseries2 = new intvec(k+1);
  for (i = k-1; i >= 0; i--)
    (*hseries2)[i] = (*work)[i];
  (*hseries2)[k] = (*work)[l];
  delete work;
  return hseries2;
}

◆ hSecondSeries0b()

bigintmat * hSecondSeries0b	(	ideal	I,
		ideal	Q,
		intvec *	wdegree,
		intvec *	shifts,
		const ring	src,
		const coeffs	biv_cf
	)

Definition at line 2041 of file hilb.cc.

{
  if (hilb_Qt==NULL) hilb_Qt=makeQt();
  poly h;
  if (id_IsModule(I,src))
    h=hFirstSeries0m(I,Q,wdegree,shifts,src,hilb_Qt);
  else
    h=hFirstSeries0p(I,Q,wdegree,src,hilb_Qt);
  int co;
  poly h2=hFirst2Second(h,hilb_Qt,co);
  p_Delete(&h,hilb_Qt);
  bigintmat *biv=hPoly2BIV(h2,hilb_Qt,biv_cf);
  p_Delete(&h2,hilb_Qt);
  return biv;
}

◆ hSeries()

static intvec * hSeries	(	ideal	S,
		intvec *	modulweight,
		intvec *	wdegree,
		ideal	Q
	)

static

Definition at line 1828 of file hilb.cc.

{
  intvec *work, *hseries1=NULL;
  int  mc;
  int64  p0;
  int  i, j, k, l, ii, mw;
  hexist = hInit(S, Q, &hNexist);
  if (hNexist==0)
  {
    hseries1=new intvec(2);
    (*hseries1)[0]=1;
    (*hseries1)[1]=0;
    return hseries1;
  }
 
  if (wdegree != NULL) hWDegree(wdegree);
 
  p0 = 1;
  hwork = (scfmon)omAlloc(hNexist * sizeof(scmon));
  hvar = (varset)omAlloc(((currRing->N) + 1) * sizeof(int));
  hpure = (scmon)omAlloc((1 + ((currRing->N) * (currRing->N))) * sizeof(int));
  stcmem = hCreate((currRing->N) - 1);
  Qpol = (int64 **)omAlloc(((currRing->N) + 1) * sizeof(int64 *));
  Ql = (int64 *)omAlloc0(((currRing->N) + 1) * sizeof(int64));
  Q0 = (int64 *)omAlloc(((currRing->N) + 1) * sizeof(int64));
  *Qpol = NULL;
  hLength = k = j = 0;
  mc = hisModule;
  if (mc!=0)
  {
    mw = hMinModulweight(modulweight);
    hstc = (scfmon)omAlloc(hNexist * sizeof(scmon));
  }
  else
  {
    mw = 0;
    hstc = hexist;
    hNstc = hNexist;
  }
  loop
  {
    if (mc!=0)
    {
      hComp(hexist, hNexist, mc, hstc, &hNstc);
      if (modulweight != NULL)
        j = (*modulweight)[mc-1]-mw;
    }
    if (hNstc!=0)
    {
      hNvar = (currRing->N);
      for (i = hNvar; i>=0; i--)
        hvar[i] = i;
      //if (notstc) // TODO: no mon divides another
        hStaircase(hstc, &hNstc, hvar, hNvar);
      hSupp(hstc, hNstc, hvar, &hNvar);
      if (hNvar!=0)
      {
        if ((hNvar > 2) && (hNstc > 10))
          hOrdSupp(hstc, hNstc, hvar, hNvar);
        hHilbEst(hstc, hNstc, hvar, hNvar);
        memset(hpure, 0, ((currRing->N) + 1) * sizeof(int));
        hPure(hstc, 0, &hNstc, hvar, hNvar, hpure, &hNpure);
        hLexS(hstc, hNstc, hvar, hNvar);
        Q0[hNvar] = 0;
        hHilbStep(hpure, hstc, hNstc, hvar, hNvar, &p0, 1);
      }
    }
    else
    {
      if(*Qpol!=NULL)
        (**Qpol)++;
      else
      {
        *Qpol = (int64 *)omAlloc(sizeof(int64));
        hLength = *Ql = **Qpol = 1;
      }
    }
    if (*Qpol!=NULL)
    {
      i = hLength;
      while ((i > 0) && ((*Qpol)[i - 1] == 0))
        i--;
      if (i > 0)
      {
        l = i + j;
        if (l > k)
        {
          work = new intvec(l);
          for (ii=0; ii<k; ii++)
            (*work)[ii] = (*hseries1)[ii];
          if (hseries1 != NULL)
            delete hseries1;
          hseries1 = work;
          k = l;
        }
        while (i > 0)
        {
          (*hseries1)[i + j - 1] += (*Qpol)[i - 1];
          (*Qpol)[i - 1] = 0;
          i--;
        }
      }
    }
    mc--;
    if (mc <= 0)
      break;
  }
  if (k==0)
  {
    hseries1=new intvec(2);
    (*hseries1)[0]=0;
    (*hseries1)[1]=0;
  }
  else
  {
    l = k+1;
    while ((*hseries1)[l-2]==0) l--;
    if (l!=k)
    {
      work = new intvec(l);
      for (ii=l-2; ii>=0; ii--)
        (*work)[ii] = (*hseries1)[ii];
      delete hseries1;
      hseries1 = work;
    }
    (*hseries1)[l-1] = mw;
  }
  for (i = 0; i <= (currRing->N); i++)
  {
    if (Ql[i]!=0)
      omFreeSize((ADDRESS)Qpol[i], Ql[i] * sizeof(int64));
  }
  omFreeSize((ADDRESS)Q0, ((currRing->N) + 1) * sizeof(int64));
  omFreeSize((ADDRESS)Ql, ((currRing->N) + 1) * sizeof(int64));
  omFreeSize((ADDRESS)Qpol, ((currRing->N) + 1) * sizeof(int64 *));
  hKill(stcmem, (currRing->N) - 1);
  omFreeSize((ADDRESS)hpure, (1 + ((currRing->N) * (currRing->N))) * sizeof(int));
  omFreeSize((ADDRESS)hvar, ((currRing->N) + 1) * sizeof(int));
  omFreeSize((ADDRESS)hwork, hNexist * sizeof(scmon));
  hDelete(hexist, hNexist);
  if (hisModule!=0)
    omFreeSize((ADDRESS)hstc, hNexist * sizeof(scmon));
  return hseries1;
}

◆ hWDegree()

static void hWDegree ( intvec * wdegree )

static

Definition at line 1601 of file hilb.cc.

{
  int i, k;
  int x;
 
  for (i=(currRing->N); i; i--)
  {
    x = (*wdegree)[i-1];
    if (x != 1)
    {
      for (k=hNexist-1; k>=0; k--)
      {
        hexist[k][i] *= x;
      }
    }
  }
}

◆ id_DelDiv_hi()

static void id_DelDiv_hi	(	ideal	id,
		BOOLEAN *	bad,
		const ring	r
	)

static

Definition at line 1205 of file hilb.cc.

{
  int k=IDELEMS(id)-1;
  while(id->m[k]==NULL) k--;
  int kk = k+1;
  long *sev=(long*)omAlloc0(kk*sizeof(long));
  BOOLEAN only_lm=r->cf->has_simple_Alloc;
  if (BIT_SIZEOF_LONG / r->N==0) // 1 bit per exp
  {
    for (int i=k; i>=0; i--)
    {
        sev[i]=p_GetShortExpVector0(id->m[i],r);
    }
  }
  else
  if (BIT_SIZEOF_LONG / r->N==1) // 1..2 bit per exp
  {
    for (int i=k; i>=0; i--)
    {
        sev[i]=p_GetShortExpVector1(id->m[i],r);
    }
  }
  else
  {
    for (int i=k; i>=0; i--)
    {
        sev[i]=p_GetShortExpVector(id->m[i],r);
    }
  }
  if (only_lm)
  {
    for (int i=0; i<k; i++)
    {
      if (bad[i] && (id->m[i] != NULL))
      {
        poly m_i=id->m[i];
        long sev_i=sev[i];
        for (int j=i+1; j<=k; j++)
        {
          if (id->m[j]!=NULL)
          {
            if (p_LmShortDivisibleBy(m_i, sev_i, id->m[j],~sev[j],r))
            {
              p_LmFree(&id->m[j],r);
            }
            else if (p_LmShortDivisibleBy(id->m[j],sev[j], m_i,~sev_i,r))
            {
              p_LmFree(&id->m[i],r);
              break;
            }
          }
        }
      }
    }
  }
  else
  {
    for (int i=0; i<k; i++)
    {
      if (bad[i] && (id->m[i] != NULL))
      {
        poly m_i=id->m[i];
        long sev_i=sev[i];
        for (int j=i+1; j<=k; j++)
        {
          if (id->m[j]!=NULL)
          {
            if (p_LmShortDivisibleBy(m_i, sev_i, id->m[j],~sev[j],r))
            {
              p_Delete(&id->m[j],r);
            }
            else if (p_LmShortDivisibleBy(id->m[j],sev[j], m_i,~sev_i,r))
            {
              p_Delete(&id->m[i],r);
              break;
            }
          }
        }
      }
    }
  }
  omFreeSize(sev,kk*sizeof(long));
}

◆ idInsertMonomial()

static void idInsertMonomial	(	ideal	I,
		poly	p
	)

static

Definition at line 252 of file hilb.cc.

{
  /*
   * It adds monomial in I and if required,
   * enlarge the size of poly-set by 16.
   * It does not make copy of  p.
   */
 
  if(I == NULL)
  {
    return;
  }
 
  int j = IDELEMS(I) - 1;
  while ((j >= 0) && (I->m[j] == NULL))
  {
    j--;
  }
  j++;
  if (j == IDELEMS(I))
  {
    pEnlargeSet(&(I->m), IDELEMS(I), 16);
    IDELEMS(I) +=16;
  }
  I->m[j] = p;
}

◆ makeQt()

static ring makeQt ( )

static

Definition at line 196 of file hilb.cc.

{
  ring Qt=(ring) omAlloc0Bin(sip_sring_bin);
  Qt->cf = nInitChar(n_Q, NULL);
  Qt->N=1;
  Qt->names=(char**)omAlloc(sizeof(char_ptr));
  Qt->names[0]=omStrDup("t");
  Qt->wvhdl=(int **)omAlloc0(3 * sizeof(int_ptr));
  Qt->order = (rRingOrder_t *) omAlloc(3 * sizeof(rRingOrder_t *));
  Qt->block0 = (int *)omAlloc0(3 * sizeof(int *));
  Qt->block1 = (int *)omAlloc0(3 * sizeof(int *));
  /* ringorder lp for the first block: var 1 */
  Qt->order[0]  = ringorder_lp;
  Qt->block0[0] = 1;
  Qt->block1[0] = 1;
  /* ringorder C for the second block: no vars */
  Qt->order[1]  = ringorder_C;
  /* the last block: everything is 0 */
  Qt->order[2]  = (rRingOrder_t)0;
  rComplete(Qt);
  return Qt;
}

◆ minimalMonomialGenSet()

static ideal minimalMonomialGenSet ( ideal I )

static

Definition at line 575 of file hilb.cc.

{
  /*
   * eliminates monomials which
   * can be generated by others in I
   */
  //first sort monomials of the ideal
 
  idSkipZeroes(I);
 
  sortMonoIdeal_pCompare(I);
 
  int i, k;
  int ICount = IDELEMS(I);
 
  for(k = ICount - 1; k >=1; k--)
  {
    for(i = 0; i < k; i++)
    {
 
      if(p_LmDivisibleBy(I->m[i], I->m[k], currRing))
      {
        pDelete(&(I->m[k]));
        break;
      }
    }
  }
 
  idSkipZeroes(I);
  return(I);
}

◆ monCompare()

static int monCompare	(	const void *	m,
		const void *	n
	)

static

Definition at line 555 of file hilb.cc.

{
  /* compares monomials */
 
 return(p_Compare(*(poly*) m, *(poly*)n, currRing));
}

◆ p_Div_hi()

static BOOLEAN p_Div_hi	(	poly	p,
		const int *	exp_q,
		const ring	src
	)

static

Definition at line 1151 of file hilb.cc.

{
  BOOLEAN bad=FALSE;
  // e=max(0,p-q) for all exps
  for(int i=src->N;i>0;i--)
  {
    int pi=p_GetExp(p,i,src)-exp_q[i];
    if (pi<0)
    {
      pi=0;
      bad=TRUE;
    }
    p_SetExp(p,i,pi,src);
  }
  #ifdef PDEBUG
  p_Setm(p,src);
  #endif
  return bad;
}

◆ positionInOrbit_FG_Case()

static int positionInOrbit_FG_Case	(	ideal	I,
		poly	,
		std::vector< ideal >	idorb,
		std::vector< poly >	,
		int	,
		int
	)

static

Definition at line 450 of file hilb.cc.

{
  /*
   * It compares the ideal I with ideals in the set 'idorb'.
   * I and ideals of 'idorb' are sorted.
   *
   * It returns 0 if I is not equal to any ideal of 'idorb'
   * else returns position of the matched ideal.
   */
  int ps = 0;
  int i, s = 0;
  int OrbCount = idorb.size();
 
  if(idIs0(I))
  {
    return(1);
  }
 
  for(i = 1; i < OrbCount; i++)
  {
    s = comapreMonoIdBases(I, idorb[i]);
    if(s)
    {
      ps = i + 1;
      break;
    }
  }
 
  return(ps);
}

◆ positionInOrbit_IG_Case()

static int positionInOrbit_IG_Case	(	ideal	I,
		poly	w,
		std::vector< ideal >	idorb,
		std::vector< poly >	polist,
		int	trInd,
		int
	)

static

Definition at line 372 of file hilb.cc.

{
  /*
   * It compares the ideal I with ideals in the set 'idorb'
   * up to total degree =
   * trInd - max(deg of w, deg of word in polist) polynomials.
   *
   * It returns 0 if I is not equal to any ideal in the
   * 'idorb' else returns position of the matched ideal.
   */
 
  int ps = 0;
  int i, s = 0;
  int orbCount = idorb.size();
 
  if(idIs0(I))
  {
    return(1);
  }
 
  int degw = p_Totaldegree(w, currRing);
  int degp;
  int dtr;
  int dtrp;
 
  dtr = trInd - degw;
  int IwCount;
 
   IwCount = CountOnIdUptoTruncationIndex(I, dtr);
 
  if(IwCount == 0)
  {
    return(1);
  }
 
  int ObCount;
 
  bool flag2 = FALSE;
 
  for(i = 1;i < orbCount; i++)
  {
    degp = p_Totaldegree(polist[i], currRing);
    if(degw > degp)
    {
      dtr = trInd - degw;
 
      ObCount = 0;
      ObCount = CountOnIdUptoTruncationIndex(idorb[i], dtr);
      if(ObCount == 0)
      {continue;}
      if(flag2)
      {
        IwCount = 0;
        IwCount = CountOnIdUptoTruncationIndex(I, dtr);
        flag2 = FALSE;
      }
    }
    else
    {
      flag2 = TRUE;
      dtrp = trInd - degp;
      ObCount = 0;
      ObCount = CountOnIdUptoTruncationIndex(idorb[i], dtrp);
      IwCount = 0;
      IwCount = CountOnIdUptoTruncationIndex(I, dtrp);
    }
 
    s = comapreMonoIdBases_IG_Case(I, IwCount, idorb[i], ObCount);
 
    if(s)
    {
      ps = i + 1;
      break;
    }
  }
  return(ps);
}

◆ positionInOrbitTruncationCase()

static int positionInOrbitTruncationCase	(	ideal	I,
		poly	w,
		std::vector< ideal >	idorb,
		std::vector< poly >	polist,
		int	,
		int	trunDegHs
	)

static

Definition at line 481 of file hilb.cc.

{
  /*
   * It compares the ideal I with ideals in the set 'idorb'.
   * I and ideals in 'idorb' are sorted.
 
   * returns 0 if I is not equal to any ideal of 'idorb'
   * else returns position of the matched ideal.
   */
 
  int ps = 0;
  int i, s = 0;
  int OrbCount = idorb.size();
  int dtr=0; int IwCount, ObCount;
  dtr = trunDegHs - 1 - p_Totaldegree(w, currRing);
 
  if(idIs0(I))
  {
    for(i = 1; i < OrbCount; i++)
    {
      if(p_Totaldegree(w, currRing) == p_Totaldegree(polist[i], currRing))
      {
        if(idIs0(idorb[i]))
          return(i+1);
        ObCount=0;
        ObCount = CountOnIdUptoTruncationIndex(idorb[i], dtr);
        if(ObCount==0)
        {
          ps = i + 1;
          break;
        }
      }
    }
 
    return(ps);
  }
 
  IwCount = CountOnIdUptoTruncationIndex(I, dtr);
 
  if(p_Totaldegree(I->m[0], currRing)==0)
  {
    for(i = 1; i < OrbCount; i++)
    {
       if(idIs0(idorb[i]))
         continue;
       if(p_Totaldegree(idorb[i]->m[0], currRing)==0)
       {
         ps = i + 1;
         break;
       }
     }
     return(ps);
  }
 
  for(i = 1; i < OrbCount; i++)
  {
    if(p_Totaldegree(w, currRing) == p_Totaldegree(polist[i], currRing))
    {
      if(idIs0(idorb[i]))
        continue;
      ObCount=0;
      ObCount = CountOnIdUptoTruncationIndex(idorb[i], dtr);
      s = comapreMonoIdBases_IG_Case(I, IwCount, idorb[i], ObCount);
      if(s)
      {
        ps = i + 1;
        break;
      }
    }
  }
 
  return(ps);
}

◆ RightColonOperation()

ideal RightColonOperation	(	ideal	S,
		poly	w,
		int	lV
	)

Definition at line 1117 of file hilb.cc.

{
  /*
   * This returns right colon ideal of a monomial two-sided ideal of
   * the free associative algebra with respect to a monomial 'w'
   * (S:_R w).
   */
  S = minimalMonomialGenSet(S);
  ideal Iw = idInit(1,1);
  Iw = colonIdeal(S, w, lV, Iw, 0);
  return (Iw);
}

◆ scDegree()

void scDegree	(	ideal	S,
		intvec *	modulweight,
		ideal	Q
	)

Definition at line 2057 of file hilb.cc.

{
  int co;
  int mu=0;
#if 0
  if (hilb_Qt==NULL) hilb_Qt=makeQt();
  poly h1;
  if (id_IsModule(S,currRing))
    h1 = hFirstSeries0p(S,Q,NULL,currRing,hilb_Qt);
  else
    h1 = hFirstSeries0m(S,Q,NULL, modulweight,currRing,hilb_Qt);
 
  poly h2=hFirst2Second(h1,hilb_Qt,co);
  int di = (currRing->N)-co;
  if (h1==NULL) di=0;
  poly p=h2;
  while(p!=NULL)
  {
    mu+=n_Int(pGetCoeff(p),hilb_Qt->cf);
    p_LmDelete(&p,hilb_Qt);
  }
#else
  bigintmat *h1=hFirstSeries0b(S,Q,NULL,modulweight,currRing,coeffs_BIGINT);
  intvec *hseries1=new intvec(1,h1->cols());
  for(int i=0;i<h1->cols();i++)
  {
    (*hseries1)[i]=n_Int(BIMATELEM(*h1,1,i+1),coeffs_BIGINT);
  }
  delete h1;
  intvec *hseries2;
  int l = hseries1->length()-1;
  if (l > 1)
    hseries2 = hSecondSeries(hseries1);
  else
    hseries2 = hseries1;
  hDegreeSeries(hseries1, hseries2, &co, &mu);
  if (l>1)
    delete hseries1;
  delete hseries2;
  if ((l == 1) &&(mu == 0))
    scPrintDegree((currRing->N)+1, 0);
  else
#endif
    scPrintDegree(co, mu);
}

◆ shiftInMon()

static poly shiftInMon	(	poly	p,
		int	i,
		int	lV,
		const ring	r
	)

static

Definition at line 607 of file hilb.cc.

{
  /*
   * shifts the variables of monomial p in the  i^th layer,
   * p remains unchanged,
   * creates new poly and returns it for the colon ideal
   */
  poly smon = p_One(r);
  int j, sh, cnt;
  cnt = r->N;
  sh = i*lV;
  int *e=(int *)omAlloc((r->N+1)*sizeof(int));
  int *s=(int *)omAlloc0((r->N+1)*sizeof(int));
  p_GetExpV(p, e, r);
 
  for(j = 1; j <= cnt; j++)
  {
    if(e[j] == 1)
    {
      s[j+sh] = e[j];
    }
  }
 
  p_SetExpV(smon, s, currRing);
  omFree(e);
  omFree(s);
 
  p_SetComp(smon, p_GetComp(p, currRing), currRing);
  p_Setm(smon, currRing);
 
  return(smon);
}

◆ sortMonoIdeal_pCompare()

static void sortMonoIdeal_pCompare ( ideal I )

static

Definition at line 562 of file hilb.cc.

{
  /*
   * sorts monomial ideal in ascending order
   * order must be a total degree
   */
 
  qsort(I->m, IDELEMS(I), sizeof(poly), monCompare);
 
}

◆ TwordMap()

static void TwordMap	(	poly	p,
		poly	w,
		int	lV,
		int	d,
		ideal	Jwi,
		bool &	flag
	)

static

Definition at line 674 of file hilb.cc.

{
  /*
   * computes T_w(p) in a new poly object and places it
   * in Jwi which stores elements of colon ideal of I,
   * p and w remain unchanged,
   * the new polys for Jwi are constructed by sub-routines
   * deleteInMon, shiftInMon, p_MDivide,
   * places the result in Jwi and deletes the new polys
   * coming in dw, smon, qmon
   */
  int i;
  poly smon, dw;
  poly qmonp = NULL;
  bool del;
 
  for(i = 0;i <= d - 1; i++)
  {
    dw = deleteInMon(w, i, lV, currRing);
    smon = shiftInMon(p, i, lV, currRing);
    del = TRUE;
 
    if(pLmDivisibleBy(smon, w))
    {
      flag = TRUE;
      del  = FALSE;
 
      pDelete(&dw);
      pDelete(&smon);
 
      //delete all monomials of Jwi
      //and make Jwi =1
 
      for(int j = 0;j < IDELEMS(Jwi); j++)
      {
        pDelete(&Jwi->m[j]);
      }
 
      idInsertMonomial(Jwi, p_One(currRing));
      break;
    }
 
    if(pLmDivisibleBy(dw, smon))
    {
      del = FALSE;
      qmonp = p_MDivide(smon, dw, currRing);
      idInsertMonomial(Jwi, shiftInMon(qmonp, -d, lV, currRing));
      pLmFree(&qmonp);
      pDelete(&dw);
      pDelete(&smon);
    }
    //in case both if are false, delete dw and smon
    if(del)
    {
      pDelete(&dw);
      pDelete(&smon);
    }
  }
 
}

Variable Documentation

◆ hilb_Qt

STATIC_VAR ring hilb_Qt =NULL

Definition at line 219 of file hilb.cc.

◆ hLength

STATIC_VAR int hLength

Definition at line 69 of file hilb.cc.

◆ Q0

STATIC_VAR int64* Q0

Definition at line 68 of file hilb.cc.

◆ Ql

STATIC_VAR int64 * Ql

Definition at line 68 of file hilb.cc.

◆ Qpol

STATIC_VAR int64** Qpol

Definition at line 67 of file hilb.cc.

Macros

Functions

Variables

Macro Definition Documentation

◆ omsai

◆ OVERFLOW_MAX

◆ OVERFLOW_MIN

Function Documentation

◆ colonIdeal()

◆ comapreMonoIdBases()

◆ comapreMonoIdBases_IG_Case()

◆ compare_rp_currRing()

◆ CountOnIdUptoTruncationIndex()

◆ deleteInMon()

◆ getModuleComp()

◆ hAddHilb()

◆ hBIV2Poly()

◆ hDegreeSeries()

◆ hFirst2Second()

◆ hFirstSeries()

◆ hFirstSeries0()

◆ hFirstSeries0b()

◆ hFirstSeries0m()

◆ hFirstSeries0p()

◆ hFirstSeries1()

◆ hHilbEst()

◆ hHilbStep()

◆ hilbert_series()

◆ HilbertSeries_OrbitData()

◆ hLastHilb()

◆ hLookSeries()

◆ hMinModulweight()

◆ hPoly2BIV()

◆ hPrintHilb()

◆ hSecondSeries()

◆ hSecondSeries0b()

◆ hSeries()

◆ hWDegree()

◆ id_DelDiv_hi()

◆ idInsertMonomial()

◆ makeQt()

◆ minimalMonomialGenSet()

◆ monCompare()

◆ p_Div_hi()

◆ positionInOrbit_FG_Case()

◆ positionInOrbit_IG_Case()

◆ positionInOrbitTruncationCase()

◆ RightColonOperation()

◆ scDegree()

◆ shiftInMon()

◆ sortMonoIdeal_pCompare()

◆ TwordMap()

Variable Documentation

◆ hilb_Qt

◆ hLength

◆ Q0

◆ Ql

◆ Qpol