1 /****************************************
2 *  Computer Algebra System SINGULAR     *
3 ****************************************/
4 /*
5 * ABSTRACT: resolutions
6 */
7 
8 #include "kernel/mod2.h"
9 
10 #include "misc/mylimits.h"
11 #include "misc/options.h"
12 #include "misc/intvec.h"
13 
14 #include "coeffs/numbers.h"
15 
16 #include "polys/monomials/ring.h"
17 #include "polys/kbuckets.h"
18 #include "polys/prCopy.h"
19 #include "polys/matpol.h"
20 
21 #include "kernel/combinatorics/stairc.h"
22 #include "kernel/combinatorics/hilb.h"
23 
24 #include "kernel/GBEngine/kstd1.h"
25 #include "kernel/GBEngine/kutil.h"
26 #include "kernel/GBEngine/syz.h"
27 
28 #include "kernel/ideals.h"
29 #include "kernel/polys.h"
30 
31 //#define SHOW_PROT
32 //#define SHOW_RED
33 //#define SHOW_Kosz
34 //#define SHOW_RESULT
35 //#define INVERT_PAIRS
36 //#define ONLY_STD
37 //#define EXPERIMENT1    //Hier stimmt was mit der Anzahl der Erzeuger in xyz11 nicht!!
38 #define EXPERIMENT2
39 #define EXPERIMENT3
40 #define WITH_BUCKET     //Use of buckets in EXPERIMENT3 (Product criterion)
41 #define WITH_SCHREYER_ORD
42 #define USE_CHAINCRIT
43 #define USE_CHAINCRIT0
44 #define USE_PROD_CRIT
45 #define USE_REGULARITY
46 #define WITH_SORT
47 //#define FULL_TOTAKE
48 VAR int discard_pairs;
49 VAR int short_pairs;
50 
51 /*3
52 * assumes the ideals old_ideal and new_ideal to be homogeneous
53 * tests wether the new_ideal is a regular extension of the old_ideal
54 */
syIsRegular(ideal old_ideal,ideal new_ideal,int deg)55 static BOOLEAN syIsRegular(ideal old_ideal,ideal new_ideal,int deg)
56 {
57   intvec * old_hilbs=hHstdSeries(old_ideal,NULL,NULL,NULL);
58   intvec * new_hilbs=hHstdSeries(new_ideal,NULL,NULL,NULL);
59   int biggest_length=si_max(old_hilbs->length()+deg,new_hilbs->length());
60   intvec * shifted_old_hilbs=new intvec(biggest_length);
61   intvec * old_hilb1=new intvec(biggest_length);
62   intvec * new_hilb1=new intvec(biggest_length);
63   int i;
64   BOOLEAN isRegular=TRUE;
65 
66   for (i=old_hilbs->length()+deg-1;i>=deg;i--)
67     (*shifted_old_hilbs)[i] = (*old_hilbs)[i-deg];
68   for (i=old_hilbs->length()-1;i>=0;i--)
69     (*old_hilb1)[i] = (*old_hilbs)[i]-(*shifted_old_hilbs)[i];
70   for (i=old_hilbs->length()+deg-1;i>=old_hilbs->length();i--)
71     (*old_hilb1)[i] = -(*shifted_old_hilbs)[i];
72   for (i=new_hilbs->length()-1;i>=0;i--)
73     (*new_hilb1)[i] = (*new_hilbs)[i];
74   i = 0;
75   while ((i<biggest_length) && isRegular)
76   {
77     isRegular = isRegular && ((*old_hilb1)[i] == (*new_hilb1)[i]);
78     i++;
79   }
80   delete old_hilbs;
81   delete new_hilbs;
82   delete old_hilb1;
83   delete new_hilb1;
84   delete shifted_old_hilbs;
85   return isRegular;
86 }
87 
88 /*3
89 * shows the resolution stored in syzstr->orderedRes
90 */
91 #if 0 /* unused*/
92 static void syShowRes(syStrategy syzstr)
93 {
94   int i=0;
95 
96   while ((i<syzstr->length) && (!idIs0(syzstr->res[i])))
97   {
98     Print("aktueller hoechster index ist: %d\n",(*syzstr->Tl)[i]);
99     Print("der %d-te modul ist:\n",i);
100     idPrint(syzstr->res[i]);
101     PrintS("Seine Darstellung:\n");
102     idPrint(syzstr->orderedRes[i]);
103     i++;
104   }
105 }
106 #endif
107 
108 /*3
109 * produces the next subresolution for a regular extension
110 */
syCreateRegularExtension(syStrategy syzstr,ideal old_ideal,ideal old_repr,int old_tl,poly next_generator,resolvente totake)111 static void syCreateRegularExtension(syStrategy syzstr,ideal old_ideal,
112             ideal old_repr,int old_tl, poly next_generator,resolvente totake)
113 {
114   int index=syzstr->length-1,i,j,start,start_ttk/*,new_tl*/;
115   poly gen=pCopy(next_generator),p;
116   poly neg_gen=pCopy(next_generator);
117   ideal current_ideal,current_repr;
118   int current_tl;
119   poly w_gen=pHead(next_generator);
120   pSetComp(w_gen,0);
121   pSetmComp(w_gen);
122 
123   //syShowRes(syzstr);
124   neg_gen = pNeg(neg_gen);
125   if (pGetComp(gen)>0)
126   {
127     p_Shift(&gen,-1,currRing);
128     p_Shift(&neg_gen,-1,currRing);
129   }
130   while (index>0)
131   {
132     if (index%2==0)
133       p = gen;
134     else
135       p = neg_gen;
136     if (index>1)
137     {
138       current_ideal = syzstr->res[index-1];
139       current_repr = syzstr->orderedRes[index-1];
140       current_tl = (*syzstr->Tl)[index-1];
141     }
142     else
143     {
144       current_ideal = old_ideal;
145       current_repr = old_repr;
146       current_tl = old_tl;
147     }
148     if (!idIs0(current_ideal))
149     {
150       if (idIs0(syzstr->res[index]))
151       {
152         syzstr->res[index] = idInit(IDELEMS(current_ideal),
153           current_ideal->rank+current_tl);
154         syzstr->orderedRes[index] = idInit(IDELEMS(current_ideal),
155           current_ideal->rank);
156         start = 0;
157       }
158       else
159       {
160         start = IDELEMS(syzstr->res[index]);
161         while ((start>0) && (syzstr->res[index]->m[start-1]==NULL)) start--;
162         if (IDELEMS(syzstr->res[index])<start+IDELEMS(current_ideal))
163         {
164           pEnlargeSet(&syzstr->res[index]->m,IDELEMS(syzstr->res[index]),
165                    IDELEMS(current_ideal));
166           IDELEMS(syzstr->res[index]) += IDELEMS(current_ideal);
167           pEnlargeSet(&syzstr->orderedRes[index]->m,IDELEMS(syzstr->orderedRes[index]),
168                    IDELEMS(current_ideal));
169           IDELEMS(syzstr->orderedRes[index]) += IDELEMS(current_ideal);
170         }
171       }
172       if (idIs0(totake[index]))
173       {
174         totake[index] = idInit(IDELEMS(current_ideal),
175           current_ideal->rank+current_tl);
176         start_ttk = 0;
177       }
178       else
179       {
180         start_ttk = IDELEMS(totake[index]);
181         while ((start_ttk>0) && (totake[index]->m[start_ttk-1]==NULL)) start_ttk--;
182         if (IDELEMS(totake[index])<start_ttk+IDELEMS(current_ideal))
183         {
184           pEnlargeSet(&totake[index]->m,IDELEMS(totake[index]),
185                    IDELEMS(current_ideal));
186           for (j=IDELEMS(totake[index]);j<IDELEMS(totake[index])+
187                                   IDELEMS(current_ideal);j++)
188             totake[index]->m[j] = NULL;
189           IDELEMS(totake[index]) += IDELEMS(current_ideal);
190         }
191       }
192       for (i=0;i<IDELEMS(current_ideal);i++)
193       {
194         if (current_ideal->m[i]!=NULL)
195         {
196           syzstr->res[index]->m[i+start] = pCopy(current_ideal->m[i]);
197           syzstr->res[index]->m[i+start] = pMult_mm(syzstr->res[index]->m[i+start],w_gen);
198           p_Shift(&syzstr->res[index]->m[i+start],current_tl,currRing);
199           syzstr->res[index]->m[i+start] = pAdd(syzstr->res[index]->m[i+start],
200             ppMult_qq(current_repr->m[i],p));
201           syzstr->orderedRes[index]->m[i+start] = pCopy(current_repr->m[i]);
202           syzstr->orderedRes[index]->m[i+start] =
203             pMult_mm(syzstr->orderedRes[index]->m[i+start],w_gen);
204           if ((*syzstr->Tl)[index]!=0)
205             p_Shift(&syzstr->orderedRes[index]->m[i+start],(*syzstr->Tl)[index],currRing);
206         }
207       }
208       for (i=0;i<IDELEMS(totake[index-1]);i++)
209       {
210         if (totake[index-1]->m[i]!=NULL)
211         {
212           if ((index==1) && ((i==IDELEMS(current_ideal) ||
213                (totake[index-1]->m[i+1]==NULL)))) break;
214           totake[index]->m[i+start_ttk] =
215             pMult_mm(pCopy(totake[index-1]->m[i]),w_gen);
216           p_Shift(&totake[index]->m[i+start_ttk],current_tl,currRing);
217 #ifdef FULL_TOTAKE
218           poly pp=pCopy(p);
219           p_Shift(&pp,i+1,currRing);
220           totake[index]->m[i+start_ttk] = pAdd(totake[index]->m[i+start_ttk],pp);
221 #endif
222         }
223       }
224       (*syzstr->Tl)[index] += current_tl;
225     }
226     index--;
227   }
228   pDelete(&gen);
229   pDelete(&neg_gen);
230   pDelete(&w_gen);
231   //syShowRes(syzstr);
232 }
233 
234 /*3
235 * proves the consistence of the pairset resPairs with the corresponding
236 * set of generators;
237 * only for tests
238 */
syTestPairs(SSet resPairs,int length,ideal old_generators)239 static void syTestPairs(SSet resPairs,int length,ideal old_generators)
240 {
241   int i=0;
242 
243   while (i<length)
244   {
245     if (resPairs[i].lcm!=NULL)
246     {
247       if (resPairs[i].p1!=NULL)
248         assume(resPairs[i].p1==old_generators->m[resPairs[i].ind1]);
249       if (resPairs[i].p2!=NULL)
250         assume(resPairs[i].p2==old_generators->m[resPairs[i].ind2]);
251     }
252     i++;
253   }
254 }
255 
256 /*3
257 * cancels the weight monomials given by the leading terms of totake
258 * from the resolution res;
259 * works in place on res, but reads only from totake
260 */
syReorder_Kosz(syStrategy syzstr)261 void syReorder_Kosz(syStrategy syzstr)
262 {
263   int length=syzstr->length;
264   int syzIndex=length-1,i,j;
265   resolvente res=syzstr->fullres;
266   poly p;
267 
268   while ((syzIndex!=0) && (res[syzIndex]==NULL)) syzIndex--;
269   while (syzIndex>0)
270   {
271     for(i=0;i<IDELEMS(res[syzIndex]);i++)
272     {
273 #ifdef USE_REGULARITY
274       if ((syzstr->regularity>0) && (res[syzIndex]->m[i]!=NULL))
275       {
276         if (p_FDeg(res[syzIndex]->m[i],currRing)>=syzstr->regularity+syzIndex)
277           pDelete(&res[syzIndex]->m[i]);
278       }
279 #endif
280       p = res[syzIndex]->m[i];
281       while (p!=NULL)
282       {
283         if (res[syzIndex-1]->m[pGetComp(p)-1]!=NULL)
284         {
285           for(j=1;j<=(currRing->N);j++)
286           {
287             pSetExp(p,j,pGetExp(p,j)
288                         -pGetExp(res[syzIndex-1]->m[pGetComp(p)-1],j));
289           }
290         }
291         else
292           PrintS("error in the resolvent\n");
293         pSetm(p);
294         pIter(p);
295       }
296     }
297     syzIndex--;
298   }
299 }
300 
301 /*3
302 * updates the pairset resPairs by generating all pairs including the
303 * new_generators in the 0-th modul;
304 * the new_generators are inserted in the old_generators;
305 * new_generators is empty after the procedure;
306 */
updatePairs(SSet * resPairs,int * l_pairs,syStrategy syzstr,int index,ideal new_generators,ideal new_repr,int crit_comp)307 static void updatePairs(SSet *resPairs,int *l_pairs,syStrategy syzstr,
308        int index,ideal new_generators,ideal new_repr,int crit_comp)
309 {
310   if (idIs0(new_generators)) return;
311   ideal old_generators=syzstr->res[index];
312   ideal old_repr=syzstr->orderedRes[index];
313   int i=0,j,k,kk,og_elem=0,og_idel=IDELEMS(old_generators),l=*l_pairs,jj,ll,j1;
314   int og_ini=0;
315   ideal pairs=idInit(og_idel+IDELEMS(new_generators),old_generators->rank);
316   polyset prs=pairs->m;
317   poly p=NULL;
318   SObject tso;
319 
320   syInitializePair(&tso);
321   while ((og_elem<og_idel) && (old_generators->m[og_elem]!=NULL))
322   {
323     if ((index>0) && (pGetComp(old_generators->m[og_elem])<=crit_comp))
324       og_ini = og_elem;
325     og_elem++;
326   }
327   while ((l>0) && ((*resPairs)[l-1].lcm==NULL)) l--;
328   while ((i<IDELEMS(new_generators)) && (new_generators->m[i]!=NULL))
329   {
330     syTestPairs(*resPairs,*l_pairs,old_generators);
331     if (IDELEMS(old_generators)==og_elem)
332     {
333       pEnlargeSet(&old_generators->m,IDELEMS(old_generators),16);
334       IDELEMS(old_generators) += 16;
335       pEnlargeSet(&old_repr->m,IDELEMS(old_repr),16);
336       IDELEMS(old_repr) += 16;
337     }
338     k = p_FDeg(new_generators->m[i],currRing);
339     kk = pGetComp(new_generators->m[i]);
340     j = og_ini;
341     while ((j<og_elem) && (old_generators->m[j]!=NULL) &&
342            (pGetComp(old_generators->m[j])<kk)) j++;
343     while ((j<og_elem) && (old_generators->m[j]!=NULL) &&
344            (p_FDeg(old_generators->m[j],currRing)<=k)) j++;
345     for (jj=og_elem;jj>j;jj--)
346     {
347       old_generators->m[jj] = old_generators->m[jj-1];
348       old_repr->m[jj] = old_repr->m[jj-1];
349     }
350     old_generators->m[j] = new_generators->m[i];
351     new_generators->m[i] = NULL;
352     old_repr->m[j] = new_repr->m[i];
353     new_repr->m[i] = NULL;
354     og_elem++;
355     for (jj=0;jj<*l_pairs;jj++)
356     {
357       if ((*resPairs)[jj].lcm!=NULL)
358       {
359         if ((*resPairs)[jj].ind1>=j) (*resPairs)[jj].ind1++;
360         if ((*resPairs)[jj].ind2>=j) (*resPairs)[jj].ind2++;
361       }
362     }
363     syTestPairs(*resPairs,*l_pairs,old_generators);
364     for (jj=og_ini;jj<og_elem;jj++)
365     {
366       if ((j!=jj) && (pGetComp(old_generators->m[jj])==pGetComp(old_generators->m[j])))
367       {
368         p = pOne();
369         pLcm(old_generators->m[jj],old_generators->m[j],p);
370         pSetComp(p,j+1);
371         pSetm(p);
372         j1 = 0;
373         while (j1<jj)
374         {
375           if (prs[j1]!=NULL)
376           {
377             if (pLmDivisibleByNoComp(prs[j1],p))
378             {
379               pDelete(&p);
380               break;
381             }
382             else if (pLmDivisibleByNoComp(p,prs[j1]))
383             {
384               pDelete(&(prs[j1]));
385             }
386 #ifdef USE_CHAINCRIT0
387             else
388             {
389               poly p1,p2;
390               int ip=(currRing->N);
391               p1 = pMDivide(p,old_generators->m[jj]);
392               p2 = pMDivide(prs[j1],old_generators->m[j1]);
393               while ((ip>0) && (pGetExp(p1,ip)*pGetExp(p2,ip)==0)) ip--;
394               if (ip==0)
395               {
396                 int ti=0;
397                 while ((ti<l) && (((*resPairs)[ti].ind1!=j1)|| ((*resPairs)[ti].ind2!=jj))) ti++;
398                 if (ti<l)
399                 {
400                   if (TEST_OPT_PROT) PrintS("cc");
401                   syDeletePair(&(*resPairs)[ti]);
402                   syCompactifyPairSet(*resPairs,*l_pairs,ti);
403                   l--;
404                 }
405               }
406               pDelete(&p1);
407               pDelete(&p2);
408             }
409 #endif
410           }
411           j1++;
412         }
413         if (p!=NULL)
414           prs[jj] = p;
415       }
416     }
417     for (jj=og_ini;jj<og_elem;jj++)
418     {
419       if (prs[jj] !=NULL)
420       {
421         if (l>=*l_pairs)
422         {
423           SSet temp = (SSet)omAlloc0((*l_pairs+16)*sizeof(SObject));
424           for (ll=0;ll<*l_pairs;ll++)
425           {
426             temp[ll].p = (*resPairs)[ll].p;
427             temp[ll].p1 = (*resPairs)[ll].p1;
428             temp[ll].p2 = (*resPairs)[ll].p2;
429             temp[ll].syz = (*resPairs)[ll].syz;
430             temp[ll].lcm = (*resPairs)[ll].lcm;
431             temp[ll].ind1 = (*resPairs)[ll].ind1;
432             temp[ll].ind2 = (*resPairs)[ll].ind2;
433             temp[ll].syzind = (*resPairs)[ll].syzind;
434             temp[ll].order = (*resPairs)[ll].order;
435             temp[ll].isNotMinimal = (*resPairs)[ll].isNotMinimal;
436           }
437           omFreeSize((ADDRESS)(*resPairs),*l_pairs*sizeof(SObject));
438           *l_pairs += 16;
439           (*resPairs) = temp;
440         }
441         tso.lcm = prs[jj];
442         prs[jj] = NULL;
443         tso.order = p_FDeg(tso.lcm,currRing);
444         tso.p1 = old_generators->m[jj];
445         tso.p2 = old_generators->m[j];
446         tso.ind1 = jj;
447         tso.ind2 = j;
448         tso.syzind = -1;
449         tso.isNotMinimal = NULL;
450         tso.p = NULL;
451         tso.syz = NULL;
452         SSet rP=*resPairs;
453 #ifdef SHOW_PROT
454 Print("erzeuge Paar im Modul %d,%d mit: \n",index,tso.order);
455 PrintS("poly1: ");pWrite(tso.p1);
456 PrintS("poly2: ");pWrite(tso.p2);
457 PrintS("syz: ");pWrite(tso.syz);
458 PrintS("sPoly: ");pWrite(tso.p);
459 PrintLn();
460 #endif
461         syEnterPair(rP,&tso,&l,index);
462         syInitializePair(&tso);
463       }
464     }
465     i++;
466   }
467   idDelete(&pairs);
468 }
469 
470 /*3
471 * performs the modification of a single reduction on the syzygy-level
472 */
sySPRedSyz_Kosz(syStrategy syzstr,poly redWith,poly syz,poly q=NULL,int l_syz=-1)473 inline void sySPRedSyz_Kosz(syStrategy syzstr,poly redWith,poly syz,poly q=NULL,int l_syz=-1)
474 {
475   poly p=pMDivide(q,redWith);
476   pSetCoeff(p,nDiv(pGetCoeff(q),pGetCoeff(redWith)));
477   kBucket_Minus_m_Mult_p(syzstr->syz_bucket,p,syz,&l_syz,NULL);
478   pDelete(&p);
479 }
480 
481 /*3
482 * normalizes the poly bucket by the ideal;
483 * stops the reduction whenever the leading component is less than the
484 * crit_comp;
485 * returns the changing status
486 */
syRedSyz(kBucket_pt bucket,ideal red,int crit_comp,int * g_l)487 static BOOLEAN syRedSyz(kBucket_pt bucket,ideal red,int crit_comp,int* g_l)
488 {
489   poly p = kBucketGetLm(bucket);
490   int j = 0,i=IDELEMS(red)-1;
491   number n;
492   BOOLEAN isChanged=FALSE;
493 
494   loop
495   {
496     if ((j>=i) || (p==NULL) || (pGetComp(p)<=crit_comp)) break;
497     if ((red->m[j]!=NULL) && (pDivisibleBy(red->m[j],p)))
498     {
499       n = kBucketPolyRed(bucket,red->m[j], g_l[j], NULL);
500       nDelete(&n);
501       p = kBucketGetLm(bucket);
502       isChanged = TRUE;
503       j = 0;
504     }
505     else
506       j++;
507   }
508   return isChanged;
509 }
510 
511 /*3
512 * a tail reduction for the syzygies yielding new generators
513 */
syRedTailSyz(poly tored,ideal red,ideal sec_red,int crit_comp,syStrategy syzstr,int * gen_length,int * secgen_length,int * tored_length)514 static poly syRedTailSyz(poly tored,ideal red,ideal sec_red,int crit_comp,syStrategy syzstr,
515             int * gen_length,int * secgen_length,int * tored_length)
516 {
517   int i=IDELEMS(red)-1,num_mon,num_tail;
518   poly h,hn;
519   // BOOLEAN dummy;
520 
521   while ((i>0) && (red->m[i-1]==NULL)) i--;
522   i--;
523   h = tored;
524   if ((h!=NULL) && (pGetComp(h)>crit_comp))
525   {
526     num_mon = 1;
527     hn = pNext(h);
528     num_tail = *tored_length-1;
529     while (hn!=NULL)
530     {
531       kBucketInit(syzstr->syz_bucket,hn,num_tail);
532       /*dummy =*/ (void) syRedSyz(syzstr->syz_bucket,red,crit_comp,gen_length);
533       kBucketClear(syzstr->syz_bucket,&hn,&num_tail);
534       pNext(h) = hn;
535       if ((hn==NULL) || (pGetComp(hn)<=crit_comp))
536         break;
537       else
538       {
539         pIter(h);
540         pIter(hn);
541         num_mon++;
542         num_tail--;
543       }
544     }
545     if (sec_red!=NULL)
546     {
547       while (hn!=NULL)
548       {
549         kBucketInit(syzstr->syz_bucket,hn,num_tail);
550         /*dummy =*/ (void) syRedSyz(syzstr->syz_bucket,sec_red,crit_comp,secgen_length);
551         kBucketClear(syzstr->syz_bucket,&hn,&num_tail);
552         pNext(h) = hn;
553         if (hn==NULL)
554           break;
555         else
556         {
557           pIter(h);
558           pIter(hn);
559           num_mon++;
560           num_tail--;
561         }
562       }
563     }
564     *tored_length = num_mon+num_tail;
565   }
566   assume(pLength(tored)==*tored_length);
567   return tored;
568 }
569 
570 #if 0
571 // unused
572 /*3
573 * the complete reduction of a single pair which is just stored
574 * in bucket and syz_bucket
575 */
576 static BOOLEAN syRedSyzPair(syStrategy syzstr,int index,int* g_l,int* orp_l)
577 {
578   kBucket_pt bucket=syzstr->bucket;
579   poly p = kBucketGetLm(bucket);
580   ideal red=syzstr->res[index],repr=syzstr->orderedRes[index];
581   int j = 0,i=IDELEMS(red)-1;
582   number n;
583   BOOLEAN isChanged=FALSE;
584 
585   loop
586   {
587     if ((j>=i) || (p==NULL)) break;
588     if ((red->m[j]!=NULL) && (pDivisibleBy(red->m[j],p)))
589     {
590       sySPRedSyz_Kosz(syzstr,red->m[j],repr->m[j],p,orp_l[j]);
591       n = kBucketPolyRed(bucket,red->m[j], g_l[j], NULL);
592       nDelete(&n);
593       p = kBucketGetLm(bucket);
594       isChanged = TRUE;
595       j = 0;
596     }
597     else
598       j++;
599   }
600   return isChanged;
601 }
602 #endif
603 
604 /*3
605 * the tailreduction for generators (which includes the correction of
606 * the corresponding representation)
607 */
608 #if 0 /*unused*/
609 static void syRedTailSyzPair(SObject tso,syStrategy syzstr,int index,
610             int * gen_length,int* orp_l,int * tored_l,int * syzred_l)
611 {
612   int num_mon,num_tail,syz_l;
613   poly h,hn;
614   BOOLEAN dummy;
615 
616   h = tso.p;
617   kBucketInit(syzstr->syz_bucket,tso.syz,*syzred_l);
618   if (h!=NULL)
619   {
620     num_mon = 1;
621     hn = pNext(h);
622     num_tail = *tored_l-1;
623     while (hn!=NULL)
624     {
625       kBucketInit(syzstr->bucket,hn,num_tail);
626       dummy = syRedSyzPair(syzstr,index,gen_length,orp_l);
627       kBucketClear(syzstr->bucket,&hn,&num_tail);
628       pNext(h) = hn;
629       if (hn==NULL)
630         break;
631       else
632       {
633         pIter(h);
634         pIter(hn);
635         num_mon++;
636         num_tail--;
637       }
638     }
639     *tored_l = num_mon+num_tail;
640   }
641   kBucketClear(syzstr->syz_bucket,&tso.syz,&syz_l);
642   assume(pLength(tso.syz)==syz_l);
643   assume(pLength(tso.p)==*tored_l);
644 }
645 #endif
646 
647 /*3
648 * the reduction of a pair in the 0-th module
649 */
redOnePair(SSet resPairs,int itso,int l,ideal syzygies,int crit_comp,syStrategy syzstr,int index,ideal new_generators,ideal new_repr,int * ogm_l,int * orp_l)650 static void redOnePair(SSet resPairs,int itso,int l, ideal syzygies,
651             int crit_comp, syStrategy syzstr,int index,ideal new_generators,
652             ideal new_repr,int * ogm_l,int * orp_l)
653 {
654   SObject tso = resPairs[itso];
655   assume (tso.lcm!=NULL);
656   ideal old_generators=syzstr->res[index];
657   ideal old_repr=syzstr->orderedRes[index];
658   int og_idel=IDELEMS(old_generators),ng_place=IDELEMS(new_generators);
659   int toReplace=0;
660   int i,j,syz_l;
661   number /*coefgcd,*/n;
662   polyset ogm=old_generators->m;
663   poly p;
664   BOOLEAN deleteP=FALSE;
665 #ifdef EXPERIMENT1
666   poly syzp;
667 #endif
668   int syz_place=IDELEMS(syzygies);
669 
670   while ((syz_place>0) && (syzygies->m[syz_place-1]==NULL)) syz_place--;
671   while ((ng_place>0) && (new_generators->m[ng_place-1]==NULL)) ng_place--;
672   while ((og_idel>0) && (old_generators->m[og_idel-1]==NULL)) og_idel--;
673   assume (tso.ind1<og_idel);
674   assume (tso.ind2<og_idel);
675   assume (tso.ind1!=tso.ind2);
676   assume (tso.p1 == old_generators->m[tso.ind1]);
677   assume (tso.p2 == old_generators->m[tso.ind2]);
678   tso.p1 = old_generators->m[tso.ind1];
679   tso.p2 = old_generators->m[tso.ind2];
680   if ((tso.p1!=NULL) && (tso.p2!=NULL))
681   {
682     if (TEST_OPT_PROT)
683       PrintS(".");
684     if (index==0)
685     {
686 /*--- tests wether a generator must be replaced (lt(f1)|lt(f2)!)--*/
687       if (p_FDeg(tso.p1,currRing)==p_FDeg(tso.lcm,currRing))
688         toReplace = tso.ind1+1;
689       else if (p_FDeg(tso.p2,currRing)==p_FDeg(tso.lcm,currRing))
690         toReplace = tso.ind2+1;
691     }
692 #ifdef EXPERIMENT3
693 /*--- tests wether the product criterion applies --------------*/
694     if ((index==0) && (old_generators->rank==1) &&
695         (p_FDeg(tso.p1,currRing)+p_FDeg(tso.p2,currRing)==tso.order))
696     {
697       tso.p = NULL;
698       p = pCopy(tso.p1);
699       p_Shift(&p,-1,currRing);
700 #ifdef WITH_BUCKET
701       poly pp;
702       pp = pMult_mm(pCopy(old_repr->m[tso.ind2]),p);
703       kBucketInit(syzstr->syz_bucket,pp,-1);
704       pLmDelete(&p);
705       p = pNeg(p);
706       pp = pCopy(old_repr->m[tso.ind2]);
707       int il=-1;
708       while (p!=NULL)
709       {
710         kBucket_Minus_m_Mult_p(syzstr->syz_bucket,p,pp,&il,NULL);
711         pLmDelete(&p);
712       }
713       pDelete(&pp);
714       p = pCopy(tso.p2);
715       p_Shift(&p,-1,currRing);
716       pp = pCopy(old_repr->m[tso.ind1]);
717       il=-1;
718       while (p!=NULL)
719       {
720         kBucket_Minus_m_Mult_p(syzstr->syz_bucket,p,pp,&il,NULL);
721         pLmDelete(&p);
722       }
723       pDelete(&pp);
724       kBucketClear(syzstr->syz_bucket,&tso.syz,&j);
725 #else
726       tso.syz = pMult(p,pCopy(old_repr->m[tso.ind2]));
727       p = pCopy(tso.p2);
728       p_Shift(&p,-1,currRing);
729       tso.syz = pSub(tso.syz,pMult(p,pCopy(old_repr->m[tso.ind1])));
730 #endif
731     }
732     else
733 #endif
734 /*--- the product criterion does not apply --------------------*/
735     {
736       tso.p = ksOldCreateSpoly(tso.p2,tso.p1);
737       number coefgcd = n_Gcd(pGetCoeff(tso.p1),pGetCoeff(tso.p2),currRing->cf);
738       assume (old_repr->m[tso.ind1]!=NULL);
739       tso.syz = pCopy(old_repr->m[tso.ind1]);
740       poly tt = pMDivide(tso.lcm,tso.p1);
741       pSetComp(tt,0);
742       pSetmComp(tt);
743       pSetCoeff(tt,nDiv(pGetCoeff(tso.p1),coefgcd));
744       tso.syz = pMult_mm(tso.syz,tt);
745       pDelete(&tt);
746       coefgcd = nInpNeg(coefgcd);
747       assume (old_repr->m[tso.ind2]!=NULL);
748       p = pCopy(old_repr->m[tso.ind2]);
749       tt = pMDivide(tso.lcm,tso.p2);
750       pSetComp(tt,0);
751       pSetmComp(tt);
752       pSetCoeff(tt,nDiv(pGetCoeff(tso.p2),coefgcd));
753       p = pMult_mm(p,tt);
754       pDelete(&tt);
755       tso.syz = pAdd(p,tso.syz);
756 #ifdef EXPERIMENT2
757       if ((tso.syz!=NULL) && (pGetComp(tso.syz)<=crit_comp))
758       {
759 /*--- breaks when the leading component is less than crit_comp ------*/
760         deleteP = TRUE;
761         discard_pairs++;
762       }
763 #endif
764       nDelete(&coefgcd);
765     }                             //End of the else-part of EXPERIMENT3
766 #ifdef SHOW_PROT
767 Print("reduziere Paar im Module %d mit: \n",index);
768 PrintS("poly1: ");pWrite(tso.p1);
769 PrintS("poly2: ");pWrite(tso.p2);
770 PrintS("syz: ");pWrite(tso.syz);
771 PrintS("sPoly: ");pWrite(tso.p);
772 #endif
773     assume(tso.syz!=NULL);
774     kBucketInit(syzstr->syz_bucket,tso.syz,-1);
775     if ((tso.p!=NULL) && (!deleteP))
776     {
777       kBucketInit(syzstr->bucket,tso.p,-1);
778       p = kBucketGetLm(syzstr->bucket);
779       j = 0;
780       loop
781       {
782         if (j>=og_idel)
783         {
784 /*--- reduction with generators computed in this procedure ---*/
785           j = 0;
786           while ((j<ng_place) && (!pDivisibleBy(new_generators->m[j],p))) j++;
787           if (j>=ng_place) break;
788           assume (new_repr->m[j]!=NULL);
789           sySPRedSyz_Kosz(syzstr,new_generators->m[j],new_repr->m[j],p);
790           n = kBucketPolyRed(syzstr->bucket,new_generators->m[j],pLength(new_generators->m[j]), NULL);
791           p = kBucketGetLm(syzstr->bucket);
792 #ifdef EXPERIMENT1
793           syzp = kBucketGetLm(syzstr->syz_bucket);
794           if ((syzp!=NULL) && (pGetComp(syzp)<=crit_comp))
795           {
796             deleteP =TRUE;
797             break;
798           }
799           //if (syzp==NULL)
800             //assume(p==NULL);
801           //else
802             //if (pGetComp(syzp)<=crit_comp) short_pairs++;
803 #endif
804           if (p==NULL) break;
805           j = 0;
806         }
807         if (pDivisibleBy(ogm[j],p))
808         {
809 /*--- reduction with general old generators ---------------------*/
810           assume (old_repr->m[j]!=NULL);
811           sySPRedSyz_Kosz(syzstr,ogm[j],old_repr->m[j],p,orp_l[j]);
812           n = kBucketPolyRed(syzstr->bucket,ogm[j],ogm_l[j], NULL);
813           p = kBucketGetLm(syzstr->bucket);
814 #ifdef EXPERIMENT1
815           syzp = kBucketGetLm(syzstr->syz_bucket);
816           if ((syzp!=NULL) && (pGetComp(syzp)<=crit_comp))
817           {
818             break;
819             deleteP =TRUE;
820           }
821           //if (syzp==NULL)
822             //assume(p==NULL);
823           //else
824             //if ((pGetComp(syzp)<=crit_comp) && (p!=NULL)) short_pairs++;
825 #endif
826           if (p==NULL) break;
827           j = 0;
828         }
829         else
830           j++;
831       }
832       kBucketClear(syzstr->bucket,&tso.p,&tso.length);
833     }
834     kBucketClear(syzstr->syz_bucket,&tso.syz,&syz_l);
835     if (deleteP)
836     {
837       pDelete(&tso.p);
838       pDelete(&tso.syz);
839     }
840   }
841   else
842   {
843     PrintS("Shit happens!\n");
844   }
845 #ifdef SHOW_PROT
846 Print("erhalte Paar im Module %d mit: \n",index);
847 PrintS("syz: ");pWrite(tso.syz);
848 PrintS("sPoly: ");pWrite(tso.p);
849 PrintLn();
850 #endif
851   if (toReplace)
852   {
853 /*-- replaces the generator if neccesary ------------------*/
854     pDelete(&old_generators->m[toReplace-1]);
855     pDelete(&old_repr->m[toReplace-1]);
856     for (i=toReplace-1;i<og_idel-1;i++)
857     {
858       old_generators->m[i] = old_generators->m[i+1];
859       old_repr->m[i] = old_repr->m[i+1];
860     }
861     old_generators->m[og_idel-1] = NULL;
862     old_repr->m[og_idel-1] = NULL;
863     for (i=itso+1;i<l;i++)
864     {
865       if (resPairs[i].lcm!=NULL)
866       {
867         if ((resPairs[i].ind1==toReplace-1)||(resPairs[i].ind2==toReplace-1))
868           syDeletePair(&resPairs[i]);
869         else
870         {
871           if (resPairs[i].ind1>=toReplace)
872             (resPairs[i].ind1)--;
873           if (resPairs[i].ind2>=toReplace)
874             (resPairs[i].ind2)--;
875         }
876       }
877     }
878     syCompactifyPairSet(resPairs,l,itso+1);
879   }
880   if (tso.p!=NULL)
881   {
882 /*-- stores the new generator ---------------------------------*/
883     //syRedTailSyzPair(tso,syzstr,index,ogm_l,orp_l,&tso.length,&syz_l);
884     if (ng_place>=IDELEMS(new_generators))
885     {
886       pEnlargeSet(&new_generators->m,IDELEMS(new_generators),16);
887       IDELEMS(new_generators) += 16;
888       pEnlargeSet(&new_repr->m,IDELEMS(new_repr),16);
889       IDELEMS(new_repr) += 16;
890     }
891     if (!nIsOne(pGetCoeff(tso.p)))
892     {
893       n=nInvers(pGetCoeff(tso.p));
894       pNorm(tso.p);
895       tso.syz=__p_Mult_nn(tso.syz,n,currRing);
896       nDelete(&n);
897     }
898     new_generators->m[ng_place] = tso.p;
899     tso.p = NULL;
900     new_repr->m[ng_place] = tso.syz;
901     tso.syz = NULL;
902   }
903   else
904   {
905 /*--- takes the syzygy as new generator of the next module ---*/
906     if (tso.syz==NULL)
907     {
908 #ifndef EXPERIMENT2
909 #ifdef EXPERIMENT3
910       short_pairs++;
911 #endif
912 #endif
913     }
914     else if (pGetComp(tso.syz)<=crit_comp)
915     {
916       pDelete(&tso.syz);
917     }
918     else
919     {
920       if (syz_place>=IDELEMS(syzygies))
921       {
922         pEnlargeSet(&syzygies->m,IDELEMS(syzygies),16);
923         IDELEMS(syzygies) += 16;
924       }
925       syzygies->m[syz_place] = tso.syz;
926       tso.syz = NULL;
927       pNorm(syzygies->m[syz_place]);
928     }
929   }
930   resPairs[itso] = tso;
931   syDeletePair(&resPairs[itso]);
932   syTestPairs(resPairs,l,old_generators);
933 }
934 
935 /*3
936 * reduction of all pairs of a fixed degree of the 0-th module
937 */
redPairs(SSet resPairs,int l_pairs,ideal syzygies,ideal new_generators,ideal new_repr,int crit_comp,syStrategy syzstr,int index)938 static BOOLEAN redPairs(SSet resPairs,int l_pairs, ideal syzygies,
939   ideal new_generators,ideal new_repr, int crit_comp,syStrategy syzstr,
940   int index)
941 {
942   if (resPairs[0].lcm==NULL) return TRUE;
943   int i,j,actdeg=resPairs[0].order;
944   int * ogm_l=(int*)omAlloc0(IDELEMS(syzstr->res[index])*sizeof(int));
945   int * orp_l=(int*)omAlloc0(IDELEMS(syzstr->orderedRes[index])*sizeof(int));
946   // int t1=IDELEMS(syzstr->res[index]),t2=IDELEMS(syzstr->orderedRes[index]);
947 
948   for (j=IDELEMS(syzstr->res[index])-1;j>=0;j--)
949   {
950     if (syzstr->res[index]->m[j]!=NULL)
951       ogm_l[j] = pLength(syzstr->res[index]->m[j]);
952   }
953   for (j=IDELEMS(syzstr->orderedRes[index])-1;j>=0;j--)
954   {
955     if (syzstr->orderedRes[index]->m[j]!=NULL)
956       orp_l[j] = pLength(syzstr->orderedRes[index]->m[j]);
957   }
958   loop
959   {
960     i = 0;
961     if (TEST_OPT_PROT)
962       Print("(%d,%d)",index,resPairs[0].order);
963     while (resPairs[i].order==actdeg)
964     {
965       syTestPairs(resPairs,l_pairs,syzstr->res[index]);
966       redOnePair(resPairs,i,l_pairs,syzygies,crit_comp,syzstr,index,
967                  new_generators, new_repr,ogm_l,orp_l);
968       i++;
969       syTestPairs(resPairs,l_pairs,syzstr->res[index]);
970     }
971     syTestPairs(resPairs,l_pairs,syzstr->res[index]);
972     syCompactifyPairSet(resPairs,l_pairs,0);
973     syTestPairs(resPairs,l_pairs,syzstr->res[index]);
974     if (!idIs0(new_generators))
975       break;
976     else if (resPairs[0].lcm==NULL)  //there are no pairs left and no new_gens
977     {
978       omFreeSize((ADDRESS)ogm_l,IDELEMS(syzstr->res[index])*sizeof(int));
979       omFreeSize((ADDRESS)orp_l,IDELEMS(syzstr->orderedRes[index])*sizeof(int));
980       return TRUE;
981     }
982     else
983       actdeg = resPairs[0].order;
984   }
985   syTestPairs(resPairs,l_pairs,syzstr->res[index]);
986   omFreeSize((ADDRESS)ogm_l,IDELEMS(syzstr->res[index])*sizeof(int));
987   omFreeSize((ADDRESS)orp_l,IDELEMS(syzstr->orderedRes[index])*sizeof(int));
988   return FALSE;
989 }
990 
991 /*3
992 * extends the standard basis old_generators with new_generators;
993 * returns the syzygies which involve the new elements;
994 * assumes that the components of the new_generators are sperated
995 * from those of old_generators, i.e. whenever the leading term
996 * of a syzygy lies in the part of the old_generators, the syzygy
997 * lie just in the module old_generators
998 * assumes that the new_generators are reduced w.r.t. old_generators
999 */
kosz_std(ideal new_generators,ideal new_repr,syStrategy syzstr,int index,int next_comp)1000 static ideal kosz_std(ideal new_generators,ideal new_repr,syStrategy syzstr,
1001                       int index,int next_comp)
1002 {
1003   int og_idel=IDELEMS(syzstr->res[index]);
1004   int l_pairs=2*og_idel;
1005   ideal syzygies=idInit(16,syzstr->res[index]->rank+1);
1006   if ((idIs0(new_generators)) || (new_generators->m[0]==NULL))
1007   {
1008     WerrorS("Hier ist was faul!\n");
1009     return NULL;
1010   }
1011   SSet resPairs=(SSet)omAlloc0(l_pairs*sizeof(SObject));
1012   loop
1013   {
1014     updatePairs(&resPairs,&l_pairs,syzstr,index,
1015                 new_generators,new_repr,next_comp);
1016     if (redPairs(resPairs,l_pairs,syzygies, new_generators,new_repr,
1017                  next_comp,syzstr,index)) break;
1018   }
1019   omFreeSize((SSet)resPairs,l_pairs*sizeof(SObject));
1020   return syzygies;
1021 }
1022 
1023 /*3
1024 * normalizes the incoming generators
1025 */
normalize(poly next_p,ideal add_generators,syStrategy syzstr,int * g_l,int * p_l,int crit_comp)1026 static poly normalize(poly next_p,ideal add_generators, syStrategy syzstr,
1027                       int * g_l,int * p_l,int crit_comp)
1028 {
1029   int j=0,i=IDELEMS(add_generators);
1030   kBucketInit(syzstr->bucket,next_p,pLength(next_p));
1031   poly p = kBucketGetLm(syzstr->bucket),result;
1032   number n;
1033 
1034   loop
1035   {
1036     if ((j>=i) || (p==NULL) || (pGetComp(p)<=crit_comp)) break;
1037     if ((add_generators->m[j]!=NULL) && (pDivisibleBy(add_generators->m[j],p)))
1038     {
1039       n = kBucketPolyRed(syzstr->bucket,add_generators->m[j], g_l[j], NULL);
1040       nDelete(&n);
1041       p = kBucketGetLm(syzstr->bucket);
1042       j = 0;
1043     }
1044     else
1045       j++;
1046   }
1047   kBucketClear(syzstr->bucket,&result,p_l);
1048   return result;
1049 }
1050 
1051 /*3
1052 * updates the pairs inthe higher modules
1053 */
updatePairsHIndex(SSet * resPairs,int * l_pairs,syStrategy,int index,ideal add_generators,ideal,ideal,ideal,int,int * first_new)1054 static void updatePairsHIndex(SSet *resPairs,int *l_pairs,syStrategy /*syzstr*/,
1055        int index,ideal add_generators,ideal /*add_repr*/,ideal /*new_generators*/,
1056        ideal /*new_repr*/,int /*crit_comp*/,int* first_new)
1057 {
1058   int i=*first_new,l=*l_pairs,j,ll,j1,add_idel=IDELEMS(add_generators);
1059   ideal pairs=idInit(add_idel,add_generators->rank);
1060   polyset prs=pairs->m;
1061   poly p=NULL;
1062   SObject tso;
1063 
1064   syInitializePair(&tso);
1065   while ((l>0) && ((*resPairs)[l-1].lcm==NULL)) l--;
1066   while ((i<add_idel) && (add_generators->m[i]!=NULL))
1067   {
1068     for (j=0;j<i;j++)
1069     {
1070       if (pGetComp(add_generators->m[j]) == pGetComp(add_generators->m[i]))
1071       {
1072         p = pOne();
1073         pLcm(add_generators->m[j],add_generators->m[i],p);
1074         pSetComp(p,i+1);
1075         pSetm(p);
1076         j1 = 0;
1077         while (j1<j)
1078         {
1079           if (prs[j1]!=NULL)
1080           {
1081             if (pLmDivisibleByNoComp(prs[j1],p))
1082             {
1083               pDelete(&p);
1084               break;
1085             }
1086             else if (pLmDivisibleByNoComp(p,prs[j1]))
1087             {
1088               pDelete(&(prs[j1]));
1089             }
1090 #ifdef USE_CHAINCRIT
1091             else
1092             {
1093               poly p1,p2;
1094               int ip=(currRing->N);
1095               p1 = pMDivide(p,add_generators->m[j]);
1096               p2 = pMDivide(prs[j1],add_generators->m[j1]);
1097               while ((ip>0) && (pGetExp(p1,ip)*pGetExp(p2,ip)==0)) ip--;
1098               if (ip==0)
1099               {
1100                 int ti=0;
1101                 while ((ti<l) && (((*resPairs)[ti].ind1!=j1)|| ((*resPairs)[ti].ind2!=j))) ti++;
1102                 if (ti<l)
1103                 {
1104                   if (TEST_OPT_PROT) PrintS("cc");
1105                   syDeletePair(&(*resPairs)[ti]);
1106                   syCompactifyPairSet(*resPairs,*l_pairs,ti);
1107                   l--;
1108                 }
1109               }
1110               pDelete(&p1);
1111               pDelete(&p2);
1112             }
1113 #endif
1114           }
1115           j1++;
1116         }
1117         if (p!=NULL)
1118           prs[j] = p;
1119       }
1120     }
1121     for (j=0;j<i;j++)
1122     {
1123       if (prs[j] !=NULL)
1124       {
1125         if (l>=*l_pairs)
1126         {
1127           SSet temp = (SSet)omAlloc0((*l_pairs+16)*sizeof(SObject));
1128           for (ll=0;ll<*l_pairs;ll++)
1129           {
1130             temp[ll].p = (*resPairs)[ll].p;
1131             temp[ll].p1 = (*resPairs)[ll].p1;
1132             temp[ll].p2 = (*resPairs)[ll].p2;
1133             temp[ll].syz = (*resPairs)[ll].syz;
1134             temp[ll].lcm = (*resPairs)[ll].lcm;
1135             temp[ll].ind1 = (*resPairs)[ll].ind1;
1136             temp[ll].ind2 = (*resPairs)[ll].ind2;
1137             temp[ll].syzind = (*resPairs)[ll].syzind;
1138             temp[ll].order = (*resPairs)[ll].order;
1139             temp[ll].isNotMinimal = (*resPairs)[ll].isNotMinimal;
1140           }
1141           omFreeSize((ADDRESS)(*resPairs),*l_pairs*sizeof(SObject));
1142           *l_pairs += 16;
1143           (*resPairs) = temp;
1144         }
1145         tso.lcm = prs[j];
1146         prs[j] = NULL;
1147         tso.order = p_FDeg(tso.lcm,currRing);
1148         tso.p1 = add_generators->m[j];
1149         tso.p2 = add_generators->m[i];
1150         tso.ind1 = j;
1151         tso.ind2 = i;
1152         tso.syzind = -1;
1153         tso.isNotMinimal = NULL;
1154         tso.p = NULL;
1155         tso.syz = NULL;
1156         SSet rP=*resPairs;
1157 #ifdef SHOW_PROT
1158 Print("erzeuge Paar im Modul %d,%d mit: \n",index,tso.order);
1159 PrintS("poly1: ");pWrite(tso.p1);
1160 PrintS("poly2: ");pWrite(tso.p2);
1161 PrintS("syz: ");pWrite(tso.syz);
1162 PrintS("sPoly: ");pWrite(tso.p);
1163 PrintLn();
1164 #endif
1165         syEnterPair(rP,&tso,&l,index);
1166         syInitializePair(&tso);
1167       }
1168     }
1169     i++;
1170   }
1171   *first_new = i;
1172   idDelete(&pairs);
1173 }
1174 
1175 /*3
1176 * reduction of a single pair in the higher moduls
1177 */
1178 #ifdef SHOW_PROT
redOnePairHIndex(SSet resPairs,int itso,int crit_comp,syStrategy syzstr,int index,ideal add_generators,ideal add_repr,ideal new_generators,ideal new_repr,int * next_place_add,int ** g_l,poly deg_soc)1179 static void redOnePairHIndex(SSet resPairs,int itso, int crit_comp,
1180             syStrategy syzstr,int index,ideal add_generators, ideal add_repr,
1181             ideal new_generators, ideal new_repr,int * next_place_add,int ** g_l,
1182             poly deg_soc)
1183 #else
1184 static void redOnePairHIndex(SSet resPairs,int itso, int crit_comp,
1185             syStrategy syzstr,int /*index*/,ideal add_generators, ideal add_repr,
1186             ideal new_generators, ideal new_repr,int * next_place_add,int ** g_l,
1187             poly deg_soc)
1188 #endif
1189 {
1190   SObject tso = resPairs[itso];
1191   assume (tso.lcm!=NULL);
1192   int ng_place=IDELEMS(new_generators);
1193   int i,j;
1194   number n;
1195   poly p;
1196 #ifdef EXPERIMENT1
1197   poly syzp;
1198 #endif
1199 
1200   assume (tso.ind1<*next_place_add);
1201   assume (tso.ind2<*next_place_add);
1202   assume (tso.ind1!=tso.ind2);
1203   assume (tso.p1 == add_generators->m[tso.ind1]);
1204   assume (tso.p2 == add_generators->m[tso.ind2]);
1205   tso.p1 = add_generators->m[tso.ind1];
1206   tso.p2 = add_generators->m[tso.ind2];
1207   if ((tso.p1!=NULL) && (tso.p2!=NULL))
1208   {
1209     if (TEST_OPT_PROT)
1210       PrintS(".");
1211 #ifdef USE_PROD_CRIT
1212     if (p_FDeg(tso.p1,currRing)+p_FDeg(tso.p2,currRing)==tso.order+p_FDeg(deg_soc,currRing))
1213     {
1214       if (TEST_OPT_PROT) PrintS("pc");
1215       int ac=pGetComp(tso.p1);
1216       assume(ac=pGetComp(tso.p2));
1217       poly p1=pCopy(tso.p1);
1218       poly p2=pCopy(tso.p2);
1219       poly pp1,pp2,tp1,tp2;
1220       poly sp1=pCopy(add_repr->m[tso.ind1]),sp2=pCopy(add_repr->m[tso.ind2]);
1221       pp1 = p1;
1222       pp2 = p2;
1223       loop
1224       {
1225         assume(pp1!=NULL);
1226         for(i=(int)(currRing->N); i; i--)
1227           pSetExp(pp1,i, pGetExp(pp1,i)- pGetExp(deg_soc,i));
1228         pSetComp(pp1, 0);
1229         pSetm(pp1);
1230         if ((pNext(pp1)!=NULL) && (pGetComp(pNext(pp1))!=ac))  break;
1231         pIter(pp1);
1232       }
1233       loop
1234       {
1235         assume(pp2!=NULL);
1236         for(i=(int)(currRing->N); i; i--)
1237           pSetExp(pp2,i, pGetExp(pp2,i)- pGetExp(deg_soc,i));
1238         pSetComp(pp2, 0);
1239         pSetm(pp2);
1240         if ((pNext(pp2)!=NULL) && (pGetComp(pNext(pp2))!=ac)) break;
1241         pIter(pp2);
1242       }
1243       tp1 = pNext(pp1);
1244       tp2 = pNext(pp2);
1245       pNext(pp1) = NULL;
1246       pNext(pp2) = NULL;
1247       //p_Shift(&p1,-ac,currRing);
1248       //p_Shift(&p2,-ac,currRing);
1249       tp1 = pMult(tp1,pCopy(p2));
1250       tp2 = pMult(tp2,pCopy(p1));
1251       sp1 = pMult(p2,sp1);
1252       sp2 = pMult(p1,sp2);
1253       tso.p = pSub(tp1,tp2);
1254       tso.syz = pSub(sp1,sp2);
1255     }
1256     else
1257 #endif
1258     {
1259       tso.p = ksOldCreateSpoly(tso.p2,tso.p1);
1260       number coefgcd = n_Gcd(pGetCoeff(tso.p1),pGetCoeff(tso.p2),currRing->cf);
1261       assume (add_repr->m[tso.ind1]!=NULL);
1262       tso.syz = pCopy(add_repr->m[tso.ind1]);
1263       poly tt = pMDivide(tso.lcm,tso.p1);
1264       pSetComp(tt,0);
1265       pSetmComp(tt);
1266       pSetCoeff(tt,nDiv(pGetCoeff(tso.p1),coefgcd));
1267       tso.syz = pMult_mm(tso.syz,tt);
1268       pDelete(&tt);
1269       coefgcd = nInpNeg(coefgcd);
1270       assume (add_repr->m[tso.ind2]!=NULL);
1271       p = pCopy(add_repr->m[tso.ind2]);
1272       tt = pMDivide(tso.lcm,tso.p2);
1273       pSetComp(tt,0);
1274       pSetmComp(tt);
1275       pSetCoeff(tt,nDiv(pGetCoeff(tso.p2),coefgcd));
1276       p = pMult_mm(p,tt);
1277       pDelete(&tt);
1278       tso.syz = pAdd(p,tso.syz);
1279       nDelete(&coefgcd);
1280     }
1281 #ifdef SHOW_PROT
1282 Print("reduziere Paar im Module %d mit: \n",index);
1283 PrintS("poly1: ");pWrite(tso.p1);
1284 PrintS("poly2: ");pWrite(tso.p2);
1285 PrintS("syz: ");pWrite(tso.syz);
1286 PrintS("sPoly: ");pWrite(tso.p);
1287 #endif
1288     assume(tso.syz!=NULL);
1289     kBucketInit(syzstr->syz_bucket,tso.syz,-1);
1290     if (tso.p!=NULL)
1291     {
1292       kBucketInit(syzstr->bucket,tso.p,-1);
1293       p = kBucketGetLm(syzstr->bucket);
1294       j = 0;
1295       loop
1296       {
1297         if (j>=*next_place_add) break;
1298         if (pDivisibleBy(add_generators->m[j],p))
1299         {
1300           assume (add_repr->m[j]!=NULL);
1301           sySPRedSyz_Kosz(syzstr,add_generators->m[j],add_repr->m[j],p);
1302           n = kBucketPolyRed(syzstr->bucket,add_generators->m[j],
1303                    pLength(add_generators->m[j]), NULL);
1304           p = kBucketGetLm(syzstr->bucket);
1305           if ((p==NULL) || (pGetComp(p)<=crit_comp)) break;
1306           j = 0;
1307         }
1308         else
1309           j++;
1310       }
1311       kBucketClear(syzstr->bucket,&tso.p,&tso.length);
1312     }
1313     kBucketClear(syzstr->syz_bucket,&tso.syz,&j);
1314   }
1315   else
1316   {
1317     PrintS("Shit happens!\n");
1318   }
1319 #ifdef SHOW_PROT
1320 Print("erhalte Paar im Module %d mit: \n",index);
1321 PrintS("syz: ");pWrite(tso.syz);
1322 PrintS("sPoly: ");pWrite(tso.p);
1323 PrintLn();
1324 #endif
1325   if (tso.p!=NULL)
1326   {
1327     if (!nIsOne(pGetCoeff(tso.p)))
1328     {
1329       n=nInvers(pGetCoeff(tso.p));
1330       pNorm(tso.p);
1331       tso.syz=__p_Mult_nn(tso.syz,n,currRing);
1332       nDelete(&n);
1333     }
1334   }
1335   if ((TEST_OPT_PROT) && (tso.syz==NULL)) PrintS("null");
1336   if ((tso.p!=NULL) && (pGetComp(tso.p)>crit_comp))
1337   {
1338     if (*next_place_add>=IDELEMS(add_generators))
1339     {
1340       pEnlargeSet(&add_generators->m,IDELEMS(add_generators),16);
1341       pEnlargeSet(&add_repr->m,IDELEMS(add_repr),16);
1342       *g_l = (int*)omRealloc0Size((ADDRESS)*g_l, IDELEMS(add_generators)*sizeof(int),
1343                             (IDELEMS(add_generators)+16)*sizeof(int));
1344       IDELEMS(add_generators) += 16;
1345       IDELEMS(add_repr) += 16;
1346     }
1347     assume(add_repr->m[*next_place_add]==NULL);
1348     add_generators->m[*next_place_add] = tso.p;
1349     add_repr->m[*next_place_add] = tso.syz;
1350     (*g_l)[*next_place_add] = tso.length;
1351     (*next_place_add)++;
1352   }
1353   else
1354   {
1355     while ((ng_place>0) && (new_generators->m[ng_place-1]==NULL) &&
1356           (new_repr->m[ng_place-1]==NULL)) ng_place--;
1357     if (ng_place>=IDELEMS(new_generators))
1358     {
1359       pEnlargeSet(&new_generators->m,IDELEMS(new_generators),16);
1360       IDELEMS(new_generators) += 16;
1361       pEnlargeSet(&new_repr->m,IDELEMS(new_repr),16);
1362       IDELEMS(new_repr) += 16;
1363     }
1364     new_generators->m[ng_place] = tso.p;
1365     new_repr->m[ng_place] = tso.syz;
1366   }
1367   tso.p = NULL;
1368   tso.syz = NULL;
1369   resPairs[itso] = tso;
1370   syDeletePair(&resPairs[itso]);
1371 }
1372 
1373 /*3
1374 * reduction of all pairs of a fixed degree of a fixed module
1375 */
reducePairsHIndex(SSet resPairs,int l_pairs,syStrategy syzstr,int index,ideal add_generators,ideal add_repr,ideal new_generators,ideal new_repr,int crit_comp,int * red_deg,int * next_place_add,int ** g_l,resolvente totake)1376 static BOOLEAN reducePairsHIndex(SSet resPairs,int l_pairs,syStrategy syzstr,
1377        int index,ideal add_generators,ideal add_repr,ideal new_generators,
1378        ideal new_repr,int crit_comp,int * red_deg,int * next_place_add,int **g_l,
1379        resolvente totake)
1380 {
1381   if (resPairs[0].lcm==NULL) return FALSE;
1382   int i=0;
1383   poly deg_soc;
1384 
1385   if (TEST_OPT_PROT)
1386     Print("(%d,%d)",index,resPairs[0].order);
1387   while ((i<l_pairs) && (resPairs[i].order==*red_deg))
1388   {
1389     assume(totake[index-1]!=NULL);
1390     assume(pGetComp(resPairs[i].p1)<=IDELEMS(totake[index-1]));
1391     assume(totake[index-1]->m[pGetComp(resPairs[i].p1)-1]!=NULL);
1392     deg_soc = totake[index-1]->m[pGetComp(resPairs[i].p1)-1];
1393     redOnePairHIndex(resPairs,i,crit_comp,syzstr,index, add_generators,add_repr,
1394                      new_generators, new_repr,next_place_add,g_l,deg_soc);
1395     i++;
1396   }
1397   syCompactifyPairSet(resPairs,l_pairs,0);
1398   if (resPairs[0].lcm==NULL)  //there are no pairs left and no new_gens
1399     return FALSE;
1400   else
1401     *red_deg = resPairs[0].order;
1402   return TRUE;
1403 }
1404 
1405 /*3
1406 * we proceed the generators of the next module;
1407 * they are stored in add_generators and add_repr;
1408 * if the normal form of a new genrators w.r.t. add_generators has
1409 * pGetComp<crit_comp it is skipped from the reduction;
1410 * new_generators and new_repr (which are empty) stores the result of the
1411 * reduction which is normalized afterwards
1412 */
procedeNextGenerators(ideal temp_generators,ideal,ideal new_generators,ideal new_repr,ideal add_generators,ideal add_repr,syStrategy syzstr,int index,int crit_comp,resolvente totake)1413 static void procedeNextGenerators(ideal temp_generators,ideal /*temp_repr*/,
1414       ideal new_generators, ideal new_repr, ideal add_generators,
1415       ideal add_repr, syStrategy syzstr,int index, int crit_comp,
1416       resolvente totake)
1417 {
1418   int i=0,j,next_new_el;
1419   int idel_temp=IDELEMS(temp_generators);
1420   int next_place_add;
1421   int p_length,red_deg,l_pairs=IDELEMS(add_generators);
1422   poly next_p;
1423   int * gen_length=(int*)omAlloc0(IDELEMS(add_generators)*sizeof(int));
1424   int * secgen_length=(int*)omAlloc0(IDELEMS(syzstr->res[index])*sizeof(int));
1425   BOOLEAN pairs_left;
1426   SSet resPairs=(SSet)omAlloc0(l_pairs*sizeof(SObject));
1427 
1428   for (j=IDELEMS(syzstr->res[index])-1;j>=0;j--)
1429   {
1430     if (syzstr->res[index]->m[j]!=NULL)
1431       secgen_length[j] = pLength(syzstr->res[index]->m[j]);
1432   }
1433   assume(idIs0(new_generators));
1434   next_place_add = IDELEMS(add_generators);
1435   while ((next_place_add>0) && (add_generators->m[next_place_add-1]==NULL))
1436     next_place_add--;
1437   int next_deg = p_FDeg(temp_generators->m[i],currRing);
1438   next_new_el = next_place_add;
1439 /*--- loop about all all elements-----------------------------------*/
1440   while ((i<idel_temp) && (temp_generators->m[i]!=NULL))
1441   {
1442 /*--- separates elements of equal degree----------------------------*/
1443 #ifdef USE_REGULARITY
1444     if (syzstr->regularity>0)
1445     {
1446       if (next_deg >= syzstr->regularity+index)
1447       {
1448         while ((i<idel_temp) && (temp_generators->m[i]!=NULL))
1449         {
1450           pDelete(&temp_generators->m[i]);
1451           i++;
1452         }
1453         break;
1454       }
1455     }
1456 #endif
1457     while ((i<idel_temp) && (p_FDeg(temp_generators->m[i],currRing)==next_deg))
1458     {
1459       next_p = temp_generators->m[i];
1460       temp_generators->m[i] = NULL;
1461       next_p = normalize(next_p,add_generators,syzstr,gen_length,&p_length,
1462                 crit_comp);
1463       if (next_p!=NULL)
1464       {
1465         if (pGetComp(next_p)<=crit_comp)
1466         {
1467           pDelete(&next_p);
1468           //if (TEST_OPT_PROT) Print("u(%d)",index);
1469         }
1470         else
1471         {
1472           next_p = syRedTailSyz(next_p,add_generators,syzstr->res[index],crit_comp,syzstr,
1473             gen_length,secgen_length,&p_length);
1474           if (!nIsOne(pGetCoeff(next_p)))
1475             pNorm(next_p);
1476           if (next_place_add>=IDELEMS(add_generators))
1477           {
1478             pEnlargeSet(&add_generators->m,IDELEMS(add_generators),16);
1479             pEnlargeSet(&add_repr->m,IDELEMS(add_repr),16);
1480             gen_length = (int*)omRealloc0Size((ADDRESS)gen_length, IDELEMS(add_generators)*sizeof(int),
1481                                         (IDELEMS(add_generators)+16)*sizeof(int));
1482             IDELEMS(add_generators) += 16;
1483             IDELEMS(add_repr) += 16;
1484           }
1485           add_generators->m[next_place_add] = next_p;
1486           if (totake[index]==NULL)
1487             totake[index] = idInit(16,new_generators->rank);
1488           if ((*syzstr->Tl)[index]==IDELEMS(totake[index]))
1489           {
1490             pEnlargeSet(&totake[index]->m,IDELEMS(totake[index]),
1491                         (*syzstr->Tl)[index]+16-IDELEMS(totake[index]));
1492             for (j=IDELEMS(totake[index]);j<(*syzstr->Tl)[index]+16;j++)
1493               totake[index]->m[j] = NULL;
1494             IDELEMS(totake[index]) = (*syzstr->Tl)[index]+16;
1495           }
1496 #ifdef FULL_TOTAKE
1497           totake[index]->m[(*syzstr->Tl)[index]] = pCopy(next_p);
1498 #else
1499           totake[index]->m[(*syzstr->Tl)[index]] = pHead(next_p);
1500 #endif
1501           assume(add_repr->m[next_place_add]==NULL);
1502 #ifdef WITH_SCHREYER_ORD
1503           add_repr->m[next_place_add] = pHead(add_generators->m[next_place_add]);
1504 #else
1505           add_repr->m[next_place_add] = pOne();
1506 #endif
1507           ((*syzstr->Tl)[index])++;
1508           pSetComp(add_repr->m[next_place_add],(*syzstr->Tl)[index]);
1509           pSetmComp(add_repr->m[next_place_add]);
1510           gen_length[next_place_add] = p_length;
1511           next_place_add++;
1512         }
1513       }
1514       i++;
1515     }                        //end inner loop
1516     red_deg = next_deg;
1517     if (i<idel_temp)
1518       next_deg = p_FDeg(temp_generators->m[i],currRing);
1519     else
1520       next_deg = -1;
1521     if ((next_place_add>next_new_el) || (next_deg<0))  //there are new generators or pairs
1522     {
1523 /*-reducing and generating pairs untill the degree of the next generators-*/
1524       pairs_left = TRUE;
1525       while (pairs_left && ((next_deg<0) || (red_deg<= next_deg)))
1526       {
1527         updatePairsHIndex(&resPairs,&l_pairs,syzstr,index,add_generators,
1528           add_repr,new_generators,new_repr,crit_comp,&next_new_el);
1529         pairs_left = reducePairsHIndex(resPairs,l_pairs,syzstr,index,add_generators,
1530            add_repr,new_generators,new_repr,crit_comp,&red_deg,&next_place_add,&gen_length,
1531            totake);
1532       }
1533     }
1534   }
1535   omFreeSize((SSet)resPairs,l_pairs*sizeof(SObject));
1536   omFreeSize((ADDRESS)gen_length,IDELEMS(add_generators)*sizeof(int));
1537   omFreeSize((ADDRESS)secgen_length,IDELEMS(syzstr->res[index])*sizeof(int));
1538 }
1539 
1540 /*3
1541 * normalizes the part of the next reduction lying within the block
1542 * of former generators (old_generators);
1543 */
normalizeOldPart(ideal new_generators,ideal new_repr,syStrategy syzstr,int index,int)1544 static ideal normalizeOldPart(ideal new_generators,ideal new_repr,
1545                       syStrategy syzstr,int index,int /*crit_comp*/)
1546 {
1547   ideal old_generators= syzstr->res[index];
1548   ideal old_repr= syzstr->orderedRes[index];
1549   int i,j=0,ii=IDELEMS(old_generators)-1,dummy;
1550   poly p;
1551   number n;
1552   int * g_l=(int*)omAlloc0(IDELEMS(old_generators)*sizeof(int));
1553 
1554   for (i=0;i<IDELEMS(old_generators);i++)
1555   {
1556     if (old_generators->m[i]!=NULL)
1557     {
1558       g_l[i] = pLength(old_generators->m[i]);
1559     }
1560   }
1561   for (i=IDELEMS(new_generators)-1;i>=0;i--)
1562   {
1563     if (new_generators->m[i]!=NULL)
1564     {
1565       kBucketInit(syzstr->bucket,new_generators->m[i],
1566                    pLength(new_generators->m[i]));
1567       kBucketInit(syzstr->syz_bucket,new_repr->m[i],
1568                    pLength(new_repr->m[i]));
1569       p = kBucketGetLm(syzstr->bucket);
1570       loop
1571       {
1572         if ((j>=ii) || (p==NULL)) break;
1573         if ((old_generators->m[j]!=NULL) &&
1574             (pDivisibleBy(old_generators->m[j],p)))
1575         {
1576           sySPRedSyz_Kosz(syzstr,old_generators->m[j],old_repr->m[j],p);
1577           n = kBucketPolyRed(syzstr->bucket,old_generators->m[j], g_l[j], NULL);
1578           nDelete(&n);
1579           p = kBucketGetLm(syzstr->bucket);
1580           j = 0;
1581         }
1582         else
1583           j++;
1584       }
1585       assume (p==NULL);
1586       kBucketClear(syzstr->bucket,&new_generators->m[i],&dummy);
1587       kBucketClear(syzstr->syz_bucket,&new_repr->m[i],&dummy);
1588     }
1589   }
1590   ideal result=idInit(IDELEMS(new_repr),new_repr->rank);
1591   for (j=IDELEMS(new_repr)-1;j>=0;j--)
1592   {
1593     result->m[j] = new_repr->m[j];
1594     if ((result->m[j]!=NULL) && (!nIsOne(pGetCoeff(result->m[j]))))
1595       pNorm(result->m[j]);
1596     new_repr->m[j] = NULL;
1597   }
1598   omFreeSize((ADDRESS)g_l,IDELEMS(old_generators)*sizeof(int));
1599   return result;
1600 }
1601 
1602 /*3
1603 * constructs the new subresolution for a nonregular extension
1604 */
kosz_ext(ideal new_generators,ideal new_repr,syStrategy syzstr,int index,int next_comp,resolvente totake)1605 static ideal kosz_ext(ideal new_generators,ideal new_repr,syStrategy syzstr,
1606                       int index,int next_comp,resolvente totake)
1607 {
1608   ideal temp_generators =idInit(IDELEMS(new_generators),new_generators->rank);
1609   ideal temp_repr=idInit(IDELEMS(new_repr),new_repr->rank);
1610   ideal add_generators =idInit(IDELEMS(new_generators),new_generators->rank);
1611   ideal add_repr=idInit(IDELEMS(new_repr),new_repr->rank);
1612   int min_deg=-1;
1613   int j,jj,k,deg_p,idel_temp=IDELEMS(temp_generators);
1614   poly p;
1615 /*--reorder w.r.t. the degree----------------------------------------*/
1616   for (j=IDELEMS(new_generators)-1;j>=0;j--)
1617   {
1618     if (new_generators->m[j]!=NULL)
1619     {
1620       p = new_generators->m[j];
1621       new_generators->m[j] = NULL;
1622       deg_p = p_FDeg(p,currRing);
1623       if (min_deg<0)
1624       {
1625         min_deg = deg_p;
1626       }
1627       else
1628       {
1629         if (deg_p<min_deg) min_deg = deg_p;
1630       }
1631       k = 0;
1632       while ((k<idel_temp) && (temp_generators->m[k]!=NULL) &&
1633              (p_FDeg(temp_generators->m[k],currRing)<=deg_p)) k++;
1634       for (jj=idel_temp-1;jj>k;jj--)
1635       {
1636         temp_generators->m[jj] = temp_generators->m[jj-1];
1637       }
1638       temp_generators->m[k] = p;
1639     }
1640   }
1641 /*--- computing the standard basis in the resolution of the extension -*/
1642   procedeNextGenerators(temp_generators,temp_repr,new_generators,new_repr,
1643     add_generators,add_repr,syzstr,index,next_comp,totake);
1644   j = IDELEMS(syzstr->res[index]);
1645   while ((j>0) && (syzstr->res[index]->m[j-1]==NULL)) j--;
1646   jj = IDELEMS(add_generators);
1647   while ((jj>0) && (add_generators->m[jj-1]==NULL)) jj--;
1648   if (j+jj>=IDELEMS(syzstr->res[index]))
1649   {
1650     pEnlargeSet(&syzstr->res[index]->m,IDELEMS(syzstr->res[index]),
1651                 j+jj+1-IDELEMS(syzstr->res[index]));
1652     IDELEMS(syzstr->res[index]) = j+jj+1;
1653     pEnlargeSet(&syzstr->orderedRes[index]->m,IDELEMS(syzstr->orderedRes[index]),
1654                 j+jj+1-IDELEMS(syzstr->orderedRes[index]));
1655     IDELEMS(syzstr->orderedRes[index]) = j+jj+1;
1656   }
1657   for (k=0;k<jj;k++)
1658   {
1659     syzstr->res[index]->m[j+k] = add_generators->m[k];
1660     syzstr->orderedRes[index]->m[j+k] = add_repr->m[k];
1661     add_generators->m[k] = NULL;
1662     add_repr->m[k] = NULL;
1663   }
1664   assume(idIs0(add_generators));
1665   assume(idIs0(add_repr));
1666   idDelete(&add_generators);
1667   idDelete(&add_repr);
1668   idDelete(&temp_generators);
1669   idDelete(&temp_repr);
1670 /*--- normalizing the rest to get the syzygies ------------------------*/
1671   return normalizeOldPart(new_generators,new_repr,syzstr,index,next_comp);
1672 }
1673 
1674 /*
1675 * this procedure assumes that the first order is C !!!
1676 * INPUT: old_generators - the generators of the actual module
1677 *                         computed so far (they are mixed vectors)
1678 *        old_repr       - the representations of the old generators
1679 *        new_generators - generators coming from reductions below,
1680 *                         they must have leading terms in new components
1681 *                         (they live only in the module part)
1682 *  (*syzstr->Tl)[index] - the last used component in the syzygy
1683 * OUTPUT: old_generators is updated
1684 *         new_generators is empty
1685 *         the return value is a set of new generators for the syzygies,
1686 */
syAppendSyz(ideal new_generators,syStrategy syzstr,int index,int crit_comp,resolvente totake)1687 static ideal syAppendSyz(ideal new_generators, syStrategy syzstr,int index,int crit_comp,
1688                          resolvente totake)
1689 {
1690   int i,j;
1691   ideal result;
1692   int rk_new_gens = id_RankFreeModule(new_generators,currRing);
1693   if (syzstr->res[index]==NULL)
1694   {
1695     syzstr->res[index] = idInit(1,si_max(rk_new_gens,1));
1696     syzstr->orderedRes[index] = idInit(1,si_max(rk_new_gens,1));
1697   }
1698   int ng_idel=IDELEMS(new_generators);
1699   ideal new_repr =idInit(ng_idel, crit_comp+ng_idel);
1700 
1701   if (index==0)
1702   {
1703     //int * og_l=(int*)omAlloc0(IDELEMS(syzstr->res[0])*sizeof(int));
1704     //for (i=IDELEMS(syzstr->res[0])-1;i>=0;i--)
1705     //{
1706       //if (syzstr->res[0]->m[i]!=NULL)
1707         //og_l[i] = pLength(syzstr->res[0]->m[i]);
1708     //}
1709     for (i=0;i<ng_idel;i++)
1710     {
1711       if (new_generators->m[i]!=NULL)
1712       {
1713         //int ng_l=pLength(new_generators->m[i]);
1714         //new_generators->m[i] = syRedTailSyz(new_generators->m[i],syzstr->res[0],NULL,0,syzstr,
1715             //og_l,NULL,&ng_l);
1716         if (totake[index]==NULL)
1717           totake[index] = idInit(16,new_generators->rank);
1718         if ((*syzstr->Tl)[index]>=IDELEMS(totake[index]))
1719         {
1720           pEnlargeSet(&totake[index]->m,IDELEMS(totake[index]),
1721                       (*syzstr->Tl)[index]+16-IDELEMS(totake[index]));
1722           for (j=IDELEMS(totake[index]);j<(*syzstr->Tl)[index]+16;j++)
1723             totake[index]->m[j] = NULL;
1724           IDELEMS(totake[index]) = (*syzstr->Tl)[index]+16;
1725         }
1726 #ifdef FULL_TOTAKE
1727         totake[index]->m[(*syzstr->Tl)[index]] = pCopy(new_generators->m[i]);
1728 #else
1729         totake[index]->m[(*syzstr->Tl)[index]] = pHead(new_generators->m[i]);
1730 #endif
1731 #ifdef WITH_SCHREYER_ORD
1732         new_repr->m[i] = pHead(new_generators->m[i]);
1733 #else
1734         new_repr->m[i] = pOne();
1735 #endif
1736         ((*syzstr->Tl)[index])++;
1737         pSetComp(new_repr->m[i],(*syzstr->Tl)[index]);
1738         pSetmComp(new_repr->m[i]);
1739       }
1740     }
1741     //omFreeSize((ADDRESS)og_l,IDELEMS(syzstr->res[0])*sizeof(int));
1742 #ifdef SHOW_PROT
1743 PrintS("Add new generators:\n");
1744 idPrint(new_generators);
1745 PrintS("with representaions:\n");
1746 idPrint(new_repr);
1747 #endif
1748     result = kosz_std(new_generators,new_repr,syzstr,index,crit_comp);
1749   }
1750   else
1751   {
1752     result = kosz_ext(new_generators,new_repr,syzstr,index,crit_comp,totake);
1753   }
1754   idSkipZeroes(result);
1755   assume(idIs0(new_repr));
1756   idDelete(&new_repr);
1757   return result;
1758 }
1759 
1760 /*
1761 * main call of the extended Koszul-resolution
1762 */
syKosz(ideal arg,int * length)1763 syStrategy syKosz(ideal arg,int * length)
1764 {
1765   int i,j,jj,k=0,index=0,rk_arg/*,next_syz=0*/;
1766   int crit_comp,t_comp,next_deg,old_tl;
1767   ideal temp=NULL,old_ideal,old_repr;
1768   ring origR = currRing;
1769   poly next_gen;
1770   BOOLEAN isRegular;
1771 
1772   discard_pairs = 0;
1773   short_pairs = 0;
1774   if (idIs0(arg)) return NULL;
1775   rk_arg = id_RankFreeModule(arg,currRing);
1776   syStrategy syzstr=(syStrategy)omAlloc0(sizeof(ssyStrategy));
1777 /*--- changes to a Cdp-ring ----------------------------*/
1778   syzstr->syRing = rAssure_C_dp(origR); rChangeCurrRing(syzstr->syRing);
1779 /*--- initializes the data structures---------------*/
1780   syzstr->length = *length = (syzstr->syRing->N)+2;
1781   syzstr->regularity = -1;
1782   if (origR!=syzstr->syRing)
1783     temp = idrCopyR(arg, origR, syzstr->syRing);
1784   else
1785     temp = idCopy(arg);
1786   if (rk_arg==0)
1787   {
1788     id_Shift(temp,1,currRing);
1789   }
1790   idSkipZeroes(temp);
1791 #ifdef WITH_SORT
1792   if (temp->m[0]!=NULL)
1793   {
1794     int md;
1795     int maxdeg=p_FDeg(temp->m[IDELEMS(temp)-1],currRing);
1796     ideal temp1=idInit(IDELEMS(temp),temp->rank);
1797     for (j=IDELEMS(temp)-2;j>=0;j--)
1798     {
1799       jj = p_FDeg(temp->m[j],currRing);
1800       if (jj>maxdeg) maxdeg = jj;
1801     }
1802     while (!idIs0(temp))
1803     {
1804       md = maxdeg;
1805       for (j=IDELEMS(temp)-1;j>=0;j--)
1806       {
1807         if (temp->m[j]!=NULL)
1808         {
1809           jj = p_FDeg(temp->m[j],currRing);
1810           if (jj<md) md = jj;
1811         }
1812       }
1813       for (j=0;j<IDELEMS(temp);j++)
1814       {
1815         if ((temp->m[j]!=NULL) && (p_FDeg(temp->m[j],currRing)==md))
1816         {
1817           temp1->m[k] = temp->m[j];
1818           temp->m[j] = NULL;
1819           k++;
1820         }
1821       }
1822     }
1823     idDelete(&temp);
1824     temp = temp1;
1825     temp1 = NULL;
1826   }
1827 #endif
1828 #ifdef USE_REGULARITY
1829   int last_generator=IDELEMS(temp)-1;
1830   while ((last_generator>=0) && (temp->m[last_generator]==NULL))
1831     last_generator--;
1832 #endif
1833   syzstr->res = (resolvente)omAlloc0((*length+1)*sizeof(ideal));
1834   syzstr->orderedRes = (resolvente)omAlloc0((*length+1)*sizeof(ideal));
1835   resolvente totake=(resolvente)omAlloc0((*length+1)*sizeof(ideal));
1836   syzstr->Tl = new intvec(*length+1);
1837   syzstr->bucket = kBucketCreate(currRing);
1838   syzstr->syz_bucket = kBucketCreate(currRing);
1839   ideal new_generators=idInit(1,si_max(rk_arg,1));
1840   ideal temp_gens,old_std;
1841   syzstr->res[0] = idInit(1,1);
1842   if (rk_arg>1) syzstr->res[0]->rank = rk_arg;
1843   syzstr->orderedRes[0] = idInit(1,1);
1844 /*--- computes the resolution ----------------------*/
1845   i = 0;
1846   while (i<IDELEMS(temp))
1847   {
1848     if (temp->m[i]!=NULL)
1849     {
1850       new_generators->m[0] = kNF(syzstr->res[0],currRing->qideal,temp->m[i]);
1851       if (!nIsOne(pGetCoeff(new_generators->m[0])))
1852         pNorm(new_generators->m[0]);
1853       next_deg = p_FDeg(new_generators->m[0],currRing);
1854       next_gen = pCopy(new_generators->m[0]);
1855     }
1856     if (!idIs0(new_generators))
1857     {
1858       index = 0;
1859       while (index<=*length)
1860       {
1861         if (index==0)
1862         {
1863           old_ideal = idCopy(syzstr->res[0]);
1864           old_repr = idCopy(syzstr->orderedRes[0]);
1865           old_tl = (*syzstr->Tl)[0];
1866           old_std = id_Head(syzstr->res[0],currRing);
1867         }
1868         t_comp = (*syzstr->Tl)[index];
1869         if (index==0) crit_comp = t_comp;
1870         temp_gens = syAppendSyz(new_generators,syzstr, index,crit_comp,totake);
1871         crit_comp = t_comp;
1872         if (index==0)
1873         {
1874           isRegular = syIsRegular(old_std,syzstr->res[0],next_deg);
1875 #ifndef ONLY_STD
1876           if (isRegular)
1877             syCreateRegularExtension(syzstr,old_ideal,old_repr,old_tl,next_gen,
1878                                      totake);
1879 #ifdef USE_REGULARITY
1880         if ((index==0) && (!isRegular) && (i==last_generator))
1881         {
1882 /*----------- we are computing the regularity -----------------------*/
1883           ideal initial=id_Head(syzstr->res[0],currRing);
1884           int len=0,reg=0;
1885           intvec *w=NULL;
1886           ring dp_C_ring = rAssure_dp_C(currRing); rChangeCurrRing(dp_C_ring);
1887           initial = idrMoveR_NoSort(initial, syzstr->syRing, dp_C_ring);
1888           resolvente res = sySchreyerResolvente(initial,-1,&len,TRUE, TRUE);
1889           intvec * dummy = syBetti(res,len,&reg, w);
1890           syzstr->regularity = reg+2;
1891           delete dummy;
1892           delete w;
1893           for (j=0;j<len;j++)
1894           {
1895             if (res[j]!=NULL) idDelete(&(res[j]));
1896           }
1897           omFreeSize((ADDRESS)res,len*sizeof(ideal));
1898           idDelete(&initial);
1899           rChangeCurrRing(syzstr->syRing);
1900           rDelete(dp_C_ring);
1901         }
1902 #endif
1903 #endif
1904           idDelete(&old_ideal);
1905           idDelete(&old_repr);
1906           idDelete(&old_std);
1907           if (TEST_OPT_PROT)
1908           {
1909             if (isRegular)
1910               PrintS("\n regular\n");
1911             else
1912               PrintS("\n not regular\n");
1913           }
1914           if (next_gen!=NULL)
1915             pDelete(&next_gen);
1916           if (isRegular)
1917           {
1918             idDelete(&temp_gens);
1919             break;
1920           }
1921         }
1922         idDelete(&new_generators);
1923         new_generators = temp_gens;
1924 #ifdef ONLY_STD
1925         break;
1926 #endif
1927         if (idIs0(new_generators)) break;
1928         index++;
1929       }
1930       if (!idIs0(new_generators))
1931       {
1932         for (j=0;j<IDELEMS(new_generators);j++)
1933         {
1934           if (new_generators->m[j]!=NULL)
1935           {
1936             pDelete(&new_generators->m[j]);
1937             new_generators->m[j] = NULL;
1938           }
1939         }
1940       }
1941     }
1942     i++;
1943   }
1944   if (idIs0(new_generators) && new_generators!=NULL) idDelete(&new_generators);
1945   if (temp!=NULL) idDelete(&temp);
1946   kBucketDestroy(&(syzstr->bucket));
1947   kBucketDestroy(&(syzstr->syz_bucket));
1948   index = 0;
1949   syzstr->fullres = syzstr->res;
1950   syzstr->res = NULL;
1951   index = 0;
1952   while ((index<=*length) && (syzstr->fullres[index]!=NULL))
1953   {
1954 #ifdef SHOW_RESULT
1955     Print("The %d-th syzygy-module is now:\n",index);
1956     ideal ttt=id_Head(syzstr->fullres[index],currRing);
1957     idShow(ttt);
1958     idDelete(&ttt);
1959     //if (index>0)
1960     //{
1961       //Print("The related module is: \n");
1962       //idPrint(totake[index-1]);
1963     //}
1964     //Print("The %d-th module of the minimal resolution is:\n",index);
1965     if (!idIs0(totake[index]))
1966       idShow(totake[index]);
1967     //Print("with standard basis:\n");
1968     //idPrint(syzstr->fullres[index]);
1969     //if ((index<*length) && (totake[index+1]!=NULL))
1970     //{
1971       //Print("The %d-th syzygy-module is now:\n",index+1);
1972       //idPrint(totake[index+1]);
1973       //matrix m1=idModule2Matrix(totake[index]);
1974       //matrix m2=idModule2Matrix(totake[index+1]);
1975       //matrix m3=mpMult(m1,m2);
1976       //idPrint((ideal)m3);
1977     //}
1978 #endif
1979     if (!idIs0(totake[index]))
1980     {
1981       for(i=0;i<IDELEMS(totake[index]);i++)
1982       {
1983         if (totake[index]->m[i]!=NULL)
1984         {
1985           j=0;
1986           while ((j<IDELEMS(syzstr->fullres[index])) &&
1987             ((syzstr->fullres[index]->m[j]==NULL) ||
1988             (!pLmEqual(syzstr->fullres[index]->m[j],totake[index]->m[i])))) j++;
1989           if (j<IDELEMS(syzstr->fullres[index]))
1990           {
1991             pDelete(&totake[index]->m[i]);
1992             totake[index]->m[i] = syzstr->fullres[index]->m[j];
1993             syzstr->fullres[index]->m[j] = NULL;
1994           }
1995           else
1996           {
1997             PrintS("Da ist was faul!!!\n");
1998             Print("Aber: Regularitaet %d, Grad %ld\n",
1999                    syzstr->regularity,p_FDeg(totake[index]->m[i],currRing));
2000           }
2001         }
2002       }
2003       idDelete(&syzstr->fullres[index]);
2004       syzstr->fullres[index] = totake[index];
2005     }
2006 #ifdef SHOW_RESULT
2007     idShow(syzstr->fullres[index]);
2008 #endif
2009     index++;
2010   }
2011   syReorder_Kosz(syzstr);
2012   index = 0;
2013   while ((index<=*length) && (syzstr->orderedRes[index]!=NULL))
2014   {
2015     idDelete(&(syzstr->orderedRes[index]));
2016     index++;
2017   }
2018   if (origR!=syzstr->syRing)
2019   {
2020     rChangeCurrRing(origR);
2021     index = 0;
2022     while ((index<=*length) && (syzstr->fullres[index]!=NULL))
2023     {
2024       syzstr->fullres[index] = idrMoveR(syzstr->fullres[index],syzstr->syRing, origR);
2025       index++;
2026     }
2027   }
2028   delete syzstr->Tl;
2029   syzstr->Tl = NULL;
2030   rDelete(syzstr->syRing);
2031   syzstr->syRing = NULL;
2032   omFreeSize((ADDRESS)totake,(*length+1)*sizeof(ideal));
2033   omFreeSize((ADDRESS)syzstr->orderedRes,(*length+1)*sizeof(ideal));
2034 //Print("Pairs to discard: %d\n",discard_pairs);
2035 //Print("Pairs shorter reduced: %d\n",short_pairs);
2036 //discard_pairs = 0;
2037 //short_pairs = 0;
2038   return syzstr;
2039 }
2040 
2041