source: git/Singular/syz3.cc @ 50cbdc

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