source: git/kernel/syz3.cc @ 7b25fe

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