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

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