source: git/kernel/GBEngine/syz3.cc @ 922890

spielwiese
Last change on this file since 922890 was 922890, checked in by Hans Schoenemann <hannes@…>, 9 years ago
moved stairc.h to kernel/combinatorics/
  • Property mode set to 100644
File size: 61.2 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: resolutions
6*/
7
8
9
10
11#include <kernel/mod2.h>
12#include <misc/mylimits.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>
18#include <kernel/combinatorics/stairc.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#include <polys/matpol.h>
29
30//#define SHOW_PROT
31//#define SHOW_RED
32//#define SHOW_Kosz
33//#define SHOW_RESULT
34//#define INVERT_PAIRS
35//#define ONLY_STD
36//#define EXPERIMENT1    //Hier stimmt was mit der Anzahl der Erzeuger in xyz11 nicht!!
37#define EXPERIMENT2
38#define EXPERIMENT3
39#define WITH_BUCKET     //Use of buckets in EXPERIMENT3 (Product criterion)
40#define WITH_SCHREYER_ORD
41#define USE_CHAINCRIT
42#define USE_CHAINCRIT0
43#define USE_PROD_CRIT
44#define USE_REGULARITY
45#define WITH_SORT
46//#define FULL_TOTAKE
47int discard_pairs;
48int short_pairs;
49
50/*3
51* assumes the ideals old_ideal and new_ideal to be homogeneous
52* tests wether the new_ideal is a regular extension of the old_ideal
53*/
54static BOOLEAN syIsRegular(ideal old_ideal,ideal new_ideal,int deg)
55{
56  intvec * old_hilbs=hHstdSeries(old_ideal,NULL,NULL,NULL);
57  intvec * new_hilbs=hHstdSeries(new_ideal,NULL,NULL,NULL);
58  int biggest_length=si_max(old_hilbs->length()+deg,new_hilbs->length());
59  intvec * shifted_old_hilbs=new intvec(biggest_length);
60  intvec * old_hilb1=new intvec(biggest_length);
61  intvec * new_hilb1=new intvec(biggest_length);
62  int i;
63  BOOLEAN isRegular=TRUE;
64
65  for (i=old_hilbs->length()+deg-1;i>=deg;i--)
66    (*shifted_old_hilbs)[i] = (*old_hilbs)[i-deg];
67  for (i=old_hilbs->length()-1;i>=0;i--)
68    (*old_hilb1)[i] = (*old_hilbs)[i]-(*shifted_old_hilbs)[i];
69  for (i=old_hilbs->length()+deg-1;i>=old_hilbs->length();i--)
70    (*old_hilb1)[i] = -(*shifted_old_hilbs)[i];
71  for (i=new_hilbs->length()-1;i>=0;i--)
72    (*new_hilb1)[i] = (*new_hilbs)[i];
73  i = 0;
74  while ((i<biggest_length) && isRegular)
75  {
76    isRegular = isRegular && ((*old_hilb1)[i] == (*new_hilb1)[i]);
77    i++;
78  }
79  delete old_hilbs;
80  delete new_hilbs;
81  delete old_hilb1;
82  delete new_hilb1;
83  delete shifted_old_hilbs;
84  return isRegular;
85}
86
87/*3
88* shows the resolution stored in syzstr->orderedRes
89*/
90#if 0 /* unused*/
91static void syShowRes(syStrategy syzstr)
92{
93  int i=0;
94
95  while ((i<syzstr->length) && (!idIs0(syzstr->res[i])))
96  {
97    Print("aktueller hoechster index ist: %d\n",(*syzstr->Tl)[i]);
98    Print("der %d-te modul ist:\n",i);
99    idPrint(syzstr->res[i]);
100    PrintS("Seine Darstellung:\n");
101    idPrint(syzstr->orderedRes[i]);
102    i++;
103  }
104}
105#endif
106
107/*3
108* produces the next subresolution for a regular extension
109*/
110static void syCreateRegularExtension(syStrategy syzstr,ideal old_ideal,
111            ideal old_repr,int old_tl, poly next_generator,resolvente totake)
112{
113  int index=syzstr->length-1,i,j,start,start_ttk/*,new_tl*/;
114  poly gen=pCopy(next_generator),p;
115  poly neg_gen=pCopy(next_generator);
116  ideal current_ideal,current_repr;
117  int current_tl;
118  poly w_gen=pHead(next_generator);
119  pSetComp(w_gen,0);
120  pSetmComp(w_gen);
121
122  //syShowRes(syzstr);
123  neg_gen = pNeg(neg_gen);
124  if (pGetComp(gen)>0)
125  {
126    p_Shift(&gen,-1,currRing);
127    p_Shift(&neg_gen,-1,currRing);
128  }
129  while (index>0)
130  {
131    if (index%2==0)
132      p = gen;
133    else
134      p = neg_gen;
135    if (index>1)
136    {
137      current_ideal = syzstr->res[index-1];
138      current_repr = syzstr->orderedRes[index-1];
139      current_tl = (*syzstr->Tl)[index-1];
140    }
141    else
142    {
143      current_ideal = old_ideal;
144      current_repr = old_repr;
145      current_tl = old_tl;
146    }
147    if (!idIs0(current_ideal))
148    {
149      if (idIs0(syzstr->res[index]))
150      {
151        syzstr->res[index] = idInit(IDELEMS(current_ideal),
152          current_ideal->rank+current_tl);
153        syzstr->orderedRes[index] = idInit(IDELEMS(current_ideal),
154          current_ideal->rank);
155        start = 0;
156      }
157      else
158      {
159        start = IDELEMS(syzstr->res[index]);
160        while ((start>0) && (syzstr->res[index]->m[start-1]==NULL)) start--;
161        if (IDELEMS(syzstr->res[index])<start+IDELEMS(current_ideal))
162        {
163          pEnlargeSet(&syzstr->res[index]->m,IDELEMS(syzstr->res[index]),
164                   IDELEMS(current_ideal));
165          IDELEMS(syzstr->res[index]) += IDELEMS(current_ideal);
166          pEnlargeSet(&syzstr->orderedRes[index]->m,IDELEMS(syzstr->orderedRes[index]),
167                   IDELEMS(current_ideal));
168          IDELEMS(syzstr->orderedRes[index]) += IDELEMS(current_ideal);
169        }
170      }
171      if (idIs0(totake[index]))
172      {
173        totake[index] = idInit(IDELEMS(current_ideal),
174          current_ideal->rank+current_tl);
175        start_ttk = 0;
176      }
177      else
178      {
179        start_ttk = IDELEMS(totake[index]);
180        while ((start_ttk>0) && (totake[index]->m[start_ttk-1]==NULL)) start_ttk--;
181        if (IDELEMS(totake[index])<start_ttk+IDELEMS(current_ideal))
182        {
183          pEnlargeSet(&totake[index]->m,IDELEMS(totake[index]),
184                   IDELEMS(current_ideal));
185          for (j=IDELEMS(totake[index]);j<IDELEMS(totake[index])+
186                                  IDELEMS(current_ideal);j++)
187            totake[index]->m[j] = NULL;
188          IDELEMS(totake[index]) += IDELEMS(current_ideal);
189        }
190      }
191      for (i=0;i<IDELEMS(current_ideal);i++)
192      {
193        if (current_ideal->m[i]!=NULL)
194        {
195          syzstr->res[index]->m[i+start] = pCopy(current_ideal->m[i]);
196          syzstr->res[index]->m[i+start] = pMult_mm(syzstr->res[index]->m[i+start],w_gen);
197          p_Shift(&syzstr->res[index]->m[i+start],current_tl,currRing);
198          syzstr->res[index]->m[i+start] = pAdd(syzstr->res[index]->m[i+start],
199            ppMult_qq(current_repr->m[i],p));
200          syzstr->orderedRes[index]->m[i+start] = pCopy(current_repr->m[i]);
201          syzstr->orderedRes[index]->m[i+start] =
202            pMult_mm(syzstr->orderedRes[index]->m[i+start],w_gen);
203          if ((*syzstr->Tl)[index]!=0)
204            p_Shift(&syzstr->orderedRes[index]->m[i+start],(*syzstr->Tl)[index],currRing);
205        }
206      }
207      for (i=0;i<IDELEMS(totake[index-1]);i++)
208      {
209        if (totake[index-1]->m[i]!=NULL)
210        {
211          if ((index==1) && ((i==IDELEMS(current_ideal) ||
212               (totake[index-1]->m[i+1]==NULL)))) break;
213          totake[index]->m[i+start_ttk] =
214            pMult_mm(pCopy(totake[index-1]->m[i]),w_gen);
215          p_Shift(&totake[index]->m[i+start_ttk],current_tl,currRing);
216#ifdef FULL_TOTAKE
217          poly pp=pCopy(p);
218          p_Shift(&pp,i+1,currRing);
219          totake[index]->m[i+start_ttk] = pAdd(totake[index]->m[i+start_ttk],pp);
220#endif
221        }
222      }
223      (*syzstr->Tl)[index] += current_tl;
224    }
225    index--;
226  }
227  pDelete(&gen);
228  pDelete(&neg_gen);
229  pDelete(&w_gen);
230  //syShowRes(syzstr);
231}
232
233/*3
234* proves the consistence of the pairset resPairs with the corresponding
235* set of generators;
236* only for tests
237*/
238#ifdef SING_NDEBUG
239static void syTestPairs(SSet resPairs,int length,ideal /*old_generators*/)
240#else
241static void syTestPairs(SSet resPairs,int length,ideal old_generators)
242#endif
243{
244  int i=0;
245
246  while (i<length)
247  {
248    if (resPairs[i].lcm!=NULL)
249    {
250      if (resPairs[i].p1!=NULL)
251        assume(resPairs[i].p1==old_generators->m[resPairs[i].ind1]);
252      if (resPairs[i].p2!=NULL)
253        assume(resPairs[i].p2==old_generators->m[resPairs[i].ind2]);
254    }
255    i++;
256  }
257}
258
259/*3
260* cancels the weight monomials given by the leading terms of totake
261* from the resolution res;
262* works in place on res, but reads only from totake
263*/
264void syReorder_Kosz(syStrategy syzstr)
265{
266  int length=syzstr->length;
267  int syzIndex=length-1,i,j;
268  resolvente res=syzstr->fullres;
269  poly p;
270
271  while ((syzIndex!=0) && (res[syzIndex]==NULL)) syzIndex--;
272  while (syzIndex>0)
273  {
274    for(i=0;i<IDELEMS(res[syzIndex]);i++)
275    {
276#ifdef USE_REGULARITY
277      if ((syzstr->regularity>0) && (res[syzIndex]->m[i]!=NULL))
278      {
279        if (p_FDeg(res[syzIndex]->m[i],currRing)>=syzstr->regularity+syzIndex)
280          pDelete(&res[syzIndex]->m[i]);
281      }
282#endif
283      p = res[syzIndex]->m[i];
284      while (p!=NULL)
285      {
286        if (res[syzIndex-1]->m[pGetComp(p)-1]!=NULL)
287        {
288          for(j=1;j<=(currRing->N);j++)
289          {
290            pSetExp(p,j,pGetExp(p,j)
291                        -pGetExp(res[syzIndex-1]->m[pGetComp(p)-1],j));
292          }
293        }
294        else
295          PrintS("error in the resolvent\n");
296        pSetm(p);
297        pIter(p);
298      }
299    }
300    syzIndex--;
301  }
302}
303
304/*3
305* updates the pairset resPairs by generating all pairs including the
306* new_generators in the 0-th modul;
307* the new_generators are inserted in the old_generators;
308* new_generators is empty after the procedure;
309*/
310static void updatePairs(SSet *resPairs,int *l_pairs,syStrategy syzstr,
311       int index,ideal new_generators,ideal new_repr,int crit_comp)
312{
313  if (idIs0(new_generators)) return;
314  ideal old_generators=syzstr->res[index];
315  ideal old_repr=syzstr->orderedRes[index];
316  int i=0,j,k,kk,og_elem=0,og_idel=IDELEMS(old_generators),l=*l_pairs,jj,ll,j1;
317  int og_ini=0;
318  ideal pairs=idInit(og_idel+IDELEMS(new_generators),old_generators->rank);
319  polyset prs=pairs->m;
320  poly p=NULL;
321  SObject tso;
322
323  syInitializePair(&tso);
324  while ((og_elem<og_idel) && (old_generators->m[og_elem]!=NULL))
325  {
326    if ((index>0) && (pGetComp(old_generators->m[og_elem])<=crit_comp))
327      og_ini = og_elem;
328    og_elem++;
329  }
330  while ((l>0) && ((*resPairs)[l-1].lcm==NULL)) l--;
331  while ((i<IDELEMS(new_generators)) && (new_generators->m[i]!=NULL))
332  {
333    syTestPairs(*resPairs,*l_pairs,old_generators);
334    if (IDELEMS(old_generators)==og_elem)
335    {
336      pEnlargeSet(&old_generators->m,IDELEMS(old_generators),16);
337      IDELEMS(old_generators) += 16;
338      pEnlargeSet(&old_repr->m,IDELEMS(old_repr),16);
339      IDELEMS(old_repr) += 16;
340    }
341    k = p_FDeg(new_generators->m[i],currRing);
342    kk = pGetComp(new_generators->m[i]);
343    j = og_ini;
344    while ((j<og_elem) && (old_generators->m[j]!=NULL) &&
345           (pGetComp(old_generators->m[j])<kk)) j++;
346    while ((j<og_elem) && (old_generators->m[j]!=NULL) &&
347           (p_FDeg(old_generators->m[j],currRing)<=k)) j++;
348    for (jj=og_elem;jj>j;jj--)
349    {
350      old_generators->m[jj] = old_generators->m[jj-1];
351      old_repr->m[jj] = old_repr->m[jj-1];
352    }
353    old_generators->m[j] = new_generators->m[i];
354    new_generators->m[i] = NULL;
355    old_repr->m[j] = new_repr->m[i];
356    new_repr->m[i] = NULL;
357    og_elem++;
358    for (jj=0;jj<*l_pairs;jj++)
359    {
360      if ((*resPairs)[jj].lcm!=NULL)
361      {
362        if ((*resPairs)[jj].ind1>=j) (*resPairs)[jj].ind1++;
363        if ((*resPairs)[jj].ind2>=j) (*resPairs)[jj].ind2++;
364      }
365    }
366    syTestPairs(*resPairs,*l_pairs,old_generators);
367    for (jj=og_ini;jj<og_elem;jj++)
368    {
369      if ((j!=jj) && (pGetComp(old_generators->m[jj])==pGetComp(old_generators->m[j])))
370      {
371        p = pOne();
372        pLcm(old_generators->m[jj],old_generators->m[j],p);
373        pSetComp(p,j+1);
374        pSetm(p);
375        j1 = 0;
376        while (j1<jj)
377        {
378          if (prs[j1]!=NULL)
379          {
380            if (pLmDivisibleByNoComp(prs[j1],p))
381            {
382              pDelete(&p);
383              break;
384            }
385            else if (pLmDivisibleByNoComp(p,prs[j1]))
386            {
387              pDelete(&(prs[j1]));
388            }
389#ifdef USE_CHAINCRIT0
390            else
391            {
392              poly p1,p2;
393              int ip=(currRing->N);
394              p1 = pDivide(p,old_generators->m[jj]);
395              p2 = pDivide(prs[j1],old_generators->m[j1]);
396              while ((ip>0) && (pGetExp(p1,ip)*pGetExp(p2,ip)==0)) ip--;
397              if (ip==0)
398              {
399                int ti=0;
400                while ((ti<l) && (((*resPairs)[ti].ind1!=j1)|| ((*resPairs)[ti].ind2!=jj))) ti++;
401                if (ti<l)
402                {
403                  if (TEST_OPT_PROT) PrintS("cc");
404                  syDeletePair(&(*resPairs)[ti]);
405                  syCompactifyPairSet(*resPairs,*l_pairs,ti);
406                  l--;
407                }
408              }
409              pDelete(&p1);
410              pDelete(&p2);
411            }
412#endif
413          }
414          j1++;
415        }
416        if (p!=NULL)
417          prs[jj] = p;
418      }
419    }
420    for (jj=og_ini;jj<og_elem;jj++)
421    {
422      if (prs[jj] !=NULL)
423      {
424        if (l>=*l_pairs)
425        {
426          SSet temp = (SSet)omAlloc0((*l_pairs+16)*sizeof(SObject));
427          for (ll=0;ll<*l_pairs;ll++)
428          {
429            temp[ll].p = (*resPairs)[ll].p;
430            temp[ll].p1 = (*resPairs)[ll].p1;
431            temp[ll].p2 = (*resPairs)[ll].p2;
432            temp[ll].syz = (*resPairs)[ll].syz;
433            temp[ll].lcm = (*resPairs)[ll].lcm;
434            temp[ll].ind1 = (*resPairs)[ll].ind1;
435            temp[ll].ind2 = (*resPairs)[ll].ind2;
436            temp[ll].syzind = (*resPairs)[ll].syzind;
437            temp[ll].order = (*resPairs)[ll].order;
438            temp[ll].isNotMinimal = (*resPairs)[ll].isNotMinimal;
439          }
440          omFreeSize((ADDRESS)(*resPairs),*l_pairs*sizeof(SObject));
441          *l_pairs += 16;
442          (*resPairs) = temp;
443        }
444        tso.lcm = prs[jj];
445        prs[jj] = NULL;
446        tso.order = p_FDeg(tso.lcm,currRing);
447        tso.p1 = old_generators->m[jj];
448        tso.p2 = old_generators->m[j];
449        tso.ind1 = jj;
450        tso.ind2 = j;
451        tso.syzind = -1;
452        tso.isNotMinimal = NULL;
453        tso.p = NULL;
454        tso.syz = NULL;
455        SSet rP=*resPairs;
456#ifdef SHOW_PROT
457Print("erzeuge Paar im Modul %d,%d mit: \n",index,tso.order);
458PrintS("poly1: ");pWrite(tso.p1);
459PrintS("poly2: ");pWrite(tso.p2);
460PrintS("syz: ");pWrite(tso.syz);
461PrintS("sPoly: ");pWrite(tso.p);
462PrintLn();
463#endif
464        syEnterPair(rP,&tso,&l,index);
465        syInitializePair(&tso);
466      }
467    }
468    i++;
469  }
470  idDelete(&pairs);
471}
472
473/*3
474* performs the modification of a single reduction on the syzygy-level
475*/
476inline void sySPRedSyz_Kosz(syStrategy syzstr,poly redWith,poly syz,poly q=NULL,int l_syz=-1)
477{
478  poly p=pDivide(q,redWith);
479  pSetCoeff(p,nDiv(pGetCoeff(q),pGetCoeff(redWith)));
480  kBucket_Minus_m_Mult_p(syzstr->syz_bucket,p,syz,&l_syz,NULL);
481  pDelete(&p);
482}
483
484/*3
485* normalizes the poly bucket by the ideal;
486* stops the reduction whenever the leading component is less than the
487* crit_comp;
488* returns the changing status
489*/
490static BOOLEAN syRedSyz(kBucket_pt bucket,ideal red,int crit_comp,int* g_l)
491{
492  poly p = kBucketGetLm(bucket);
493  int j = 0,i=IDELEMS(red)-1;
494  number n;
495  BOOLEAN isChanged=FALSE;
496
497  loop
498  {
499    if ((j>=i) || (p==NULL) || (pGetComp(p)<=crit_comp)) break;
500    if ((red->m[j]!=NULL) && (pDivisibleBy(red->m[j],p)))
501    {
502      n = kBucketPolyRed(bucket,red->m[j], g_l[j], NULL);
503      nDelete(&n);
504      p = kBucketGetLm(bucket);
505      isChanged = TRUE;
506      j = 0;
507    }
508    else
509      j++;
510  }
511  return isChanged;
512}
513
514/*3
515* a tail reduction for the syzygies yielding new generators
516*/
517static poly syRedTailSyz(poly tored,ideal red,ideal sec_red,int crit_comp,syStrategy syzstr,
518            int * gen_length,int * secgen_length,int * tored_length)
519{
520  int i=IDELEMS(red)-1,num_mon,num_tail;
521  poly h,hn;
522  // BOOLEAN dummy;
523
524  while ((i>0) && (red->m[i-1]==NULL)) i--;
525  i--;
526  h = tored;
527  if ((h!=NULL) && (pGetComp(h)>crit_comp))
528  {
529    num_mon = 1;
530    hn = pNext(h);
531    num_tail = *tored_length-1;
532    while (hn!=NULL)
533    {
534      kBucketInit(syzstr->syz_bucket,hn,num_tail);
535      /*dummy =*/ (void) syRedSyz(syzstr->syz_bucket,red,crit_comp,gen_length);
536      kBucketClear(syzstr->syz_bucket,&hn,&num_tail);
537      pNext(h) = hn;
538      if ((hn==NULL) || (pGetComp(hn)<=crit_comp))
539        break;
540      else
541      {
542        pIter(h);
543        pIter(hn);
544        num_mon++;
545        num_tail--;
546      }
547    }
548    if (sec_red!=NULL)
549    {
550      while (hn!=NULL)
551      {
552        kBucketInit(syzstr->syz_bucket,hn,num_tail);
553        /*dummy =*/ (void) syRedSyz(syzstr->syz_bucket,sec_red,crit_comp,secgen_length);
554        kBucketClear(syzstr->syz_bucket,&hn,&num_tail);
555        pNext(h) = hn;
556        if (hn==NULL)
557          break;
558        else
559        {
560          pIter(h);
561          pIter(hn);
562          num_mon++;
563          num_tail--;
564        }
565      }
566    }
567    *tored_length = num_mon+num_tail;
568  }
569  assume(pLength(tored)==*tored_length);
570  return tored;
571}
572
573/*3
574* the complete reduction of a single pair which is just stored
575* in bucket and syz_bucket
576*/
577static BOOLEAN syRedSyzPair(syStrategy syzstr,int index,int* g_l,int* orp_l)
578{
579  kBucket_pt bucket=syzstr->bucket;
580  poly p = kBucketGetLm(bucket);
581  ideal red=syzstr->res[index],repr=syzstr->orderedRes[index];
582  int j = 0,i=IDELEMS(red)-1;
583  number n;
584  BOOLEAN isChanged=FALSE;
585
586  loop
587  {
588    if ((j>=i) || (p==NULL)) break;
589    if ((red->m[j]!=NULL) && (pDivisibleBy(red->m[j],p)))
590    {
591      sySPRedSyz_Kosz(syzstr,red->m[j],repr->m[j],p,orp_l[j]);
592      n = kBucketPolyRed(bucket,red->m[j], g_l[j], NULL);
593      nDelete(&n);
594      p = kBucketGetLm(bucket);
595      isChanged = TRUE;
596      j = 0;
597    }
598    else
599      j++;
600  }
601  return isChanged;
602}
603
604/*3
605* the tailreduction for generators (which includes the correction of
606* the corresponding representation)
607*/
608#if 0 /*unused*/
609static 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*/
650static 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 = pDivide(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 = pDivide(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
767Print("reduziere Paar im Module %d mit: \n",index);
768PrintS("poly1: ");pWrite(tso.p1);
769PrintS("poly2: ");pWrite(tso.p2);
770PrintS("syz: ");pWrite(tso.syz);
771PrintS("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
846Print("erhalte Paar im Module %d mit: \n",index);
847PrintS("syz: ");pWrite(tso.syz);
848PrintS("sPoly: ");pWrite(tso.p);
849PrintLn();
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      pMult_nn(tso.syz,n);
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*/
938static 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*/
1000static 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    Werror("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*/
1026static 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*/
1054static 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 = pDivide(p,add_generators->m[j]);
1096              p2 = pDivide(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
1158Print("erzeuge Paar im Modul %d,%d mit: \n",index,tso.order);
1159PrintS("poly1: ");pWrite(tso.p1);
1160PrintS("poly2: ");pWrite(tso.p2);
1161PrintS("syz: ");pWrite(tso.syz);
1162PrintS("sPoly: ");pWrite(tso.p);
1163PrintLn();
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
1179static 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
1184static 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 = pDivide(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 = pDivide(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
1282Print("reduziere Paar im Module %d mit: \n",index);
1283PrintS("poly1: ");pWrite(tso.p1);
1284PrintS("poly2: ");pWrite(tso.p2);
1285PrintS("syz: ");pWrite(tso.syz);
1286PrintS("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
1320Print("erhalte Paar im Module %d mit: \n",index);
1321PrintS("syz: ");pWrite(tso.syz);
1322PrintS("sPoly: ");pWrite(tso.p);
1323PrintLn();
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      pMult_nn(tso.syz,n);
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*/
1376static 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*/
1413static 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)<=(unsigned)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*/
1544static 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*/
1605static 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*/
1687static 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
1743PrintS("Add new generators:\n");
1744idPrint(new_generators);
1745PrintS("with representaions:\n");
1746idPrint(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*/
1763syStrategy 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    for (j=0;j<IDELEMS(temp);j++)
1789    {
1790      if (temp->m[j]!=NULL)
1791        p_Shift(&temp->m[j],1,currRing);
1792    }
1793  }
1794  idSkipZeroes(temp);
1795#ifdef WITH_SORT
1796  if (temp->m[0]!=NULL)
1797  {
1798    int md;
1799    int maxdeg=p_FDeg(temp->m[IDELEMS(temp)-1],currRing);
1800    ideal temp1=idInit(IDELEMS(temp),temp->rank);
1801    for (j=IDELEMS(temp)-2;j>=0;j--)
1802    {
1803      jj = p_FDeg(temp->m[j],currRing);
1804      if (jj>maxdeg) maxdeg = jj;
1805    }
1806    while (!idIs0(temp))
1807    {
1808      md = maxdeg;
1809      for (j=IDELEMS(temp)-1;j>=0;j--)
1810      {
1811        if (temp->m[j]!=NULL)
1812        {
1813          jj = p_FDeg(temp->m[j],currRing);
1814          if (jj<md) md = jj;
1815        }
1816      }
1817      for (j=0;j<IDELEMS(temp);j++)
1818      {
1819        if ((temp->m[j]!=NULL) && (p_FDeg(temp->m[j],currRing)==md))
1820        {
1821          temp1->m[k] = temp->m[j];
1822          temp->m[j] = NULL;
1823          k++;
1824        }
1825      }
1826    }
1827    idDelete(&temp);
1828    temp = temp1;
1829    temp1 = NULL;
1830  }
1831#endif
1832#ifdef USE_REGULARITY
1833  int last_generator=IDELEMS(temp)-1;
1834  while ((last_generator>=0) && (temp->m[last_generator]==NULL))
1835    last_generator--;
1836#endif
1837  syzstr->res = (resolvente)omAlloc0((*length+1)*sizeof(ideal));
1838  syzstr->orderedRes = (resolvente)omAlloc0((*length+1)*sizeof(ideal));
1839  resolvente totake=(resolvente)omAlloc0((*length+1)*sizeof(ideal));
1840  syzstr->Tl = new intvec(*length+1);
1841  syzstr->bucket = kBucketCreate(currRing);
1842  syzstr->syz_bucket = kBucketCreate(currRing);
1843  ideal new_generators=idInit(1,si_max(rk_arg,1));
1844  ideal temp_gens,old_std;
1845  syzstr->res[0] = idInit(1,1);
1846  if (rk_arg>1) syzstr->res[0]->rank = rk_arg;
1847  syzstr->orderedRes[0] = idInit(1,1);
1848/*--- computes the resolution ----------------------*/
1849  i = 0;
1850  while (i<IDELEMS(temp))
1851  {
1852    if (temp->m[i]!=NULL)
1853    {
1854      new_generators->m[0] = kNF(syzstr->res[0],currRing->qideal,temp->m[i]);
1855      if (!nIsOne(pGetCoeff(new_generators->m[0])))
1856        pNorm(new_generators->m[0]);
1857      next_deg = p_FDeg(new_generators->m[0],currRing);
1858      next_gen = pCopy(new_generators->m[0]);
1859    }
1860    if (!idIs0(new_generators))
1861    {
1862      index = 0;
1863      while (index<=*length)
1864      {
1865        if (index==0)
1866        {
1867          old_ideal = idCopy(syzstr->res[0]);
1868          old_repr = idCopy(syzstr->orderedRes[0]);
1869          old_tl = (*syzstr->Tl)[0];
1870          old_std = id_Head(syzstr->res[0],currRing);
1871        }
1872        t_comp = (*syzstr->Tl)[index];
1873        if (index==0) crit_comp = t_comp;
1874        temp_gens = syAppendSyz(new_generators,syzstr, index,crit_comp,totake);
1875        crit_comp = t_comp;
1876        if (index==0)
1877        {
1878          isRegular = syIsRegular(old_std,syzstr->res[0],next_deg);
1879#ifndef ONLY_STD
1880          if (isRegular)
1881            syCreateRegularExtension(syzstr,old_ideal,old_repr,old_tl,next_gen,
1882                                     totake);
1883#ifdef USE_REGULARITY
1884        if ((index==0) && (!isRegular) && (i==last_generator))
1885        {
1886/*----------- we are computing the regularity -----------------------*/
1887          ideal initial=id_Head(syzstr->res[0],currRing);
1888          int len=0,reg=0;
1889          intvec *w=NULL;
1890          ring dp_C_ring = rAssure_dp_C(currRing); rChangeCurrRing(dp_C_ring); 
1891          initial = idrMoveR_NoSort(initial, syzstr->syRing, dp_C_ring);
1892          resolvente res = sySchreyerResolvente(initial,-1,&len,TRUE, TRUE);
1893          intvec * dummy = syBetti(res,len,&reg, w);
1894          syzstr->regularity = reg+2;
1895          delete dummy;
1896          delete w;
1897          for (j=0;j<len;j++)
1898          {
1899            if (res[j]!=NULL) idDelete(&(res[j]));
1900          }
1901          omFreeSize((ADDRESS)res,len*sizeof(ideal));
1902          idDelete(&initial);
1903          rChangeCurrRing(syzstr->syRing);
1904          rDelete(dp_C_ring);
1905        }
1906#endif
1907#endif
1908          idDelete(&old_ideal);
1909          idDelete(&old_repr);
1910          idDelete(&old_std);
1911          if (TEST_OPT_PROT)
1912          {
1913            if (isRegular)
1914              PrintS("\n regular\n");
1915            else
1916              PrintS("\n not regular\n");
1917          }
1918          if (next_gen!=NULL)
1919            pDelete(&next_gen);
1920          if (isRegular)
1921          {
1922            idDelete(&temp_gens);
1923            break;
1924          }
1925        }
1926        idDelete(&new_generators);
1927        new_generators = temp_gens;
1928#ifdef ONLY_STD
1929        break;
1930#endif
1931        if (idIs0(new_generators)) break;
1932        index++;
1933      }
1934      if (!idIs0(new_generators))
1935      {
1936        for (j=0;j<IDELEMS(new_generators);j++)
1937        {
1938          if (new_generators->m[j]!=NULL)
1939          {
1940            pDelete(&new_generators->m[j]);
1941            new_generators->m[j] = NULL;
1942          }
1943        }
1944      }
1945    }
1946    i++;
1947  }
1948  if (idIs0(new_generators) && new_generators!=NULL) idDelete(&new_generators);
1949  if (temp!=NULL) idDelete(&temp);
1950  kBucketDestroy(&(syzstr->bucket));
1951  kBucketDestroy(&(syzstr->syz_bucket));
1952  index = 0;
1953  syzstr->fullres = syzstr->res;
1954  syzstr->res = NULL;
1955  index = 0;
1956  while ((index<=*length) && (syzstr->fullres[index]!=NULL))
1957  {
1958#ifdef SHOW_RESULT
1959    Print("The %d-th syzygy-module is now:\n",index);
1960    ideal ttt=id_Head(syzstr->fullres[index],currRing);
1961    idShow(ttt);
1962    idDelete(&ttt);
1963    //if (index>0)
1964    //{
1965      //Print("The related module is: \n");
1966      //idPrint(totake[index-1]);
1967    //}
1968    //Print("The %d-th module of the minimal resolution is:\n",index);
1969    if (!idIs0(totake[index]))
1970      idShow(totake[index]);
1971    //Print("with standard basis:\n");
1972    //idPrint(syzstr->fullres[index]);
1973    //if ((index<*length) && (totake[index+1]!=NULL))
1974    //{
1975      //Print("The %d-th syzygy-module is now:\n",index+1);
1976      //idPrint(totake[index+1]);
1977      //matrix m1=idModule2Matrix(totake[index]);
1978      //matrix m2=idModule2Matrix(totake[index+1]);
1979      //matrix m3=mpMult(m1,m2);
1980      //idPrint((ideal)m3);
1981    //}
1982#endif
1983    if (!idIs0(totake[index]))
1984    {
1985      for(i=0;i<IDELEMS(totake[index]);i++)
1986      {
1987        if (totake[index]->m[i]!=NULL)
1988        {
1989          j=0;
1990          while ((j<IDELEMS(syzstr->fullres[index])) &&
1991            ((syzstr->fullres[index]->m[j]==NULL) ||
1992            (!pLmEqual(syzstr->fullres[index]->m[j],totake[index]->m[i])))) j++;
1993          if (j<IDELEMS(syzstr->fullres[index]))
1994          {
1995            pDelete(&totake[index]->m[i]);
1996            totake[index]->m[i] = syzstr->fullres[index]->m[j];
1997            syzstr->fullres[index]->m[j] = NULL;
1998          }
1999          else
2000          {
2001            PrintS("Da ist was faul!!!\n");
2002            Print("Aber: Regularitaet %d, Grad %ld\n",
2003                   syzstr->regularity,p_FDeg(totake[index]->m[i],currRing));
2004          }
2005        }
2006      }
2007      idDelete(&syzstr->fullres[index]);
2008      syzstr->fullres[index] = totake[index];
2009    }
2010#ifdef SHOW_RESULT
2011    idShow(syzstr->fullres[index]);
2012#endif
2013    index++;
2014  }
2015  syReorder_Kosz(syzstr);
2016  index = 0;
2017  while ((index<=*length) && (syzstr->orderedRes[index]!=NULL))
2018  {
2019    idDelete(&(syzstr->orderedRes[index]));
2020    index++;
2021  }
2022  if (origR!=syzstr->syRing)
2023  {
2024    rChangeCurrRing(origR);
2025    index = 0;
2026    while ((index<=*length) && (syzstr->fullres[index]!=NULL))
2027    {
2028      syzstr->fullres[index] = idrMoveR(syzstr->fullres[index],syzstr->syRing, origR);
2029      index++;
2030    }
2031  }
2032  delete syzstr->Tl;
2033  syzstr->Tl = NULL;
2034  rDelete(syzstr->syRing);
2035  syzstr->syRing = NULL;
2036  omFreeSize((ADDRESS)totake,(*length+1)*sizeof(ideal));
2037  omFreeSize((ADDRESS)syzstr->orderedRes,(*length+1)*sizeof(ideal));
2038//Print("Pairs to discard: %d\n",discard_pairs);
2039//Print("Pairs shorter reduced: %d\n",short_pairs);
2040//discard_pairs = 0;
2041//short_pairs = 0;
2042  return syzstr;
2043}
2044
Note: See TracBrowser for help on using the repository browser.