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