Changeset eb836c in git

Timestamp:

Jul 21, 2014, 11:41:11 AM (9 years ago)

Author:

Yue Ren <ren@…>

Branches:

(u'spielwiese', '8e0ad00ce244dfd0756200662572aef8402f13d5')

Children:

d4049cf8c7b6c2765919c8fcc08a3b038d0e1c29

Parents:

e744d9067cafc8137765bab643ab9ebf8120f339

git-author:

Yue Ren <ren@mathematik.uni-kl.de>2014-07-21 11:41:11+02:00

git-committer:

Yue Ren <ren@mathematik.uni-kl.de>2015-02-06 13:47:03+01:00

Message:

chg: first working version

Files:

• Singular/dyn_modules/gfanlib/adjustWeights.cc (modified) (1 diff)
• Singular/dyn_modules/gfanlib/groebnerCone.cc (modified) (5 diffs)
• Singular/dyn_modules/gfanlib/groebnerCone.h (modified) (1 diff)
• Singular/dyn_modules/gfanlib/tropicalCurves.cc (modified) (9 diffs)
• Singular/dyn_modules/gfanlib/tropicalCurves.h (modified) (1 diff)
• dyn_modules/callgfanlib/std_wrapper.cc (modified) (1 diff)
• dyn_modules/callgfanlib/tropicalTraversal.cc (modified) (1 diff)

Singular/dyn_modules/gfanlib/adjustWeights.cc

@@ -64 +64 @@
   /* find k such that e+k*w is strictly positive,
    * i.e. k such that e[i]+k*w[i] is strictly positive
-   * for all i=0,...,n */
+   * for all i=0,...,n
+   * note that the division is rounded towards zero,
+   * hence we increment the value by 1 */
   gfan::Integer k((long)0);
   if (e[0].sign()<=0)
-    k = gfan::Integer((long)1)-(gfan::Integer((long)e[0].sign())*e[0])/w[0];
+    k = gfan::Integer((long)1)-(e[0]/w[0]);
   for (unsigned i=1; i<e.size(); i++)
   {
-    gfan::Integer kk = gfan::Integer((long)1)-(gfan::Integer((long)e[i].sign())*e[i])/w[i];
-    if (k<kk)
-      k = kk;
+    if (e[i].sign()<=0)
+    {
+      gfan::Integer kk = gfan::Integer((long)1)-(e[i]/w[i]);
+      if (k<kk)
+        k = kk;
+    }
   }
   /* compute e+k*w, check it for correctness and return it */

Singular/dyn_modules/gfanlib/groebnerCone.cc

@@ -41 +41 @@
                 << "weight: " << std::endl
                 << v << std::endl;
-    }
-    return zc.contains(v);
+      return false;
+    }
+    return true;
   }
   return (zc.dimension()==0);
@@ -144 +145 @@
 }
 
+/***
+ * Computes the groebner cone of I around u+e*w for e>0 sufficiently small.
+ * Assumes that this cone is a face of the maximal Groenbner cone given by the ordering of r.
+ **/
+groebnerCone::groebnerCone(const ideal I, const ring r, const gfan::ZVector& u, const gfan::ZVector& w, const tropicalStrategy& currentCase):
+  polynomialIdeal(NULL),
+  polynomialRing(NULL),
+  currentStrategy(&currentCase)
+{
+  assume(checkPolynomialInput(I,r));
+  if (I) polynomialIdeal = id_Copy(I,r);
+  if (r) polynomialRing = rCopy(r);
+
+  int n = rVar(r);
+  gfan::ZMatrix inequalities = gfan::ZMatrix(0,n);
+  gfan::ZMatrix equations = gfan::ZMatrix(0,n);
+  int* expv = (int*) omAlloc((n+1)*sizeof(int));
+  for (int i=0; i<idSize(polynomialIdeal); i++)
+  {
+    poly g = polynomialIdeal->m[i];
+    p_GetExpV(g,expv,polynomialRing);
+    gfan::ZVector leadexpv = intStar2ZVector(n,expv);
+    long d1 = wDeg(g,polynomialRing,u);
+    long d2 = wDeg(g,polynomialRing,w);
+    for (pIter(g); g; pIter(g))
+    {
+      p_GetExpV(g,expv,polynomialRing);
+      gfan::ZVector tailexpv = intStar2ZVector(n,expv);
+      if (wDeg(g,polynomialRing,u)==d1 && wDeg(g,polynomialRing,w)==d2)
+        equations.appendRow(leadexpv-tailexpv);
+      else
+      {
+        assume(wDeg(g,polynomialRing,u)<d1 || wDeg(g,polynomialRing,w)<d2);
+        inequalities.appendRow(leadexpv-tailexpv);
+      }
+    }
+  }
+  omFreeSize(expv,(n+1)*sizeof(int));
+
+  polyhedralCone = gfan::ZCone(inequalities,equations);
+  interiorPoint = polyhedralCone.getRelativeInteriorPoint();
+  assume(checkOrderingAndCone(polynomialRing,polyhedralCone));
+}
+
+
 groebnerCone::groebnerCone(const ideal I, const ideal inI, const ring r, const tropicalStrategy& currentCase):
   polynomialIdeal(id_Copy(I,r)),
