shiftgb.cc
Go to the documentation of this file.
1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /*
5 * ABSTRACT: kernel: utils for shift GB and free GB
6 */
7 
8 
9 
10 
11 #include <kernel/mod2.h>
12 
13 #ifdef HAVE_SHIFTBBA
14 #include <polys/monomials/ring.h>
15 #include <kernel/polys.h>
16 #include <coeffs/numbers.h>
17 #include <kernel/ideals.h>
18 #include <polys/matpol.h>
19 #include <polys/kbuckets.h>
20 #include <kernel/GBEngine/kstd1.h>
21 #include <polys/sbuckets.h>
23 #include <kernel/GBEngine/kutil.h>
24 #include <kernel/structs.h>
25 #include <omalloc/omalloc.h>
26 #include <kernel/GBEngine/khstd.h>
27 #include <polys/kbuckets.h>
28 #include <polys/weight.h>
29 #include <misc/intvec.h>
30 #include <kernel/structs.h>
33 #include <polys/weight.h>
34 #include <misc/intvec.h>
36 #include <polys/nc/sca.h>
37 
38 
39 #define freeT(A,v) omFreeSize((ADDRESS)A,(v+1)*sizeof(int))
40 
41 
42 /* TODO: write p* stuff as instances of p_* for all the functions */
43 /* p_* functions are new, p* are old */
44 
45 poly p_LPshiftT(poly p, int sh, int uptodeg, int lV, kStrategy strat, const ring r)
46 {
47  /* assume shift takes place, shifts the poly p by sh */
48  /* p is like TObject: lm in currRing = r, tail in tailRing */
49  /* copies p */
50 
51  if (p==NULL) return(p);
52 
55 
56  /* assume sh and uptodeg agree TODO check */
57 
58  if (sh == 0) return(p); /* the zero shift */
59 
60  poly q = NULL;
61  poly s = p_mLPshift(p_Head(p,r), sh, uptodeg, lV, r); // lm in currRing
62  /* pNext(s) will be fixed below */
63  poly pp = pNext(p);
64 
65  while (pp != NULL)
66  {
67  poly h=p_mLPshift(p_Head(pp,strat->tailRing),sh,uptodeg,lV,strat->tailRing);
68  pIter(pp);
69 
70  q = p_Add_q(q, h,strat->tailRing);
71  }
72  pNext(s) = q;
73  /* int version: returns TRUE if it was successful */
74  return(s);
75 }
76 
77 poly p_LPshift(poly p, int sh, int uptodeg, int lV, const ring r)
78 {
79  /* assume shift takes place */
80  /* shifts the poly p from the ring r by sh */
81 
82  /* assume sh and uptodeg agree TODO check */
83  assume(sh>=0);
84 
85  if (sh == 0) return(p); /* the zero shift */
86 
87  poly q = NULL;
88  poly pp = p;
89  while (pp!=NULL)
90  {
91  poly h=pp;
92  pIter(pp);
93  pNext(h)=NULL;
94  h=p_mLPshift(h,sh,uptodeg,lV,r);
95  q = p_Add_q(q, h,r);
96  }
97  return(q);
98 }
99 
100 poly p_mLPshift(poly p, int sh, int uptodeg, int lV, const ring r)
101 {
102  /* p is a monomial from the ring r */
103 
104  if (sh == 0) return(p); /* the zero shift */
105 
106  assume(sh>=0);
107  int L = p_mLastVblock(p,lV,r);
108  assume(L+sh-1<=uptodeg);
109 
110  int *e=(int *)omAlloc0((r->N+1)*sizeof(int));
111  int *s=(int *)omAlloc0((r->N+1)*sizeof(int));
112  p_GetExpV(p,e,r);
113 
114  int j;
115  // for (j=1; j<=r->N; j++)
116  // L*lV gives the last position of the last block
117  for (j=1; j<= L*lV ; j++)
118  {
119  assume(e[j]<=1);
120  if (e[j]==1)
121  {
122  assume(j + (sh*lV)<=r->N);
123  s[j + (sh*lV)] = e[j]; /* actually 1 */
124  }
125  }
126  p_SetExpV(p,s,r);
127  freeT(e, r->N);
128  freeT(s, r->N);
129  /* pSetm(m); */ /* done in the pSetExpV */
130  /* think on the component and coefficient */
131  // number c = pGetCoeff(p);
132  // p_SetCoeff0(m,p_GetCoeff(p,r),r);
133  return(p);
134 }
135 
136 int p_LastVblockT(poly p, int lV, kStrategy strat, const ring r)
137 {
138  /* returns the number of maximal block */
139  /* appearing among the monomials of p */
140  /* the 0th block is the 1st one */
141 
142  /* p is like TObject: lm in currRing = r, tail in tailRing */
145 
146  int ans = p_mLastVblock(p, lV, r); // Block of LM
147  poly q = pNext(p);
148  int ansnew = 0;
149  while (q != NULL)
150  {
151  ansnew = p_mLastVblock(q, lV, strat->tailRing);
152  ans = si_max(ans,ansnew);
153  pIter(q);
154  }
155  /* do not need to delete q */
156  return(ans);
157 }
158 
159 int p_LastVblock(poly p, int lV, const ring r)
160 {
161  /* returns the number of maximal block */
162  /* appearing among the monomials of p */
163  /* the 0th block is the 1st one */
164  poly q = p;
165  int ans = 0;
166  int ansnew = 0;
167  while (q!=NULL)
168  {
169  ansnew = p_mLastVblock(q, lV, r);
170  ans = si_max(ans,ansnew);
171  pIter(q);
172  }
173  return(ans);
174 }
175 
176 int p_mLastVblock(poly p, int lV, const ring r)
177 {
178  /* for a monomial p, returns the number of the last block */
179  /* where a nonzero exponent is sitting */
180  if (p_LmIsConstant(p,r))
181  {
182  return(0);
183  }
184  int *e=(int *)omAlloc0((r->N+1)*sizeof(int));
185  p_GetExpV(p,e,r);
186  int j,b;
187  j = r->N;
188  while ( (!e[j]) && (j>=1) ) j--;
189  freeT(e, r->N);
190  assume(j>0);
191  b = (int)((j+lV-1)/lV); /* the number of the block, >=1 */
192  return (b);
193 }
194 
195 int pFirstVblock(poly p, int lV)
196 {
197  /* returns the number of maximal block */
198  /* appearing among the monomials of p */
199  /* the 0th block is the 1st one */
200  poly q = p; //p_Copy(p,currRing); /* need it ? */
201  int ans = 0;
202  int ansnew = 0;
203  while (q!=NULL)
204  {
205  ansnew = pmFirstVblock(q,lV);
206  ans = si_min(ans,ansnew);
207  pIter(q);
208  }
209  /* do not need to delete q */
210  return(ans);
211 }
212 
213 int pmFirstVblock(poly p, int lV)
214 {
215  if (pIsConstantPoly(p))
216  {
217  return(int(0));
218  }
219  /* for a monomial p, returns the number of the first block */
220  /* where a nonzero exponent is sitting */
221  int *e=(int *)omAlloc0((currRing->N+1)*sizeof(int));
222  pGetExpV(p,e);
223  int j,b;
224  j = 1;
225  while ( (!e[j]) && (j<=currRing->N-1) ) j++;
226  if (j==currRing->N + 1)
227  {
228 #ifdef PDEBUG
229  PrintS("pmFirstVblock: unexpected zero exponent vector\n");
230 #endif
231  return(j);
232  }
233  b = (int)(j/lV)+1; /* the number of the block, 1<= N <= currRing->N */
234  return (b);
235 }
236 
237  /* there should be two routines: */
238  /* 1. test place-squarefreeness: in homog this suffices: isInV */
239  /* 2. test the presence of a hole -> in the tail??? */
240 
241 int isInV(poly p, int lV)
242 {
243  /* investigate only the leading monomial of p in currRing */
244  if ( pTotaldegree(p)==0 ) return(1);
245  if (lV <= 0) return(0);
246  /* returns 1 iff p is in V */
247  /* that is in each block up to a certain one there is only one nonzero exponent */
248  /* lV = the length of V = the number of orig vars */
249  int *e = (int *)omAlloc0((currRing->N+1)*sizeof(int));
250  int b = (int)((currRing->N +lV-1)/lV); /* the number of blocks */
251  //int b = (int)(currRing->N)/lV;
252  int *B = (int *)omAlloc0((b+1)*sizeof(int)); /* the num of elements in a block */
253  pGetExpV(p,e);
254  int i,j;
255  for (j=1; j<=b; j++)
256  {
257  /* we go through all the vars */
258  /* by blocks in lV vars */
259  for (i=(j-1)*lV + 1; i<= j*lV; i++)
260  {
261  if (e[i]) B[j] = B[j]+1;
262  }
263  }
264  // j = b;
265  // while ( (!B[j]) && (j>=1)) j--;
266  for (j=b; j>=1; j--)
267  {
268  if (B[j]!=0) break;
269  }
270  /* do not need e anymore */
271  freeT(e, currRing->N);
272 
273  if (j==0) goto ret_true;
274 // {
275 // /* it is a zero exp vector, which is in V */
276 // freeT(B, b);
277 // return(1);
278 // }
279  /* now B[j] != 0 and we test place-squarefreeness */
280  for (; j>=1; j--)
281  {
282  if (B[j]!=1)
283  {
284  freeT(B, b);
285  return(0);
286  }
287  }
288  ret_true:
289  freeT(B, b);
290  return(1);
291 }
292 
293 int poly_isInV(poly p, int lV)
294 {
295  /* tests whether the whole polynomial p in in V */
296  poly q = p;
297  while (q!=NULL)
298  {
299  if ( !isInV(q,lV) )
300  {
301  return(0);
302  }
303  q = pNext(q);
304  }
305  return(1);
306 }
307 
308 int ideal_isInV(ideal I, int lV)
309 {
310  /* tests whether each polynomial of an ideal I lies in in V */
311  int i;
312  int s = IDELEMS(I)-1;
313  for(i = 0; i <= s; i++)
314  {
315  if ( !poly_isInV(I->m[i],lV) )
316  {
317  return(0);
318  }
319  }
320  return(1);
321 }
322 
323 
324 int itoInsert(poly p, int uptodeg, int lV, const ring r)
325 {
326  /* for poly in lmCR/tailTR presentation */
327  /* the below situation (commented out) might happen! */
328 // if (r == currRing)
329 // {
330 // "Current ring is not expected in toInsert";
331 // return(0);
332 // }
333  /* compute the number of insertions */
334  int i = p_mLastVblock(p, lV, currRing);
335  if (pNext(p) != NULL)
336  {
337  i = si_max(i, p_LastVblock(pNext(p), lV, r) );
338  }
339  // i = uptodeg - i +1;
340  i = uptodeg - i;
341  // p_wrp(p,currRing,r); Print("----i:%d",i); PrintLn();
342  return(i);
343 }
344 
345 poly p_ShrinkT(poly p, int lV, kStrategy strat, const ring r)
346 //poly p_Shrink(poly p, int uptodeg, int lV, kStrategy strat, const ring r)
347 {
348  /* p is like TObject: lm in currRing = r, tail in tailRing */
349  /* proc shrinks the poly p in ring r */
350  /* lV = the length of V = the number of orig vars */
351  /* check assumes/exceptions */
352  /* r->N is a multiple of lV */
353 
354  if (p==NULL) return(p);
355 
358 
359  poly q = NULL;
360  poly s = p_mShrink(p, lV, r); // lm in currRing
361  poly pp = pNext(p);
362 
363  while (pp != NULL)
364  {
365  // q = p_Add_q(q, p_mShrink(pp,uptodeg,lV,strat->tailRing),strat->tailRing);
366  q = p_Add_q(q, p_mShrink(pp,lV,strat->tailRing),strat->tailRing);
367  pIter(pp);
368  }
369  pNext(s) = q;
370  return(s);
371 }
372 
373 poly p_Shrink(poly p, int lV, const ring r)
374 {
375  /* proc shrinks the poly p in ring r */
376  /* lV = the length of V = the number of orig vars */
377  /* check assumes/exceptions */
378  /* r->N is a multiple of lV */
379 
380  if (p==NULL) return(p);
382  poly q = NULL;
383  poly pp = p;
384 
385  while (pp != NULL)
386  {
387  q = p_Add_q(q, p_mShrink(pp,lV,r),r);
388  pIter(pp);
389  }
390  return(q);
391 }
392 
393 poly p_mShrink(poly p, int lV, const ring r)
394 {
395  /* shrinks the monomial p in ring r */
396  /* lV = the length of V = the number of orig vars */
397 
398  /* check assumes/exceptions */
399  /* r->N is a multiple of lV */
400 
401  int *e = (int *)omAlloc0((r->N+1)*sizeof(int));
402  int b = (int)((r->N +lV-1)/lV); /* the number of blocks */
403  // int *B = (int *)omAlloc0((b+1)*sizeof(int)); /* the num of elements in a block */
404  int *S = (int *)omAlloc0((r->N+1)*sizeof(int)); /* the shrinked exponent */
405  p_GetExpV(p,e,r);
406  int i,j; int cnt = 1; //counter for blocks in S
407  for (j=1; j<=b; j++)
408  {
409  /* we go through all the vars */
410  /* by blocks in lV vars */
411  for (i=(j-1)*lV + 1; i<= j*lV; i++)
412  {
413  if (e[i]==1)
414  {
415  // B[j] = B[j]+1; // for control in V?
416  S[(cnt-1)*lV + (i - (j-1)*lV)] = e[i];
417  /* assuming we are in V, can interrupt here */
418  cnt++;
419  // break; //results in incomplete shrink!
420  i = j*lV; // manual break under assumption p is in V
421  }
422  }
423  }
424 #ifdef PDEBUG
425  // Print("p_mShrink: cnt = [%d], b = %d\n",cnt,b);
426 #endif
427  // cnt -1 <= b must hold!
428  // freeT(B, b);
429  poly s = p_One(r);
430  p_SetExpV(s,S,r);
431  freeT(e, r->N);
432  freeT(S, r->N);
433  /* p_Setm(s,r); // done by p_SetExpV */
434  p_SetComp(s,p_GetComp(p,r),r); // component is preserved
435  p_SetCoeff(s,n_Copy(p_GetCoeff(p,r),r->cf),r); // coeff is preserved
436 #ifdef PDEBUG
437  // Print("p_mShrink: from "); p_wrp(p,r); Print(" to "); p_wrp(s,r); PrintLn();
438 #endif
439  return(s);
440 }
441 
442 /* shiftgb stuff */
443 
444 
445 /*2
446  *if the leading term of p
447  *divides the leading term of some T[i] it will be canceled
448  */
449 // static inline void clearSShift (poly p, unsigned long p_sev,int l, int* at, int* k,
450 // kStrategy strat)
451 // {
452 // assume(p_sev == pGetShortExpVector(p));
453 // if (!pLmShortDivisibleBy(p,p_sev, strat->T[*at].p, ~ strat->sevT[*at])) return;
454 // // if (l>=strat->lenS[*at]) return;
455 // if (TEST_OPT_PROT)
456 // PrintS("!");
457 // mflush();
458 // //pDelete(&strat->S[*at]);
459 // deleteInS((*at),strat);
460 // (*at)--;
461 // (*k)--;
462 // // assume(lenS_correct(strat));
463 // }
464 
465 /* remarks: cleanT : just deletion
466 enlargeT: just reallocation */
467 
468 #endif
int p_mLastVblock(poly p, int lV, const ring r)
Definition: shiftgb.cc:176
const CanonicalForm int s
Definition: facAbsFact.cc:55
int itoInsert(poly p, int uptodeg, int lV, const ring r)
Definition: shiftgb.cc:324
BOOLEAN p_LmCheckIsFromRing(poly p, ring r)
Definition: pDebug.cc:71
static int si_min(const int a, const int b)
Definition: auxiliary.h:121
Compatiblity layer for legacy polynomial operations (over currRing)
return P p
Definition: myNF.cc:203
poly p_ShrinkT(poly p, int lV, kStrategy strat, const ring r)
Definition: shiftgb.cc:345
int poly_isInV(poly p, int lV)
Definition: shiftgb.cc:293
int ideal_isInV(ideal I, int lV)
Definition: shiftgb.cc:308
static unsigned long p_SetComp(poly p, unsigned long c, ring r)
Definition: p_polys.h:242
#define p_GetComp(p, r)
Definition: monomials.h:72
static void p_GetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1443
poly p_mShrink(poly p, int lV, const ring r)
Definition: shiftgb.cc:393
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:407
static void p_SetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1451
poly pp
Definition: myNF.cc:296
static BOOLEAN p_LmIsConstant(const poly p, const ring r)
Definition: p_polys.h:949
poly p_LPshift(poly p, int sh, int uptodeg, int lV, const ring r)
Definition: shiftgb.cc:77
#define pIter(p)
Definition: monomials.h:44
poly p_Shrink(poly p, int lV, const ring r)
Definition: shiftgb.cc:373
int p_LastVblockT(poly p, int lV, kStrategy strat, const ring r)
Definition: shiftgb.cc:136
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
static poly p_Head(poly p, const ring r)
Definition: p_polys.h:812
const ring r
Definition: syzextra.cc:208
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:49
poly p_One(const ring r)
Definition: p_polys.cc:1314
int j
Definition: myNF.cc:70
static long pTotaldegree(poly p)
Definition: polys.h:264
#define assume(x)
Definition: mod2.h:394
static int si_max(const int a, const int b)
Definition: auxiliary.h:120
int i
Definition: cfEzgcd.cc:123
void PrintS(const char *s)
Definition: reporter.cc:284
int isInV(poly p, int lV)
Definition: shiftgb.cc:241
#define freeT(A, v)
Definition: shiftgb.cc:39
#define IDELEMS(i)
Definition: simpleideals.h:24
int int kStrategy strat
Definition: myNF.cc:68
int pFirstVblock(poly p, int lV)
Definition: shiftgb.cc:195
poly p_LPshiftT(poly p, int sh, int uptodeg, int lV, kStrategy strat, const ring r)
Definition: shiftgb.cc:45
BOOLEAN p_CheckIsFromRing(poly p, ring r)
Definition: pDebug.cc:101
#define NULL
Definition: omList.c:10
static FORCE_INLINE number n_Copy(number n, const coeffs r)
return a copy of &#39;n&#39;
Definition: coeffs.h:455
ring tailRing
Definition: kutil.h:331
b *CanonicalForm B
Definition: facBivar.cc:51
#define pGetExpV(p, e)
Gets a copy of (resp. set) the exponent vector, where e is assumed to point to (r->N +1)*sizeof(long)...
Definition: polys.h:96
poly p_mLPshift(poly p, int sh, int uptodeg, int lV, const ring r)
Definition: shiftgb.cc:100
#define pNext(p)
Definition: monomials.h:43
#define p_GetCoeff(p, r)
Definition: monomials.h:57
#define pIsConstantPoly(p)
return TRUE if all monomials of p are constant
Definition: polys.h:229
int p_LastVblock(poly p, int lV, const ring r)
Definition: shiftgb.cc:159
polyrec * poly
Definition: hilb.h:10
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:877
static Poly * h
Definition: janet.cc:978
const poly b
Definition: syzextra.cc:213
int pmFirstVblock(poly p, int lV)
Definition: shiftgb.cc:213
#define omAlloc0(size)
Definition: omAllocDecl.h:211