source: git/kernel/syz3.cc @ 3a67ea7

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