Changeset 1d85871 in git


Ignore:
Timestamp:
Feb 18, 2015, 9:41:53 AM (9 years ago)
Author:
Yue Ren <ren@…>
Branches:
(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', '38dfc5131670d387a89455159ed1e071997eec94')
Children:
9f7a9cd448041b035d092b0dcf11b4d5f14ac002
Parents:
4a6847d88eb8e6f833551e2a34563ba917cb5c52
Message:
new: functions for groebnerFans, groebnerComplexes
Location:
Singular/dyn_modules/gfanlib
Files:
16 edited

Legend:

Unmodified
Added
Removed

  • Singular/dyn_modules/gfanlib/Makefile.am

    r4a6847 r1d85871  
    11ACLOCAL_AMFLAGS = -I ../../m4
    22
    3 SOURCES = singularWishlist.h gfanlib_exceptions.h callgfanlib_conversion.cc callgfanlib_conversion.h bbcone.cc bbcone.h bbfan.cc bbfan.h bbpolytope.cc bbpolytope.h gfan.h gitfan.cc gitfan.h std_wrapper.cc std_wrapper.h tropicalVarietyOfPolynomials.h tropicalVarietyOfPolynomials.cc ppinitialReduction.cc ppinitialReduction.h containsMonomial.cc containsMonomial.h adjustWeights.cc adjustWeights.h tropicalStrategy.cc tropicalStrategy.h initial.cc initial.h witness.cc witness.h lift.cc lift.h flip.cc flip.h tropicalCurves.cc tropicalCurves.h groebnerCone.cc groebnerCone.h startingCone.cc startingCone.h tropicalTraversal.cc tropicalTraversal.h tropicalVarietyOfIdeals.cc tropicalVarietyOfIdeals.h tropicalVariety.cc tropicalVariety.h tropical.cc tropical.h gfanlib.cc
     3SOURCES = singularWishlist.h gfanlib_exceptions.h callgfanlib_conversion.cc callgfanlib_conversion.h bbcone.cc bbcone.h bbfan.cc bbfan.h bbpolytope.cc bbpolytope.h gfan.h gitfan.cc gitfan.h std_wrapper.cc std_wrapper.h tropicalVarietyOfPolynomials.h tropicalVarietyOfPolynomials.cc ppinitialReduction.cc ppinitialReduction.h containsMonomial.cc containsMonomial.h adjustWeights.cc adjustWeights.h tropicalStrategy.cc tropicalStrategy.h initial.cc initial.h witness.cc witness.h lift.cc lift.h flip.cc flip.h tropicalCurves.cc tropicalCurves.h groebnerCone.cc groebnerCone.h startingCone.cc startingCone.h tropicalTraversal.cc tropicalTraversal.h tropicalVarietyOfIdeals.cc tropicalVarietyOfIdeals.h tropicalVariety.cc tropicalVariety.h groebnerFan.cc groebnerFan.h groebnerComplex.cc groebnerComplex.h tropical.cc tropical.h gfanlib.cc
    44
    55MY_CPPFLAGS = -I${srcdir} -I${top_srcdir} -I${top_builddir} \

  • Singular/dyn_modules/gfanlib/bbcone.cc

    r4a6847 r1d85871  
    16661666
    16671667  return relativeInteriorPoints;
     1668}
     1669
     1670
     1671/***
     1672 * Computes a relative interior point and an outer normal vector for each facet of zc
     1673 **/
     1674std::pair<gfan::ZMatrix,gfan::ZMatrix> interiorPointsAndNormalsOfFacets(const gfan::ZCone zc, const std::set<gfan::ZVector> &exceptThesePoints, const bool onlyLowerHalfSpace)
     1675{
     1676  gfan::ZMatrix inequalities = zc.getFacets();
     1677  gfan::ZMatrix equations = zc.getImpliedEquations();
     1678  int r = inequalities.getHeight();
     1679  int c = inequalities.getWidth();
     1680
     1681  /* our cone has r facets, if r==0 return empty matrices */
     1682  gfan::ZMatrix relativeInteriorPoints = gfan::ZMatrix(0,c);
     1683  gfan::ZMatrix outerFacetNormals = gfan::ZMatrix(0,c);
     1684  if (r==0)
     1685    return std::make_pair(relativeInteriorPoints,outerFacetNormals);
     1686
     1687  /* next we iterate over each of the r facets,
     1688   * build the respective cone and add it to the list
     1689   * this is the i=0 case */
     1690  gfan::ZMatrix newInequalities = inequalities.submatrix(1,0,r,c);
     1691  gfan::ZMatrix newEquations = equations;
     1692  newEquations.appendRow(inequalities[0]);
     1693  gfan::ZCone facet = gfan::ZCone(newInequalities,newEquations);
     1694  gfan::ZVector interiorPoint = facet.getRelativeInteriorPoint();
     1695  if (onlyLowerHalfSpace==false || interiorPoint[0].sign()<0)
     1696  {
     1697    if (exceptThesePoints.count(interiorPoint)==0)
     1698    {
     1699      relativeInteriorPoints.appendRow(interiorPoint);
     1700      outerFacetNormals.appendRow(-inequalities[0]);
     1701    }
     1702  }
     1703
     1704  /* these are the cases i=1,...,r-2 */
     1705  for (int i=1; i<r-1; i++)
     1706  {
     1707    newInequalities = inequalities.submatrix(0,0,i,c);
     1708    newInequalities.append(inequalities.submatrix(i+1,0,r,c));
     1709    newEquations = equations;
     1710    newEquations.appendRow(inequalities[i]);
     1711    facet = gfan::ZCone(newInequalities,newEquations);
     1712    interiorPoint = facet.getRelativeInteriorPoint();
     1713    if (onlyLowerHalfSpace==false || interiorPoint[0].sign()<0)
     1714    {
     1715      if (exceptThesePoints.count(interiorPoint)==0)
     1716      {
     1717        relativeInteriorPoints.appendRow(interiorPoint);
     1718        outerFacetNormals.appendRow(-inequalities[i]);
     1719      }
     1720    }
     1721  }
     1722
     1723  /* this is the i=r-1 case */
     1724  newInequalities = inequalities.submatrix(0,0,r-1,c);
     1725  newEquations = equations;
     1726  newEquations.appendRow(inequalities[r-1]);
     1727  facet = gfan::ZCone(newInequalities,newEquations);
     1728  interiorPoint = facet.getRelativeInteriorPoint();
     1729  if (onlyLowerHalfSpace==false || interiorPoint[0].sign()<0)
     1730  {
     1731    if (exceptThesePoints.count(interiorPoint)==0)
     1732    {
     1733      relativeInteriorPoints.appendRow(interiorPoint);
     1734      outerFacetNormals.appendRow(-inequalities[r-1]);
     1735    }
     1736  }
     1737
     1738  return std::make_pair(relativeInteriorPoints,outerFacetNormals);
    16681739}
    16691740

  • Singular/dyn_modules/gfanlib/bbcone.h

    r4a6847 r1d85871  
    2121gfan::ZCone liftUp(const gfan::ZCone &zc);
    2222gfan::ZMatrix interiorPointsOfFacets(const gfan::ZCone &zc, const std::set<gfan::ZVector> &exceptThese=std::set<gfan::ZVector>());
     23std::pair<gfan::ZMatrix,gfan::ZMatrix> interiorPointsAndNormalsOfFacets(const gfan::ZCone zc, const std::set<gfan::ZVector> &exceptThesePoints=std::set<gfan::ZVector>(), const bool onlyLowerHalfSpace=false);
    2324
    2425#endif

  • Singular/dyn_modules/gfanlib/groebnerCone.cc

    r4a6847 r1d85871  
    444444
    445445/***
    446  * Computes a relative interior point and an outer normal vector for each facet of zc
    447  **/
    448 static std::pair<gfan::ZMatrix,gfan::ZMatrix> interiorPointsAndNormalsOfFacets(const gfan::ZCone zc)
    449 {
    450   gfan::ZMatrix inequalities = zc.getFacets();
    451   gfan::ZMatrix equations = zc.getImpliedEquations();
    452   int r = inequalities.getHeight();
    453   int c = inequalities.getWidth();
    454 
    455   /* our cone has r facets, if r==0 return empty matrices */
    456   gfan::ZMatrix relativeInteriorPoints = gfan::ZMatrix(0,c);
    457   gfan::ZMatrix outerFacetNormals = gfan::ZMatrix(0,c);
    458   if (r==0)
    459     return std::make_pair(relativeInteriorPoints,outerFacetNormals);
    460 
    461   /* next we iterate over each of the r facets,
    462    * build the respective cone and add it to the list
    463    * this is the i=0 case */
    464   gfan::ZMatrix newInequalities = inequalities.submatrix(1,0,r,c);
    465   gfan::ZMatrix newEquations = equations;
    466   newEquations.appendRow(inequalities[0]);
    467   gfan::ZCone facet = gfan::ZCone(newInequalities,newEquations);
    468   relativeInteriorPoints.appendRow(facet.getRelativeInteriorPoint());
    469   assume(zc.contains(relativeInteriorPoints[0]) && !zc.containsRelatively(relativeInteriorPoints[0]));
    470   outerFacetNormals.appendRow(-inequalities[0]);
    471 
    472   /* these are the cases i=1,...,r-2 */
    473   for (int i=1; i<r-1; i++)
    474   {
    475     newInequalities = inequalities.submatrix(0,0,i,c);
    476     newInequalities.append(inequalities.submatrix(i+1,0,r,c));
    477     newEquations = equations;
    478     newEquations.appendRow(inequalities[i]);
    479     facet = gfan::ZCone(newInequalities,newEquations);
    480     relativeInteriorPoints.appendRow(facet.getRelativeInteriorPoint());
    481     assume(zc.contains(relativeInteriorPoints[i]) && !zc.containsRelatively(relativeInteriorPoints[i]));
    482     outerFacetNormals.appendRow(-inequalities[i]);
    483   }
    484 
    485   /* this is the i=r-1 case */
    486   newInequalities = inequalities.submatrix(0,0,r-1,c);
    487   newEquations = equations;
    488   newEquations.appendRow(inequalities[r-1]);
    489   facet = gfan::ZCone(newInequalities,newEquations);
    490   relativeInteriorPoints.appendRow(facet.getRelativeInteriorPoint());
    491   assume(zc.contains(relativeInteriorPoints[r-1]) && !zc.containsRelatively(relativeInteriorPoints[r-1]));
    492   outerFacetNormals.appendRow(-inequalities[r-1]);
    493 
    494   return std::make_pair(relativeInteriorPoints,outerFacetNormals);
    495 }
    496 
    497 
    498 /***
    499446 * Returns a complete list of neighboring Groebner cones.
    500447 **/
    … …  
    565512    {
    566513      gfan::ZCone zc = sigma->getPolyhedralCone();
    567       assume(isCompatible(zf,&zc));
     514      // assume(isCompatible(zf,&zc));
    568515      zf->insert(zc);
    569516    }

  • Singular/dyn_modules/gfanlib/initial.cc

    r4a6847 r1d85871  
    1212  {
    1313    if (!w[i].fitsInInt())
     14    {
     15      WerrorS("wDeg: overflow in weight vector");
    1416      throw 0; // weightOverflow;
     17    }
    1518    d += p_GetExp(p,i+1,r)*w[i].toInt();
    1619  }

  • Singular/dyn_modules/gfanlib/ppinitialReduction.cc

    r4a6847 r1d85871  
    121121}
    122122
     123bool p_xLeadmonomDivisibleBy(const poly g, const poly f, const ring r)
     124{
     125  poly gx = p_Head(g,r);
     126  poly fx = p_Head(f,r);
     127  p_SetExp(gx,1,0,r);
     128  p_SetExp(fx,1,0,r);
     129  p_Setm(gx,r);
     130  p_Setm(fx,r);
     131  bool b = p_LeadmonomDivisibleBy(gx,fx,r);
     132  p_Delete(&gx,r);
     133  p_Delete(&fx,r);
     134  return b;
     135}
     136
     137void pReduceInhomogeneous(poly &g, const number p, const ring r)
     138{
     139  if (g==NULL)
     140    return;
     141  p_Test(g,r);
     142
     143  poly toBeChecked = pNext(g);
     144  pNext(g) = NULL; poly gEnd = g;
     145  poly gCache;
     146
     147  number coeff, pPower; int power; poly subst;
     148  while(toBeChecked)
     149  {
     150    for (gCache = g; gCache; pIter(gCache))
     151      if (p_xLeadmonomDivisibleBy(gCache,toBeChecked,r)) break;
     152    if (gCache)
     153    {
     154      n_Power(p,p_GetExp(toBeChecked,1,r)-p_GetExp(gCache,1,r),&pPower,r->cf);
     155      coeff = n_Mult(p_GetCoeff(toBeChecked,r),pPower,r->cf);
     156      p_SetCoeff(gCache,n_Add(p_GetCoeff(gCache,r),coeff,r->cf),r);
     157      n_Delete(&pPower,r->cf); n_Delete(&coeff,r->cf);
     158      toBeChecked=p_LmDeleteAndNext(toBeChecked,r);
     159    }
     160    else
     161    {
     162      if (n_DivBy(p_GetCoeff(toBeChecked,r),p,r->cf))
     163      {
     164        power=1;
     165        coeff=n_Div(p_GetCoeff(toBeChecked,r),p,r->cf);
     166        while (n_DivBy(coeff,p,r->cf))
     167        {
     168          power++;
     169          number coeff0 = n_Div(coeff,p,r->cf);
     170          n_Delete(&coeff,r->cf);
     171          coeff = coeff0;
     172          coeff0 = NULL;
     173          if (power<1)
     174          {
     175            WerrorS("pReduce: overflow in exponent");
     176            throw 0;
     177          }
     178        }
     179        subst=p_LmInit(toBeChecked,r);
     180        p_AddExp(subst,1,power,r);
     181        p_SetCoeff(subst,coeff,r);
     182        p_Setm(subst,r); p_Test(subst,r);
     183        toBeChecked=p_LmDeleteAndNext(toBeChecked,r);
     184        toBeChecked=p_Add_q(toBeChecked,subst,r);
     185        p_Test(toBeChecked,r);
     186      }
     187      else
     188      {
     189        pNext(gEnd)=toBeChecked;
     190        pIter(gEnd); pIter(toBeChecked);
     191        pNext(gEnd)=NULL;
     192        p_Test(g,r);
     193      }
     194    }
     195  }
     196  p_Test(g,r);
     197  divideByCommonGcd(g,r);
     198  return;
     199}
     200
    123201void ptNormalize(poly* gStar, const number p, const ring r)
    124202{

  • Singular/dyn_modules/gfanlib/ppinitialReduction.h

    r4a6847 r1d85871  
    1010bool isOrderingLocalInT(const ring r);
    1111void pReduce(ideal &I, const number p, const ring r);
     12void pReduceInhomogeneous(poly &g, const number p, const ring r);
    1213bool ppreduceInitially(ideal I, const ring r, const number p);
    1314/* BOOLEAN ppreduceInitially(leftv res, leftv args); */

  • Singular/dyn_modules/gfanlib/startingCone.cc

    r4a6847 r1d85871  
    55#include <lift.h>
    66#include <groebnerCone.h>
    7 // #include <neighbours.h>
    87#include <tropicalStrategy.h>
    98#include <tropicalCurves.h>
    … …  
    1413
    1514
    16 /***
    17  * checks whether sigma is contained in the tropical variety
    18  * by testing whether the initial Ideal with respect to the interior point
    19  * is monomial free.
    20  **/
    21 // static bool checkContainmentInTropicalVariety(const groebnerCone sigma)
    22 // {
    23 //   ideal I = sigma.getPolynomialIdeal();
    24 //   ring r = sigma.getPolynomialRing();
    25 //   const tropicalStrategy* currentStrategy = sigma.getTropicalStrategy();
    26 
    27 //   gfan::ZCone zc = sigma.getPolyhedralCone();
    28 //   gfan::ZMatrix zm = zc.extremeRays();
    29 //   for (int i=0; i<zm.getHeight(); i++)
    30 //   {
    31 //     gfan::ZVector w = zm[i];
    32 //     if (currentStrategy->isValuationNonTrivial() && w[0].sign()==0)
    33 //       continue;
    34 //     poly s = currentStrategy->checkInitialIdealForMonomial(I,r,w);
    35 //     if (s)
    36 //     {
    37 //       p_Delete(&s,r);
    38 //       return false;
    39 //     }
    40 //   }
    41 
    42 //   zm = zc.generatorsOfLinealitySpace();
    43 //   for (int i=0; i<zm.getHeight(); i++)
    44 //   {
    45 //     gfan::ZVector w = zm[i];
    46 //     if (currentStrategy->isValuationNonTrivial() && w[0].sign()==0)
    47 //       continue;
    48 //     poly s = currentStrategy->checkInitialIdealForMonomial(I,r,w);
    49 //     if (s)
    50 //     {
    51 //       p_Delete(&s,r);
    52 //       return false;
    53 //     }
    54 //   }
    55 
    56 //   return true;
    57 // }
    58 
    59 
    60 // static bool checkOneCodimensionalLinealitySpace(const groebnerCone sigma)
    61 // {
    62 //   gfan::ZCone zc = sigma.getPolyhedralCone();
    63 //   int linDim = zc.dimensionOfLinealitySpace();
    64 //   int dim = zc.dimension();
    65 //   return (linDim+1)==dim;
    66 // }
     15groebnerCone groebnerStartingCone(const tropicalStrategy& currentStrategy)
     16{
     17  groebnerCone sigma(currentStrategy.getStartingIdeal(), currentStrategy.getStartingRing(), currentStrategy);
     18  return sigma;
     19}
    6720
    6821

  • Singular/dyn_modules/gfanlib/startingCone.h

    r4a6847 r1d85871  
    55#include <tropicalStrategy.h>
    66
     7groebnerCone groebnerStartingCone(const tropicalStrategy& currentStrategy);
    78groebnerCone tropicalStartingCone(const tropicalStrategy& currentCase);
    89

  • Singular/dyn_modules/gfanlib/tropical.cc

    r4a6847 r1d85871  
    55#include <bbcone.h>
    66#include <ppinitialReduction.h>
    7 // #include <ttinitialReduction.h>
    87#include <containsMonomial.h>
    98#include <initial.h>
    109#include <witness.h>
    1110#include <tropicalCurves.h>
    12 // #include <neighbours.h>
    1311#include <tropicalStrategy.h>
    1412#include <startingCone.h>
     13#include <groebnerFan.h>
     14#include <groebnerComplex.h>
    1515#include <tropicalVariety.h>
    1616
     17int tropicalVerboseLevel = 0;
    1718
    1819gfan::ZCone homogeneitySpace(ideal I, ring r)
    … …  
    5152{
    5253  leftv u = args;
     54  if ((u != NULL) && (u->Typ() == POLY_CMD))
     55  {
     56    leftv v = u->next;
     57    if (v == NULL)
     58    {
     59      poly g = (poly) u->Data();
     60      ideal I = idInit(1);
     61      I->m[0] = g;
     62      res->rtyp = coneID;
     63      res->data = (void*) new gfan::ZCone(homogeneitySpace(I,currRing));
     64      I->m[0] = NULL;
     65      id_Delete(&I,currRing);
     66      return FALSE;
     67    }
     68  }
    5369  if ((u != NULL) && (u->Typ() == IDEAL_CMD))
    5470  {
    … …  
    107123{
    108124  leftv u = args;
     125  if ((u != NULL) && (u->Typ() == POLY_CMD))
     126  {
     127    leftv v = u->next;
     128    if (v == NULL)
     129    {
     130      poly g = (poly) u->Data();
     131      ideal I = idInit(1);
     132      I->m[0] = g;
     133      res->rtyp = coneID;
     134      res->data = (void*) new gfan::ZCone(lowerHomogeneitySpace(I,currRing));
     135      I->m[0] = NULL;
     136      id_Delete(&I,currRing);
     137      return FALSE;
     138    }
     139  }
    109140  if ((u != NULL) && (u->Typ() == IDEAL_CMD))
    110141  {
    … …  
    181212{
    182213  leftv u = args;
     214  if ((u != NULL) && (u->Typ() == POLY_CMD))
     215  {
     216    leftv v = u->next;
     217    if ((v !=NULL) && ((v->Typ() == BIGINTMAT_CMD) || (v->Typ() == INTVEC_CMD)))
     218    {
     219      try
     220      {
     221        poly g = (poly) u->Data();
     222        ideal I = idInit(1);
     223        I->m[0] = g;
     224        gfan::ZVector* weightVector;
     225        if (v->Typ() == INTVEC_CMD)
     226        {
     227          intvec* w0 = (intvec*) v->Data();
     228          bigintmat* w1 = iv2bim(w0,coeffs_BIGINT);
     229          w1->inpTranspose();
     230          weightVector = bigintmatToZVector(*w1);
     231          delete w1;
     232        }
     233        else
     234        {
     235          bigintmat* w1 = (bigintmat*) v->Data();
     236          weightVector = bigintmatToZVector(*w1);
     237        }
     238        res->rtyp = coneID;
     239        res->data = (void*) new gfan::ZCone(groebnerCone(I,currRing,*weightVector));
     240        delete weightVector;
     241        I->m[0] = NULL;
     242        id_Delete(&I,currRing);
     243        return FALSE;
     244      }
     245      catch (const std::exception& ex)
     246      {
     247        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     248        return TRUE;
     249      }
     250    }
     251  }
    183252  if ((u != NULL) && (u->Typ() == IDEAL_CMD))
    184253  {
    … …  
    186255    if ((v !=NULL) && ((v->Typ() == BIGINTMAT_CMD) || (v->Typ() == INTVEC_CMD)))
    187256    {
    188       ideal I = (ideal) u->Data();
    189       gfan::ZVector* weightVector;
    190       if (v->Typ() == INTVEC_CMD)
    191       {
    192         intvec* w0 = (intvec*) v->Data();
    193         bigintmat* w1 = iv2bim(w0,coeffs_BIGINT);
    194         w1->inpTranspose();
    195         weightVector = bigintmatToZVector(*w1);
    196         delete w1;
    197       }
    198       else
    199       {
    200         bigintmat* w1 = (bigintmat*) v->Data();
    201         weightVector = bigintmatToZVector(*w1);
    202       }
    203       res->rtyp = coneID;
    204       res->data = (void*) new gfan::ZCone(groebnerCone(I,currRing,*weightVector));
    205       delete weightVector;
    206       return FALSE;
     257      try
     258      {
     259        ideal I = (ideal) u->Data();
     260        gfan::ZVector* weightVector;
     261        if (v->Typ() == INTVEC_CMD)
     262        {
     263          intvec* w0 = (intvec*) v->Data();
     264          bigintmat* w1 = iv2bim(w0,coeffs_BIGINT);
     265          w1->inpTranspose();
     266          weightVector = bigintmatToZVector(*w1);
     267          delete w1;
     268        }
     269        else
     270        {
     271          bigintmat* w1 = (bigintmat*) v->Data();
     272          weightVector = bigintmatToZVector(*w1);
     273        }
     274        res->rtyp = coneID;
     275        res->data = (void*) new gfan::ZCone(groebnerCone(I,currRing,*weightVector));
     276        delete weightVector;
     277        return FALSE;
     278      }
     279      catch (const std::exception& ex)
     280      {
     281        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     282        return TRUE;
     283      }
    207284    }
    208285  }
    … …  
    243320{
    244321  leftv u = args;
     322  if ((u != NULL) && (u->Typ() == POLY_CMD))
     323  {
     324    leftv v = u->next;
     325    if (v == NULL)
     326    {
     327      try
     328      {
     329        poly g = (poly) u->Data();
     330        ideal I = idInit(1);
     331        I->m[0] = g;
     332        res->rtyp = coneID;
     333        res->data = (void*) new gfan::ZCone(maximalGroebnerCone(I,currRing));
     334        I->m[0] = NULL;
     335        id_Delete(&I,currRing);
     336        return FALSE;
     337      }
     338      catch (const std::exception& ex)
     339      {
     340        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     341        return TRUE;
     342      }
     343    }
     344  }
    245345  if ((u != NULL) && (u->Typ() == IDEAL_CMD))
    246346  {
    … …  
    248348    if (v == NULL)
    249349    {
    250       ideal I = (ideal) u->Data();
    251       res->rtyp = coneID;
    252       res->data = (void*) new gfan::ZCone(maximalGroebnerCone(I,currRing));
    253       return FALSE;
     350      try
     351      {
     352        ideal I = (ideal) u->Data();
     353        res->rtyp = coneID;
     354        res->data = (void*) new gfan::ZCone(maximalGroebnerCone(I,currRing));
     355        return FALSE;
     356      }
     357      catch (const std::exception& ex)
     358      {
     359        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     360        return TRUE;
     361      }
    254362    }
    255363  }
    … …  
    293401    if ((v !=NULL) && ((v->Typ() == BIGINTMAT_CMD) || (v->Typ() == INTVEC_CMD)))
    294402    {
    295       ideal I = (ideal) u->Data();
    296       gfan::ZVector* weightVector;
    297       if (v->Typ() == INTVEC_CMD)
    298       {
    299         intvec* w0 = (intvec*) v->Data();
    300         bigintmat* w1 = iv2bim(w0,coeffs_BIGINT);
    301         w1->inpTranspose();
    302         weightVector = bigintmatToZVector(*w1);
    303         delete w1;
    304       }
    305       else
    306       {
    307         bigintmat* w1 = (bigintmat*) v->Data();
    308         weightVector = bigintmatToZVector(*w1);
    309       }
    310       res->rtyp = IDEAL_CMD;
    311       res->data = (void*) initial(I, currRing, *weightVector);
    312       delete weightVector;
    313       return FALSE;
     403      try
     404      {
     405        ideal I = (ideal) u->Data();
     406        gfan::ZVector* weightVector;
     407        if (v->Typ() == INTVEC_CMD)
     408        {
     409          intvec* w0 = (intvec*) v->Data();
     410          bigintmat* w1 = iv2bim(w0,coeffs_BIGINT);
     411          w1->inpTranspose();
     412          weightVector = bigintmatToZVector(*w1);
     413          delete w1;
     414        }
     415        else
     416        {
     417          bigintmat* w1 = (bigintmat*) v->Data();
     418          weightVector = bigintmatToZVector(*w1);
     419        }
     420        res->rtyp = IDEAL_CMD;
     421        res->data = (void*) initial(I, currRing, *weightVector);
     422        delete weightVector;
     423        return FALSE;
     424      }
     425      catch (const std::exception& ex)
     426      {
     427        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     428        return TRUE;
     429      }
    314430    }
    315431  }
    … …  
    326442  p->iiAddCproc("","lowerHomogeneitySpace",FALSE,lowerHomogeneitySpace);
    327443  p->iiAddCproc("","initial",FALSE,initial);
     444  p->iiAddCproc("","tropicalVariety",FALSE,tropicalVariety);
     445  p->iiAddCproc("","groebnerFan",FALSE,groebnerFan);
     446  p->iiAddCproc("","groebnerComplex",FALSE,groebnerComplex);
     447#ifndef NDEBUG
    328448  // p->iiAddCproc("","tropicalNeighbours",FALSE,tropicalNeighbours);
    329 #ifndef NDEBUG
    330449  // p->iiAddCproc("","initial0",FALSE,initial0);
    331450  p->iiAddCproc("","pReduceDebug",FALSE,pReduceDebug);
    … …  
    362481  p->iiAddCproc("","nonPositiveTropicalStartingPoint",FALSE,nonPositiveTropicalStartingPoint);
    363482  p->iiAddCproc("","tropicalStartingCone",FALSE,tropicalStartingCone);
    364   p->iiAddCproc("","tropicalVariety",FALSE,tropicalVariety);
    365483#endif //NDEBUG
    366484  // p->iiAddCproc("","ppreduceInitially",FALSE,ppreduceInitially);

  • Singular/dyn_modules/gfanlib/tropical.h

    r4a6847 r1d85871  
    44#include "Singular/ipid.h"
    55
     6
     7extern int tropicalVerboseLevel;
     8
     9
    610gfan::ZCone homogeneitySpace(ideal I, ring r);
    711void tropical_setup(SModulFunctions* p);

  • Singular/dyn_modules/gfanlib/tropicalCurves.cc

    r4a6847 r1d85871  
    7070
    7171  /* construct a new ordering as describe above */
    72   bool overflow;
     72  bool overflow = false;
    7373  s->order[0] = ringorder_a;
    7474  s->block0[0] = 1;
    … …  
    9797
    9898  if (overflow)
    99     throw 0; //todo: use a proper custom exception
     99  {
     100    WerrorS("genericlyWeightedOrdering: overflow in weight vector");
     101    throw 0; // weightOverflow;
     102  }
    100103
    101104  /* complete the ring and return it */

  • Singular/dyn_modules/gfanlib/tropicalStrategy.cc

    r4a6847 r1d85871  
    134134  extraReductionAlgorithm(noExtraReduction)
    135135{
    136   assume(rField_is_Q(r) || rField_is_Zp(r));
     136  assume(rField_is_Q(r) || rField_is_Zp(r) || rField_is_Ring_Z(r));
    137137  if (!completelyHomogeneous)
    138138  {
    … …  
    225225  omFreeSize(shiftByOne,(n+1)*sizeof(int));
    226226
    227   // ring origin = currRing;
    228   // rChangeCurrRing(startingRing);
    229   // ideal startingIdeal = kNF(pt,startingRing->qideal,J);
    230   // rChangeCurrRing(origin);
    231   ideal startingIdeal = J; J = NULL;
     227  ring origin = currRing;
     228  rChangeCurrRing(startingRing);
     229  ideal startingIdeal = kNF(pt,startingRing->qideal,J); // mathematically redundant,
     230  rChangeCurrRing(origin);                              // but helps with upcoming std computation
     231  // ideal startingIdeal = J; J = NULL;
    232232  assume(startingIdeal->m[k]==NULL);
    233233  startingIdeal->m[k] = pt->m[0];
    … …  
    271271  /* map the input ideal into the new polynomial ring */
    272272  startingIdeal = constructStartingIdeal(J,s,uniformizingParameter,startingRing);
     273  reduce(startingIdeal,startingRing);
    273274
    274275  linealitySpace = homogeneitySpace(startingIdeal,startingRing);

  • Singular/dyn_modules/gfanlib/tropicalTraversal.cc

    r4a6847 r1d85871  
    33#include <tropicalCurves.h>
    44
    5 std::vector<bool> checkNecessaryFlips(const groebnerCones &tropicalVariety, const groebnerCones &workingList,
    6                                       const gfan::ZVector &interiorPoint, const gfan::ZMatrix &normalVectors)
     5std::vector<bool> checkNecessaryTropicalFlips(const groebnerCones &tropicalVariety, const groebnerCones &workingList,
     6                                              const gfan::ZVector &interiorPoint, const gfan::ZMatrix &normalVectors)
    77{
    88  int k = normalVectors.getHeight();
    … …  
    7272        id_Delete(&inI,sigma.getPolynomialRing());
    7373
    74         std::vector<bool> needToFlip = checkNecessaryFlips(tropicalVariety,workingList,interiorPoint,normalVectors);
     74        std::vector<bool> needToFlip = checkNecessaryTropicalFlips(tropicalVariety,workingList,interiorPoint,normalVectors);
    7575        for (int j=0; j<normalVectors.getHeight(); j++)
    7676        {
    … …  
    8888    workingList.erase(sigma);
    8989    tropicalVariety.insert(sigma);
    90     // std::cout << "tropicalVariety.size():" << tropicalVariety.size() << std::endl;
    91     // std::cout << "workingList.size():" << workingList.size() << std::endl;
     90    if (printlevel > 0)
     91      Print("cones finished: %lu   cones in working list: %lu\n", tropicalVariety.size(), workingList.size());
    9292  }
    9393  return tropicalVariety;
    9494}
    9595
    96 groebnerCones tropicalTraversal(const groebnerCone startingCone)
     96
     97std::vector<bool> checkNecessaryGroebnerFlips(const groebnerCones &groebnerFan, const groebnerCones &workingList,
     98                                              const gfan::ZMatrix &interiorPoints)
    9799{
    98   groebnerCones tropicalVariety;
     100  int k = interiorPoints.getHeight();
     101  std::vector<bool> needToFlip(k,true);
     102
     103  for (groebnerCones::iterator sigma = groebnerFan.begin(); sigma!=groebnerFan.end(); sigma++)
     104  {
     105    for (int i=0; i<k; i++)
     106    {
     107      if (needToFlip[i] && sigma->contains(interiorPoints[i]))
     108        needToFlip[i] = false;
     109    }
     110  }
     111
     112  for (groebnerCones::iterator sigma = workingList.begin(); sigma!=workingList.end(); sigma++)
     113  {
     114    for (int i=0; i<k; i++)
     115    {
     116      if (needToFlip[i] && sigma->contains(interiorPoints[i]))
     117        needToFlip[i] = false;
     118    }
     119  }
     120
     121  return needToFlip;
     122}
     123
     124
     125groebnerCones groebnerTraversal(const groebnerCone startingCone)
     126{
     127  const tropicalStrategy* currentStrategy = startingCone.getTropicalStrategy();
     128
     129  groebnerCones groebnerFan;
    99130  groebnerCones workingList;
    100131  workingList.insert(startingCone);
    101   // std::cout << "starting traversal" << std::endl;
     132  std::set<gfan::ZVector> finishedInteriorPoints;
     133  bool onlyLowerHalfSpace = !currentStrategy->isValuationTrivial();
     134
    102135  while(!workingList.empty())
    103136  {
    104     const groebnerCone sigma=*(workingList.begin());
    105     const groebnerCones neighbours = sigma.tropicalNeighbours();
    106     for (groebnerCones::iterator tau = neighbours.begin(); tau!=neighbours.end(); tau++)
     137    /**
     138     * Pick a maximal Groebner cone from the working list
     139     * and compute interior points on its facets as well as outer facet normals
     140     */
     141    groebnerCone sigma=*(workingList.begin());
     142    workingList.erase(workingList.begin());
     143
     144    std::pair<gfan::ZMatrix,gfan::ZMatrix> interiorPointsAndOuterFacetNormals = interiorPointsAndNormalsOfFacets(sigma.getPolyhedralCone(), finishedInteriorPoints, onlyLowerHalfSpace);
     145    gfan::ZMatrix interiorPoints = interiorPointsAndOuterFacetNormals.first;
     146    gfan::ZMatrix outerFacetNormals = interiorPointsAndOuterFacetNormals.second;
     147    std::vector<bool> needToFlip = checkNecessaryGroebnerFlips(groebnerFan,workingList, interiorPoints);
     148
     149    for (int i=0; i<interiorPoints.getHeight(); i++)
    107150    {
    108       if (tropicalVariety.count(*tau)==0)
    109         workingList.insert(*tau);
     151      gfan::ZVector interiorPoint = interiorPoints[i];
     152
     153      if (needToFlip[i]==true)
     154      {
     155        groebnerCone neighbour = sigma.flipCone(interiorPoint,outerFacetNormals[i]);
     156        workingList.insert(neighbour);
     157      }
     158      finishedInteriorPoints.insert(interiorPoints[i]);
    110159    }
    111     tropicalVariety.insert(sigma);
    112     workingList.erase(sigma);
    113     // std::cout << "tropicalVariety.size():" << tropicalVariety.size() << std::endl;
    114     // std::cout << "workingList.size():" << workingList.size() << std::endl;
     160
     161    sigma.deletePolynomialData();
     162    groebnerFan.insert(sigma);
     163    if (printlevel > 0)
     164      Print("cones finished: %lu   cones in working list: %lu\n", groebnerFan.size(), workingList.size());
    115165  }
    116   return tropicalVariety;
     166  return groebnerFan;
    117167}

  • Singular/dyn_modules/gfanlib/tropicalTraversal.h

    r4a6847 r1d85871  
    22#define TROPICAL_TRAVERSAL_H
    33
     4#include <groebnerCone.h>
     5
    46groebnerCones tropicalTraversalMinimizingFlips(const groebnerCone startingCone);
    57groebnerCones tropicalTraversal(const groebnerCone startingCone);
     8groebnerCones groebnerTraversal(const groebnerCone startingCone);
    69
    710#endif

  • Singular/dyn_modules/gfanlib/tropicalVariety.cc

    r4a6847 r1d85871  
    1414BITSET bitsetSave1, bitsetSave2;
    1515
     16
    1617/***
    1718 * sets option(redSB)
    … …  
    2223  si_opt_1|=Sy_bit(OPT_REDSB);
    2324}
     25
    2426
    2527/***
    … …  
    4042}
    4143
     44
    4245BOOLEAN tropicalVariety(leftv res, leftv args)
    4346{
    4447  leftv u = args;
     48  if ((u!=NULL) && (u->Typ()==POLY_CMD))
     49  {
     50    poly g = (poly) u->Data();
     51    leftv v = u->next;
     52    if (v==NULL)
     53    {
     54      try
     55      {
     56        ideal I = idInit(1);
     57        I->m[0] = g;
     58        tropicalStrategy currentStrategy(I,currRing);
     59        std::set<gfan::ZCone> maxCones = tropicalVariety(g,currRing,&currentStrategy);
     60        res->rtyp = fanID;
     61        res->data = (char*) toZFan(maxCones);
     62        I->m[0] = NULL;
     63        id_Delete(&I,currRing);
     64        return FALSE;
     65      }
     66      catch (const std::exception& ex)
     67      {
     68        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     69        return TRUE;
     70      }
     71    }
     72    if ((v!=NULL) && (v->Typ()==NUMBER_CMD))
     73    {
     74      try
     75      {
     76        ideal I = idInit(1);
     77        I->m[0] = g;
     78        number p = (number) v->Data();
     79        tropicalStrategy currentStrategy(I,p,currRing);
     80        ideal startingIdeal = currentStrategy.getStartingIdeal();
     81        ring startingRing = currentStrategy.getStartingRing();
     82        poly gStart = startingIdeal->m[0];
     83        std::set<gfan::ZCone> maxCones = tropicalVariety(gStart,startingRing,&currentStrategy);
     84        res->rtyp = fanID;
     85        res->data = (char*) toZFan(maxCones);
     86        I->m[0] = NULL;
     87        id_Delete(&I,currRing);
     88        return FALSE;
     89      }
     90      catch (const std::exception& ex)
     91      {
     92        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     93        return TRUE;
     94      }
     95
     96    }
     97  }
    4598  if ((u!=NULL) && (u->Typ()==IDEAL_CMD))
    4699  {
    … …  
    53106      if (v==NULL)
    54107      {
    55         tropicalStrategy currentStrategy(I,currRing);
    56         std::set<gfan::ZCone> maxCones = tropicalVariety(g,currRing,&currentStrategy);
    57         // gfan::ZFan* zf = toZFan(maxCones);
    58         // delete zf;
    59         // res->rtyp = NONE;
    60         // res->data = NULL;
    61         // while (1)
    62         // {
    63         //   omUpdateInfo();
    64         //   Print("usedBytesAfter=%ld\n",om_Info.UsedBytes);
    65         //   tropicalStrategy debugTest(I,currRing);
    66         //   maxCones = tropicalVariety(g,currRing,debugTest);
    67         //   gfan::ZFan* zf = toZFan(maxCones);
    68         //   delete zf;
    69         // }
    70         res->rtyp = fanID;
    71         res->data = (char*) toZFan(maxCones);
    72         return FALSE;
     108        try
     109        {
     110          tropicalStrategy currentStrategy(I,currRing);
     111          std::set<gfan::ZCone> maxCones = tropicalVariety(g,currRing,&currentStrategy);
     112          res->rtyp = fanID;
     113          res->data = (char*) toZFan(maxCones);
     114          return FALSE;
     115        }
     116        catch (const std::exception& ex)
     117        {
     118          WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     119          return TRUE;
     120        }
    73121      }
    74122      if ((v!=NULL) && (v->Typ()==NUMBER_CMD))
    75123      {
    76         number p = (number) v->Data();
    77         tropicalStrategy currentStrategy(I,p,currRing);
    78         ideal startingIdeal = currentStrategy.getStartingIdeal();
    79         ring startingRing = currentStrategy.getStartingRing();
    80         poly gStart = startingIdeal->m[0];
    81         std::set<gfan::ZCone> maxCones = tropicalVariety(gStart,startingRing,&currentStrategy);
    82         res->rtyp = fanID;
    83         res->data = (char*) toZFan(maxCones);
    84         return FALSE;
     124        try
     125        {
     126          number p = (number) v->Data();
     127          tropicalStrategy currentStrategy(I,p,currRing);
     128          ideal startingIdeal = currentStrategy.getStartingIdeal();
     129          ring startingRing = currentStrategy.getStartingRing();
     130          poly gStart = startingIdeal->m[0];
     131          std::set<gfan::ZCone> maxCones = tropicalVariety(gStart,startingRing,&currentStrategy);
     132          res->rtyp = fanID;
     133          res->data = (char*) toZFan(maxCones);
     134          return FALSE;
     135        }
     136        catch (const std::exception& ex)
     137        {
     138          WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     139          return TRUE;
     140        }
    85141      }
    86142    }
    … …  
    88144    if (v==NULL)
    89145    {
    90       setOptionRedSB();
    91       // ideal stdI;
    92       // if (!hasFlag(u,FLAG_STD))
    93       //   stdI = gfanlib_kStd_wrapper(I,currRing);
    94       // else
    95       //   stdI = id_Copy(I,currRing);
    96       tropicalStrategy currentStrategy(I,currRing);
    97       gfan::ZFan* tropI = tropicalVariety(currentStrategy);
    98       res->rtyp = fanID;
    99       res->data = (char*) tropI;
    100       undoSetOptionRedSB();
    101       // id_Delete(&stdI,currRing);
    102       return FALSE;
     146      try
     147      {
     148        setOptionRedSB();
     149        ideal stdI;
     150        if (!hasFlag(u,FLAG_STD))
     151          stdI = gfanlib_kStd_wrapper(I,currRing);
     152        else
     153          stdI = id_Copy(I,currRing);
     154        tropicalStrategy currentStrategy(stdI,currRing);
     155        gfan::ZFan* tropI = tropicalVariety(currentStrategy);
     156        res->rtyp = fanID;
     157        res->data = (char*) tropI;
     158        undoSetOptionRedSB();
     159        id_Delete(&stdI,currRing);
     160        return FALSE;
     161      }
     162      catch (const std::exception& ex)
     163      {
     164        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     165        return TRUE;
     166      }
    103167    }
    104168    if ((v!=NULL) && (v->Typ()==NUMBER_CMD))
    105169    {
    106       number p = (number) v->Data();
    107       // ideal stdI;
    108       // if (!hasFlag(u,FLAG_STD))
    109       //   stdI = gfanlib_kStd_wrapper(I,currRing);
    110       // else
    111       //   stdI = id_Copy(I,currRing);
    112       // tropicalStrategy currentStrategy(stdI,p,currRing);
    113       tropicalStrategy currentStrategy(I,p,currRing);
    114       gfan::ZFan* tropI = tropicalVariety(currentStrategy);
    115       res->rtyp = fanID;
    116       res->data = (char*) tropI;
    117       // id_Delete(&stdI,currRing);
    118       return FALSE;
     170      try
     171      {
     172        number p = (number) v->Data();
     173        ideal stdI;
     174        if (!hasFlag(u,FLAG_STD))
     175          stdI = gfanlib_kStd_wrapper(I,currRing);
     176        else
     177          stdI = id_Copy(I,currRing);
     178        tropicalStrategy currentStrategy(stdI,p,currRing);
     179        gfan::ZFan* tropI = tropicalVariety(currentStrategy);
     180        res->rtyp = fanID;
     181        res->data = (char*) tropI;
     182        id_Delete(&stdI,currRing);
     183        return FALSE;
     184      }
     185      catch (const std::exception& ex)
     186      {
     187        WerrorS("ERROR: "); WerrorS(ex.what()); WerrorS("\n");
     188        return TRUE;
     189      }
    119190    }
    120191    return FALSE;
Note: See TracChangeset for help on using the changeset viewer.