source: git/Singular/syz1.cc @ 0f5510

spielwiese
Last change on this file since 0f5510 was 0f5510, checked in by Thomas Siebert <siebert@…>, 26 years ago
*** empty log message *** git-svn-id: file:///usr/local/Singular/svn/trunk@1534 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 65.7 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id: syz1.cc,v 1.26 1998-04-29 07:20:43 siebert Exp $ */
5/*
6* ABSTRACT: resolutions
7*/
8
9
10#include "mod2.h"
11#include "tok.h"
12#include "mmemory.h"
13#include "polys.h"
14#include "febase.h"
15#include "kstd1.h"
16#include "kutil.h"
17#include "spolys.h"
18#include "spolys0.h"
19#include "stairc.h"
20#include "ipid.h"
21#include "cntrlc.h"
22#include "ipid.h"
23#include "intvec.h"
24#include "ipshell.h"
25#include "limits.h"
26#include "numbers.h"
27#include "modulop.h"
28#include "ideals.h"
29#include "intvec.h"
30#include "ring.h"
31#include "lists.h"
32#include "syz.h"
33
34
35/*--------------static variables------------------------*/
36/*---contains the real components wrt. cdp--------------*/
37//static int ** truecomponents=NULL;
38//static int ** backcomponents=NULL;
39//static int ** Howmuch=NULL;
40//static int ** Firstelem=NULL;
41/*---points to the real components of the actual module-*/
42static int *  currcomponents=NULL;
43/*---head-term-polynomials for the reduction------------*/
44static poly redpol=NULL;
45/*---counts number of applications of GM-criteria-------*/
46//static int crit;
47//static int zeroRed;
48//static int dsim;
49//static int simple;
50static int euler;
51/*---controls the ordering------------------------------*/
52static intvec * orderedcomp;
53static int *binomials;
54static int highdeg;
55static int highdeg_1;
56
57// uncomment to play around with new spolys for syzygies
58// #include "syzvector.cc"
59
60/*3
61* deletes all entres of a pair
62*/
63static void syDeletePair(SObject * so)
64{
65  pDelete(&(*so).p);
66  pDelete(&(*so).lcm);
67  pDelete(&(*so).syz);
68  (*so).p1 = NULL;
69  (*so).p2 = NULL;
70  (*so).ind1 = 0;
71  (*so).ind2 = 0;
72  (*so).syzind = -1;
73  (*so).order = 0;
74  (*so).isNotMinimal = NULL;
75}
76
77/*3
78* puts all entres of a pair to another
79*/
80static void syCopyPair(SObject * argso, SObject * imso)
81{
82  *imso=*argso;
83  //(*imso).p = (*argso).p;
84  //(*imso).p1 = (*argso).p1;
85  //(*imso).p2 = (*argso).p2;
86  //(*imso).lcm = (*argso).lcm;
87  //(*imso).syz = (*argso).syz;
88  //(*imso).ind1 = (*argso).ind1;
89  //(*imso).ind2 = (*argso).ind2;
90  //(*imso).syzind = (*argso).syzind;
91  //(*imso).order = (*argso).order;
92  //(*imso).isNotMinimal = (*argso).isNotMinimal;
93  (*argso).p = NULL;
94  (*argso).p1 = NULL;
95  (*argso).p2 = NULL;
96  (*argso).lcm = NULL;
97  (*argso).syz = NULL;
98  (*argso).ind1 = 0;
99  (*argso).ind2 = 0;
100  (*argso).syzind = -1;
101  (*argso).order = 0;
102  (*argso).isNotMinimal = NULL;
103}
104
105/*3
106* deletes empty objects from a pair set beginning with
107* pair first
108* assumes a pair to be empty if .lcm does so
109*/
110static void syCompactifyPairSet(SSet sPairs, int sPlength, int first)
111{
112  int k=first,kk=0;
113
114  while (k+kk<sPlength)
115  {
116    if (sPairs[k+kk].lcm!=NULL)
117    {
118      if (kk>0) syCopyPair(&sPairs[k+kk],&sPairs[k]);
119      k++;
120    }
121    else
122    {
123      kk++;
124    }
125  }
126}
127
128/*3
129* deletes empty objects from a pair set beginning with
130* pair first
131* assumes a pair to be empty if .lcm does so
132*/
133static void syCompactify1(SSet sPairs, int* sPlength, int first)
134{
135  int k=first,kk=0;
136
137  while (k+kk<=*sPlength)
138  {
139    if (sPairs[k+kk].lcm!=NULL)
140    {
141      if (kk>0) syCopyPair(&sPairs[k+kk],&sPairs[k]);
142      k++;
143    }
144    else
145    {
146      kk++;
147    }
148  }
149  *sPlength -= kk;
150}
151
152/*3
153* replaces comp1dpc during homogeneous syzygy-computations
154* compares with components of currcomponents instead of the
155* exp[0]
156*/
157static int syzcomp1dpc(poly p1, poly p2)
158{
159  /*4 compare monomials by order then revlex*/
160    int i = pVariables;
161    if ((pGetExp(p1,i) == pGetExp(p2,i)))
162    {
163      do
164      {
165        i--;
166        if (i <= 1)
167        {
168          //(*orderingdepth)[pVariables-i]++;
169           /*4 handle module case:*/
170           if (pGetComp(p1)==pGetComp(p2)) return 0;
171           else if 
172              (currcomponents[pGetComp(p1)]>currcomponents[pGetComp(p2)]) 
173                return 1;
174           else return -1;
175        }
176      } while ((pGetExp(p1,i) == pGetExp(p2,i)));
177    }
178    //(*orderingdepth)[pVariables-i]++;
179    if (pGetExp(p1,i) < pGetExp(p2,i)) return 1;
180    return -1;
181}
182
183/*3
184* replaces comp1dpc during homogeneous syzygy-computations
185* compares with components of currcomponents instead of the
186* exp[0]
187*/
188
189#ifdef SY_VEC_DEBUG
190static int normal_syzcomp2dpc(poly p1, poly p2);
191
192static int syzcomp2dpc(poly p1, poly p2)
193{
194  int n1 = normal_syzcomp2dpc(p1, p2);
195  int n2 = SyVecSyzCmp(p1, p2);
196 
197  if (n1 != n2)
198  {
199    PrintS("Error in syzcomp\n");
200    SyVecSyzCmp(p1, p2);
201    normal_syzcomp2dpc(p1, p2);
202  }
203  return n1;
204}
205
206static int normal_syzcomp2dpc(poly p1, poly p2)
207#else
208#ifdef HAVE_SY_VECTOR
209static int dummy(poly p1, poly p2)
210#else
211static int syzcomp2dpc(poly p1, poly p2)
212#endif
213#endif
214{
215  int o1=pGetOrder(p1), o2=pGetOrder(p2);
216  if (o1 > o2) return 1;
217  if (o1 < o2) return -1;
218
219  if (o1>0)
220  {
221    int i = pVariables;
222    while ((i>1) && (pGetExp(p1,i)==pGetExp(p2,i)))
223      i--;
224    //(*orderingdepth)[pVariables-i]++;
225    if (i>1)
226    {
227      if (pGetExp(p1,i) < pGetExp(p2,i)) return 1;
228      return -1;
229    }
230  }
231  o1=pGetComp(p1);
232  o2=pGetComp(p2);
233  if (o1==o2/*pGetComp(p1)==pGetComp(p2)*/) return 0;
234  if (currcomponents[o1]>currcomponents[o2]) return 1;
235  return -1;
236}
237
238/*3
239* compares only the monomial without component
240*/
241static int syzcompmonomdp(poly p1, poly p2)
242{
243  int i;
244
245  /*4 compare monomials by order then revlex*/
246  if (pGetOrder(p1) == pGetOrder(p2))
247  {
248    i = pVariables;
249    if ((pGetExp(p1,i) == pGetExp(p2,i)))
250    {
251      do
252      {
253        i--;
254        if (i <= 1)
255          return 0;
256      } while ((pGetExp(p1,i) == pGetExp(p2,i)));
257    }
258    if (pGetExp(p1,i) < pGetExp(p2,i)) return 1;
259    return -1;
260  }
261  else if (pGetOrder(p1) > pGetOrder(p2)) return 1;
262  return -1;
263}
264
265/*
266* (i,j) in binomials[j-1,i-j+1]
267* table size is: pVariables * (highdeg+1)
268* highdeg_1==highdeg+1
269*/
270static void syBinomSet()
271{
272  highdeg_1=highdeg+1;
273  for(int j=1;j<=highdeg;j++)
274  {
275    binomials[j/*0,j*/] = j;
276    for (int i=1;i<pVariables;i++)
277    {
278      binomials[i*(highdeg_1)+j/*i,j*/]
279      = binomials[(i-1)*(highdeg_1)+j/*i-1,j*/]*(j+i)/(i+1);
280    }
281  }
282  for (int i=0;i<pVariables;i++)
283  { 
284    binomials[i*(highdeg_1)/*i,0*/]=0;
285  }
286}
287
288static inline int syBinom(int i,int j)
289{
290  return binomials[(j-1)*(highdeg_1)+i-j+1/*j-1,i-j*/];
291}
292
293#if defined(HAVE_SY_VECTOR) && ! defined(SY_VEC_DEBUG)
294#define syzSetm(p) pSetm(p)
295#else
296static void syzSetm(poly p)
297{
298  int i = 0,j;
299  for (j=pVariables;j>0;j--)
300    i += pGetExp(p,j);
301
302  if (i<=highdeg)
303  {
304    i=1;
305    j = -INT_MAX;
306    int *ip=binomials+pGetExp(p,1);
307    loop
308    {
309      j += (*ip);
310      if (i==pVariables) break;
311      i++;
312      ip+=highdeg_1+pGetExp(p,i);
313    }
314    pGetOrder(p) = j;
315  }
316  else
317    pGetOrder(p)=i;
318}
319#endif
320
321static inline poly syMultT(poly p,poly m)
322{
323  poly q,result=q=pNew();
324  int j;
325 
326  if (pGetOrder(p)>0)
327  {
328    loop
329    {
330      spMemcpy(q,p);
331      for (j=pVariables;j>0;j--)
332        pAddExp(q,j, pGetExp(m,j));
333      pSetCoeff0(q,nCopy(pGetCoeff(p)));
334      syzSetm(q);
335      pIter(p);
336      if (p!=NULL)
337      {
338        pNext(q) = pNew();
339        pIter(q);
340      }
341      else break;
342    }
343  }
344  else
345  {
346    poly lastmon=NULL;
347    int i=0;
348    loop
349    {
350      if (pGetOrder(p)!=i)
351      {
352        spMemcpy(q,p);
353        for (j=pVariables;j>0;j--)
354          pAddExp(q,j,  pGetExp(m,j));
355        syzSetm(q);
356        lastmon = q;
357        i = pGetOrder(p);
358      }
359      else
360      {
361        spMemcpy(q,lastmon);
362        pSetComp(q,pGetComp(p));
363      }
364      pSetCoeff0(q,nCopy(pGetCoeff(p)));
365      pIter(p);
366      if (p!=NULL)
367      {
368        pNext(q) = pNew();
369        pIter(q);
370      }
371      else break;
372    }
373  }
374  pNext(q) = NULL;
375  return result;
376}
377
378static inline poly syMultTNeg(poly p,poly m)
379{
380  poly q,result=q=pNew();
381  int j;
382 
383  if (pGetOrder(p)>0)
384  {
385    loop
386    {
387      spMemcpy(q,p);
388      for (j=pVariables;j>0;j--)
389        pAddExp(q,j, pGetExp(m,j));
390      pSetCoeff0(q,nCopy(pGetCoeff(p)));
391      pGetCoeff(q) = nNeg(pGetCoeff(q));
392      syzSetm(q);
393      pIter(p);
394      if (p!=NULL)
395      {
396        pNext(q) = pNew();
397        pIter(q);
398      }
399      else break;
400    }
401  }
402  else
403  {
404    poly lastmon=NULL;
405    int i=0;
406    loop
407    {
408      if (pGetOrder(p)!=i)
409      {
410        spMemcpy(q,p);
411        for (j=pVariables;j>0;j--)
412          pAddExp(q,j, pGetExp(m,j));
413        syzSetm(q);
414        lastmon = q;
415        i = pGetOrder(p);
416      }
417      else
418      {
419        spMemcpy(q,lastmon);
420        pSetComp(q,pGetComp(p));
421      }
422      pSetCoeff0(q,nCopy(pGetCoeff(p)));
423      pGetCoeff(q) = nNeg(pGetCoeff(q));
424      pIter(p);
425      if (p!=NULL)
426      {
427        pNext(q) = pNew();
428        pIter(q);
429      }
430      else break;
431    }
432  }
433  pNext(q) = NULL;
434  return result;
435}
436
437static poly syMultT1(poly p,poly m)
438{
439  poly q,result=q=pNew();
440  int j;
441 
442  if (pGetOrder(p)>0)
443  {
444    loop
445    {
446      spMemcpy(q,p);
447      for (j=pVariables;j>0;j--)
448        pAddExp(q,j, pGetExp(m,j));
449      pSetCoeff0(q,nMult(pGetCoeff(p),pGetCoeff(m)));
450      syzSetm(q);
451      pIter(p);
452      if (p!=NULL)
453      {
454        pNext(q) = pNew();
455        pIter(q);
456      }
457      else break;
458    }
459  }
460  else
461  {
462    poly lastmon=NULL;
463    int i=0;
464    loop
465    {
466      if (pGetOrder(p)!=i)
467      {
468        spMemcpy(q,p);
469        for (j=pVariables;j>0;j--)
470          pAddExp(q,j, pGetExp(m,j));
471        syzSetm(q);
472        lastmon = q;
473        i = pGetOrder(p);
474      }
475      else
476      {
477        spMemcpy(q,lastmon);
478        pSetComp(q,pGetComp(p));
479      }
480      pSetCoeff0(q,nMult(pGetCoeff(p),pGetCoeff(m)));
481      pIter(p);
482      if (p!=NULL)
483      {
484        pNext(q) = pNew();
485        pIter(q);
486      }
487      else break;
488    }
489  }
490  pNext(q) = NULL;
491  return result;
492}
493
494static inline poly syAdd(poly m1,poly m2)
495{
496  if (m1==NULL)
497    return m2;
498  if (m2==NULL)
499    return m1;
500  if ((pGetOrder(m1)>0) && (pGetOrder(m2)>0))
501  {
502    return pAdd(m1,m2);
503  }
504  else
505  {
506    poly p,result=p=pNew();
507    number n;
508    loop
509    {
510      if (pGetOrder(m1)>pGetOrder(m2))
511      {
512        pNext(p) = m1;
513        pIter(p);
514        pIter(m1);
515      }
516      else if (pGetOrder(m1)<pGetOrder(m2))
517      {
518        pNext(p) = m2;
519        pIter(p);
520        pIter(m2);
521      }
522      else if (currcomponents[pGetComp(m1)]>currcomponents[pGetComp(m2)])
523      {
524        pNext(p) = m1;
525        pIter(p);
526        pIter(m1);
527      }
528      else if (currcomponents[pGetComp(m1)]<currcomponents[pGetComp(m2)])
529      {
530        pNext(p) = m2;
531        pIter(p);
532        pIter(m2);
533      }
534      else
535      {
536        n = nAdd(pGetCoeff(m1),pGetCoeff(m2));
537        if (nIsZero(n))
538        {
539          pDelete1(&m1);
540        }
541        else
542        {
543          pNext(p) = m1;
544          pIter(p);
545          pIter(m1);
546          pSetCoeff(p,n);
547        }
548        pDelete1(&m2);
549      }
550      if (m1==NULL)
551      {
552        pNext(p) = m2;
553        break;
554      }
555      if (m2==NULL)
556      {
557        pNext(p) = m1;
558        break;
559      }
560    }
561    pDelete1(&result);
562    return result;
563  }
564}
565
566poly syOrdPolySchreyer(poly p)
567{
568  return pOrdPolyMerge(p);
569}
570
571#ifdef SY_VEC_DEBUG
572static poly normal_sySPAdd(poly m1,poly m2,poly m);
573static poly sySPAdd(poly m1, poly m2, poly m)
574{
575  poly temp1 = pCopy(m1);
576 
577  poly p1 = normal_sySPAdd(m1, m2, m);
578  poly p2 = syVecSpoly(m2, temp1, m);
579 
580  pTest(p1);
581  pTest(p2);
582  if (p1 != NULL) 
583  {
584    if (p2 == NULL || ! pEqual(p1, p2))
585    {
586      PrintS("Error in SySpoly\n");
587//      poly temp5 = pCopy(temp3);
588//      poly temp6 = pCopy(temp4);
589//      p1 = normal_sySPAdd(temp3, temp4, m);
590//      p2 = syVecSpoly(temp6, temp5, m);
591    }
592  }
593  else if (p2 != NULL) 
594    PrintS("Error in SySpoly\n");
595
596  return p1;
597}
598
599static poly normal_sySPAdd(poly m1,poly m2,poly m)
600#else
601#ifdef HAVE_SY_VECTOR
602#define sySPAdd(p1, p2, m) syVecSpoly(p2, p1, m)
603static poly dummy(poly m1,poly m2,poly m)
604#else
605static poly sySPAdd(poly m1,poly m2,poly m)
606#endif
607#endif
608{
609  int i2=-INT_MAX;
610  int i1;
611  int j;
612  //register poly m1=m11;
613  poly p=redpol;
614  poly tm2=pNew();
615  number multwith=pGetCoeff(m);
616
617  pDelete1(&m1);
618  i1 = pGetOrder(m1);
619  pIter(m2);
620  spMemcpy(tm2,m2);
621  j = 1;
622  pAddExp(tm2,1, pGetExp(m,1));
623  int *ip=binomials+pGetExp(tm2,1);
624  loop
625  {
626    i2 += (*ip);
627    if (j==pVariables) break;
628    j++;
629    pAddExp(tm2,j,pGetExp(m,j));
630    ip+=highdeg_1+pGetExp(tm2,j);
631  }
632  pGetOrder(tm2) = i2;
633  pSetCoeff0(tm2,nMult(pGetCoeff(m2),multwith));
634  loop
635  {
636    j = i1-i2;
637    if (j>0/*m1->order>tm2->order*/)
638    {
639      p = pNext(p) = m1;
640      pIter(m1);
641      if (m1!=NULL)
642        i1 = pGetOrder(m1);
643      else
644        break;
645    }
646    else if (j<0/*m1->order<tm2->order*/)
647    {
648      p = pNext(p) = tm2;
649      j = pGetOrder(m2);
650      pIter(m2);
651      if (m2!=NULL)
652      {
653        j -=pGetOrder(m2);
654        tm2 = pNew();
655        if (j!=0)
656        {
657          spMemcpy(tm2,m2);
658          j = 1;
659          pAddExp(tm2,1, pGetExp(m,1));
660          int *ip=binomials+pGetExp(tm2,1);
661          i2=-INT_MAX;
662          loop
663          {
664            i2 += (*ip);
665            if (j==pVariables) break;
666            j++;
667            pAddExp(tm2,j, pGetExp(m,j));
668            ip+=highdeg_1+pGetExp(tm2,j);
669          }
670          pGetOrder(tm2) = i2;
671          //simple++;
672        }
673        else
674        {
675          spMemcpy(tm2,p);
676          pSetComp(tm2,pGetComp(m2));
677          //dsim++;
678        }
679        //pNext(tm2) = NULL;
680        pSetCoeff0(tm2,nMult(pGetCoeff(m2),multwith));
681      }
682      else
683        break;
684    }
685    else 
686    {
687      j = currcomponents[pGetComp(m1)]-currcomponents[pGetComp(m2)];
688      if (j>0/*currcomponents[pGetComp(m1)]>currcomponents[pGetComp(m2)]*/)
689      {
690        p = pNext(p) = m1;
691        pIter(m1);
692        if (m1!=NULL)
693          i1 = pGetOrder(m1);
694        else
695          break;
696      }
697      else 
698      {
699        if (j<0/*currcomponents[pGetComp(m1)]<currcomponents[pGetComp(m2)]*/)
700        {
701          p = pNext(p) = tm2;
702          j = pGetOrder(m2);
703        }
704        else
705        {
706          number n = nAdd(pGetCoeff(m1),pGetCoeff(tm2));
707          if (nIsZero(n))
708          {
709            pDelete1(&tm2);
710            nDelete(&n);
711            j = 0;
712          }
713          else
714          {
715            p = pNext(p) = tm2;
716            pSetCoeff(p,n);
717            j = pGetOrder(m2);
718          }
719          pDelete1(&m1);
720          if (m1==NULL)
721          {
722            pIter(m2);
723            break;
724          }
725          else
726            i1 = pGetOrder(m1);
727        }
728        pIter(m2);
729        if (m2!=NULL)
730        {
731          j -=pGetOrder(m2);
732          tm2 = pNew();
733          if (j!=0)
734          {
735            spMemcpy(tm2,m2);
736            j = 1;
737            pAddExp(tm2,1, pGetExp(m,1));
738            int *ip=binomials+pGetExp(tm2,1);
739            i2=-INT_MAX;
740            loop
741            {
742              i2 += (*ip);
743              if (j==pVariables) break;
744              j++;
745              pAddExp(tm2,j, pGetExp(m,j));
746              ip+=highdeg_1+pGetExp(tm2,j);
747            }
748            pGetOrder(tm2) = i2;
749            //simple++;
750          }
751          else
752          {
753            spMemcpy(tm2,p);
754            pSetComp(tm2,pGetComp(m2));
755            //dsim++;
756          }
757          //pNext(tm2) = NULL;
758          pSetCoeff0(tm2,nMult(pGetCoeff(m2),multwith));
759        }
760        else
761          break;
762      }
763    }
764  }
765  if (m2==NULL)
766  {
767    pNext(p) = m1;
768  }
769  else
770  {
771    pNext(p) = syMultT1(m2,m);
772  }
773  return pNext(redpol);
774}
775
776static inline poly sySPoly(poly m1,poly m2,poly lcm)
777{
778  poly a=pDivide(lcm,m1);
779  poly b1=NULL,b2=NULL;
780  if (pNext(m1)!=NULL)
781    b1 = syMultT(pNext(m1),a);
782  pDelete(&a);
783  a=pDivide(lcm,m2);
784  if (pNext(m2)!=NULL)
785    b2 = syMultTNeg(pNext(m2),a);
786  pDelete(&a);
787  return syAdd(b1,b2);
788}
789
790static poly sySPolyRed(poly m1,poly m2)
791{
792  if (pNext(m2)==NULL)
793  {
794    pDelete1(&m1);
795    return m1;
796  }
797  poly a=pDivide(m1,m2);
798  pSetCoeff0(a,nCopy(pGetCoeff(m1)));
799  pGetCoeff(a) = nNeg(pGetCoeff(a));
800  //return spSpolyRed(m2,m1,NULL);
801  if (pNext(m2)==NULL)
802  {
803    pDelete1(&m1);
804    return m1;
805  }
806  if (pNext(m1)==NULL)
807    return syMultT1(pNext(m2),a);
808  poly res;
809#if defined(HAVE_SY_VECTOR) && ! defined (SY_VEC_DEBUG)
810    res = sySpolyProc(m2, m1,a);
811#else   
812  if (pGetOrder(m1)>0)
813  {
814    // TBC: initialize spSpolyLoop where ordering is set!!!
815    res = spSpolyRed(m2,m1,NULL, NULL);
816  }
817  else
818  {
819    res = sySPAdd(m1,m2,a);
820  }
821#endif
822  pDelete1(&a);
823  return res;
824}
825
826static int syLength(poly p)
827{
828  int result=0;
829
830  while (p!=NULL)
831  {
832    result++;
833    pIter(p);
834  }
835  return result;
836}
837
838poly syRedtail (poly p, syStrategy syzstr, int index)
839{
840  poly h, hn;
841  int j,pos;
842  ideal redWith=syzstr->orderedRes[index];
843
844  h = p;
845  hn = pNext(h);
846  while(hn != NULL)
847  {
848    j = syzstr->Firstelem[index-1][pGetComp(hn)]-1;
849    if (j>=0)
850    {
851      pos = j+syzstr->Howmuch[index-1][pGetComp(hn)];
852      while (j < pos)
853      {
854        if (pDivisibleBy2(redWith->m[j], hn))
855        {
856           //int syL=syLength(redWith->m[j]);
857            //Print("r");
858//for(int jj=j+1;jj<pos;jj++)
859//{
860  //if (syDivisibleBy1(redWith->m[jj],hn))
861  //{
862    //int syL1=syLength(redWith->m[jj]);
863    //if (syL1<syL)
864    //{
865      //j = jj;
866      //syL = syL1;
867      //Print("take poly %d with %d monomials\n",j,syL);
868      //Print("have poly %d with %d monomials\n",jj,syL1);
869    //}
870  //}
871//}
872          //if (pGetComp(redWith->m[j])!=pGetComp(hn))
873          //{
874            //Print("Hilfe!!!!!!!\n");
875            //Print("Fehler in Modul %d bei Elem %d\n",index,j);
876            //Print("Poly p: ");pWrite(hn);
877            //Print("Poly redWith: ");pWrite(redWith->m[j]);
878          //}
879          hn = sySPolyRed(hn,redWith->m[j]);
880          if (hn == NULL)
881          {
882            pNext(h) = NULL;
883            return p;
884          }
885          j = syzstr->Firstelem[index-1][pGetComp(hn)]-1;
886          pos = j+syzstr->Howmuch[index-1][pGetComp(hn)];
887        }
888        else
889        {
890          j++;
891        }
892      }
893    }
894    h = pNext(h) = hn;
895    hn = pNext(h);
896  }
897  return p;
898}
899
900/*3
901* local procedure for of syInitRes for the module case
902*/
903static int syChMin(intvec * iv)
904{
905  int i,j=-1,r=-1;
906
907  for (i=iv->length()-1;i>=0;i--)
908  {
909    if ((*iv)[i]>=0)
910    {
911      if ((j<0) || ((*iv)[i]<j))
912      {
913        j = (*iv)[i];
914        r = i;
915      }
916    }
917  }
918  return r;
919}
920
921/*3
922* initialize the resolution and puts in the argument as
923* zeroth entre, length must be > 0
924* assumes that the basering is degree-compatible
925*/
926static SRes syInitRes(ideal arg,int * length, intvec * Tl, intvec * cw=NULL)
927{
928  if (idIs0(arg)) return NULL;
929  SRes resPairs = (SRes)Alloc0(*length*sizeof(SSet));
930  resPairs[0] = (SSet)Alloc0(IDELEMS(arg)*sizeof(SObject));
931  intvec * iv=NULL;
932  int i,j;
933
934  if (idRankFreeModule(arg)==0)
935  {
936    iv = idSort(arg);
937    for (i=0;i<IDELEMS(arg);i++)
938    {
939      (resPairs[0])[i].syz = /*pCopy*/(arg->m[(*iv)[i]-1]);
940      arg->m[(*iv)[i]-1] = NULL;
941      (resPairs[0])[i].order = pTotaldegree((resPairs[0])[i].syz);
942    } 
943  }
944  else
945  {
946    iv = new intvec(IDELEMS(arg),1,-1);
947    for (i=0;i<IDELEMS(arg);i++)
948    {
949      (*iv)[i] = pTotaldegree(arg->m[i])+(*cw)[pGetComp(arg->m[i])-1];
950    }
951    for (i=0;i<IDELEMS(arg);i++)
952    {
953      j = syChMin(iv);
954      if (j<0) break;
955      (resPairs[0])[i].syz = arg->m[j];
956      arg->m[j] = NULL;
957      (resPairs[0])[i].order = (*iv)[j];
958      (*iv)[j] = -1;
959    }
960  }
961  if (iv!=NULL)  delete iv;
962  (*Tl)[0] = IDELEMS(arg);
963  return resPairs;
964}
965
966/*3
967* determines the place of a polynomial in the right ordered resolution
968* set the vectors of truecomponents
969*/
970static void syOrder(poly p,syStrategy syzstr,int index,
971                    int realcomp)
972{
973  int i=IDELEMS(syzstr->res[index-1])+1,j=0,k,tc,orc,ie=realcomp-1;
974  int *trind1=syzstr->truecomponents[index-1];
975  int *trind=syzstr->truecomponents[index];
976  int *bc=syzstr->backcomponents[index];
977  int *F1=syzstr->Firstelem[index-1];
978  int *H1=syzstr->Howmuch[index-1];
979  poly pp;
980  polyset o_r=syzstr->orderedRes[index]->m;
981  polyset or1=syzstr->orderedRes[index-1]->m;
982
983  if (p==NULL) return;
984  if (realcomp==0) realcomp=1;
985 
986  if (index>1)
987    tc = trind1[pGetComp(p)]-1;
988  else
989    tc = pGetComp(p)-1;
990  loop         //while ((j<ie) && (trind1[orc]<=tc+1))
991  {
992    if (j>=ie)
993      break;
994    else
995    {
996      orc = pGetComp(o_r[j]);
997      if (trind1[orc]>tc+1) break;
998      j += H1[orc];
999    }
1000  }
1001  if (j>ie)
1002  {
1003    WerrorS("orderedRes to small");
1004    return;
1005  }
1006  ie++;
1007  if (o_r[j]!=NULL)
1008  {
1009    for (k=ie-1;k>j;k--)
1010    {
1011      o_r[k] = o_r[k-1];
1012      bc[k] = bc[k-1];
1013    }
1014  }
1015  o_r[j] = p;
1016  bc[j] = realcomp-1;
1017  (H1[pGetComp(p)])++;
1018  for (k=0;k<i;k++)
1019  {
1020    if (F1[k]>j)
1021      (F1[k])++;
1022  }
1023  if (F1[pGetComp(p)]==0)
1024    F1[pGetComp(p)]=j+1;
1025//Print("write in sort %d till %d\n",index-1,i-1);
1026//Print("poly: ");pWrite(p);
1027//Print("in module %d as %d -th element\n",index,j);
1028  for (k=0;k<IDELEMS((syzstr->res)[index]);k++)
1029  {
1030    if (trind[k]>j)
1031      trind[k] += 1;
1032  }
1033  for (k=IDELEMS((syzstr->res)[index])-1;k>realcomp;k--)
1034    trind[k] = trind[k-1];
1035  trind[realcomp] = j+1;
1036}
1037
1038/*3
1039* reduces all pairs of degree deg in the module index
1040* put the reduced generators to the resolvente which contains
1041* the truncated std
1042*/
1043static void syRedNextPairs(SSet nextPairs, syStrategy syzstr,
1044               int howmuch, int index)
1045{
1046  int i=howmuch-1,j,k=IDELEMS(syzstr->res[index]);
1047  int ks=IDELEMS(syzstr->res[index+1]),kk,l,ll;
1048  int ** Fin=syzstr->Firstelem;
1049  int ** Hin=syzstr->Howmuch;
1050  int ** bin=syzstr->backcomponents;
1051  number coefgcd,n;
1052  ideal o_r=syzstr->orderedRes[index+1];
1053  polyset redset=syzstr->orderedRes[index]->m;
1054  poly p=NULL,q;
1055  BOOLEAN isDivisible;
1056  SObject tso;
1057
1058  if ((nextPairs==NULL) || (howmuch==0)) return;
1059  while ((k>0) && (syzstr->res[index]->m[k-1]==NULL)) k--;
1060  while ((ks>0) && (syzstr->res[index+1]->m[ks-1]==NULL)) ks--;
1061  while (i>=0)
1062  {
1063    tso = nextPairs[i];
1064    isDivisible = FALSE;
1065    if (syzstr->res[index+1]!=NULL)
1066    {
1067      l = Fin[index][pGetComp(tso.lcm)]-1;
1068      if (l>=0)
1069      {
1070        ll = l+Hin[index][pGetComp(tso.lcm)];
1071        while ((l<ll) && (!isDivisible))
1072        {
1073          if (o_r->m[l]!=NULL)
1074          {
1075            isDivisible = isDivisible || 
1076              pDivisibleBy2(o_r->m[l],tso.lcm);
1077          }
1078          l++;
1079        }
1080      }
1081    }
1082    if (!isDivisible)
1083    {
1084      tso.p = 
1085        sySPoly(tso.p1, tso.p2,tso.lcm);
1086      coefgcd = 
1087        nGcd(pGetCoeff(tso.p1),pGetCoeff(tso.p1));
1088      tso.syz = pHead(tso.lcm);
1089      pSetm(tso.syz);
1090      p = tso.syz;
1091      pSetCoeff(p,nDiv(pGetCoeff(tso.p1),coefgcd));
1092      pGetCoeff(p) = nNeg(pGetCoeff(p));
1093      pSetComp(p,tso.ind2+1);
1094      pNext(p) = pHead(tso.lcm);
1095      pIter(p);
1096      pSetm(p);
1097      pSetComp(p,tso.ind1+1);
1098      pSetCoeff(p,nDiv(pGetCoeff(tso.p2),coefgcd));
1099      nDelete(&coefgcd);
1100      if (tso.p != NULL)
1101      {
1102        q = tso.p;
1103        j = Fin[index-1][pGetComp(q)]-1;
1104        int pos = j+Hin[index-1][pGetComp(q)];
1105        loop
1106        {
1107          if (j<0) break;
1108          if (pDivisibleBy2(redset[j],q))
1109          {
1110            //int syL=syLength(redset[j]);
1111            //Print("r");
1112//for(int jj=j+1;jj<k;jj++)
1113//{
1114  //if (syDivisibleBy(redset[jj],q))
1115  //{
1116    //int syL1=syLength(redset[jj]);
1117    //if (syL1<syL)
1118    //{
1119      //j = jj;
1120      //syL = syL1;
1121      //Print("take poly %d with %d monomials\n",j,syL);
1122      //Print("have poly %d with %d monomials\n",jj,syL1);
1123    //}
1124  //}
1125//}
1126            pNext(p) = pHead(q);
1127            pIter(p);
1128            pSetComp(p,bin[index][j]+1);
1129            pGetCoeff(p) = nNeg(pGetCoeff(p));
1130            q = sySPolyRed(q,redset[j]);
1131            if (q==NULL) break;
1132            j = Fin[index-1][pGetComp(q)]-1;
1133            pos = j+Hin[index-1][pGetComp(q)];
1134          }
1135          else
1136          {
1137            j++;
1138            if (j==pos) break;
1139          }
1140        }
1141        tso.p = q;
1142      }
1143      if (tso.p != NULL)
1144      {
1145        if (TEST_OPT_PROT) Print("g");
1146        if (k==IDELEMS((syzstr->res)[index]))
1147        {
1148          pEnlargeSet(&(syzstr->res[index]->m),IDELEMS(syzstr->res[index]),16);
1149          syzstr->truecomponents[index]=(int*)ReAlloc((ADDRESS)syzstr->truecomponents[index],
1150                                       (IDELEMS(syzstr->res[index])+1)*sizeof(int),
1151                                       (IDELEMS(syzstr->res[index])+17)*sizeof(int));
1152          syzstr->backcomponents[index]=(int*)ReAlloc((ADDRESS)syzstr->backcomponents[index],
1153                                       (IDELEMS(syzstr->res[index])+1)*sizeof(int),
1154                                       (IDELEMS(syzstr->res[index])+17)*sizeof(int));
1155          syzstr->Howmuch[index]=(int*)ReAlloc((ADDRESS)syzstr->Howmuch[index],
1156                                       (IDELEMS(syzstr->res[index])+1)*sizeof(int),
1157                                       (IDELEMS(syzstr->res[index])+17)*sizeof(int));
1158          syzstr->Firstelem[index]=(int*)ReAlloc((ADDRESS)syzstr->Firstelem[index],
1159                                       (IDELEMS(syzstr->res[index])+1)*sizeof(int),
1160                                       (IDELEMS(syzstr->res[index])+17)*sizeof(int));
1161        int mk;
1162        for (mk=IDELEMS(syzstr->res[index])+1;mk<IDELEMS(syzstr->res[index])+17;mk++)
1163        {
1164          syzstr->Howmuch[index][mk] = 0;
1165          syzstr->Firstelem[index][mk] = 0;
1166        }
1167//Print("sort %d has now size %d\n",index,IDELEMS(res[index])+17);
1168          IDELEMS(syzstr->res[index]) += 16;
1169          pEnlargeSet(&(syzstr->orderedRes[index]->m),IDELEMS(syzstr->orderedRes[index]),16);
1170          IDELEMS(syzstr->orderedRes[index]) += 16;
1171          redset=syzstr->orderedRes[index]->m;
1172        }
1173        pNext(p) = pHead(tso.p);
1174        pIter(p);
1175        if (p!=NULL)
1176        {
1177          k++;
1178          pSetComp(p,k);
1179          pGetCoeff(p) = nNeg(pGetCoeff(p));
1180        }
1181        if (tso.p!=NULL)
1182        {
1183          syzstr->res[index]->m[k-1] = tso.p;
1184          pNorm(syzstr->res[index]->m[k-1]);
1185          syOrder(syzstr->res[index]->m[k-1],syzstr,index,k);
1186        }
1187        tso.isNotMinimal = p;
1188        tso.p = NULL;
1189      }
1190      else
1191      {
1192        if (TEST_OPT_PROT) Print(".");
1193        if (index % 2==0)
1194          euler++;
1195        else 
1196          euler--;
1197      }
1198      if (ks==IDELEMS(syzstr->res[index+1]))
1199      {
1200        pEnlargeSet(&(syzstr->res[index+1]->m),IDELEMS(syzstr->res[index+1]),16);
1201        syzstr->truecomponents[index+1]=(int*)ReAlloc((ADDRESS)syzstr->truecomponents[index+1],
1202                                       (IDELEMS(syzstr->res[index+1])+1)*sizeof(int),
1203                                       (IDELEMS(syzstr->res[index+1])+17)*sizeof(int));
1204        syzstr->backcomponents[index+1]=(int*)ReAlloc((ADDRESS)syzstr->backcomponents[index+1],
1205                                       (IDELEMS(syzstr->res[index+1])+1)*sizeof(int),
1206                                       (IDELEMS(syzstr->res[index+1])+17)*sizeof(int));
1207        syzstr->Howmuch[index+1]=(int*)ReAlloc((ADDRESS)syzstr->Howmuch[index+1],
1208                                       (IDELEMS(syzstr->res[index+1])+1)*sizeof(int),
1209                                       (IDELEMS(syzstr->res[index+1])+17)*sizeof(int));
1210        syzstr->Firstelem[index+1]=(int*)ReAlloc((ADDRESS)syzstr->Firstelem[index+1],
1211                                       (IDELEMS(syzstr->res[index+1])+1)*sizeof(int),
1212                                       (IDELEMS(syzstr->res[index+1])+17)*sizeof(int));
1213        int mk;
1214        for (mk=IDELEMS(syzstr->res[index+1])+1;mk<IDELEMS(syzstr->res[index+1])+17;mk++)
1215        {
1216          syzstr->Howmuch[index+1][mk] = 0;
1217          syzstr->Firstelem[index+1][mk] = 0;
1218        }
1219//Print("sort %d has now size %d\n",index+1,IDELEMS(res[index+1])+17);
1220        IDELEMS(syzstr->res[index+1]) += 16;
1221        pEnlargeSet(&(syzstr->orderedRes[index+1]->m),IDELEMS(syzstr->orderedRes[index+1]),16);
1222        IDELEMS(syzstr->orderedRes[index+1]) += 16;
1223      }
1224      currcomponents = syzstr->truecomponents[index];
1225      syzstr->res[index+1]->m[ks] = syRedtail(tso.syz,syzstr,index+1);
1226      currcomponents = syzstr->truecomponents[index-1];
1227      syzstr->res[index+1]->m[ks] = tso.syz;
1228      pNorm(syzstr->res[index+1]->m[ks]);
1229      tso.syz =NULL;
1230      tso.syzind = ks;
1231      syOrder(syzstr->res[index+1]->m[ks],syzstr,index+1,ks+1);
1232      ks++;
1233      p = NULL;
1234      nextPairs[i] = tso;
1235    }
1236    else
1237    {
1238      syDeletePair(&nextPairs[i]);
1239      //crit++;
1240    }
1241    i--;
1242  }
1243} 
1244
1245/*3
1246* reduces the generators of the module index in degree deg
1247* (which are actual syzygies of the module index-1)
1248* wrt. the ideal generated by elements of lower degrees
1249*/
1250static void syRedGenerOfCurrDeg(syStrategy syzstr, int deg, int index)
1251{
1252  ideal res=syzstr->res[index];
1253  int i=0,j,k=IDELEMS(res),kk;
1254  SSet sPairs=syzstr->resPairs[index-1];
1255
1256  while ((k>0) && (res->m[k-1]==NULL)) k--;
1257  while ((i<(*syzstr->Tl)[index-1]) && (((sPairs)[i].syz==NULL) || 
1258          ((sPairs)[i].order<deg)))
1259    i++;
1260  if ((i>=(*syzstr->Tl)[index-1]) || ((sPairs)[i].order>deg)) return;
1261  while ((i<(*syzstr->Tl)[index-1]) && (((sPairs)[i].syz==NULL) ||
1262         ((sPairs)[i].order==deg))) 
1263  {
1264    if ((sPairs)[i].syz!=NULL)
1265    {
1266      j = k-1;
1267      while ((j>=0) && (res->m[j]!=NULL) && 
1268             ((sPairs)[i].syz!=NULL))
1269      {
1270        if (pDivisibleBy1(res->m[j],(sPairs)[i].syz))
1271        {
1272          //Print("r");
1273          (sPairs)[i].syz = 
1274            //spSpolyRed(res->m[j],(sPairs)[i].syz,NULL);
1275            sySPolyRed((sPairs)[i].syz,res->m[j]);
1276          j = k-1;
1277        }
1278        else
1279        {
1280          j--;
1281        }
1282      }
1283      if ((sPairs)[i].syz != NULL)
1284      {
1285        if (k==IDELEMS(res))
1286        {
1287          pEnlargeSet(&(res->m),IDELEMS(res),16);
1288          syzstr->truecomponents[index]=(int*)ReAlloc((ADDRESS)syzstr->truecomponents[index],
1289                                       (IDELEMS(res)+1)*sizeof(int),
1290                                       (IDELEMS(res)+17)*sizeof(int));
1291          syzstr->backcomponents[index]=(int*)ReAlloc((ADDRESS)syzstr->backcomponents[index],
1292                                       (IDELEMS(res)+1)*sizeof(int),
1293                                       (IDELEMS(res)+17)*sizeof(int));
1294          syzstr->Howmuch[index]=(int*)ReAlloc((ADDRESS)syzstr->Howmuch[index],
1295                                       (IDELEMS(res)+1)*sizeof(int),
1296                                       (IDELEMS(res)+17)*sizeof(int));
1297          syzstr->Firstelem[index]=(int*)ReAlloc((ADDRESS)syzstr->Firstelem[index],
1298                                       (IDELEMS(res)+1)*sizeof(int),
1299                                       (IDELEMS(res)+17)*sizeof(int));
1300        int mk;
1301        for (mk=IDELEMS(res)+1;mk<IDELEMS(res)+17;mk++)
1302        {
1303          syzstr->Howmuch[index][mk] = 0;
1304          syzstr->Firstelem[index][mk] = 0;
1305        }
1306//Print("sort %d has now size %d\n",index,IDELEMS(res[index])+17);
1307          IDELEMS(res) += 16;
1308          pEnlargeSet(&(syzstr->orderedRes[index]->m),IDELEMS(syzstr->orderedRes[index]),16);
1309          IDELEMS(syzstr->orderedRes[index]) += 16;
1310        }
1311        if (BTEST1(6))
1312        {
1313          if ((sPairs)[i].isNotMinimal==NULL)
1314          {
1315            PrintLn();
1316            Print("minimal generator: ");pWrite((syzstr->resPairs[index-1])[i].syz);
1317            Print("comes from: ");pWrite((syzstr->resPairs[index-1])[i].p1);
1318            Print("and: ");pWrite((syzstr->resPairs[index-1])[i].p2);
1319          }
1320        }
1321        //res->m[k] = (sPairs)[i].syz;
1322        res->m[k] = syRedtail((sPairs)[i].syz,syzstr,index);
1323        (sPairs)[i].syzind = k;
1324        pNorm(res->m[k]);
1325  //      (sPairs)[i].syz = NULL;
1326        k++;
1327        syOrder(res->m[k-1],syzstr,index,k);
1328        euler++;
1329      }
1330      else
1331        (sPairs)[i].syzind = -1;
1332    }
1333    i++;
1334  }
1335}
1336
1337/*3
1338* puts a pair into the right place in resPairs
1339*/
1340static void syEnterPair(SSet sPairs, SObject * so, int * sPlength,int index)
1341{
1342  int ll,k,no=(*so).order,sP=*sPlength,i;
1343  poly p=(*so).lcm;
1344
1345  if ((sP==0) || (sPairs[sP-1].order<=no))
1346    ll = sP;
1347  else if (sP==1)
1348    ll = 0;
1349  else
1350  {
1351    int an=0,en=sP-1;
1352    loop
1353    {
1354      if (an>=en-1)
1355      {
1356        if ((sPairs[an].order<=no) && (sPairs[an+1].order>no))
1357        {
1358          ll = an+1;
1359          break;
1360        }
1361        else if ((sPairs[en].order<=no) && (sPairs[en+1].order>no))
1362        {
1363          ll = en+1;
1364          break;
1365        }
1366        else if (sPairs[an].order>no)
1367        {
1368          ll = an;
1369          break;
1370        }
1371        else
1372        {
1373          Print("Hier ist was faul!\n");
1374          break;
1375        }
1376      }
1377      i=(an+en) / 2;
1378      if (sPairs[i].order <= no)
1379        an=i;
1380      else
1381        en=i;
1382    }
1383  }
1384  for (k=(*sPlength);k>ll;k--)
1385  {
1386    syCopyPair(&sPairs[k-1],&sPairs[k]);
1387  }
1388  syCopyPair(so,&sPairs[ll]);
1389  (*sPlength)++;
1390}
1391
1392/*3
1393* computes pairs from the new elements (beginning with the element newEl)
1394* in the module index
1395*/
1396static void syCreateNewPairs(syStrategy syzstr, int index, int newEl)
1397{
1398  SSet temp;
1399  SObject tso;
1400  int i,ii,j,k=IDELEMS(syzstr->res[index]),l=(*syzstr->Tl)[index],ll;
1401  int qc,first,pos,jj,j1;
1402  int * bci=syzstr->backcomponents[index];
1403  poly p,q;
1404  polyset rs=syzstr->res[index]->m,nPm;
1405
1406  while ((k>0) && (rs[k-1]==NULL)) k--;
1407  if (newEl>=k) return;
1408  ideal nP=idInit(k,syzstr->res[index]->rank);
1409  nPm=nP->m;
1410  while ((l>0) && ((syzstr->resPairs[index])[l-1].p1==NULL)) l--;
1411//Print("new pairs in module %d\n",index);
1412  for (j=newEl;j<k;j++)
1413  {
1414    q = rs[j];
1415    qc = pGetComp(q);
1416    first = syzstr->Firstelem[index-1][pGetComp(q)]-1;
1417    pos = first+syzstr->Howmuch[index-1][pGetComp(q)];
1418    for (i=first;i<pos;i++)
1419    {
1420      jj = bci[i];
1421      if (jj>=j) break;
1422      p = pOne();
1423      pLcm(rs[jj],q,p);
1424      pSetComp(p,j+1);
1425      ii = first;
1426      loop
1427      {
1428        j1 = bci[ii];
1429        if (nPm[j1]!=NULL)
1430        {
1431          if (pDivisibleBy2(nPm[j1],p))
1432          {
1433            pDelete(&p);
1434            break;
1435          }
1436          else if (pDivisibleBy2(p,nPm[j1]))
1437          {
1438            pDelete(&(nPm[j1]));
1439            break;
1440          }
1441        }
1442        ii++;
1443        if (ii>=pos) break;
1444      }
1445      if (p!=NULL)
1446      {
1447        nPm[jj] = p;
1448      }
1449    }
1450    //for (i=0;i<j;i++)
1451    for (i=first;i<pos;i++)
1452    {
1453      ii = bci[i];
1454      if (nPm[ii]!=NULL)
1455      {
1456        if (l>=(*syzstr->Tl)[index])
1457        {
1458          temp = (SSet)Alloc0(((*syzstr->Tl)[index]+16)*sizeof(SObject));
1459          for (ll=0;ll<(*syzstr->Tl)[index];ll++)
1460          {
1461            temp[ll].p = (syzstr->resPairs[index])[ll].p;
1462            temp[ll].p1 = (syzstr->resPairs[index])[ll].p1;
1463            temp[ll].p2 = (syzstr->resPairs[index])[ll].p2;
1464            temp[ll].syz = (syzstr->resPairs[index])[ll].syz;
1465            temp[ll].lcm = (syzstr->resPairs[index])[ll].lcm;
1466            temp[ll].ind1 = (syzstr->resPairs[index])[ll].ind1;
1467            temp[ll].ind2 = (syzstr->resPairs[index])[ll].ind2;
1468            temp[ll].syzind = (syzstr->resPairs[index])[ll].syzind;
1469            temp[ll].order = (syzstr->resPairs[index])[ll].order;
1470            temp[ll].isNotMinimal = (syzstr->resPairs[index])[ll].isNotMinimal;
1471          }
1472          Free((ADDRESS)syzstr->resPairs[index],(*syzstr->Tl)[index]*sizeof(SObject));
1473          (*syzstr->Tl)[index] += 16;
1474          syzstr->resPairs[index] = temp;
1475        }
1476        tso.lcm = p = nPm[ii];
1477        nPm[ii] = NULL;
1478        tso.order = pGetOrder(p) = pTotaldegree(p);
1479        if ((syzstr->cw!=NULL) && (index>0) && (pGetComp(q)>0))
1480        {
1481          int ii=index-1,jj=pGetComp(q);
1482          while (ii>0) 
1483          {
1484            jj = pGetComp(syzstr->res[ii]->m[jj-1]);
1485            ii--;
1486          }
1487          tso.order += (*syzstr->cw)[jj-1];
1488        }
1489        tso.p1 = rs[ii];
1490        tso.p2 = q;
1491        tso.ind1 = ii;
1492        tso.ind2 = j;
1493        tso.syzind = -1;
1494        tso.isNotMinimal = NULL;
1495        tso.p = NULL;
1496        tso.syz = NULL;
1497        syEnterPair(syzstr->resPairs[index],&tso,&l,index);
1498      }
1499    }
1500  }
1501  idDelete(&nP);
1502}
1503
1504static void sySyzTail(syStrategy syzstr, int index, int newEl)
1505{
1506  int j,ll,k=IDELEMS((syzstr->res)[index]);
1507  while ((k>0) && ((syzstr->res)[index]->m[k-1]==NULL)) k--;
1508  for (j=newEl;j<k;j++)
1509  {
1510    ll = 0;
1511    while ((ll<(*syzstr->Tl)[index]) && (syzstr->resPairs[index][ll].p1!=NULL) &&
1512           (syzstr->resPairs[index][ll].ind1!=j) && (syzstr->resPairs[index][ll].ind2!=j))
1513      ll++;
1514    if ((ll<(*syzstr->Tl)[index]) && (syzstr->resPairs[index][ll].p1!=NULL))
1515      syzstr->res[index]->m[j] = syRedtail(syzstr->res[index]->m[j],syzstr,index);
1516  }
1517}
1518
1519static SSet syChosePairsPutIn(syStrategy syzstr, int *index,
1520               int *howmuch, int * actdeg, int an, int en)
1521{
1522  int newdeg=*actdeg,newindex=-1,i,t,sldeg;
1523  poly p;
1524  SSet result;
1525  SRes resPairs=syzstr->resPairs;
1526 
1527  if (an>syzstr->length) return NULL;
1528  if (en>syzstr->length) en=syzstr->length;
1529  while (*index<en)
1530  {
1531    if (resPairs[*index]!=NULL)
1532    {
1533      sldeg = (*actdeg)+*index;
1534      i = 0;
1535      if (*index!=0)
1536      {
1537        while ((i<(*syzstr->Tl)[*index]))
1538        {
1539          if ((resPairs[*index])[i].lcm!=NULL)
1540          {
1541            if ((resPairs[*index])[i].order == sldeg)
1542            {
1543              result = &(resPairs[*index])[i];
1544              *howmuch =1;
1545              i++;
1546              while ((i<(*syzstr->Tl)[*index]) && ((resPairs[*index])[i].lcm!=NULL)
1547                      && ((resPairs[*index])[i].order == sldeg))
1548              {
1549                i++;
1550                (*howmuch)++;
1551              }
1552              return result;
1553            }
1554          }
1555          i++;
1556        }
1557      }
1558      else
1559      {
1560        while ((i<(*syzstr->Tl)[*index]))
1561        {
1562          if ((resPairs[*index])[i].syz!=NULL)
1563          {
1564            if ((resPairs[*index])[i].order == sldeg)
1565            {
1566              result = &(resPairs[*index])[i];
1567              (*howmuch) =1;
1568              i++;
1569              while ((i<(*syzstr->Tl)[*index]) && ((resPairs[*index])[i].syz!=NULL)
1570                      && ((resPairs[*index])[i].order == *actdeg))
1571              {
1572                i++;
1573                (*howmuch)++;
1574              }
1575              return result;
1576            }
1577          }
1578          i++;
1579        }
1580      }
1581    }
1582    (*index)++;
1583  }
1584  *index = an;
1585  //if (TEST_OPT_PROT) Print("(Euler:%d)",euler);
1586  while (*index<en)
1587  {
1588    if (resPairs[*index]!=NULL)
1589    {
1590      i = 0;
1591      while ((i<(*syzstr->Tl)[*index]))
1592      {
1593        t = *actdeg+*index;
1594        if (((resPairs[*index])[i].lcm!=NULL) || 
1595              ((resPairs[*index])[i].syz!=NULL))
1596        {
1597          if ((resPairs[*index])[i].order > t)
1598            t = (resPairs[*index])[i].order;
1599        }
1600        if ((t>*actdeg+*index) && ((newdeg==*actdeg) || (t<newdeg+*index)))
1601        {
1602          newdeg = t-*index;
1603          newindex = *index;
1604          break;
1605        }
1606        i++;
1607      } 
1608    }
1609    (*index)++;
1610  }
1611  if (newdeg>*actdeg)
1612  {
1613    *actdeg = newdeg;
1614    *index = newindex;
1615    return syChosePairsPutIn(syzstr,index,howmuch,actdeg,an,en);
1616  }
1617  else return NULL;
1618}
1619
1620/*3
1621* FOR THE HOMOGENEOUS CASE ONLY!
1622* looks through the pair set and the given module for
1623* remaining pairs or generators to consider
1624* returns a pointer to the first pair and the number of them in the given module
1625* works with slanted degree (i.e. deg=realdeg-index)
1626*/
1627static SSet syChosePairs(syStrategy syzstr, int *index,
1628               int *howmuch, int * actdeg)
1629{
1630  return syChosePairsPutIn(syzsts,index,howmuch,actdeg,0,syzstr->length);
1631}
1632
1633/*3
1634* FOR THE INHOMOGENEOUS CASE ONLY!
1635* looks through the pair set and the given module for
1636* remaining pairs or generators to consider
1637* returns a pointer to the first pair and the number of them in the given module
1638* works with slanted degree (i.e. deg=realdeg-index)
1639* looks first through the 0 and 1 module then through the other
1640*/
1641static SSet syChosePairsIH(syStrategy syzstr, int *index,
1642               int *howmuch, int * actdeg, int mindeg)
1643{
1644  SSet result=NULL;
1645
1646  result = syChosePairsPutIn(syzstr,index,howmuch,actdeg,0,2);
1647  if (result == NULL)
1648  {
1649    *actdeg = mindeg;
1650    result = syChosePairsPutIn(syzstr,index,howmuch,actdeg,2,syzstr->length);
1651  }
1652  return result;
1653}
1654
1655/*3
1656* looks through the pair set and the given module for
1657* remaining pairs or generators to consider
1658* returns a pointer to the first pair and the number of them in the given module
1659* works deg by deg
1660*/
1661/*
1662*static SSet syChosePairs1(SRes resPairs,intvec * Tl, int *index, int *howmuch,
1663*                   int length,int * actdeg)
1664*{
1665*  int newdeg=*actdeg,newindex=-1,i,t;
1666*  SSet result;
1667* 
1668*  while (*index>=0)
1669*  {
1670*    if (resPairs[*index]!=NULL)
1671*    {
1672*      i = 0;
1673*      if (*index!=0)
1674*      {
1675*        while ((i<(*Tl)[*index]))
1676*        {
1677*          if ((resPairs[*index])[i].lcm!=NULL)
1678*          {
1679*            if (pGetOrder((resPairs[*index])[i].lcm) == *actdeg)
1680*            {
1681*              result = &(resPairs[*index])[i];
1682*              *howmuch =1;
1683*              i++;
1684*              while ((i<(*Tl)[*index]) && ((resPairs[*index])[i].lcm!=NULL)
1685*                      && (pGetOrder((resPairs[*index])[i].lcm) == *actdeg))
1686*              {
1687*                i++;
1688*                (*howmuch)++;
1689*              }
1690*              return result;
1691*            }
1692*          }
1693*          i++;
1694*        }
1695*      }
1696*      else
1697*      {
1698*        while ((i<(*Tl)[*index]))
1699*        {
1700*          if ((resPairs[*index])[i].syz!=NULL)
1701*          {
1702*            if ((resPairs[*index])[i].order == *actdeg)
1703*            {
1704*              result = &(resPairs[*index])[i];
1705*              (*howmuch) =1;
1706*              i++;
1707*              while ((i<(*Tl)[*index]) && ((resPairs[*index])[i].syz!=NULL)
1708*                      && ((resPairs[*index])[i].order == *actdeg))
1709*              {
1710*                i++;
1711*                (*howmuch)++;
1712*              }
1713*              return result;
1714*            }
1715*          }
1716*          i++;
1717*        }
1718*      }
1719*    }
1720*    (*index)--;
1721*  }
1722*  *index = length-1;
1723*  while (*index>=0)
1724*  {
1725*    if (resPairs[*index]!=NULL)
1726*    {
1727*      i = 0;
1728*      while ((i<(*Tl)[*index]))
1729*      {
1730*        t = *actdeg;
1731*        if ((resPairs[*index])[i].lcm!=NULL)
1732*        {
1733*          if (pGetOrder((resPairs[*index])[i].lcm) > *actdeg)
1734*            t = pGetOrder((resPairs[*index])[i].lcm);
1735*        }
1736*        else if ((resPairs[*index])[i].syz!=NULL)
1737*        {
1738*          if ((resPairs[*index])[i].order > *actdeg)
1739*            t = (resPairs[*index])[i].order;
1740*        }
1741*        if ((t>*actdeg) && ((newdeg==*actdeg) || (t<newdeg)))
1742*        {
1743*          newdeg = t;
1744*          newindex = *index;
1745*          break;
1746*        }
1747*        i++;
1748*      }
1749*    }
1750*    (*index)--;
1751*  }
1752*  if (newdeg>*actdeg)
1753*  {
1754*    *actdeg = newdeg;
1755*    *index = newindex;
1756*    return syChosePairs1(resPairs,Tl,index,howmuch,length,actdeg);
1757*  }
1758*  else return NULL;
1759*}
1760*/
1761/*3
1762* statistics of the resolution
1763*/
1764static void syStatistics(resolvente res,int length)
1765{
1766  int i,j=1,k,deg=0;
1767
1768  PrintLn();
1769  while ((j<length) && (res[j]!=NULL))
1770  {
1771    Print("In module %d: \n",j);
1772    k = 0;
1773    while ((k<IDELEMS(res[j])) && (res[j]->m[k]!=NULL))
1774    {
1775      i = 1;
1776      deg = pGetOrder(res[j]->m[k]);
1777      k++;
1778      while ((k<IDELEMS(res[j])) && (res[j]->m[k]!=NULL) &&
1779              (pGetOrder(res[j]->m[k])==deg))
1780      {
1781        i++;
1782        k++;
1783      }
1784      Print("%d elements of degree %d\n",i,deg);
1785    }
1786    j++;
1787  }
1788}
1789
1790/*3
1791* initialize a module
1792*/
1793static int syInitSyzMod(syStrategy syzstr, int index, int init=17)
1794{
1795  int result;
1796
1797  if (syzstr->res[index]==NULL) 
1798  {
1799    syzstr->res[index] = idInit(init-1,1);
1800    syzstr->truecomponents[index] = (int*)Alloc0(init*sizeof(int));
1801    if (index==0)
1802    {
1803      for (int i=0;i<init;i++)
1804        syzstr->truecomponents[0][i] = i;
1805    }
1806    syzstr->backcomponents[index] = (int*)Alloc0(init*sizeof(int));
1807    syzstr->Howmuch[index] = (int*)Alloc0(init*sizeof(int));
1808    syzstr->Firstelem[index] = (int*)Alloc0(init*sizeof(int));
1809//Print("sort %d has now size %d\n",index,init);
1810    syzstr->orderedRes[index] = idInit(init-1,1);
1811    result = 0;
1812  }
1813  else
1814  {
1815    result = IDELEMS(syzstr->res[index]);
1816    while ((result>0) && (syzstr->res[index]->m[result-1]==NULL)) result--;
1817  }
1818  return result;
1819}
1820
1821void syKillComputation(syStrategy syzstr)
1822{
1823//Print("ref: %d\n",syzstr->references);
1824  if (syzstr->references>0)
1825  {
1826    (syzstr->references)--;
1827  }
1828  else
1829  {
1830    int i,j;
1831 
1832    if (syzstr->resPairs!=NULL)
1833    {
1834      for (i=0;i<syzstr->length;i++)
1835      {
1836        for (j=0;j<(*syzstr->Tl)[i];j++)
1837        {
1838          if ((syzstr->resPairs[i])[j].lcm!=NULL)
1839            pDelete(&((syzstr->resPairs[i])[j].lcm));
1840        }
1841        if (syzstr->orderedRes[i]!=NULL)
1842        {
1843          for (j=0;j<IDELEMS(syzstr->orderedRes[i]);j++)
1844          {
1845            syzstr->orderedRes[i]->m[j] = NULL;
1846          }
1847        }
1848        idDelete(&(syzstr->orderedRes[i]));
1849        if (syzstr->truecomponents[i]!=NULL)
1850        {
1851          Free((ADDRESS)syzstr->truecomponents[i],(IDELEMS(syzstr->res[i])+1)*sizeof(int));
1852          syzstr->truecomponents[i]=NULL;
1853        }
1854        if (syzstr->backcomponents[i]!=NULL)
1855        {
1856          Free((ADDRESS)syzstr->backcomponents[i],(IDELEMS(syzstr->res[i])+1)*sizeof(int));
1857          syzstr->backcomponents[i]=NULL;
1858        }
1859        if (syzstr->Howmuch[i]!=NULL)
1860        {
1861          Free((ADDRESS)syzstr->Howmuch[i],(IDELEMS(syzstr->res[i])+1)*sizeof(int));
1862          syzstr->Howmuch[i]=NULL;
1863        }
1864        if (syzstr->Firstelem[i]!=NULL)
1865        {
1866          Free((ADDRESS)syzstr->Firstelem[i],(IDELEMS(syzstr->res[i])+1)*sizeof(int));
1867          syzstr->Firstelem[i]=NULL;
1868        }
1869        if (syzstr->res[i]!=NULL)
1870        {
1871          for (j=0;j<IDELEMS(syzstr->res[i]);j++)
1872          {
1873            if (syzstr->res[i]->m[j]!=NULL)
1874              pDelete(&(syzstr->res[i]->m[j]));
1875          }
1876        }
1877        if (syzstr->weights!=0)
1878        {
1879          if (syzstr->weights[i]!=NULL)
1880          {
1881            delete syzstr->weights[i];
1882          }
1883          Free((ADDRESS)syzstr->weights,syzstr->length*sizeof(intvec*));
1884        }
1885        idDelete(&(syzstr->res[i]));
1886        Free((ADDRESS)syzstr->resPairs[i],(*syzstr->Tl)[i]*sizeof(SObject));
1887      }
1888      Free((ADDRESS)syzstr->resPairs,syzstr->length*sizeof(SObject*));
1889      Free((ADDRESS)syzstr->res,(syzstr->length+1)*sizeof(ideal));
1890      Free((ADDRESS)syzstr->orderedRes,(syzstr->length+1)*sizeof(ideal));
1891      Free((ADDRESS)syzstr->truecomponents,(syzstr->length+1)*sizeof(int*));
1892      Free((ADDRESS)syzstr->backcomponents,(syzstr->length+1)*sizeof(int*));
1893      Free((ADDRESS)syzstr->Howmuch,(syzstr->length+1)*sizeof(int*));
1894      Free((ADDRESS)syzstr->Firstelem,(syzstr->length+1)*sizeof(int*));
1895      Free((ADDRESS)syzstr->binom,pVariables*(syzstr->highdeg_1)*sizeof(int));
1896    }
1897    if (syzstr->minres!=NULL)
1898    {
1899      for (i=0;i<syzstr->length;i++)
1900      {
1901        if (syzstr->minres[i]!=NULL)
1902        {
1903          for (j=0;j<IDELEMS(syzstr->minres[i]);j++)
1904          {
1905            if (syzstr->minres[i]->m[j]!=NULL)
1906              pDelete(&(syzstr->minres[i]->m[j]));
1907          }
1908        }
1909        idDelete(&(syzstr->minres[i]));
1910      }
1911      Free((ADDRESS)syzstr->minres,(syzstr->length+1)*sizeof(ideal));
1912    }
1913    if (syzstr->fullres!=NULL)
1914    {
1915      for (i=0;i<syzstr->length;i++)
1916      {
1917        if (syzstr->fullres[i]!=NULL)
1918        {
1919          for (j=0;j<IDELEMS(syzstr->fullres[i]);j++)
1920          {
1921            if (syzstr->fullres[i]->m[j]!=NULL)
1922              pDelete(&(syzstr->fullres[i]->m[j]));
1923          }
1924        }
1925        idDelete(&(syzstr->fullres[i]));
1926      }
1927      Free((ADDRESS)syzstr->fullres,(syzstr->length+1)*sizeof(ideal));
1928    }
1929    if (syzstr->cw!=NULL)
1930      delete syzstr->cw;
1931    if (syzstr->resolution!=NULL)
1932      delete syzstr->resolution;
1933    Free((ADDRESS)syzstr,sizeof(ssyStrategy));
1934  }
1935}
1936
1937intvec * syBettiOfComputation(syStrategy syzstr)
1938{
1939  int dummy;
1940  if (syzstr->resPairs!=NULL)
1941  {
1942    int i,j=-1,jj=-1,l;
1943    SRes rP=syzstr->resPairs;
1944   
1945    l = syzstr->length;
1946    while ((l>0) && (rP[l-1]==NULL)) l--;
1947    if (l==0) return NULL;
1948    l--;
1949    while (l>=0)
1950    {
1951      i = 0;
1952      while ((i<(*syzstr->Tl)[l]) && 
1953        ((rP[l][i].lcm!=NULL) || (rP[l][i].syz!=NULL)))
1954      {
1955        if (rP[l][i].isNotMinimal==NULL) 
1956        {
1957          if (j<rP[l][i].order-l)
1958            j = rP[l][i].order-l;
1959          if (jj<l)
1960            jj = l;
1961        }
1962        i++;
1963      }
1964      l--;
1965    }
1966    jj=jj+2;
1967    intvec *result=new intvec(j,jj,0);
1968    IMATELEM(*result,1,1) = max(1,idRankFreeModule(syzstr->res[1]));
1969    for (i=0;i<jj;i++)
1970    {
1971      j = 0;
1972      while ((j<(*syzstr->Tl)[i]) && 
1973        ((rP[i][j].lcm!=NULL) || (rP[i][j].syz!=NULL)))
1974      {
1975        if (rP[i][j].isNotMinimal==NULL)
1976          IMATELEM(*result,rP[i][j].order-i,i+2)++;
1977        j++;
1978      }
1979    }
1980    return result;
1981  }
1982  else if (syzstr->fullres!=NULL)
1983    return syBetti(syzstr->fullres,syzstr->length,&dummy);
1984  else
1985    return syBetti(syzstr->minres,syzstr->length,&dummy);
1986}
1987
1988int syLength(syStrategy syzstr)
1989{
1990  return syzstr->length;
1991}
1992
1993int sySize(syStrategy syzstr)
1994{
1995  resolvente r=syzstr->res;
1996  if (r==NULL)
1997    r = syzstr->fullres;
1998  if (r==NULL)
1999    r = syzstr->minres;
2000  if (r==NULL)
2001  {
2002    WerrorS("No resolution found");
2003    return 0;
2004  }
2005  int i=syzstr->length;
2006  while ((i>0) && (r[i-1]==NULL)) i--;
2007  return i;
2008}
2009
2010int syDim(syStrategy syzstr)
2011{
2012  int i,j=-1,l;
2013  if (syzstr->resPairs!=NULL)
2014  {
2015    SRes rP=syzstr->resPairs;
2016   
2017    l = syzstr->length;
2018    while ((l>0) && (rP[l-1]==NULL)) l--;
2019    if (l==0) return -1;
2020    l--;
2021    while (l>=0)
2022    {
2023      i = 0;
2024      while ((i<(*syzstr->Tl)[l]) && 
2025        ((rP[l][i].lcm!=NULL) || (rP[l][i].syz!=NULL)) &&
2026        (rP[l][i].isNotMinimal!=NULL))
2027      {
2028        i++;
2029      }
2030      if ((i<(*syzstr->Tl)[l]) &&
2031        ((rP[l][i].lcm!=NULL) || (rP[l][i].syz!=NULL)) &&
2032        (rP[l][i].isNotMinimal==NULL))
2033        return l;
2034      l--;
2035    }
2036    return l;
2037  }
2038  else 
2039    return sySize(syzstr);
2040}
2041
2042syStrategy syCopy(syStrategy syzstr)
2043{
2044  syStrategy result=syzstr;
2045  (result->references)++;
2046  return result;
2047}
2048
2049static void syPrintEmptySpaces(int i)
2050{
2051  if (i!=0)
2052  {
2053    Print(" ");
2054    syPrintEmptySpaces(i/10);
2055  }
2056}
2057
2058static void syPrintEmptySpaces1(int i)
2059{
2060  if (i!=0)
2061  {
2062    Print(" ");
2063    syPrintEmptySpaces1(i-1);
2064  }
2065}
2066
2067static int syLengthInt(int i)
2068{
2069  int j=0;
2070
2071  if (i==0) return 1;
2072  while (i!=0)
2073  {
2074    j++;
2075    i = i/10;
2076  }
2077  return j;
2078}
2079
2080void syPrint(syStrategy syzstr)
2081{
2082  if ((syzstr->resPairs==NULL) && (syzstr->fullres==NULL) 
2083     && (syzstr->minres==NULL))
2084  {
2085    Print("No resolution defined\n");
2086    return;
2087  }
2088  int l=0;
2089  if (syzstr->resolution==NULL)
2090  {
2091    int j;
2092    if (syzstr->resPairs!=NULL)
2093    {
2094      syzstr->resolution = new intvec(syzstr->length+1);
2095      SRes rP=syzstr->resPairs;
2096      (*syzstr->resolution)[0] = max(1,idRankFreeModule(syzstr->res[1]));
2097      while ((l<syzstr->length) && (rP[l]!=NULL))
2098      {
2099        j=0;
2100        while ((j<(*syzstr->Tl)[l]) && 
2101          ((rP[l][j].lcm!=NULL) || (rP[l][j].syz!=NULL)))
2102        {
2103          if (rP[l][j].isNotMinimal==NULL)
2104            ((*syzstr->resolution)[l+1])++;
2105          j++;
2106        }
2107        l++;
2108      }
2109    }
2110    else
2111    {
2112      resolvente rr;
2113      syzstr->resolution = new intvec(syzstr->length+2);
2114      if (syzstr->minres!=NULL)
2115        rr = syzstr->minres;
2116      else
2117        rr = syzstr->fullres;
2118      (*syzstr->resolution)[0] = max(1,idRankFreeModule(rr[0]));
2119      while ((l<syzstr->length) && (rr[l]!=NULL))
2120      {
2121        j = IDELEMS(rr[l]);
2122        while ((l>0) && (rr[l]->m[j-1]==NULL)) j--;
2123        ((*syzstr->resolution)[l+1]) = j;
2124        l++;
2125      }
2126    }
2127  }
2128  char *sn=currRingHdl->id;
2129  int sl=strlen(sn);
2130  syPrintEmptySpaces1(sl);
2131  l = 0;
2132  loop
2133  {
2134    Print("%d",(*syzstr->resolution)[l]);
2135    if ((l>=syzstr->resolution->length()) || ((*syzstr->resolution)[l]==0))
2136      break;
2137    syPrintEmptySpaces1(sl+5);
2138    l++;
2139  }
2140  PrintLn();
2141  l = 0;
2142  loop
2143  {
2144    Print(sn);
2145    if ((l>=syzstr->resolution->length()) || ((*syzstr->resolution)[l]==0))
2146      break;
2147    syPrintEmptySpaces((*syzstr->resolution)[l]);
2148    Print(" <-- ");
2149    l++;
2150  }
2151  PrintLn();
2152  PrintLn();
2153  l = 0;
2154  loop
2155  {
2156    Print("%d",l);
2157    if ((l>=syzstr->resolution->length()) || ((*syzstr->resolution)[l]==0))
2158      break;
2159    syPrintEmptySpaces1(sl+5+syLengthInt((*syzstr->resolution)[l])-
2160                         syLengthInt(l));
2161    l++;
2162  }
2163  PrintLn();
2164  if (syzstr->minres==NULL)
2165  {
2166    Print("resolution not minimized yet");
2167    PrintLn();
2168  }
2169}
2170
2171static resolvente syReorder(resolvente res,int length,
2172        syStrategy syzstr,BOOLEAN toCopy=TRUE,resolvente totake=NULL)
2173{
2174  int i,j,l;
2175  poly p,q,tq;
2176  polyset ri1;
2177  resolvente fullres;
2178  fullres = (resolvente)Alloc0((length+1)*sizeof(ideal));
2179  if (totake==NULL) 
2180    totake = res;
2181  for (i=length-1;i>0;i--)
2182  {
2183    if (res[i]!=NULL)
2184    {
2185      if (i>1)
2186      {
2187        j = IDELEMS(res[i-1]);
2188        while ((j>0) && (res[i-1]->m[j-1]==NULL)) j--;
2189        fullres[i-1] = idInit(IDELEMS(res[i]),j);
2190        ri1 = totake[i-1]->m;
2191        for (j=IDELEMS(res[i])-1;j>=0;j--)
2192        {
2193          p = res[i]->m[j];
2194          q = NULL;
2195          while (p!=NULL)
2196          {
2197            if (toCopy)
2198            {
2199              tq = pHead(p);
2200              pIter(p);
2201            }
2202            else
2203            {
2204              tq = p;
2205              pIter(p);
2206              pNext(tq) = NULL;
2207            }
2208            for (l=pVariables;l>0;l--)
2209              pSubExp(tq,l, pGetExp(ri1[pGetComp(tq)-1],l));
2210            pSetm(tq);
2211            q = pAdd(q,tq);
2212          }
2213          fullres[i-1]->m[j] = q;
2214        }
2215      }
2216      else
2217      {
2218        if (toCopy)
2219          fullres[i-1] = idCopy(res[i]);
2220        else
2221          fullres[i-1] = res[i];
2222        for (j=IDELEMS(res[i])-1;j>=0;j--)
2223          fullres[i-1]->m[j] = pOrdPolyInsertSetm(fullres[i-1]->m[j]);
2224      }
2225    }
2226  }
2227  if (!toCopy)
2228    Free((ADDRESS)res,(length+1)*sizeof(ideal));
2229  //syzstr->length = length;
2230  return fullres;
2231}
2232
2233lists syConvRes(syStrategy syzstr)
2234{
2235  if ((syzstr->fullres==NULL) && (syzstr->minres==NULL))
2236  {
2237    syzstr->fullres = syReorder(syzstr->res,syzstr->length,syzstr);
2238  }
2239  resolvente tr;
2240  int typ0=IDEAL_CMD;
2241  if (syzstr->minres!=NULL)
2242    tr = syzstr->minres;
2243  else
2244    tr = syzstr->fullres;
2245  resolvente trueres=NULL;
2246  intvec ** w=NULL;
2247  if (syzstr->length>0)
2248  {
2249    trueres=(resolvente)Alloc0((syzstr->length)*sizeof(ideal));
2250    for (int i=(syzstr->length)-1;i>=0;i--)
2251    {
2252      if (tr[i]!=NULL)
2253      {
2254        trueres[i] = idCopy(tr[i]);
2255      }
2256    }
2257    if (idRankFreeModule(trueres[0]) > 0)
2258      typ0 = MODUL_CMD;
2259    if (syzstr->weights!=NULL)
2260    {
2261      w = (intvec**)Alloc0((syzstr->length)*sizeof(intvec*));
2262      for (int i=(syzstr->length)-1;i>=0;i--)
2263      {
2264        if (syzstr->weights[i]!=NULL) w[i] = ivCopy(syzstr->weights[i]);
2265      }
2266    }
2267  }
2268  return liMakeResolv(trueres,syzstr->length,-1,typ0,w);
2269}
2270
2271syStrategy syConvList(lists li)
2272{
2273  int typ0;
2274  syStrategy result=(syStrategy)Alloc0(sizeof(ssyStrategy));
2275
2276  resolvente fr = liFindRes(li,&(result->length),&typ0,&(result->weights));
2277  result->fullres = (resolvente)Alloc0((result->length+1)*sizeof(ideal));
2278  for (int i=result->length-1;i>=0;i--)
2279  {
2280    if (fr[i]!=NULL)
2281      result->fullres[i] = idCopy(fr[i]);
2282  }
2283  Free((ADDRESS)fr,(result->length)*sizeof(ideal));
2284  return result;
2285}
2286
2287syStrategy syForceMin(lists li)
2288{
2289  int typ0;
2290  syStrategy result=(syStrategy)Alloc0(sizeof(ssyStrategy));
2291
2292  resolvente fr = liFindRes(li,&(result->length),&typ0);
2293  result->minres = (resolvente)Alloc0((result->length+1)*sizeof(ideal));
2294  for (int i=result->length-1;i>=0;i--)
2295  {
2296    if (fr[i]!=NULL)
2297      result->minres[i] = idCopy(fr[i]);
2298  }
2299  Free((ADDRESS)fr,(result->length)*sizeof(ideal));
2300  return result;
2301}
2302
2303static poly syStripOut(poly p,intvec * toStrip)
2304{
2305  if (toStrip==NULL) return p;
2306  poly pp=p;
2307 
2308  while ((pp!=NULL) && ((*toStrip)[pGetComp(pp)]!=0))
2309    pDelete1(&pp);
2310  p = pp;
2311  if (pp!=NULL) 
2312  {
2313    while (pNext(pp)!=NULL)
2314    {
2315      if ((*toStrip)[pGetComp(pNext(pp))]!=0)
2316        pDelete1(&pNext(pp));
2317      else
2318        pIter(pp);
2319    }
2320  }
2321  return p;
2322}
2323
2324static poly syMinimizeP(poly toMin,syStrategy syzstr,int pNum,int index,
2325                        intvec * toStrip)
2326{
2327  int i,j,tc,lastin,newin=pNum-1;
2328  poly p,pp=pCopy(toMin),q=NULL,tq,pisN;
2329  SSet sPairs=syzstr->resPairs[index];
2330
2331  pp = syStripOut(pp,toStrip);
2332  while (newin>=0)
2333  {
2334    lastin = newin;
2335    while ((newin>=0) && (sPairs[lastin].order==sPairs[newin].order))
2336      newin--;
2337//Print("Hier lastin:%d newin:%d\n",lastin,newin);
2338    for (i=newin+1;i<=lastin;i++)
2339    {
2340      if (sPairs[i].isNotMinimal!=NULL)
2341      {
2342        pisN = sPairs[i].isNotMinimal;
2343        tc = pGetComp(pisN);
2344        p = pp;
2345        while (p!=NULL)
2346        {
2347          if (pGetComp(p)==tc)
2348          {
2349            tq = pInit();
2350            for(j=pVariables; j>0; j--)
2351              pSetExp(tq,j, pGetExp(p,j)-pGetExp(pisN,j));
2352            pSetComp(tq, 0);
2353            pSetCoeff0(tq,nDiv(pGetCoeff(p),pGetCoeff(pisN)));
2354            pGetCoeff(tq) = nNeg(pGetCoeff(tq));
2355            q = syAdd(q,syStripOut(syMultT1(syzstr->res[index+1]->m[sPairs[i].syzind],tq),toStrip));
2356          } 
2357          pIter(p);
2358        }
2359        if (q!=NULL)
2360        {
2361          pp = syAdd(pp,q);
2362          q = NULL;
2363        }
2364      }
2365    }
2366  }
2367  return pp;
2368}
2369
2370static void syKillEmptyEntres(resolvente res,int length)
2371{
2372  int i,j,jj,k,rj;
2373  intvec * changes;
2374  poly p;
2375  ideal ri;
2376
2377  for (i=0;i<length;i++)
2378  {
2379    ri = res[i];
2380    if (ri!=NULL)
2381    {
2382      rj = IDELEMS(ri);
2383      changes = new intvec(rj+1,1,-1);
2384      while ((rj>0) && (ri->m[rj-1]==NULL)) rj--;
2385      j = k = 0;
2386      while (j+k<rj)
2387      {
2388        if (ri->m[j+k]!=NULL)
2389        {
2390          ri->m[j] = ri->m[j+k];
2391          (*changes)[j+k+1] = j+1; 
2392          j++;
2393        }
2394        else
2395        {
2396          k++;
2397        }
2398      }
2399      for (jj=j;jj<rj;jj++)
2400        ri->m[jj] = NULL;
2401      if (res[i+1]!=NULL)
2402      {
2403        ri = res[i+1];
2404        for (j=IDELEMS(ri)-1;j>=0;j--)
2405        {
2406          p = ri->m[j];
2407          while (p!=NULL)
2408          {
2409            pSetComp(p,(*changes)[pGetComp(p)]);
2410            pIter(p);
2411          }
2412        }
2413      }
2414    }
2415  }
2416}
2417
2418static intvec * syToStrip(syStrategy syzstr, int index)
2419{
2420  intvec * result=NULL;
2421
2422  if ((syzstr->resPairs[index-1]!=NULL) && (!idIs0(syzstr->res[index])))
2423  {
2424    result=new intvec(IDELEMS(syzstr->res[index])+1);
2425    for (int i=(*syzstr->Tl)[index-1]-1;i>=0;i--)
2426    {
2427      if (syzstr->resPairs[index-1][i].isNotMinimal!=NULL)
2428      {
2429        (*result)[syzstr->resPairs[index-1][i].syzind+1] = 1;
2430      }
2431    }
2432  }
2433  return result;
2434}
2435
2436static resolvente syReadOutMinimalRes(syStrategy syzstr, 
2437           BOOLEAN computeStd=FALSE)
2438{
2439  intvec * Strip;
2440  resolvente tres=(resolvente)Alloc0((syzstr->length+1)*sizeof(ideal));
2441  sip_sring tmpR;
2442  ring origR = currRing;
2443
2444  if (computeStd)
2445  {
2446    tres[0] = syzstr->res[1];
2447    syzstr->res[1] = idInit(IDELEMS(tres[0]),tres[0]->rank);
2448    return tres;
2449  }
2450  int i,j,l,index,o1,ii,i1;
2451  poly p;
2452  ideal rs;
2453  SSet sPairs;
2454  int * ord,*b0,*b1;
2455  pSetmProc oldSetm=pSetm;
2456  pCompProc oldComp0=pComp0;
2457 
2458  if ((currRing->order[0]==ringorder_dp)
2459  &&  (currRing->order[1]==ringorder_C)
2460  &&  (currRing->order[2]==0))
2461  {
2462    ord=NULL;
2463  }
2464/*--- changes to a dpC-ring with special comp0------*/
2465  else
2466  {
2467    ord = (int*)Alloc0(3*sizeof(int));
2468    b0 = (int*)Alloc0(3*sizeof(int));
2469    b1 = (int*)Alloc0(3*sizeof(int));
2470    ord[0] = ringorder_dp;
2471    ord[1] = ringorder_C;
2472    b0[1] = 1;
2473    b1[0] = pVariables;
2474    tmpR  = *origR;
2475    tmpR.order = ord;
2476    tmpR.block0 = b0;
2477    tmpR.block1 = b1;
2478    rComplete(&tmpR);
2479    //pChangeRing(pVariables,1,ord,b0,b1,currRing->wvhdl);
2480    rChangeCurrRing(&tmpR, TRUE);
2481  }
2482#if ! defined(HAVE_SY_VECTOR) || defined(SY_VEC_DEBUG) 
2483  pSetm =syzSetm;
2484#endif
2485  pComp0 = syzcomp2dpc;
2486  binomials = syzstr->binom;
2487  highdeg_1 = syzstr->highdeg_1;
2488  for (index=syzstr->length-1;index>0;index--)
2489  {
2490    if (syzstr->resPairs[index]!=NULL)
2491    {
2492      currcomponents = syzstr->truecomponents[index];
2493      sPairs = syzstr->resPairs[index];
2494      Strip = syToStrip(syzstr,index);
2495      tres[index+1] = idInit(IDELEMS(syzstr->res[index+1]),syzstr->res[index+1]->rank);
2496      i1 = (*syzstr->Tl)[index];
2497      for (i=0;i<i1;i++)
2498      {
2499        if ((sPairs[i].isNotMinimal==NULL) && (sPairs[i].lcm!=NULL))
2500        {
2501          o1 = sPairs[i].order;
2502          ii = i;
2503          while ((ii<i1) && (sPairs[ii].lcm!=NULL) && (sPairs[ii].order==o1))
2504            ii++;
2505          l = sPairs[i].syzind;
2506          tres[index+1]->m[l] =
2507            syMinimizeP(syzstr->res[index+1]->m[l],syzstr,ii,index,Strip);
2508        }
2509      }
2510      delete Strip;
2511      Strip = NULL;
2512    }
2513  }
2514  currcomponents = syzstr->truecomponents[0];
2515  tres[1] = idInit(IDELEMS(syzstr->res[1]),syzstr->res[1]->rank);
2516  sPairs = syzstr->resPairs[0];
2517  for (i=(*syzstr->Tl)[0]-1;i>=0;i--)
2518  {
2519    if (sPairs[i].syzind>=0)
2520    {
2521      tres[1]->m[sPairs[i].syzind] = pCopy(syzstr->res[1]->m[sPairs[i].syzind]);
2522    }
2523  }
2524/*--- changes to the original ring------------------*/
2525  if (ord!=NULL)
2526  {
2527    //pChangeRing(pVariables,currRing->OrdSgn,currRing->order,
2528    //currRing->block0,currRing->block1,currRing->wvhdl);
2529    rChangeCurrRing(origR,TRUE);
2530  }
2531  else
2532  {
2533    pSetm=oldSetm;
2534    pComp0=oldComp0;
2535  }
2536  tres = syReorder(tres,syzstr->length,syzstr,FALSE,syzstr->res);
2537  syKillEmptyEntres(tres,syzstr->length);
2538  idSkipZeroes(tres[0]);
2539  return tres;
2540}
2541
2542syStrategy syMinimize(syStrategy syzstr)
2543{
2544  if (syzstr->minres==NULL)
2545  {
2546    if (syzstr->resPairs!=NULL)
2547    {
2548      syzstr->minres = syReadOutMinimalRes(syzstr);
2549    }
2550    else if (syzstr->fullres!=NULL)
2551    {
2552      syMinimizeResolvente(syzstr->fullres,syzstr->length,1);
2553      syzstr->minres = syzstr->fullres;
2554      syzstr->fullres = NULL;
2555    }
2556  }
2557  (syzstr->references)++;
2558  return syzstr;
2559}
2560
2561static void sySetHighdeg()
2562{
2563  long long t=1, h_d=1;
2564  long long h_n=1+pVariables;
2565  while ((t=((t*h_n)/h_d))<INT_MAX)
2566  {
2567    h_d++;
2568    h_n++;
2569  }
2570  h_d--;
2571  highdeg = h_d;
2572  //Print("max deg=%d\n",highdeg);
2573}
2574
2575/*2
2576* implementation of LaScala's algorithm
2577* assumes that the given module is homogeneous
2578* works with slanted degree, uses syChosePairs
2579*/
2580syStrategy syLaScala3(ideal arg,int * length)
2581{
2582  BOOLEAN noPair=FALSE;
2583  int i,j,* ord,*b0,*b1,actdeg=32000,index=0,reg=-1;
2584  int startdeg,howmuch;
2585  poly p;
2586  ideal temp;
2587  SSet nextPairs;
2588  syStrategy syzstr=(syStrategy)Alloc0(sizeof(ssyStrategy));
2589  sip_sring tmpR;
2590  ring origR = currRing;
2591
2592  if ((idIs0(arg)) || 
2593      ((idRankFreeModule(arg)>0) && (!idHomModule(arg,NULL,&(syzstr->cw)))))
2594  {
2595    syzstr->minres = (resolvente)Alloc0(sizeof(ideal));
2596    syzstr->length = 1;
2597    syzstr->minres[0] = idInit(1,arg->rank);
2598    return syzstr;
2599  }
2600  pSetmProc oldSetm=pSetm;
2601  pCompProc oldComp0=pComp0;
2602
2603  //crit = 0;
2604  //zeroRed = 0;
2605  //simple = 0;
2606  //dsim = 0;
2607  euler = -1;
2608  redpol = pInit();
2609  //orderingdepth = new intvec(pVariables+1);
2610  syzstr->length = *length = pVariables+2;
2611  if (idIs0(arg)) return NULL;
2612  if ((currRing->order[0]==ringorder_dp)
2613  &&  (currRing->order[1]==ringorder_C)
2614  &&  (currRing->order[2]==0))
2615  {
2616    ord=NULL;
2617  } 
2618/*--- changes to a dpC-ring with special comp0------*/
2619  else
2620  {
2621    ord = (int*)Alloc0(3*sizeof(int));
2622    b0 = (int*)Alloc0(3*sizeof(int));
2623    b1 = (int*)Alloc0(3*sizeof(int));
2624    ord[0] = ringorder_dp;
2625    ord[1] = ringorder_C;
2626    b0[1] = 1;
2627    b1[0] = pVariables;
2628    tmpR  = *origR;
2629    tmpR.order = ord;
2630    tmpR.block0 = b0;
2631    tmpR.block1 = b1;
2632    rComplete(&tmpR);
2633    //pChangeRing(pVariables,1,ord,b0,b1,currRing->wvhdl);
2634    rChangeCurrRing(&tmpR, TRUE);
2635  } 
2636/*--- initializes the data structures---------------*/
2637  syzstr->Tl = new intvec(*length);
2638  temp = idInit(IDELEMS(arg),arg->rank);
2639  for (i=0;i<IDELEMS(arg);i++)
2640  {
2641    temp->m[i] = pOrdPolyInsertSetm(pCopy(arg->m[i]));
2642    if (temp->m[i]!=NULL) 
2643    {
2644      j = pTotaldegree(temp->m[i]);
2645      if (j<actdeg) actdeg = j;
2646    }
2647  }
2648  idTest(temp);
2649  idSkipZeroes(temp);
2650  idTest(temp);
2651  sySetHighdeg();
2652  binomials = (int*)Alloc(pVariables*(highdeg+1)*sizeof(int));
2653  syBinomSet();
2654#if ! defined(HAVE_SY_VECTOR) || defined(SY_VEC_DEBUG) 
2655  pSetm =syzSetm;
2656  for (i=0;i<IDELEMS(temp);i++)
2657  {
2658    p = temp->m[i];
2659    while (p!=NULL)
2660    {
2661      pSetm(p);
2662      pIter(p);
2663    }
2664  }
2665#else
2666  SetSpolyProc();
2667#endif
2668  pComp0 = syzcomp2dpc;
2669  currcomponents = (int*)Alloc0((arg->rank+1)*sizeof(int));
2670  for (i=0;i<=arg->rank;i++)
2671    currcomponents[i] = i;
2672  syzstr->resPairs = syInitRes(temp,length,syzstr->Tl,syzstr->cw);
2673  Free((ADDRESS)currcomponents,(arg->rank+1)*sizeof(int));
2674  syzstr->res = (resolvente)Alloc0((*length+1)*sizeof(ideal));
2675  syzstr->orderedRes = (resolvente)Alloc0((*length+1)*sizeof(ideal));
2676  syzstr->truecomponents = (int**)Alloc0((*length+1)*sizeof(int*));
2677  syzstr->backcomponents = (int**)Alloc0((*length+1)*sizeof(int*));
2678  syzstr->Howmuch = (int**)Alloc0((*length+1)*sizeof(int*));
2679  syzstr->Firstelem = (int**)Alloc0((*length+1)*sizeof(int*));
2680  int len0=idRankFreeModule(arg)+1;
2681  //syzstr->truecomponents[0] = (int*)Alloc(len0*sizeof(int));
2682  //syzstr->backcomponents[0] = (int*)Alloc(len0*sizeof(int));
2683  //syzstr->Howmuch[0] = (int*)Alloc(len0*sizeof(int));
2684  //syzstr->Firstelem[0] = (int*)Alloc(len0*sizeof(int));
2685//Print("sort %d has now size %d\n",0,len0);
2686  startdeg = actdeg;
2687  nextPairs = syChosePairs(syzstr,&index,&howmuch,&actdeg);
2688  //if (TEST_OPT_PROT) Print("(%d,%d)",howmuch,index);
2689/*--- computes the resolution ----------------------*/ 
2690  while (nextPairs!=NULL)
2691  {
2692    if (TEST_OPT_PROT) Print("%d",actdeg);
2693    if (TEST_OPT_PROT) Print("(m%d)",index);
2694    if (index==0)
2695      i = syInitSyzMod(syzstr,index,len0);
2696    else
2697      i = syInitSyzMod(syzstr,index);
2698    currcomponents = syzstr->truecomponents[max(index-1,0)];
2699    j = syInitSyzMod(syzstr,index+1);
2700    if (index>0)
2701    {
2702      syRedNextPairs(nextPairs,syzstr,howmuch,index);
2703      syCompactifyPairSet(syzstr->resPairs[index],(*syzstr->Tl)[index],0);
2704    }
2705    else
2706      syRedGenerOfCurrDeg(syzstr,actdeg,index+1);
2707/*--- creates new pairs -----------------------------*/     
2708    syCreateNewPairs(syzstr,index,i);
2709    if (index<(*length)-1) 
2710    {
2711      syCreateNewPairs(syzstr,index+1,j);
2712      //currcomponents = syzstr->truecomponents[index];
2713      //sySyzTail(syzstr,index+1,j);
2714      //currcomponents = syzstr->truecomponents[index-1];
2715    }
2716    index++;
2717    nextPairs = syChosePairs(syzstr,&index,&howmuch,&actdeg);
2718    //if (TEST_OPT_PROT) Print("(%d,%d)",howmuch,index);
2719  }
2720  if (temp!=NULL) idDelete(&temp);
2721  if (ord!=NULL)
2722  {
2723    //pChangeRing(pVariables,currRing->OrdSgn,currRing->order,
2724    //currRing->block0,currRing->block1,currRing->wvhdl);
2725    rChangeCurrRing(origR,TRUE);
2726  }
2727  else
2728  {
2729    pSetm=oldSetm;
2730    pComp0=oldComp0;
2731  }
2732  if (ord!=NULL)
2733  {
2734    Free((ADDRESS)ord,3*sizeof(int));
2735    Free((ADDRESS)b0,3*sizeof(int));
2736    Free((ADDRESS)b1,3*sizeof(int));
2737  } 
2738  syzstr->binom = binomials;
2739  syzstr->highdeg_1 = highdeg_1;
2740  pDelete1(&redpol);
2741  if (TEST_OPT_PROT) PrintLn();
2742  return syzstr;
2743}
Note: See TracBrowser for help on using the repository browser.