@@ -272 +318 @@
   }
   /* check whether facet normal points outwards */
-  gfan::ZMatrix halfSpaceInequality(0,facetNormal.size());
-  halfSpaceInequality.appendRow(facetNormal);
-  gfan::ZCone halfSpace = gfan::ZCone(halfSpaceInequality,gfan::ZMatrix(0,facetNormal.size()));
-  hopefullyAFacet = intersection(polyhedralCone, halfSpace);
-  if (hopefullyAFacet.dimension()!=(polyhedralCone.dimension()-1))
+  gfan::ZCone dual = polyhedralCone.dualCone();
+  if(dual.contains(facetNormal))
   {
     std::cout << "ERROR: facetNormal is not pointing outwards!" << std::endl
@@ -302 +345 @@
   std::pair<ideal,ring> flipped = flip(polynomialIdeal,polynomialRing,interiorPoint,facetNormal,*currentStrategy);
   assume(checkPolynomialInput(flipped.first,flipped.second));
-  groebnerCone flippedCone(flipped.first, flipped.second, facetNormal, *currentStrategy);
+  groebnerCone flippedCone(flipped.first, flipped.second, interiorPoint, facetNormal, *currentStrategy);
   return flippedCone;
 }
@@ -432 +475 @@
   {
     ideal initialIdeal = initial(polynomialIdeal,polynomialRing,interiorPoints[i]);
-    std::set<gfan::ZVector> rays = raysOfTropicalCurve(initialIdeal,polynomialRing,*currentStrategy);
-    groebnerCones neighbours;
+    std::set<gfan::ZVector> rays = raysOfTropicalStar(initialIdeal,polynomialRing,interiorPoints[i],*currentStrategy);
     for (std::set<gfan::ZVector>::iterator ray = rays.begin(); ray!=rays.end(); ray++)
       neighbours.insert(this->flipCone(interiorPoints[i],*ray));

Singular/dyn_modules/gfanlib/groebnerCone.h

@@ -32 +32 @@
   groebnerCone(const ideal I, const ring r, const tropicalStrategy& currentCase);
   groebnerCone(const ideal I, const ring r, const gfan::ZVector& w, const tropicalStrategy& currentCase);
+  groebnerCone(const ideal I, const ring r, const gfan::ZVector& u, const gfan::ZVector& w, const tropicalStrategy& currentCase);
   groebnerCone(const ideal I, const ideal inI, const ring r, const tropicalStrategy& currentCase);
   groebnerCone(const groebnerCone& sigma);

Singular/dyn_modules/gfanlib/tropicalCurves.cc

@@ -4 +4 @@
 #include <callgfanlib_conversion.h>
 #include <gfanlib_exceptions.h>
+#include <std_wrapper.h>
 #include <containsMonomial.h>
 #include <initial.h>
@@ -54 +55 @@
 {
   int n = rVar(r);
-  int h = E.getHeight();
+  int h = W.getHeight();
 
   /* create a copy s of r and delete its ordering */
@@ -78 +79 @@
   s->block0[1] = 1;
   s->block1[1] = n;
-  gfan::ZVector wAdjusted = currentStrategy.adjustWeightUnterHomogeneity(w,uAdjusted);
+  gfan::ZVector wAdjusted = currentStrategy.adjustWeightUnderHomogeneity(w,uAdjusted);
   s->wvhdl[1] = ZVectorToIntStar(wAdjusted,overflow);
-  for (int j=0; j<h; j++)
+  for (int j=0; j<h-1; j++)
   {
     s->order[j+2] = ringorder_a;
@@ -88 +89 @@
     s->wvhdl[j+2] = ZVectorToIntStar(wAdjusted,overflow);
   }
-  s->order[h+2] = ringorder_wp;
-  s->block0[h+2] = 1;
-  s->block1[h+2] = n;
-  wAdjusted = currentStrategy.adjustWeightUnderHomogeneity(W[j],uAdjusted);
-  s->wvhdl[h+2] = ZVectorToIntStar(wAdjusted,overflow);
-  s->order[h+3] = ringorder_C;
+  s->order[h+1] = ringorder_wp;
+  s->block0[h+1] = 1;
+  s->block1[h+1] = n;
+  wAdjusted = currentStrategy.adjustWeightUnderHomogeneity(W[h-1],uAdjusted);
+  s->wvhdl[h+1] = ZVectorToIntStar(wAdjusted,overflow);
+  s->order[h+2] = ringorder_C;
 
   if (overflow)
@@ -130 +131 @@
   for (std::set<gfan::ZCone>::iterator zc=C.begin(); zc!=C.end(); zc++)
   {
-    gfan::ZVector v = zc->getRelativeInteriorPoint();
+    gfan::ZVector w = zc->getRelativeInteriorPoint();
     gfan::ZMatrix W = zc->generatorsOfSpan();
 
-    ring s = genericlyWeightedOrdering(r,v,W,currentStrategy);
+    ring s = genericlyWeightedOrdering(r,u,w,W,currentStrategy);
     nMapFunc identity = n_SetMap(r->cf,s->cf);
     ideal inIs = idInit(k);
@@ -140 +141 @@
 
     inIs = gfanlib_kStd_wrapper(inIs,s,isHomog);
-    ideal ininIs = initial(inIs,s,E[E.getHeight()-1]);
+    ideal ininIs = initial(inIs,s,W[W.getHeight()-1]);
 
     poly mons = checkForMonomialViaSuddenSaturation(inIs,s);
@@ -149 +150 @@
       nMapFunc mMap = n_SetMap(s->cf,r->cf);
       poly gr = p_PermPoly(gs,NULL,s,r,mMap,NULL,0);
-      idInsertPoly(I,gr);
+      idInsertPoly(inI,gr);
       k++;
       p_Delete(&mons,s);
@@ -155 +156 @@
       zc = C.begin();
     }
-    id_Delete(&Is,s);
     id_Delete(&inIs,s);
+    id_Delete(&ininIs,s);
     rDelete(s);
   }
@@ -163 +164 @@
 
 
-std::set<gfan::ZVector> raysOfTropicalCurve(ideal I, const ring r, const tropicalStrategy& currentStrategy)
-{
-  std::set<gfan::ZCone> C = tropicalCurve(I,r,currentStrategy);
+std::set<gfan::ZVector> raysOfTropicalStar(ideal I, const ring r, const gfan::ZVector u, const tropicalStrategy& currentStrategy)
+{
+  std::set<gfan::ZCone> C = tropicalStar(I,r,u,currentStrategy);
   std::set<gfan::ZVector> raysOfC;
   for (std::set<gfan::ZCone>::iterator zc=C.begin(); zc!=C.end(); zc++)

Singular/dyn_modules/gfanlib/tropicalCurves.h

@@ -7 +7 @@
 #include <tropicalStrategy.h>
 
-std::set<gfan::ZCone> tropicalCurve(const ideal I, const ring r, const tropicalStrategy& currentStrategy);
-std::set<gfan::ZVector> raysOfTropicalCurve(ideal I, const ring r, const tropicalStrategy& currentStrategy);
+std::set<gfan::ZCone> tropicalStar(const ideal I, const ring r, const gfan::ZVector u,
+                                   const tropicalStrategy& currentStrategy);
+std::set<gfan::ZVector> raysOfTropicalStar(ideal I, const ring r, const gfan::ZVector u,
+                                           const tropicalStrategy& currentStrategy);
 
 #ifndef NDEBUG

dyn_modules/callgfanlib/std_wrapper.cc

@@ -1 +1 @@
 #include <kernel/kstd1.h>
 #include <kernel/polys.h>
+#include <kernel/ideals.h>
 
 ideal gfanlib_kStd_wrapper(ideal I, ring r, tHomog h=testHomog)

dyn_modules/callgfanlib/tropicalTraversal.cc

@@ -10 +10 @@
     const groebnerCone sigma=*(workingList.begin());
     const groebnerCones neighbours = sigma.tropicalNeighbours();
-    workingList.insert(neighbours.begin(),neighbours.end());
+    for (groebnerCones::iterator tau = neighbours.begin(); tau!=neighbours.end(); tau++)
+    {
+      if (tropicalVariety.count(*tau)==0)
+        workingList.insert(*tau);
+    }
     tropicalVariety.insert(sigma);
     workingList.erase(sigma);
+    std::cout << "tropicalVariety.size():" << tropicalVariety.size() << std::endl;
+    std::cout << "workingList.size():" << workingList.size() << std::endl;
   }
   return tropicalVariety;