syz0.cc
Go to the documentation of this file.
1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /*
5 * ABSTRACT: resolutions
6 */
7 
8 
9 
10 
11 
12 #include <kernel/mod2.h>
13 #include <misc/options.h>
14 #include <omalloc/omalloc.h>
15 #include <kernel/polys.h>
16 #include <kernel/GBEngine/kstd1.h>
17 #include <kernel/GBEngine/kutil.h>
19 //#include "cntrlc.h"
20 #include <misc/intvec.h>
21 #include <coeffs/numbers.h>
22 #include <kernel/ideals.h>
23 #include <misc/intvec.h>
24 #include <polys/monomials/ring.h>
25 #include <kernel/GBEngine/syz.h>
26 #include <polys/kbuckets.h>
27 #include <polys/prCopy.h>
28 
29 static void syInitSort(ideal arg,intvec **modcomp)
30 {
31  int i,j,k,kk,kkk,jj;
32  idSkipZeroes(arg);
33  polyset F,oldF=arg->m;
34  int Fl=IDELEMS(arg);
35  int rkF=id_RankFreeModule(arg,currRing);
36  int syComponentOrder=currRing->ComponentOrder;
37 
38  while ((Fl!=0) && (oldF[Fl-1]==NULL)) Fl--;
39  if (*modcomp!=NULL) delete modcomp;
40  *modcomp = new intvec(rkF+2);
41  F=(polyset)omAlloc0(IDELEMS(arg)*sizeof(poly));
42  j=0;
43  for(i=0;i<=rkF;i++)
44  {
45  k=0;
46  jj = j;
47  (**modcomp)[i] = j;
48  while (k<Fl)
49  {
50  while ((k<Fl) && (pGetComp(oldF[k]) != i)) k++;
51  if (k<Fl)
52  {
53  kk=jj;
54  while ((kk<Fl) && (F[kk]) && (pLmCmp(oldF[k],F[kk])!=syComponentOrder))
55  {
56  kk++;
57  }
58  for (kkk=j;kkk>kk;kkk--)
59  {
60  F[kkk] = F[kkk-1];
61  }
62  F[kk] = oldF[k];
63 //Print("Element %d: ",kk);pWrite(F[kk]);
64  j++;
65  k++;
66  }
67  }
68  }
69  (**modcomp)[rkF+1] = Fl;
70  arg->m = F;
71  omFreeSize((ADDRESS)oldF,IDELEMS(arg)*sizeof(poly));
72 }
73 
74 static void syCreatePairs(polyset F,int lini,int wend,int k,int j,int i,
75  polyset pairs,int regularPairs=0,ideal mW=NULL)
76 {
77  int l,ii=0,jj;
78  poly p,q;
79 
80  while (((k<wend) && (pGetComp(F[k]) == i)) ||
81  ((currRing->qideal!=NULL) && (k<regularPairs+IDELEMS(currRing->qideal))))
82  {
83  p = pOne();
84  if ((k<wend) && (pGetComp(F[k]) == i) && (k!=j))
85  pLcm(F[j],F[k],p);
86  else if (ii<IDELEMS(currRing->qideal))
87  {
88  q = pHead(F[j]);
89  if (mW!=NULL)
90  {
91  for(jj=1;jj<=(currRing->N);jj++)
92  pSetExp(q,jj,pGetExp(q,jj) -pGetExp(mW->m[pGetComp(q)-1],jj));
93  pSetm(q);
94  }
95  pLcm(q,currRing->qideal->m[ii],p);
96  if (mW!=NULL)
97  {
98  for(jj=1;jj<=(currRing->N);jj++)
99  pSetExp(p,jj,pGetExp(p,jj) +pGetExp(mW->m[pGetComp(p)-1],jj));
100  pSetm(p);
101  }
102  pDelete(&q);
103  k = regularPairs+ii;
104  ii++;
105  }
106  l=lini;
107  while ((l<k) && ((pairs[l]==NULL) || (!pDivisibleBy(pairs[l],p))))
108  {
109  if ((pairs[l]!=NULL) && (pDivisibleBy(p,pairs[l])))
110  pDelete(&(pairs[l]));
111  l++;
112  }
113  if (l==k)
114  {
115  pSetm(p);
116  pairs[l] = p;
117  }
118  else
119  pDelete(&p);
120  k++;
121  }
122 }
123 
124 static poly syRedtail2(poly p, polyset redWith, intvec *modcomp)
125 {
126  poly h, hn;
127  int hncomp,nxt;
128  int j;
129 
130  h = p;
131  hn = pNext(h);
132  while(hn != NULL)
133  {
134  hncomp = pGetComp(hn);
135  j = (*modcomp)[hncomp];
136  nxt = (*modcomp)[hncomp+1];
137  while (j < nxt)
138  {
139  if (pLmDivisibleByNoComp(redWith[j], hn))
140  {
141  //if (TEST_OPT_PROT) PrintS("r");
142  hn = ksOldSpolyRed(redWith[j],hn);
143  if (hn == NULL)
144  {
145  pNext(h) = NULL;
146  return p;
147  }
148  hncomp = pGetComp(hn);
149  j = (*modcomp)[hncomp];
150  nxt = (*modcomp)[hncomp+1];
151  }
152  else
153  {
154  j++;
155  }
156  }
157  h = pNext(h) = hn;
158  hn = pNext(h);
159  }
160  return p;
161 }
162 
163 /*2
164 * computes the Schreyer syzygies in the local case
165 * input: arg (only allocated: Shdl, Smax)
166 * output: Shdl, Smax
167 */
168 static ideal sySchreyersSyzygiesFM(ideal arg,intvec ** modcomp)
169 {
170  int Fl=IDELEMS(arg);
171  while ((Fl!=0) && (arg->m[Fl-1]==NULL)) Fl--;
172  ideal result=idInit(16,arg->rank+Fl);
173  polyset F=arg->m,*Shdl=&(result->m);
174  if (Fl==0) return result;
175 
176  int i,j,l,k,totalToRed,ecartToRed,kk;
177  int bestEcart,totalmax,rkF,Sl=0,smax,tmax,tl;
178  int *ecartS, *ecartT, *totalS,
179  *totalT=NULL, *temp=NULL;
180  polyset pairs,S,T,ST/*,oldF*/;
181  poly p,q,toRed;
182  BOOLEAN notFound = FALSE;
183  intvec * newmodcomp = new intvec(Fl+2);
184  intvec * tempcomp;
185 
186 //Print("Naechster Modul\n");
187 //idPrint(arg);
188 /*-------------initializing the sets--------------------*/
189  ST=(polyset)omAlloc0(Fl*sizeof(poly));
190  S=(polyset)omAlloc0(Fl*sizeof(poly));
191  ecartS=(int*)omAlloc(Fl*sizeof(int));
192  totalS=(int*)omAlloc(Fl*sizeof(int));
193  T=(polyset)omAlloc0(2*Fl*sizeof(poly));
194  ecartT=(int*)omAlloc(2*Fl*sizeof(int));
195  totalT=(int*)omAlloc(2*Fl*sizeof(int));
196  pairs=(polyset)omAlloc0(Fl*sizeof(poly));
197 
198  smax = Fl;
199  tmax = 2*Fl;
200  for (j=1;j<IDELEMS(arg);j++)
201  {
202  if (arg->m[j] != NULL)
203  {
204  assume (arg->m[j-1] != NULL);
205  assume (pGetComp(arg->m[j-1])-pGetComp(arg->m[j])<=0);
206  }
207  }
208  rkF=id_RankFreeModule(arg,currRing);
209 /*----------------construction of the new ordering----------*/
210  if (rkF>0)
211  rSetSyzComp(rkF, currRing);
212  else
213  rSetSyzComp(1, currRing);
214 /*----------------creating S--------------------------------*/
215  for(j=0;j<Fl;j++)
216  {
217  S[j] = pCopy(F[j]);
218  totalS[j] = p_LDeg(S[j],&k,currRing);
219  ecartS[j] = totalS[j]-p_FDeg(S[j],currRing);
220 //Print("%d", pGetComp(S[j]));PrintS(" ");
221  p = S[j];
222  if (rkF==0) pSetCompP(p,1);
223  while (pNext(p)!=NULL) pIter(p);
224  pNext(p) = pHead(F[j]);
225  pIter(p);
226  if (rkF==0)
227  pSetComp(p,j+2);
228  else
229  pSetComp(p,rkF+j+1);
230  pSetmComp(p);
231  }
232 //PrintLn();
233  if (rkF==0) rkF = 1;
234 /*---------------creating the initial for T----------------*/
235  j=0;
236  l=-1;
237  totalmax=-1;
238  for (k=0;k<smax;k++)
239  if (totalS[k]>totalmax) totalmax=totalS[k];
240  for (kk=1;kk<=rkF;kk++)
241  {
242  for (k=0;k<=totalmax;k++)
243  {
244  for (l=0;l<smax;l++)
245  {
246  if ((pGetComp(S[l])==kk) && (totalS[l]==k))
247  {
248  ST[j] = S[l];
249  totalT[j] = totalS[l];
250  ecartT[j] = ecartS[l];
251 //Print("%d", totalS[l]);PrintS(" ");
252  j++;
253  }
254  }
255  }
256  }
257 //PrintLn();
258  for (j=0;j<smax;j++)
259  {
260  totalS[j] = totalT[j];
261  ecartS[j] = ecartT[j];
262  }
263 
264 /*---------------computing---------------------------------*/
265  for(j=0;j<smax;j++)
266  {
267  (*newmodcomp)[j+1] = Sl;
268  i = pGetComp(S[j]);
269  int syComponentOrder= currRing->ComponentOrder;
270  int lini,wend;
271  if (syComponentOrder==1)
272  {
273  lini=k=j+1;
274  wend=Fl;
275  }
276  else
277  {
278  lini=k=0;
279  while ((k<j) && (pGetComp(S[k]) != i)) k++;
280  wend=j;
281  }
282  if (TEST_OPT_PROT)
283  {
284  Print("(%d)",Fl-j);
285  mflush();
286  }
287  syCreatePairs(S,lini,wend,k,j,i,pairs);
288  for (k=lini;k<wend;k++)
289  {
290  if (pairs[k]!=NULL)
291  {
292 /*--------------creating T----------------------------------*/
293  for (l=0;l<smax;l++)
294  {
295  ecartT[l] = ecartS[l];
296  totalT[l] = totalS[l];
297  T[l] = ST[l];
298  }
299  tl = smax;
300  tempcomp = ivCopy(*modcomp);
301 /*--------------begin to reduce-----------------------------*/
302  toRed = ksOldCreateSpoly(S[j],S[k]);
303  ecartToRed = 1;
304  bestEcart = 1;
305  if (TEST_OPT_DEBUG)
306  {
307  PrintS("pair: ");pWrite0(S[j]);PrintS(" ");pWrite(S[k]);
308  }
309  if (TEST_OPT_PROT)
310  {
311  PrintS(".");
312  mflush();
313  }
314 //Print("Reduziere Paar %d,%d (ecart %d): \n",j,k,ecartToRed);
315 //Print("Poly %d: ",j);pWrite(S[j]);
316 //Print("Poly %d: ",k);pWrite(S[k]);
317 //Print("Spoly: ");pWrite(toRed);
318  while (pGetComp(toRed)<=rkF)
319  {
320  if (TEST_OPT_DEBUG)
321  {
322  PrintS("toRed: ");pWrite(toRed);
323  }
324 /*
325 * if ((bestEcart) || (ecartToRed!=0))
326 * {
327 */
328  totalToRed = p_LDeg(toRed,&kk,currRing);
329  ecartToRed = totalToRed-p_FDeg(toRed,currRing);
330 /*
331 * }
332 */
333 //Print("toRed now (neuer ecart %d): ",ecartToRed);pWrite(toRed);
334  notFound = TRUE;
335  bestEcart = 32000; //a very large integer
336  p = NULL;
337  int l=0;
338 #define OLD_SEARCH
339 #ifdef OLD_SEARCH
340  while ((l<tl) && (pGetComp(T[l])<pGetComp(toRed))) l++;
341  //assume (l==(**modcomp)[pGetComp(toRed)]);
342  while ((l<tl) && (notFound))
343 #else
344  l = (**modcomp)[pGetComp(toRed)];
345  int kkk = (**modcomp)[pGetComp(toRed)+1];
346  while ((l<kkk) && (notFound))
347 #endif
348  {
349  if ((ecartT[l]<bestEcart) && (pDivisibleBy(T[l],toRed)))
350  {
351  if (ecartT[l]<=ecartToRed) notFound = FALSE;
352  p = T[l];
353  bestEcart = ecartT[l];
354  }
355  l++;
356  }
357  if (p==NULL)
358  {
359  pDelete(&toRed);
360  for(k=j;k<Fl;k++) pDelete(&(pairs[k]));
361  omFreeSize((ADDRESS)pairs,Fl*sizeof(poly));
362  omFreeSize((ADDRESS)ST,Fl*sizeof(poly));
363  omFreeSize((ADDRESS)S,Fl*sizeof(poly));
364  omFreeSize((ADDRESS)T,tmax*sizeof(poly));
365  omFreeSize((ADDRESS)ecartT,tmax*sizeof(int));
366  omFreeSize((ADDRESS)totalT,tmax*sizeof(int));
367  omFreeSize((ADDRESS)ecartS,Fl*sizeof(int));
368  omFreeSize((ADDRESS)totalS,Fl*sizeof(int));
369  for(k=0;k<IDELEMS(result);k++) pDelete(&((*Shdl)[k]));
370  WerrorS("ideal not a standard basis");//no polynom for reduction
371  return result;
372  }
373  else
374  {
375 //Print("reduced with (ecart %d): ",bestEcart);wrp(p);PrintLn();
376  if (notFound)
377  {
378  if (tl>=tmax)
379  {
380  pEnlargeSet(&T,tmax,16);
381  tmax += 16;
382  temp = (int*)omAlloc((tmax+16)*sizeof(int));
383  for(l=0;l<tmax;l++) temp[l]=totalT[l];
384  totalT = temp;
385  temp = (int*)omAlloc((tmax+16)*sizeof(int));
386  for(l=0;l<tmax;l++) temp[l]=ecartT[l];
387  ecartT = temp;
388  }
389 //PrintS("t");
390  int comptR=pGetComp(toRed);
391  for (l=tempcomp->length()-1;l>comptR;l--)
392  {
393  if ((*tempcomp)[l]>0)
394  (*tempcomp)[l]++;
395  }
396  l=0;
397  while ((l<tl) && (comptR>pGetComp(T[l]))) l++;
398  while ((l<tl) && (totalT[l]<=totalToRed)) l++;
399  for (kk=tl;kk>l;kk--)
400  {
401  T[kk]=T[kk-1];
402  totalT[kk]=totalT[kk-1];
403  ecartT[kk]=ecartT[kk-1];
404  }
405  q = pCopy(toRed);
406  pNorm(q);
407  T[l] = q;
408  totalT[l] = totalToRed;
409  ecartT[l] = ecartToRed;
410  tl++;
411  }
412  toRed = ksOldSpolyRed(p,toRed);
413  }
414  }
415 //Print("toRed finally (neuer ecart %d): ",ecartToRed);pWrite(toRed);
416 //PrintS("s");
417  if (pGetComp(toRed)>rkF)
418  {
419  if (Sl>=IDELEMS(result))
420  {
421  pEnlargeSet(Shdl,IDELEMS(result),16);
422  IDELEMS(result) += 16;
423  }
424  //p_Shift(&toRed,-rkF,currRing);
425  pNorm(toRed);
426  (*Shdl)[Sl] = toRed;
427  Sl++;
428  }
429 /*----------------deleting all polys not from ST--------------*/
430  for(l=0;l<tl;l++)
431  {
432  kk=0;
433  while ((kk<smax) && (T[l] != S[kk])) kk++;
434  if (kk>=smax)
435  {
436  pDelete(&T[l]);
437 //Print ("#");
438  }
439  }
440  delete tempcomp;
441  }
442  }
443  for(k=lini;k<wend;k++) pDelete(&(pairs[k]));
444  }
445  (*newmodcomp)[Fl+1] = Sl;
446  omFreeSize((ADDRESS)pairs,Fl*sizeof(poly));
447  omFreeSize((ADDRESS)ST,Fl*sizeof(poly));
448  omFreeSize((ADDRESS)S,Fl*sizeof(poly));
449  omFreeSize((ADDRESS)T,tmax*sizeof(poly));
450  omFreeSize((ADDRESS)ecartT,tmax*sizeof(int));
451  omFreeSize((ADDRESS)totalT,tmax*sizeof(int));
452  omFreeSize((ADDRESS)ecartS,Fl*sizeof(int));
453  omFreeSize((ADDRESS)totalS,Fl*sizeof(int));
454  delete *modcomp;
455  *modcomp = newmodcomp;
456  return result;
457 }
458 
459 /*3
460 *special Normalform for Schreyer in factor rings
461 */
462 poly sySpecNormalize(poly toNorm,ideal mW=NULL)
463 {
464  int j,i=0;
465  poly p;
466 
467  if (toNorm==NULL) return NULL;
468  p = pHead(toNorm);
469  if (mW!=NULL)
470  {
471  for(j=1;j<=(currRing->N);j++)
472  pSetExp(p,j,pGetExp(p,j) -pGetExp(mW->m[pGetComp(p)-1],j));
473  }
474  while ((p!=NULL) && (i<IDELEMS(currRing->qideal)))
475  {
476  if (pDivisibleBy(currRing->qideal->m[i],p))
477  {
478  //pNorm(toNorm);
479  toNorm = ksOldSpolyRed(currRing->qideal->m[i],toNorm);
480  pDelete(&p);
481  if (toNorm==NULL) return NULL;
482  p = pHead(toNorm);
483  if (mW!=NULL)
484  {
485  for(j=1;j<=(currRing->N);j++)
486  pSetExp(p,j,pGetExp(p,j) -pGetExp(mW->m[pGetComp(p)-1],j));
487  }
488  i = 0;
489  }
490  else
491  {
492  i++;
493  }
494  }
495  pDelete(&p);
496  return toNorm;
497 }
498 
499 /*2
500 * computes the Schreyer syzygies in the global case
501 * input: F
502 * output: Shdl, Smax
503 * modcomp, length stores the start position of the module comp. in arg
504 */
505 static ideal sySchreyersSyzygiesFB(ideal arg,intvec ** modcomp,ideal mW,BOOLEAN redTail=TRUE)
506 {
507  kBucket_pt sy0buck = kBucketCreate(currRing);
508 
509  int Fl=IDELEMS(arg);
510  while ((Fl!=0) && (arg->m[Fl-1]==NULL)) Fl--;
511  ideal result=idInit(16,Fl);
512  int i,j,l,k,kkk,/*rkF,*/Sl=0,syComponentOrder=currRing->ComponentOrder;
513  int /*fstart,*/wend,lini,ltR,gencQ=0;
514  intvec *newmodcomp;
515  int *Flength;
516  polyset pairs,F=arg->m,*Shdl=&(result->m);
517  poly /*p,*/q,toRed,syz,lastmonom,multWith;
518  BOOLEAN isNotReduced=TRUE;
519 
520 //#define WRITE_BUCKETS
521 #ifdef WRITE_BUCKETS
522  PrintS("Input: \n");
523  ideal twr=idHead(arg);
524  idPrint(arg);
525  idDelete(&twr);
526  if (modcomp!=NULL) (*modcomp)->show(0,0);
527 #endif
528 
529  newmodcomp = new intvec(Fl+2);
530  //for (j=0;j<Fl;j++) pWrite(F[j]);
531  //PrintLn();
532  if (currRing->qideal==NULL)
533  pairs=(polyset)omAlloc0(Fl*sizeof(poly));
534  else
535  {
536  gencQ = IDELEMS(currRing->qideal);
537  pairs=(polyset)omAlloc0((Fl+gencQ)*sizeof(poly));
538  }
539  // rkF=id_RankFreeModule(arg,currRing);
540  Flength = (int*)omAlloc0(Fl*sizeof(int));
541  for(j=0;j<Fl;j++)
542  {
543  Flength[j] = pLength(F[j]);
544  }
545  for(j=0;j<Fl;j++)
546  {
547  (*newmodcomp)[j+1] = Sl;
548  if (TEST_OPT_PROT)
549  {
550  Print("(%d)",Fl-j);
551  mflush();
552  }
553  i = pGetComp(F[j]);
554  if (syComponentOrder==1)
555  {
556  lini=k=j+1;
557  wend=Fl;
558  }
559  else
560  {
561  lini=k=0;
562  while ((k<j) && (pGetComp(F[k]) != i)) k++;
563  wend=j;
564  }
565  syCreatePairs(F,lini,wend,k,j,i,pairs,Fl,mW);
566  if (currRing->qideal!=NULL) wend = Fl+gencQ;
567  for (k=lini;k<wend;k++)
568  {
569  if (pairs[k]!=NULL)
570  {
571  if (TEST_OPT_PROT)
572  {
573  PrintS(".");
574  mflush();
575  }
576  //begins to construct the syzygy
577  if (k<Fl)
578  {
579  number an=nCopy(pGetCoeff(F[k])),bn=nCopy(pGetCoeff(F[j]));
580  /*int ct =*/ (void) ksCheckCoeff(&an, &bn, currRing->cf);
581  syz = pCopy(pairs[k]);
582  //syz->coef = nCopy(F[k]->coef);
583  syz->coef = an;
584  //syz->coef = nInpNeg(syz->coef);
585  pNext(syz) = pairs[k];
586  lastmonom = pNext(syz);
587  //lastmonom->coef = nCopy(F[j]->coef);
588  lastmonom->coef = bn;
589  lastmonom->coef = nInpNeg(lastmonom->coef);
590  pSetComp(lastmonom,k+1);
591  }
592  else
593  {
594  syz = pairs[k];
595  syz->coef = nCopy(currRing->qideal->m[k-Fl]->coef);
596  syz->coef = nInpNeg(syz->coef);
597  lastmonom = syz;
598  multWith = pDivide(syz,F[j]);
599  multWith->coef = nCopy(currRing->qideal->m[k-Fl]->coef);
600  }
601  pSetComp(syz,j+1);
602  pairs[k] = NULL;
603  //the next term of the syzygy
604  //constructs the spoly
605  if (TEST_OPT_DEBUG)
606  {
607  if (k<Fl)
608  {
609  PrintS("pair: ");pWrite0(F[j]);PrintS(" ");pWrite(F[k]);
610  }
611  else
612  {
613  PrintS("pair: ");pWrite0(F[j]);PrintS(" ");pWrite(currRing->qideal->m[k-Fl]);
614  }
615  }
616  if (k<Fl)
617  toRed = ksOldCreateSpoly(F[j],F[k]);
618  else
619  {
620  q = pMult_mm(pCopy(F[j]),multWith);
621  toRed = sySpecNormalize(q,mW);
622  pDelete(&multWith);
623  }
624  kBucketInit(sy0buck,toRed,-1);
625  toRed = kBucketGetLm(sy0buck);
626  isNotReduced = TRUE;
627  while (toRed!=NULL)
628  {
629  if (TEST_OPT_DEBUG)
630  {
631  PrintS("toRed: ");pWrite(toRed);
632  }
633 // l=0;
634 // while ((l<Fl) && (!pDivisibleBy(F[l],toRed))) l++;
635 // if (l>=Fl)
636  l = (**modcomp)[pGetComp(toRed)+1]-1;
637  kkk = (**modcomp)[pGetComp(toRed)];
638  while ((l>=kkk) && (!pDivisibleBy(F[l],toRed))) l--;
639 #ifdef WRITE_BUCKETS
640  kBucketClear(sy0buck,&toRed,&ltR);
641  printf("toRed in Pair[%d, %d]:", j, k);
642  pWrite(toRed);
643  kBucketInit(sy0buck,toRed,-1);
644 #endif
645 
646  if (l<kkk)
647  {
648  if ((currRing->qideal!=NULL) && (isNotReduced))
649  {
650  kBucketClear(sy0buck,&toRed,&ltR);
651  toRed = sySpecNormalize(toRed,mW);
652 #ifdef WRITE_BUCKETS
653  printf("toRed in Pair[%d, %d]:", j, k);
654  pWrite(toRed);
655 #endif
656  kBucketInit(sy0buck,toRed,-1);
657  toRed = kBucketGetLm(sy0buck);
658  isNotReduced = FALSE;
659  }
660  else
661  {
662  pDelete(&toRed);
663 
664  pDelete(&syz);
665  for(k=j;k<Fl;k++) pDelete(&(pairs[k]));
666  omFreeSize((ADDRESS)pairs,(Fl + gencQ)*sizeof(poly));
667 
668 
669  for(k=0;k<IDELEMS(result);k++) pDelete(&((*Shdl)[k]));
670 
671  kBucketDestroy(&(sy0buck));
672 
673  //no polynom for reduction
674  WerrorS("ideal not a standard basis");
675 
676  return result;
677  }
678  }
679  else
680  {
681  //the next monom of the syzygy
682  isNotReduced = TRUE;
683  if (TEST_OPT_DEBUG)
684  {
685  PrintS("reduced with: ");pWrite(F[l]);
686  }
687  pNext(lastmonom) = pHead(toRed);
688  pIter(lastmonom);
689  lastmonom->coef = nDiv(lastmonom->coef,F[l]->coef);
690  //lastmonom->coef = nInpNeg(lastmonom->coef);
691  pSetComp(lastmonom,l+1);
692  //computes the new toRed
693  number up = kBucketPolyRed(sy0buck,F[l],Flength[l],NULL);
694  if (! nIsOne(up))
695  {
696  // Thomas: Now do whatever you need to do
697 #ifdef WRITE_BUCKETS
698  PrintS("multiplied with: ");nWrite(up);PrintLn();
699 #endif
700  pMult_nn(syz,up);
701  }
702  nDelete(&up);
703 
704  toRed = kBucketGetLm(sy0buck);
705  //the module component of the new monom
706 //pWrite(toRed);
707  }
708  }
709  kBucketClear(sy0buck,&toRed,&ltR); //Zur Sichereheit
710 //PrintLn();
711  if (syz!=NULL)
712  {
713  if (Sl>=IDELEMS(result))
714  {
715  pEnlargeSet(Shdl,IDELEMS(result),16);
716  IDELEMS(result) += 16;
717  }
718  pNorm(syz);
719  if (BTEST1(OPT_REDTAIL) && redTail)
720  {
721  (*newmodcomp)[j+2] = Sl;
722  (*Shdl)[Sl] = syRedtail2(syz,*Shdl,newmodcomp);
723  (*newmodcomp)[j+2] = 0;
724  }
725  else
726  (*Shdl)[Sl] = syz;
727  Sl++;
728  }
729  }
730  }
731 // for(k=j;k<Fl;k++) pDelete(&(pairs[k]));
732  }
733  (*newmodcomp)[Fl+1] = Sl;
734  if (currRing->qideal==NULL)
735  omFreeSize((ADDRESS)pairs,Fl*sizeof(poly));
736  else
737  omFreeSize((ADDRESS)pairs,(Fl+IDELEMS(currRing->qideal))*sizeof(poly));
738  omFreeSize((ADDRESS)Flength,Fl*sizeof(int));
739  delete *modcomp;
740  *modcomp = newmodcomp;
741 
742  kBucketDestroy(&(sy0buck));
743  return result;
744 }
745 
747 {
748  int syzIndex=length-1,i,j;
749  poly p;
750 
751  while ((syzIndex!=0) && (res[syzIndex]==NULL)) syzIndex--;
752  while (syzIndex>=initial)
753  {
754  for(i=0;i<IDELEMS(res[syzIndex]);i++)
755  {
756  p = res[syzIndex]->m[i];
757 
758  while (p!=NULL)
759  {
760  if (res[syzIndex-1]->m[pGetComp(p)-1]!=NULL)
761  {
762  for(j=1;j<=(currRing->N);j++)
763  {
764  pSetExp(p,j,pGetExp(p,j)
765  -pGetExp(res[syzIndex-1]->m[pGetComp(p)-1],j));
766  }
767  }
768  else
769  PrintS("error in the resolvent\n");
770  pSetm(p);
771  pIter(p);
772  }
773  }
774  syzIndex--;
775  }
776 }
777 
778 #if 0
779 static void syMergeSortResolventFB(resolvente res,int length, int initial=1)
780 {
781  int syzIndex=length-1,i,j;
782  poly qq,pp,result=NULL;
783  poly p;
784 
785  while ((syzIndex!=0) && (res[syzIndex]==NULL)) syzIndex--;
786  while (syzIndex>=initial)
787  {
788  for(i=0;i<IDELEMS(res[syzIndex]);i++)
789  {
790  p = res[syzIndex]->m[i];
791  if (p != NULL)
792  {
793  for (;;)
794  {
795  qq = p;
796  for(j=1;j<=(currRing->N);j++)
797  {
798  pSetExp(p,j,pGetExp(p,j)
799  -pGetExp(res[syzIndex-1]->m[pGetComp(p)-1],j));
800  }
801  pSetm(p);
802  for (;;)
803  {
804  if (pNext(p) == NULL)
805  {
806  pAdd(result, qq);
807  break;
808  }
809  pp = pNext(p);
810  for(j=1;j<=(currRing->N);j++)
811  {
813  -pGetExp(res[syzIndex-1]->m[pGetComp(pp)-1],j));
814  }
815  pSetm(pp);
816  if (pCmp(p,pNext(p)) != 1)
817  {
818  pp = p;
819  pIter(p);
820  pNext(pp) = NULL;
821  result = pAdd(result, qq);
822  break;
823  }
824  pIter(p);
825  }
826  }
827  }
828  res[syzIndex]->m[i] = p;
829  }
830  syzIndex--;
831  }
832 }
833 #endif
834 
836 {
838  if (i == 0) return FALSE;
839  int j=0;
840 
841  while ((currRing->order[j]!=ringorder_c) && (currRing->order[j]!=ringorder_C))
842  j++;
843  if (currRing->order[j+1]!=0)
844  return TRUE;
845  return FALSE;
846 }
847 
848 #if 0 /*debug only */
849 static void syPrintResolution(resolvente res,int start,int length)
850 {
851  while ((start < length) && (res[start]))
852  {
853  Print("Syz(%d): \n",start);
854  idTest(res[start]);
855  //idPrint(res[start]);
856  start++;
857  }
858 }
859 #endif
860 
861 resolvente sySchreyerResolvente(ideal arg, int maxlength, int * length,
862  BOOLEAN isMonomial, BOOLEAN /*notReplace*/)
863 {
864  ideal mW=NULL;
865  int i,syzIndex = 0,j=0;
866  intvec * modcomp=NULL,*w=NULL;
867  // int ** wv=NULL;
868  tHomog hom=(tHomog)idHomModule(arg,NULL,&w);
869  ring origR = currRing;
870  ring syRing = NULL;
871 
872  if ((!isMonomial) && syTestOrder(arg))
873  {
874  WerrorS("sres only implemented for modules with ordering ..,c or ..,C");
875  return NULL;
876  }
877  *length = 4;
878  resolvente res = (resolvente)omAlloc0(4*sizeof(ideal)),newres;
879  res[0] = idCopy(arg);
880 
881  while ((!idIs0(res[syzIndex])) && ((maxlength==-1) || (syzIndex<maxlength)))
882  {
883  i = IDELEMS(res[syzIndex]);
884  //while ((i!=0) && (!res[syzIndex]->m[i-1])) i--;
885  if (syzIndex+1==*length)
886  {
887  newres = (resolvente)omAlloc0((*length+4)*sizeof(ideal));
888  // for (j=0;j<*length+4;j++) newres[j] = NULL;
889  for (j=0;j<*length;j++) newres[j] = res[j];
890  omFreeSize((ADDRESS)res,*length*sizeof(ideal));
891  *length += 4;
892  res=newres;
893  }
894 
895  if ((hom==isHomog)|| (rHasGlobalOrdering(origR)))
896  {
897  if (syzIndex==0) syInitSort(res[0],&modcomp);
898 
899  if ((syzIndex==0) && !rRing_has_CompLastBlock(currRing))
900  res[syzIndex+1] = sySchreyersSyzygiesFB(res[syzIndex],&modcomp,mW,FALSE);
901  else
902  res[syzIndex+1] = sySchreyersSyzygiesFB(res[syzIndex],&modcomp,mW);
903 
904  if (errorreported)
905  {
906  for (j=0;j<*length;j++) idDelete( &res[j] );
907  omFreeSize((ADDRESS)res,*length*sizeof(ideal));
908  return NULL;
909  }
910 
911  mW = res[syzIndex];
912  }
913 //idPrint(res[syzIndex+1]);
914 
915  if ( /*(*/ syzIndex==0 /*)*/ )
916  {
917  if ((hom==isHomog)|| (rHasGlobalOrdering(origR)))
918  {
919  syRing = rAssure_CompLastBlock(origR, TRUE);
920  if (syRing != origR)
921  {
922  rChangeCurrRing(syRing);
923  for (i=0; i<IDELEMS(res[1]); i++)
924  {
925  res[1]->m[i] = prMoveR( res[1]->m[i], origR, syRing);
926  }
927  }
928  idTest(res[1]);
929  }
930  else
931  {
932  syRing = rAssure_SyzComp_CompLastBlock(origR);
933  if (syRing != origR)
934  {
935  rChangeCurrRing(syRing);
936  for (i=0; i<IDELEMS(res[0]); i++)
937  {
938  res[0]->m[i] = prMoveR( res[0]->m[i], origR, syRing);
939  }
940  }
941  idTest(res[0]);
942  }
943  }
944  if ((hom!=isHomog) && (rHasLocalOrMixedOrdering(origR)))
945  {
946  if (syzIndex==0) syInitSort(res[0],&modcomp);
947  res[syzIndex+1] = sySchreyersSyzygiesFM(res[syzIndex],&modcomp);
948  if (errorreported)
949  {
950  for (j=0;j<*length;j++) idDelete( &res[j] );
951  omFreeSize((ADDRESS)res,*length*sizeof(ideal));
952  return NULL;
953  }
954  }
955  syzIndex++;
956  if (TEST_OPT_PROT) Print("[%d]\n",syzIndex);
957  }
958  //syPrintResolution(res,1,*length);
959  if ((hom!=isHomog) && (rHasLocalOrMixedOrdering(origR)))
960  {
961  syzIndex = 1;
962  while ((syzIndex < *length) && (!idIs0(res[syzIndex])))
963  {
964  id_Shift(res[syzIndex],-rGetMaxSyzComp(syzIndex, currRing),currRing);
965  syzIndex++;
966  }
967  }
968  if ((hom==isHomog) || (rHasGlobalOrdering(origR)))
969  syzIndex = 1;
970  else
971  syzIndex = 0;
972  syReOrderResolventFB(res,*length,syzIndex+1);
973  if (/*ringOrderChanged:*/ origR!=syRing && syRing != NULL)
974  {
975  rChangeCurrRing(origR);
976  // Thomas: Here I assume that all (!) polys of res live in tmpR
977  while ((syzIndex < *length) && (res[syzIndex]))
978  {
979  for (i=0;i<IDELEMS(res[syzIndex]);i++)
980  {
981  if (res[syzIndex]->m[i])
982  {
983  res[syzIndex]->m[i] = prMoveR( res[syzIndex]->m[i], syRing, origR);
984  }
985  }
986  syzIndex++;
987  }
988 // j = 0; while (currRing->order[j]!=0) j++; // What was this for???!
989  rDelete(syRing);
990  }
991  else
992  {
993  // Thomas -- are you sure that you have to "reorder" here?
994  while ((syzIndex < *length) && (res[syzIndex]))
995  {
996  for (i=0;i<IDELEMS(res[syzIndex]);i++)
997  {
998  if (res[syzIndex]->m[i])
999  res[syzIndex]->m[i] = pSortCompCorrect(res[syzIndex]->m[i]);
1000  }
1001  syzIndex++;
1002  }
1003  }
1004  if ((hom==isHomog) || (rHasGlobalOrdering(origR)))
1005  {
1006  if (res[1]!=NULL)
1007  {
1009  for (i=0;i<IDELEMS(res[1]);i++)
1010  {
1011  if (res[1]->m[i])
1012  res[1]->m[i] = pSort(res[1]->m[i]);
1013  }
1014  }
1015  }
1016  //syPrintResolution(res,0,*length);
1017 
1018  //syMergeSortResolventFB(res,*length);
1019  if (modcomp!=NULL) delete modcomp;
1020  if (w!=NULL) delete w;
1021  return res;
1022 }
1023 
1024 syStrategy sySchreyer(ideal arg, int maxlength)
1025 {
1026  int rl;
1027  resolvente fr = sySchreyerResolvente(arg,maxlength,&(rl));
1028  if (fr==NULL) return NULL;
1029 
1030  // int typ0;
1032  result->length=rl;
1033  result->fullres = (resolvente)omAlloc0((rl /*result->length*/+1)*sizeof(ideal));
1034  for (int i=rl /*result->length*/-1;i>=0;i--)
1035  {
1036  if (fr[i]!=NULL)
1037  {
1038  result->fullres[i] = fr[i];
1039  fr[i] = NULL;
1040  }
1041  }
1042  if (currRing->qideal!=NULL)
1043  {
1044  for (int i=0; i<rl; i++)
1045  {
1046  if (result->fullres[i]!=NULL)
1047  {
1048  ideal t=kNF(currRing->qideal,NULL,result->fullres[i]);
1049  idDelete(&result->fullres[i]);
1050  result->fullres[i]=t;
1051  if (i<rl-1)
1052  {
1053  for(int j=IDELEMS(t)-1;j>=0; j--)
1054  {
1055  if ((t->m[j]==NULL) && (result->fullres[i+1]!=NULL))
1056  {
1057  for(int k=IDELEMS(result->fullres[i+1])-1;k>=0; k--)
1058  {
1059  if (result->fullres[i+1]->m[k]!=NULL)
1060  {
1061  pDeleteComp(&(result->fullres[i+1]->m[k]),j+1);
1062  }
1063  }
1064  }
1065  }
1066  }
1067  idSkipZeroes(result->fullres[i]);
1068  }
1069  }
1070  if ((rl>maxlength) && (result->fullres[rl-1]!=NULL))
1071  {
1072  idDelete(&result->fullres[rl-1]);
1073  }
1074  }
1075  omFreeSize((ADDRESS)fr,(rl /*result->length*/)*sizeof(ideal));
1076  return result;
1077 }
1078 
BOOLEAN rHasLocalOrMixedOrdering(const ring r)
Definition: ring.h:752
void syReOrderResolventFB(resolvente res, int length, int initial)
Definition: syz0.cc:746
#define pSetmComp(p)
TODO:
Definition: polys.h:255
void kBucketClear(kBucket_pt bucket, poly *p, int *length)
Definition: kbuckets.cc:495
resolvente sySchreyerResolvente(ideal arg, int maxlength, int *length, BOOLEAN isMonomial, BOOLEAN)
Definition: syz0.cc:861
#define pDivide(a, b)
Definition: polys.h:275
#define nWrite(n)
Definition: numbers.h:29
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
Definition: kstd1.cc:2971
#define pSetm(p)
Definition: polys.h:253
void kBucketInit(kBucket_pt bucket, poly lm, int length)
Definition: kbuckets.cc:467
void PrintLn()
Definition: reporter.cc:310
#define Print
Definition: emacs.cc:83
#define pAdd(p, q)
Definition: polys.h:186
#define idDelete(H)
delete an ideal
Definition: ideals.h:29
template CanonicalForm tmax(const CanonicalForm &, const CanonicalForm &)
#define TEST_OPT_PROT
Definition: options.h:98
static poly syRedtail2(poly p, polyset redWith, intvec *modcomp)
Definition: syz0.cc:124
#define pSetExp(p, i, v)
Definition: polys.h:42
#define FALSE
Definition: auxiliary.h:94
Compatiblity layer for legacy polynomial operations (over currRing)
return P p
Definition: myNF.cc:203
number kBucketPolyRed(kBucket_pt bucket, poly p1, int l1, poly spNoether)
Definition: kbuckets.cc:1053
syStrategy sySchreyer(ideal arg, int maxlength)
Definition: syz0.cc:1024
#define pLmCmp(p, q)
returns 0|1|-1 if p=q|p>q|p<q w.r.t monomial ordering
Definition: polys.h:105
poly prMoveR(poly &p, ring src_r, ring dest_r)
Definition: prCopy.cc:91
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
const poly kBucketGetLm(kBucket_pt bucket)
Definition: kbuckets.cc:480
int ksCheckCoeff(number *a, number *b)
intvec * ivCopy(const intvec *o)
Definition: intvec.h:126
#define pSort(p)
Definition: polys.h:201
KINLINE poly ksOldCreateSpoly(poly p1, poly p2, poly spNoether, ring r)
Definition: kInline.h:1073
#define BTEST1(a)
Definition: options.h:32
#define pCmp(p1, p2)
pCmp: args may be NULL returns: (p2==NULL ? 1 : (p1 == NULL ? -1 : p_LmCmp(p1, p2))) ...
Definition: polys.h:115
#define TRUE
Definition: auxiliary.h:98
#define nIsOne(n)
Definition: numbers.h:25
#define pLcm(a, b, m)
Definition: polys.h:277
void * ADDRESS
Definition: auxiliary.h:115
void pWrite(poly p)
Definition: polys.h:290
void WerrorS(const char *s)
Definition: feFopen.cc:24
int k
Definition: cfEzgcd.cc:93
#define TEST_OPT_DEBUG
Definition: options.h:103
static ideal sySchreyersSyzygiesFM(ideal arg, intvec **modcomp)
Definition: syz0.cc:168
#define pMult_mm(p, m)
Definition: polys.h:185
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
Definition: monomials.h:51
#define omAlloc(size)
Definition: omAllocDecl.h:210
void pWrite0(poly p)
Definition: polys.h:291
ring rAssure_SyzComp_CompLastBlock(const ring r)
makes sure that c/C ordering is last ordering and SyzIndex is first
Definition: ring.cc:4667
#define pGetComp(p)
Component.
Definition: polys.h:37
poly pp
Definition: myNF.cc:296
ideal idHead(ideal h)
#define mflush()
Definition: reporter.h:57
#define pIter(p)
Definition: monomials.h:44
poly res
Definition: myNF.cc:322
static BOOLEAN idHomModule(ideal m, ideal Q, intvec **w)
Definition: ideals.h:96
#define M
Definition: sirandom.c:24
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
#define pGetExp(p, i)
Exponent.
Definition: polys.h:41
poly initial(const poly p, const ring r, const gfan::ZVector &w)
Returns the initial form of p with respect to w.
Definition: initial.cc:32
void id_Shift(ideal M, int s, const ring r)
#define idPrint(id)
Definition: ideals.h:46
void kBucketDestroy(kBucket_pt *bucket_pt)
Definition: kbuckets.cc:200
#define pSortCompCorrect(p)
Assume: If considerd only as poly in any component of p (say, monomials of other components of p are ...
Definition: polys.h:210
Definition: intvec.h:14
long id_RankFreeModule(ideal s, ring lmRing, ring tailRing)
return the maximal component number found in any polynomial in s
#define OPT_REDTAIL
Definition: options.h:86
tHomog
Definition: structs.h:37
int j
Definition: myNF.cc:70
void pairs()
#define pLmDivisibleByNoComp(a, b)
like pLmDivisibleBy, does not check components
Definition: polys.h:142
#define pSetCompP(a, i)
Definition: polys.h:285
#define assume(x)
Definition: mod2.h:394
int rGetMaxSyzComp(int i, const ring r)
return the max-comonent wchich has syzIndex i Assume: i<= syzIndex_limit
Definition: ring.cc:5061
#define nInpNeg(n)
Definition: numbers.h:21
static void syInitSort(ideal arg, intvec **modcomp)
Definition: syz0.cc:29
static long p_FDeg(const poly p, const ring r)
Definition: p_polys.h:375
#define pSetComp(p, v)
Definition: polys.h:38
int m
Definition: cfEzgcd.cc:119
#define pMult_nn(p, n)
Definition: polys.h:183
int i
Definition: cfEzgcd.cc:123
void PrintS(const char *s)
Definition: reporter.cc:284
ring rAssure_CompLastBlock(ring r, BOOLEAN complete)
makes sure that c/C ordering is last ordering
Definition: ring.cc:4612
static long p_LDeg(const poly p, int *l, const ring r)
Definition: p_polys.h:376
#define pOne()
Definition: polys.h:297
static unsigned pLength(poly a)
Definition: p_polys.h:189
#define pHead(p)
returns newly allocated copy of Lm(p), coef is copied, next=NULL, p might be NULL ...
Definition: polys.h:67
#define IDELEMS(i)
Definition: simpleideals.h:24
void idSkipZeroes(ideal ide)
gives an ideal/module the minimal possible size
BOOLEAN syTestOrder(ideal M)
Definition: syz0.cc:835
#define nDelete(n)
Definition: numbers.h:16
ideal idCopy(ideal A)
Definition: ideals.h:60
short errorreported
Definition: feFopen.cc:23
void rSetSyzComp(int k, const ring r)
Definition: ring.cc:4989
poly sySpecNormalize(poly toNorm, ideal mW=NULL)
Definition: syz0.cc:462
void rChangeCurrRing(ring r)
Definition: polys.cc:12
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:38
#define nDiv(a, b)
Definition: numbers.h:32
#define NULL
Definition: omList.c:10
poly * polyset
Definition: hutil.h:15
#define pDivisibleBy(a, b)
returns TRUE, if leading monom of a divides leading monom of b i.e., if there exists a expvector c > ...
Definition: polys.h:138
void pEnlargeSet(poly **p, int l, int increment)
Definition: p_polys.cc:3602
int length() const
Definition: intvec.h:86
BOOLEAN rHasGlobalOrdering(const ring r)
Definition: ring.h:751
void rDelete(ring r)
unconditionally deletes fields in r
Definition: ring.cc:448
void pNorm(poly p, const ring R=currRing)
Definition: polys.h:345
const CanonicalForm & w
Definition: facAbsFact.cc:55
#define pDelete(p_ptr)
Definition: polys.h:169
#define nCopy(n)
Definition: numbers.h:15
static void syCreatePairs(polyset F, int lini, int wend, int k, int j, int i, polyset pairs, int regularPairs=0, ideal mW=NULL)
Definition: syz0.cc:74
#define pNext(p)
Definition: monomials.h:43
ideal * resolvente
Definition: ideals.h:18
KINLINE poly ksOldSpolyRed(poly p1, poly p2, poly spNoether)
Definition: kInline.h:1053
static jList * T
Definition: janet.cc:37
polyrec * poly
Definition: hilb.h:10
#define pDeleteComp(p, k)
Definition: polys.h:343
kBucket_pt kBucketCreate(const ring bucket_ring)
Creation/Destruction of buckets
Definition: kbuckets.cc:193
static Poly * h
Definition: janet.cc:978
int BOOLEAN
Definition: auxiliary.h:85
BOOLEAN idIs0(ideal h)
returns true if h is the zero ideal
BOOLEAN rRing_has_CompLastBlock(ring r)
Definition: ring.cc:5107
#define omAlloc0(size)
Definition: omAllocDecl.h:211
return result
Definition: facAbsBiFact.cc:76
int l
Definition: cfEzgcd.cc:94
#define pCopy(p)
return a copy of the poly
Definition: polys.h:168
ssyStrategy * syStrategy
Definition: syz.h:35
#define idTest(id)
Definition: ideals.h:47
static ideal sySchreyersSyzygiesFB(ideal arg, intvec **modcomp, ideal mW, BOOLEAN redTail=TRUE)
Definition: syz0.cc:505