Changeset 598f45b in git


Ignore:
Timestamp:
Jan 3, 2018, 10:49:57 AM (5 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'a657104b677b4c461d018cbf3204d72d34ad66a9')
Children:
fcbac7b9322bbd4185f67888844e361a5d59b607
Parents:
0dcc85c3cf9920d43b003b7f46b04bf39068415e
Message:
fix: system("hilbroune",I): memory leaks, format
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • Singular/extra.cc

    r0dcc85c r598f45b  
    30103010       if (strcmp(sys_cmd, "hilbroune") == 0)
    30113011       {
    3012          ideal I;
    30133012         if ((h!=NULL) && (h->Typ()==IDEAL_CMD))
    30143013         {
    3015            I=(ideal)h->CopyD();
    3016            slicehilb(I);
     3014           slicehilb((ideal)h->Data());
    30173015         }
    30183016         else return TRUE;
  • kernel/combinatorics/hilb.cc

    r0dcc85c r598f45b  
    265265//adds the new polynomial at the coresponding position
    266266//and simplifies the ideal
    267 static ideal SortByDeg_p(ideal I, poly p)
    268 {
    269     int i,j;
    270     if((I == NULL) || (idIs0(I)))
    271     {
    272         ideal res = idInit(1,1);
    273         res->m[0] = p;
    274         return(res);
    275     }
     267static void SortByDeg_p(ideal I, poly p)
     268{
     269  int i,j;
     270  if(idIs0(I))
     271  {
     272    I->m[0] = p;
     273    return;
     274  }
     275  idSkipZeroes(I);
     276  #if 1
     277  for(i = 0; (i<IDELEMS(I)) && (p_Totaldegree(I->m[i],currRing)<=p_Totaldegree(p,currRing)); i++)
     278  {
     279    if(p_DivisibleBy( I->m[i],p, currRing))
     280    {
     281      return;
     282    }
     283  }
     284  for(i = IDELEMS(I)-1; (i>=0) && (p_Totaldegree(I->m[i],currRing)>=p_Totaldegree(p,currRing)); i--)
     285  {
     286    if(p_DivisibleBy(p,I->m[i], currRing))
     287    {
     288      I->m[i] = NULL;
     289    }
     290  }
     291  if(idIs0(I))
     292  {
    276293    idSkipZeroes(I);
    277     #if 1
    278     for(i = 0; (i<IDELEMS(I)) && (p_Totaldegree(I->m[i],currRing)<=p_Totaldegree(p,currRing)); i++)
    279     {
    280         if(p_DivisibleBy( I->m[i],p, currRing))
    281         {
    282             return(I);
    283         }
    284     }
    285     for(i = IDELEMS(I)-1; (i>=0) && (p_Totaldegree(I->m[i],currRing)>=p_Totaldegree(p,currRing)); i--)
    286     {
    287         if(p_DivisibleBy(p,I->m[i], currRing))
    288         {
    289             I->m[i] = NULL;
    290         }
    291     }
    292     if(idIs0(I))
    293     {
    294         idSkipZeroes(I);
    295         I->m[0] = p;
    296         return(I);
    297     }
    298     #endif
     294    I->m[0] = p;
     295    return;
     296  }
     297  #endif
     298  idSkipZeroes(I);
     299  //First I take the case when all generators have the same degree
     300  if(p_Totaldegree(I->m[0],currRing) == p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
     301  {
     302    if(p_Totaldegree(p,currRing)<p_Totaldegree(I->m[0],currRing))
     303    {
     304      idInsertPoly(I,p);
     305      idSkipZeroes(I);
     306      for(i=IDELEMS(I)-1;i>=1; i--)
     307      {
     308        I->m[i] = I->m[i-1];
     309      }
     310      I->m[0] = p;
     311      return;
     312    }
     313    if(p_Totaldegree(p,currRing)>=p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
     314    {
     315      idInsertPoly(I,p);
     316      idSkipZeroes(I);
     317      return;
     318    }
     319  }
     320  if(p_Totaldegree(p,currRing)<=p_Totaldegree(I->m[0],currRing))
     321  {
     322    idInsertPoly(I,p);
    299323    idSkipZeroes(I);
    300     //First I take the case when all generators have the same degree
    301     if(p_Totaldegree(I->m[0],currRing) == p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
    302     {
    303         if(p_Totaldegree(p,currRing)<p_Totaldegree(I->m[0],currRing))
    304         {
    305             idInsertPoly(I,p);
    306             idSkipZeroes(I);
    307             for(i=IDELEMS(I)-1;i>=1; i--)
    308             {
    309                 I->m[i] = I->m[i-1];
    310             }
    311             I->m[0] = p;
    312             return(I);
    313         }
    314         if(p_Totaldegree(p,currRing)>=p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
    315         {
    316             idInsertPoly(I,p);
    317             idSkipZeroes(I);
    318             return(I);
    319         }
    320     }
    321     if(p_Totaldegree(p,currRing)<=p_Totaldegree(I->m[0],currRing))
    322     {
    323         idInsertPoly(I,p);
    324         idSkipZeroes(I);
    325         for(i=IDELEMS(I)-1;i>=1; i--)
    326         {
    327             I->m[i] = I->m[i-1];
    328         }
    329         I->m[0] = p;
    330         return(I);
    331     }
    332     if(p_Totaldegree(p,currRing)>=p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
    333     {
    334         idInsertPoly(I,p);
    335         idSkipZeroes(I);
    336         return(I);
    337     }
    338     for(i = IDELEMS(I)-2; ;)
    339     {
    340         if(p_Totaldegree(p,currRing)==p_Totaldegree(I->m[i],currRing))
    341         {
    342             idInsertPoly(I,p);
    343             idSkipZeroes(I);
    344             for(j = IDELEMS(I)-1; j>=i+1;j--)
    345             {
    346                 I->m[j] = I->m[j-1];
    347             }
    348             I->m[i] = p;
    349             return(I);
    350         }
    351         if(p_Totaldegree(p,currRing)>p_Totaldegree(I->m[i],currRing))
    352         {
    353             idInsertPoly(I,p);
    354             idSkipZeroes(I);
    355             for(j = IDELEMS(I)-1; j>=i+2;j--)
    356             {
    357                 I->m[j] = I->m[j-1];
    358             }
    359             I->m[i+1] = p;
    360             return(I);
    361         }
    362         i--;
    363     }
     324    for(i=IDELEMS(I)-1;i>=1; i--)
     325    {
     326      I->m[i] = I->m[i-1];
     327    }
     328    I->m[0] = p;
     329    return;
     330  }
     331  if(p_Totaldegree(p,currRing)>=p_Totaldegree(I->m[IDELEMS(I)-1],currRing))
     332  {
     333    idInsertPoly(I,p);
     334    idSkipZeroes(I);
     335    return;
     336  }
     337  for(i = IDELEMS(I)-2; ;)
     338  {
     339    if(p_Totaldegree(p,currRing)==p_Totaldegree(I->m[i],currRing))
     340    {
     341      idInsertPoly(I,p);
     342      idSkipZeroes(I);
     343      for(j = IDELEMS(I)-1; j>=i+1;j--)
     344      {
     345        I->m[j] = I->m[j-1];
     346      }
     347      I->m[i] = p;
     348      return;
     349    }
     350    if(p_Totaldegree(p,currRing)>p_Totaldegree(I->m[i],currRing))
     351    {
     352      idInsertPoly(I,p);
     353      idSkipZeroes(I);
     354      for(j = IDELEMS(I)-1; j>=i+2;j--)
     355      {
     356        I->m[j] = I->m[j-1];
     357      }
     358      I->m[i+1] = p;
     359      return;
     360    }
     361    i--;
     362  }
    364363}
    365364
     
    367366static ideal SortByDeg(ideal I)
    368367{
    369     if(idIs0(I))
    370     {
    371         return(I);
    372     }
    373     int i;
    374     ideal res;
    375     idSkipZeroes(I);
    376     res = idInit(1,1);
    377     res->m[0] = poly(0);
    378     for(i = 0; i<=IDELEMS(I)-1;i++)
    379     {
    380         res = SortByDeg_p(res, I->m[i]);
    381     }
    382     idSkipZeroes(res);
    383     //idDegSortTest(res);
    384     return(res);
     368  if(idIs0(I))
     369  {
     370    return id_Copy(I,currRing);
     371  }
     372  int i;
     373  ideal res;
     374  idSkipZeroes(I);
     375  res = idInit(1,1);
     376  res->m[0] = poly(0);
     377  for(i = 0; i<=IDELEMS(I)-1;i++)
     378  {
     379    SortByDeg_p(res, I->m[i]);
     380  }
     381  idSkipZeroes(res);
     382  //idDegSortTest(res);
     383  return(res);
    385384}
    386385
     
    390389    if(idIs0(Iorig))
    391390    {
    392         ideal res = idInit(1,1);
    393         res->m[0] = poly(0);
    394         return(res);
     391      ideal res = idInit(1,1);
     392      res->m[0] = poly(0);
     393      return(res);
    395394    }
    396395    if(idIs0(p))
    397396    {
    398         ideal res = idInit(1,1);
    399         res->m[0] = pOne();
    400         return(res);
    401     }
    402     ideal I = idCopy(Iorig);
     397      ideal res = idInit(1,1);
     398      res->m[0] = pOne();
     399      return(res);
     400    }
     401    ideal I = id_Head(Iorig,currRing);
    403402    ideal res = idInit(IDELEMS(I),1);
    404403    int i,j;
     
    406405    for(i = 0; i<IDELEMS(I); i++)
    407406    {
    408         res->m[i] = p_Copy(I->m[i], currRing);
    409         for(j = 1; (j<=currRing->N) ; j++)
    410         {
    411             dummy = p_GetExp(p->m[0], j, currRing);
    412             if(dummy > 0)
    413             {
    414                 if(p_GetExp(I->m[i], j, currRing) < dummy)
    415                 {
    416                     p_SetExp(res->m[i], j, 0, currRing);
    417                 }
    418                 else
    419                 {
    420                     p_SetExp(res->m[i], j, p_GetExp(I->m[i], j, currRing) - dummy, currRing);
    421                 }
    422             }
    423         }
    424         p_Setm(res->m[i], currRing);
    425         if(p_Totaldegree(res->m[i],currRing) == p_Totaldegree(I->m[i],currRing))
    426         {
    427             res->m[i] = NULL; // pDelete
    428         }
    429         else
    430         {
    431             I->m[i] = NULL; // pDelete
    432         }
     407      res->m[i] = p_Head(I->m[i], currRing);
     408      for(j = 1; (j<=currRing->N) ; j++)
     409      {
     410        dummy = p_GetExp(p->m[0], j, currRing);
     411        if(dummy > 0)
     412        {
     413          if(p_GetExp(I->m[i], j, currRing) < dummy)
     414          {
     415            p_SetExp(res->m[i], j, 0, currRing);
     416          }
     417          else
     418          {
     419            p_SetExp(res->m[i], j, p_GetExp(I->m[i], j, currRing) - dummy, currRing);
     420          }
     421        }
     422      }
     423      p_Setm(res->m[i], currRing);
     424      if(p_Totaldegree(res->m[i],currRing) == p_Totaldegree(I->m[i],currRing))
     425      {
     426        p_Delete(&res->m[i],currRing);
     427      }
     428      else
     429      {
     430        p_Delete(&I->m[i],currRing);
     431      }
    433432    }
    434433    idSkipZeroes(res);
     
    438437      for(i = 0; i<=IDELEMS(res)-1; i++)
    439438      {
    440         I = SortByDeg_p(I,res->m[i]);
    441       }
    442     }
     439        SortByDeg_p(I,res->m[i]);
     440      }
     441    }
     442    id_Delete(&res,currRing);
    443443    //idDegSortTest(I);
    444444    return(I);
     
    446446
    447447//id_Add for monomials
    448 static ideal idAddMon(ideal I, ideal p)
    449 {
    450     #if 1
    451     I = SortByDeg_p(I,p->m[0]);
    452     #else
    453     I = id_Add(I,p,currRing);
    454     #endif
    455     //idSkipZeroes(I);
    456     return(I);
     448static void idAddMon(ideal I, ideal p)
     449{
     450  SortByDeg_p(I,p->m[0]);
     451  //idSkipZeroes(I);
    457452}
    458453
     
    686681static poly ChoosePJL(ideal I)
    687682{
     683  int i,j,dummy;
     684  bool flag = TRUE;
     685  poly m = p_ISet(1,currRing);
     686  for(i = IDELEMS(I)-1;(i>=0) && (flag);i--)
     687  {
     688    flag = TRUE;
     689    for(j=1;(j<=currRing->N) && (flag);j++)
     690    {
     691      dummy = p_GetExp(I->m[i],j,currRing);
     692      if(dummy >= 2)
     693      {
     694        p_SetExp(m,j,dummy-1,currRing);
     695        p_Setm(m,currRing);
     696        flag = FALSE;
     697      }
     698    }
     699    if(!p_IsOne(m, currRing))
     700    {
     701      return(m);
     702    }
     703  }
     704  p_Delete(&m,currRing);
     705  m = ChoosePVar(I);
     706  return(m);
     707}
     708
     709#if 0
     710//choice JF: last entry just variable with power -1 (xy10z15 -> y9)
     711static poly ChoosePJF(ideal I)
     712{
    688713    int i,j,dummy;
    689714    bool flag = TRUE;
    690715    poly m = p_ISet(1,currRing);
    691     for(i = IDELEMS(I)-1;(i>=0) && (flag);i--)
     716    for(i = 0;(i<=IDELEMS(I)-1) && (flag);i++)
    692717    {
    693718        flag = TRUE;
     
    710735    return(m);
    711736}
    712 
    713 #if 0
    714 //choice JF: last entry just variable with power -1 (xy10z15 -> y9)
    715 static poly ChoosePJF(ideal I)
    716 {
    717     int i,j,dummy;
    718     bool flag = TRUE;
    719     poly m = p_ISet(1,currRing);
    720     for(i = 0;(i<=IDELEMS(I)-1) && (flag);i++)
    721     {
    722         flag = TRUE;
    723         for(j=1;(j<=currRing->N) && (flag);j++)
    724         {
    725             dummy = p_GetExp(I->m[i],j,currRing);
    726             if(dummy >= 2)
    727             {
    728                 p_SetExp(m,j,dummy-1,currRing);
    729                 p_Setm(m,currRing);
    730                 flag = FALSE;
    731             }
    732         }
    733         if(!p_IsOne(m, currRing))
    734         {
    735             return(m);
    736         }
    737     }
    738     m = ChoosePVar(I);
    739     return(m);
    740 }
    741737#endif
    742738
     
    744740static poly ChooseP(ideal I)
    745741{
    746     poly m;
    747     //  TEST TO SEE WHICH ONE IS BETTER
    748     //m = ChoosePXL(I);
    749     //m = ChoosePXF(I);
    750     //m = ChoosePOL(I);
    751     //m = ChoosePOF(I);
    752     //m = ChoosePVL(I);
    753     //m = ChoosePVF(I);
    754     m = ChoosePJL(I);
    755     //m = ChoosePJF(I);
    756     return(m);
     742  poly m;
     743  //  TEST TO SEE WHICH ONE IS BETTER
     744  //m = ChoosePXL(I);
     745  //m = ChoosePXF(I);
     746  //m = ChoosePOL(I);
     747  //m = ChoosePOF(I);
     748  //m = ChoosePVL(I);
     749  //m = ChoosePVF(I);
     750  m = ChoosePJL(I);
     751  //m = ChoosePJF(I);
     752  return(m);
    757753}
    758754
     
    822818    if(JustVar(I) == TRUE)
    823819    {
    824         if(IDELEMS(I) == variables)
    825         {
    826             mpz_init(dummy);
    827             if((variables % 2) == 0)
    828                 {mpz_set_ui(dummy, 1);}
    829             else
    830                 {mpz_set_si(dummy, -1);}
    831             mpz_add(ec, ec, dummy);
    832         }
    833         //mpz_clear(dummy);
    834         return;
     820      if(IDELEMS(I) == variables)
     821      {
     822        mpz_init(dummy);
     823        if((variables % 2) == 0)
     824          mpz_set_ui(dummy, 1);
     825        else
     826          mpz_set_si(dummy, -1);
     827        mpz_add(ec, ec, dummy);
     828        mpz_clear(dummy);
     829      }
     830      return;
    835831    }
    836832    ideal p = idInit(1,1);
     
    845841    for(i = 1;i<=currRing->N;i++)
    846842    {
    847         if(p_GetExp(p->m[0],i,currRing)>0)
    848         {
    849             howmanyvarinp++;
    850         }
     843      if(p_GetExp(p->m[0],i,currRing)>0)
     844      {
     845        howmanyvarinp++;
     846      }
    851847    }
    852848    eulerchar(Ip, variables-howmanyvarinp, ec);
    853849    id_Delete(&Ip, currRing);
    854     I = idAddMon(I,p);
     850    idAddMon(I,p);
     851    id_Delete(&p, currRing);
    855852  }
    856853}
     
    888885static bool IsIn(poly p, ideal I)
    889886{
    890     //assumes that I is ordered by degree
    891     if(idIs0(I))
    892     {
    893         if(p==poly(0))
    894         {
    895             return(TRUE);
    896         }
    897         else
    898         {
    899             return(FALSE);
    900         }
    901     }
     887  //assumes that I is ordered by degree
     888  if(idIs0(I))
     889  {
    902890    if(p==poly(0))
    903891    {
    904         return(FALSE);
    905     }
    906     int i,j;
    907     bool flag;
    908     for(i = 0;i<IDELEMS(I);i++)
    909     {
    910         flag = TRUE;
    911         for(j = 1;(j<=currRing->N) &&(flag);j++)
    912         {
    913             if(p_GetExp(p, j, currRing)<p_GetExp(I->m[i], j, currRing))
    914             {
    915                 flag = FALSE;
    916             }
    917         }
    918         if(flag)
    919         {
    920             return(TRUE);
    921         }
    922     }
     892      return(TRUE);
     893    }
     894    else
     895    {
     896      return(FALSE);
     897    }
     898  }
     899  if(p==poly(0))
     900  {
    923901    return(FALSE);
     902  }
     903  int i,j;
     904  bool flag;
     905  for(i = 0;i<IDELEMS(I);i++)
     906  {
     907    flag = TRUE;
     908    for(j = 1;(j<=currRing->N) &&(flag);j++)
     909    {
     910      if(p_GetExp(p, j, currRing)<p_GetExp(I->m[i], j, currRing))
     911      {
     912        flag = FALSE;
     913      }
     914    }
     915    if(flag)
     916    {
     917      return(TRUE);
     918    }
     919  }
     920  return(FALSE);
    924921}
    925922
     
    927924static poly LCMmon(ideal I)
    928925{
    929     if(idIs0(I))
    930     {
    931         return(NULL);
    932     }
    933     poly m;
    934     int dummy,i,j;
    935     m = p_ISet(1,currRing);
    936     for(i=1;i<=currRing->N;i++)
    937     {
    938         dummy=0;
    939         for(j=IDELEMS(I)-1;j>=0;j--)
    940         {
    941             if(p_GetExp(I->m[j],i,currRing) > dummy)
    942             {
    943                 dummy = p_GetExp(I->m[j],i,currRing);
    944             }
    945         }
    946         p_SetExp(m,i,dummy,currRing);
    947     }
    948     p_Setm(m,currRing);
    949     return(m);
     926  if(idIs0(I))
     927  {
     928    return(NULL);
     929  }
     930  poly m;
     931  int dummy,i,j;
     932  m = p_ISet(1,currRing);
     933  for(i=1;i<=currRing->N;i++)
     934  {
     935    dummy=0;
     936    for(j=IDELEMS(I)-1;j>=0;j--)
     937    {
     938      if(p_GetExp(I->m[j],i,currRing) > dummy)
     939      {
     940        dummy = p_GetExp(I->m[j],i,currRing);
     941      }
     942    }
     943    p_SetExp(m,i,dummy,currRing);
     944  }
     945  p_Setm(m,currRing);
     946  return(m);
    950947}
    951948
     
    964961    for(i=IDELEMS(S)-1;i>=0;i--)
    965962    {
    966         if(IsIn(S->m[i],I))
    967         {
    968             S->m[i]=NULL;
    969             prune++;
    970         }
     963      if(IsIn(S->m[i],I))
     964      {
     965        p_Delete(&S->m[i],currRing);
     966        prune++;
     967      }
    971968    }
    972969    idSkipZeroes(S);
     
    974971    for(i=IDELEMS(I)-1;i>=0;i--)
    975972    {
    976         m = p_Copy(I->m[i],currRing);
    977         for(j=1;j<=currRing->N;j++)
    978         {
    979             dummy = p_GetExp(m,j,currRing);
    980             if(dummy > 0)
    981             {
    982                 p_SetExp(m,j,dummy-1,currRing);
    983             }
    984         }
    985         p_Setm(m, currRing);
    986         if(IsIn(m,S))
    987         {
    988             I->m[i]=NULL;
    989             //printf("\n Deleted, since pi(m) is in S\n");pWrite(m);
    990         }
     973      m = p_Head(I->m[i],currRing);
     974      for(j=1;j<=currRing->N;j++)
     975      {
     976        dummy = p_GetExp(m,j,currRing);
     977        if(dummy > 0)
     978        {
     979          p_SetExp(m,j,dummy-1,currRing);
     980        }
     981      }
     982      p_Setm(m, currRing);
     983      if(IsIn(m,S))
     984      {
     985        p_Delete(&I->m[i],currRing);
     986        //printf("\n Deleted, since pi(m) is in S\n");pWrite(m);
     987      }
     988      p_Delete(&m,currRing);
    991989    }
    992990    idSkipZeroes(I);
     
    995993    if(m != NULL)
    996994    {
    997         for(i=0;i<IDELEMS(S);i++)
    998         {
    999             if(!(p_DivisibleBy(S->m[i], m, currRing)))
    1000             {
    1001                 S->m[i] = NULL;
    1002                 j++;
    1003                 moreprune++;
    1004             }
    1005             else
    1006             {
    1007                 if(pLmEqual(S->m[i],m))
    1008                 {
    1009                     S->m[i] = NULL;
    1010                     moreprune++;
    1011                 }
    1012             }
    1013         }
    1014     idSkipZeroes(S);
    1015     }
     995      for(i=0;i<IDELEMS(S);i++)
     996      {
     997        if(!(p_DivisibleBy(S->m[i], m, currRing)))
     998        {
     999          S->m[i] = NULL;
     1000          j++;
     1001          moreprune++;
     1002        }
     1003        else
     1004        {
     1005          if(pLmEqual(S->m[i],m))
     1006          {
     1007            S->m[i] = NULL;
     1008            moreprune++;
     1009          }
     1010        }
     1011      }
     1012      idSkipZeroes(S);
     1013    }
     1014    p_Delete(&m,currRing);
    10161015    /*printf("\n---------------------------\n");
    10171016    printf("\n      I\n");idPrint(I);
     
    10221021    if(idIs0(I))
    10231022    {
    1024         id_Delete(&I, currRing);
    1025         id_Delete(&S, currRing);
    1026         p_Delete(&m, currRing);
    1027         break;
     1023      id_Delete(&I, currRing);
     1024      id_Delete(&S, currRing);
     1025      break;
    10281026    }
    10291027    m = LCMmon(I);
    10301028    if(!p_DivisibleBy(x,m, currRing))
    10311029    {
    1032         //printf("\nx does not divide lcm(I)");
    1033         //printf("\nEmpty set");pWrite(q);
    1034         id_Delete(&I, currRing);
    1035         id_Delete(&S, currRing);
    1036         p_Delete(&m, currRing);
    1037         break;
    1038     }
     1030      //printf("\nx does not divide lcm(I)");
     1031      //printf("\nEmpty set");pWrite(q);
     1032      id_Delete(&I, currRing);
     1033      id_Delete(&S, currRing);
     1034      p_Delete(&m, currRing);
     1035      break;
     1036    }
     1037    p_Delete(&m, currRing);
    10391038    m = SqFree(I);
    10401039    if(m==NULL)
    10411040    {
    1042         //printf("\n      Corner: ");
    1043         //pWrite(q);
    1044         //printf("\n      With the facets of the dual simplex:\n");
    1045         //idPrint(I);
    1046         mpz_t ec;
    1047         mpz_init(ec);
    1048         mpz_ptr ec_ptr = ec;
    1049         eulerchar(I, currRing->N, ec_ptr);
    1050         bool flag = FALSE;
    1051         if(NNN==0)
    1052             {
    1053                 hilbertcoef = (mpz_ptr)omAlloc((NNN+1)*sizeof(mpz_t));
    1054                 hilbpower = (int*)omAlloc((NNN+1)*sizeof(int));
    1055                 mpz_init( &hilbertcoef[NNN]);
    1056                 mpz_set(  &hilbertcoef[NNN], ec);
    1057                 mpz_clear(ec);
    1058                 hilbpower[NNN] = p_Totaldegree(q,currRing);
    1059                 NNN++;
    1060             }
    1061         else
    1062         {
    1063             //I look if the power appears already
    1064             for(i = 0;(i<NNN)&&(flag == FALSE)&&(p_Totaldegree(q,currRing)>=hilbpower[i]);i++)
    1065             {
    1066                 if((hilbpower[i]) == (p_Totaldegree(q,currRing)))
    1067                 {
    1068                     flag = TRUE;
    1069                     mpz_add(&hilbertcoef[i],&hilbertcoef[i],ec_ptr);
    1070                 }
    1071             }
    1072             if(flag == FALSE)
    1073             {
    1074                 hilbertcoef = (mpz_ptr)omRealloc(hilbertcoef, (NNN+1)*sizeof(mpz_t));
    1075                 hilbpower = (int*)omRealloc(hilbpower, (NNN+1)*sizeof(int));
    1076                 mpz_init(&hilbertcoef[NNN]);
    1077                 for(j = NNN; j>i; j--)
    1078                 {
    1079                     mpz_set(&hilbertcoef[j],&hilbertcoef[j-1]);
    1080                     hilbpower[j] = hilbpower[j-1];
    1081                 }
    1082                 mpz_set(  &hilbertcoef[i], ec);
    1083                 mpz_clear(ec);
    1084                 hilbpower[i] = p_Totaldegree(q,currRing);
    1085                 NNN++;
    1086             }
    1087         }
    1088         break;
    1089     }
     1041      //printf("\n      Corner: ");
     1042      //pWrite(q);
     1043      //printf("\n      With the facets of the dual simplex:\n");
     1044      //idPrint(I);
     1045      mpz_t ec;
     1046      mpz_init(ec);
     1047      mpz_ptr ec_ptr = ec;
     1048      eulerchar(I, currRing->N, ec_ptr);
     1049      bool flag = FALSE;
     1050      if(NNN==0)
     1051      {
     1052        hilbertcoef = (mpz_ptr)omAlloc((NNN+1)*sizeof(mpz_t));
     1053        hilbpower = (int*)omAlloc((NNN+1)*sizeof(int));
     1054        mpz_init_set( &hilbertcoef[NNN], ec);
     1055        hilbpower[NNN] = p_Totaldegree(q,currRing);
     1056        NNN++;
     1057      }
     1058      else
     1059      {
     1060        //I look if the power appears already
     1061        for(i = 0;(i<NNN)&&(flag == FALSE)&&(p_Totaldegree(q,currRing)>=hilbpower[i]);i++)
     1062        {
     1063          if((hilbpower[i]) == (p_Totaldegree(q,currRing)))
     1064          {
     1065            flag = TRUE;
     1066            mpz_add(&hilbertcoef[i],&hilbertcoef[i],ec_ptr);
     1067          }
     1068        }
     1069        if(flag == FALSE)
     1070        {
     1071          hilbertcoef = (mpz_ptr)omRealloc(hilbertcoef, (NNN+1)*sizeof(mpz_t));
     1072          hilbpower = (int*)omRealloc(hilbpower, (NNN+1)*sizeof(int));
     1073          mpz_init(&hilbertcoef[NNN]);
     1074          for(j = NNN; j>i; j--)
     1075          {
     1076            mpz_set(&hilbertcoef[j],&hilbertcoef[j-1]);
     1077            hilbpower[j] = hilbpower[j-1];
     1078          }
     1079          mpz_set(  &hilbertcoef[i], ec);
     1080          hilbpower[i] = p_Totaldegree(q,currRing);
     1081          NNN++;
     1082        }
     1083      }
     1084      mpz_clear(ec);
     1085      id_Delete(&I, currRing);
     1086      id_Delete(&S, currRing);
     1087      break;
     1088    }
     1089    else
     1090      p_Delete(&m, currRing);
    10901091    m = ChooseP(I);
    10911092    p = idInit(1,1);
     
    10951096    poly pq = pp_Mult_mm(q,m,currRing);
    10961097    rouneslice(Ip, Sp, pq, x, prune, moreprune, steps, NNN, hilbertcoef,hilbpower);
    1097     //id_Delete(&Ip, currRing);
    1098     //id_Delete(&Sp, currRing);
    1099     S = idAddMon(S,p);
     1098    idAddMon(S,p);
    11001099    p->m[0]=NULL;
    11011100    id_Delete(&p, currRing); // p->m[0] was also in S
     
    11131112    int *hilbpower;
    11141113    ideal S = idInit(1,1);
    1115     poly q = p_ISet(1,currRing);
     1114    poly q = p_One(currRing);
    11161115    ideal X = idInit(1,1);
    11171116    X->m[0]=p_One(currRing);
    11181117    for(i=1;i<=currRing->N;i++)
    11191118    {
    1120             p_SetExp(X->m[0],i,1,currRing);
     1119      p_SetExp(X->m[0],i,1,currRing);
    11211120    }
    11221121    p_Setm(X->m[0],currRing);
    11231122    I = id_Mult(I,X,currRing);
    1124     I = SortByDeg(I);
     1123    ideal Itmp = SortByDeg(I);
     1124    id_Delete(&I,currRing);
     1125    I = Itmp;
    11251126    //printf("\n-------------RouneSlice--------------\n");
    11261127    rouneslice(I,S,q,X->m[0],prune, moreprune, steps, NNN, hilbertcoef, hilbpower);
     1128    id_Delete(&X,currRing);
     1129    p_Delete(&q,currRing);
    11271130    //printf("\nIn total Prune got rid of %i elements\n",prune);
    11281131    //printf("\nIn total More Prune got rid of %i elements\n",moreprune);
    11291132    //printf("\nSteps of rouneslice: %i\n\n", steps);
    1130     mpz_t coefhilb;
    1131     mpz_t dummy;
    1132     mpz_init(coefhilb);
    1133     mpz_init(dummy);
    11341133    printf("\n//  %8d t^0",1);
    11351134    for(i = 0; i<NNN; i++)
    11361135    {
    1137         if(mpz_sgn(&hilbertcoef[i])!=0)
    1138         {
    1139             gmp_printf("\n//  %8Zd t^%d",&hilbertcoef[i],hilbpower[i]);
    1140         }
    1141     }
     1136      if(mpz_sgn(&hilbertcoef[i])!=0)
     1137      {
     1138        gmp_printf("\n//  %8Zd t^%d",&hilbertcoef[i],hilbpower[i]);
     1139      }
     1140    }
     1141    PrintLn();
    11421142    omFreeSize(hilbertcoef, (NNN)*sizeof(mpz_t));
    11431143    omFreeSize(hilbpower, (NNN)*sizeof(int));
Note: See TracChangeset for help on using the changeset viewer.