source: git/kernel/tgb.cc @ e9f4c9

spielwiese
Last change on this file since e9f4c9 was e9f4c9, checked in by Michael Brickenstein <bricken@…>, 18 years ago
*bricken: slim_nsize exported git-svn-id: file:///usr/local/Singular/svn/trunk@8990 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 74.6 KB
Line 
1//! \file tgb.cc
2//       multiple rings
3//       shorten_tails und dessen Aufrufe pruefen wlength!!!
4/****************************************
5*  Computer Algebra System SINGULAR     *
6****************************************/
7/* $Id: tgb.cc,v 1.77 2006-02-28 12:02:52 bricken Exp $ */
8/*
9* ABSTRACT: slimgb and F4 implementation
10*/
11//#include <vector>
12//using namespace std;
13#include "mod2.h"
14#include "tgb.h"
15#include "tgb_internal.h"
16#include "tgbgauss.h"
17#include "F4.h"
18#include "digitech.h"
19#include "gring.h"
20
21#include "longrat.h"
22#define SR_HDL(A) ((long)(A))
23static const int bundle_size=100;
24static const int delay_factor=3;
25int QlogSize(number n);
26#if 1
27static omBin lm_bin=NULL;
28
29static void simplify_poly(poly p, ring r) {
30     assume(r==currRing);
31     if (!rField_is_Zp(r)){ 
32        pCleardenom(p);
33        pContent(p); //is a duplicate call, but belongs here
34       
35      }
36      else                     
37        pNorm(p);
38}
39//static const BOOLEAN up_to_radical=TRUE;
40
41int slim_nsize(number n, ring r) {
42    if (rField_is_Zp(r)){
43        return 1;
44    }
45 if (rField_is_Q(r)){
46      return QlogSize(n);
47
48  }
49    else{
50      return n_Size(n,r);
51     
52    }
53}
54static BOOLEAN monomial_root(poly m, ring r){
55    BOOLEAN changed=FALSE;
56    int i;
57    for(i=1;i<=rVar(r);i++){
58        int e=p_GetExp(m,i,r);
59        if (e>1){
60            p_SetExp(m,i,1,r);
61            changed=TRUE;
62        }
63    }
64    if (changed) {
65        p_Setm(m,r);
66    }
67    return changed;
68}
69static BOOLEAN polynomial_root(poly h, ring r){
70  poly got=gcd_of_terms(h,r);
71  BOOLEAN changed=FALSE;
72  if((got!=NULL) &&(TEST_V_UPTORADICAL)) {
73    poly copy=p_Copy(got,r);
74    //p_wrp(got,c->r);
75    changed=monomial_root(got,r);
76    if (changed)
77    {
78         poly div_by=pDivide(copy, got);
79         poly iter=h;
80         while(iter){
81            pExpVectorSub(iter,div_by);
82            pIter(iter);
83         }
84         p_Delete(&div_by, r);
85         if (TEST_OPT_PROT)
86             PrintS("U");
87    }
88    p_Delete(&copy,r);
89  }
90  p_Delete(&got,r);
91  return changed;
92}
93static inline poly p_Init_Special(const ring r)
94{
95  return p_Init(r,lm_bin);
96}
97static inline poly pOne_Special(const ring r=currRing)
98{
99  poly rc = p_Init_Special(r);
100  pSetCoeff0(rc,r->cf->nInit(1));
101  return rc;
102}
103// zum Initialiseren: in t_rep_gb plazieren:
104
105#endif
106#define LEN_VAR3
107#define degbound(p) assume(pTotaldegree(p)<10)
108//#define inDebug(p) assume((debug_Ideal==NULL)||(kNF(debug_Ideal,NULL,p,0,0)==0))
109
110//die meisten Varianten stossen sich an coef_buckets
111
112
113
114#ifdef LEN_VAR1
115// erste Variante: Laenge: Anzahl der Monome
116inline int pSLength(poly p, int l) {
117  return l; }
118inline int kSBucketLength(kBucket* bucket, poly lm) {return bucket_guess(bucket);}
119#endif
120
121#ifdef LEN_VAR2
122// 2. Variante: Laenge: Platz fuer die Koeff.
123int pSLength(poly p,int l)
124{
125  int s=0;
126  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
127  return s;
128}
129int kSBucketLength(kBucket* b, poly lm)
130{
131  int s=0;
132  int i;
133  for (i=MAX_BUCKET;i>=0;i--)
134  {
135    s+=pSLength(b->buckets[i],0);
136  }
137  return s;
138}
139#endif
140
141
142// 3.Variante: Laenge: Platz fuer Leitk * Monomanzahl
143
144
145 
146 
147int QlogSize(number n){
148   
149    if (SR_HDL(n) & SR_INT){
150       long i=SR_TO_INT(n);
151       if (i==0) return 0;
152       
153       unsigned long v;
154       v=(i>=0)?i:-i;
155       int r = 0;
156
157       while (v >>= 1)
158       {
159        r++;
160       }
161       return r+1;
162    }
163    //assume denominator is 0
164    return mpz_sizeinbase(&n->z,2);
165}
166
167
168#ifdef LEN_VAR3
169inline int pSLength(poly p,int l)
170{
171  int c;
172  number coef=pGetCoeff(p);
173  if (rField_is_Q(currRing)){
174    c=QlogSize(coef);
175  }
176  else
177    c=nSize(coef);
178 
179  return c*l /*pLength(p)*/;
180}
181//! TODO CoefBuckets berücksichtigen
182int kSBucketLength(kBucket* b, poly lm=NULL)
183{
184  int s=0;
185  int c;
186  number coef;
187  if(lm==NULL)
188    coef=pGetCoeff(kBucketGetLm(b));
189    //c=nSize(pGetCoeff(kBucketGetLm(b)));
190  else
191    coef=pGetCoeff(lm);
192    //c=nSize(pGetCoeff(lm));
193  if (rField_is_Q(currRing)){
194    c=QlogSize(coef);
195  }
196  else
197    c=nSize(coef);
198 
199  int i;
200  for (i=b->buckets_used;i>=0;i--)
201  {
202    assume((b->buckets_length[i]==0)||(b->buckets[i]!=NULL));
203    s+=b->buckets_length[i] /*pLength(b->buckets[i])*/;
204  }
205  #ifdef HAVE_COEF_BUCKETS
206  assume(b->buckets[0]==kBucketGetLm(b));
207  if (b->coef[0]!=NULL){
208   
209    if (rField_is_Q(currRing)){
210      int modifier=QlogSize(pGetCoeff(b->coef[0]));
211      c+=modifier;
212  }
213    else{
214      int modifier=nSize(pGetCoeff(b->coef[0]));
215      c*=modifier;}
216    }
217  #endif
218  return s*c;
219}
220#endif
221#ifdef LEN_VAR5
222inline wlen_type pSLength(poly p,int l)
223{
224  int c;
225  number coef=pGetCoeff(p);
226  if (rField_is_Q(currRing)){
227    c=QlogSize(coef);
228  }
229  else
230    c=nSize(coef);
231  wlen_type erg=l;
232  erg*=c;
233  erg*=c;
234  //PrintS("lenvar 5");
235  assume(erg>=0);
236  return erg; /*pLength(p)*/;
237}
238//! TODO CoefBuckets berücksichtigen
239wlen_type kSBucketLength(kBucket* b, poly lm=NULL)
240{
241  wlen_type s=0;
242  int c;
243  number coef;
244  if(lm==NULL)
245    coef=pGetCoeff(kBucketGetLm(b));
246    //c=nSize(pGetCoeff(kBucketGetLm(b)));
247  else
248    coef=pGetCoeff(lm);
249    //c=nSize(pGetCoeff(lm));
250  if (rField_is_Q(currRing)){
251    c=QlogSize(coef);
252  }
253  else
254    c=nSize(coef);
255 
256  int i;
257  for (i=b->buckets_used;i>=0;i--)
258  {
259    assume((b->buckets_length[i]==0)||(b->buckets[i]!=NULL));
260    s+=b->buckets_length[i] /*pLength(b->buckets[i])*/;
261  }
262  #ifdef HAVE_COEF_BUCKETS
263  assume(b->buckets[0]==kBucketGetLm(b));
264  if (b->coef[0]!=NULL){
265   
266    if (rField_is_Q(currRing)){
267      int modifier=QlogSize(pGetCoeff(b->coef[0]));
268      c+=modifier;
269  }
270    else{
271      int modifier=nSize(pGetCoeff(b->coef[0]));
272      c*=modifier;}
273    }
274  #endif
275  wlen_type erg=s;
276  erg*=c;
277  erg*=c;
278  return erg;
279}
280#endif
281
282#ifdef LEN_VAR4
283// 4.Variante: Laenge: Platz fuer Leitk * (1+Platz fuer andere Koeff.)
284int pSLength(poly p, int l)
285{
286  int s=1;
287  int c=nSize(pGetCoeff(p));
288  pIter(p);
289  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
290  return s*c;
291}
292int kSBucketLength(kBucket* b)
293{
294  int s=1;
295  int c=nSize(pGetCoeff(kBucketGetLm(b)));
296  int i;
297  for (i=MAX_BUCKET;i>0;i--)
298  {
299    if(b->buckets[i]==NULL) continue;
300    s+=pSLength(b->buckets[i],0);
301  }
302  return s*c;
303}
304#endif
305static int do_pELength(poly p, slimgb_alg* c, int dlm=-1){
306  if(p==NULL) return 0;
307  int s=0;
308  poly pi=p;
309  if(dlm<0){
310    dlm=pTotaldegree(p,c->r);
311    s=1;
312    pi=p->next;
313  }
314 
315  while(pi){
316    int d=pTotaldegree(pi,c->r);
317    if(d>dlm)
318      s+=1+d-dlm;
319    else
320      ++s;
321    pi=pi->next;
322  }
323  return s;
324}
325int kEBucketLength(kBucket* b, poly lm,slimgb_alg* ca)
326{
327  int s=0;
328  if(lm==NULL){
329    lm=kBucketGetLm(b);
330  }
331  if(lm==NULL) return 0;
332  int d=pTotaldegree(lm,ca->r);
333  int i;
334  for (i=b->buckets_used;i>=0;i--)
335  {
336    if(b->buckets[i]==NULL) continue;
337    s+=do_pELength(b->buckets[i],ca,d);
338  }
339  return s;
340}
341
342static inline int pELength(poly p, slimgb_alg* c,int l){
343  if (p==NULL) return 0;
344  return do_pELength(p,c);
345}
346
347
348
349
350static inline wlen_type pQuality(poly p, slimgb_alg* c, int l=-1){
351 
352  if(l<0)
353    l=pLength(p);
354  if(c->is_char0) {
355    if((pLexOrder) &&(!c->is_homog)){
356      int cs;
357      number coef=pGetCoeff(p);
358      if (rField_is_Q(currRing)){
359         cs=QlogSize(coef);
360      }
361      else
362  cs=nSize(coef);
363     wlen_type erg=cs;
364     //erg*=cs;//for quadratic
365     erg*=pELength(p,c,l);
366    //FIXME: not quadratic coeff size
367      //return cs*pELength(p,c,l);
368      return erg;
369    }
370    //Print("I am here");
371    wlen_type r=pSLength(p,l);
372    assume(r>=0);
373    return r;
374  }
375  if((pLexOrder) &&(!c->is_homog)) return pELength(p,c,l);
376  return l;
377}
378
379static inline int pTotaldegree_full(poly p){
380  int r=0;
381  while(p){
382    int d=pTotaldegree(p);
383    r=si_max(r,d);
384    pIter(p);
385  }
386  return r;
387}
388
389wlen_type red_object::guess_quality(slimgb_alg* c){
390    //works at the moment only for lenvar 1, because in different
391    //case, you have to look on coefs
392    wlen_type s=0;
393    if (c->is_char0){
394      //s=kSBucketLength(bucket,this->p);
395      if((pLexOrder) &&(!c->is_homog)){
396    int cs;
397    number coef;
398   
399    coef=pGetCoeff(kBucketGetLm(bucket));
400    //c=nSize(pGetCoeff(kBucketGetLm(b)));
401   
402    //c=nSize(pGetCoeff(lm));
403    if (rField_is_Q(currRing)){
404      cs=QlogSize(coef);
405    }
406    else
407      cs=nSize(coef);
408    #ifdef HAVE_COEF_BUCKETS
409    if (bucket->coef[0]!=NULL){
410      if (rField_is_Q(currRing)){
411        int modifier=QlogSize(pGetCoeff(bucket->coef[0]));
412        cs+=modifier;
413      }
414      else{
415        int modifier=nSize(pGetCoeff(bucket->coef[0]));
416        cs*=modifier;}
417    }
418    #endif
419    //FIXME:not quadratic
420    wlen_type erg=kEBucketLength(this->bucket,this->p,c);
421    //erg*=cs;//for quadratic
422    erg*=cs;
423    //return cs*kEBucketLength(this->bucket,this->p,c);
424    return erg;
425      }
426      s=kSBucketLength(bucket,NULL);
427    }
428    else 
429    {
430      if((pLexOrder) &&(!c->is_homog))
431  //if (false)
432  s=kEBucketLength(this->bucket,this->p,c);
433      else s=bucket_guess(bucket);
434    }
435 
436    return s;
437}
438
439
440
441static void finalize_reduction_step(reduction_step* r){
442  delete r;
443}
444static int LObject_better_gen(const void* ap, const void* bp)
445{
446  LObject* a=*(LObject**)ap;
447  LObject* b=*(LObject**)bp;
448  return(pLmCmp(a->p,b->p));
449}
450static int red_object_better_gen(const void* ap, const void* bp)
451{
452
453
454  return(pLmCmp(((red_object*) ap)->p,((red_object*) bp)->p));
455}
456
457
458static int pLmCmp_func_inverted(const void* ap1, const void* ap2){
459    poly p1,p2;
460  p1=*((poly*) ap1);
461  p2=*((poly*)ap2);
462
463  return -pLmCmp(p1,p2);
464}
465
466int tgb_pair_better_gen2(const void* ap,const void* bp){
467  return(-tgb_pair_better_gen(ap,bp));
468}
469int kFindDivisibleByInS_easy(kStrategy strat,const red_object & obj){
470  int i;
471  long not_sev=~obj.sev;
472  poly p=obj.p;
473  for(i=0;i<=strat->sl;i++){
474    if (pLmShortDivisibleBy(strat->S[i],strat->sevS[i],p,not_sev))
475      return i;
476  }
477  return -1;
478}
479int kFindDivisibleByInS_easy(kStrategy strat,poly p, long sev){
480  int i;
481  long not_sev=~sev;
482  for(i=0;i<=strat->sl;i++){
483    if (pLmShortDivisibleBy(strat->S[i],strat->sevS[i],p,not_sev))
484      return i;
485  }
486  return -1;
487}
488static int posInPairs (sorted_pair_node**  p, int pn, sorted_pair_node* qe,slimgb_alg* c,int an=0)
489{
490  if(pn==0) return 0;
491
492  int length=pn-1;
493  int i;
494  //int an = 0;
495  int en= length;
496
497  if (pair_better(qe,p[en],c))
498    return length+1;
499
500  while(1)
501    {
502      //if (an >= en-1)
503      if(en-1<=an)
504      {
505        if (pair_better(p[an],qe,c)) return an;
506        return en;
507      }
508      i=(an+en) / 2;
509        if (pair_better(p[i],qe,c))
510          en=i;
511      else an=i;
512    }
513}
514
515static BOOLEAN  ascending(int* i,int top){
516  if(top<1) return TRUE;
517  if(i[top]<i[top-1]) return FALSE;
518  return ascending(i,top-1);
519}
520
521sorted_pair_node**  spn_merge(sorted_pair_node** p, int pn,sorted_pair_node **q, int qn,slimgb_alg* c){
522  int i;
523  int* a= (int*) omalloc(qn*sizeof(int));
524//   int mc;
525//   PrintS("Debug\n");
526//   for(mc=0;mc<qn;mc++)
527// {
528
529//     wrp(q[mc]->lcm_of_lm);
530//     PrintS("\n");
531// }
532//    PrintS("Debug they are in\n");
533//   for(mc=0;mc<pn;mc++)
534// {
535
536//     wrp(p[mc]->lcm_of_lm);
537//     PrintS("\n");
538// }
539  int lastpos=0;
540  for(i=0;i<qn;i++){
541    lastpos=posInPairs(p,pn,q[i],c, si_max(lastpos-1,0));
542    //   cout<<lastpos<<"\n";
543    a[i]=lastpos;
544
545  }
546  if((pn+qn)>c->max_pairs){
547    p=(sorted_pair_node**) omrealloc(p,2*(pn+qn)*sizeof(sorted_pair_node*));
548    c->max_pairs=2*(pn+qn);
549  }
550  for(i=qn-1;i>=0;i--){
551    size_t size;
552    if(qn-1>i)
553      size=(a[i+1]-a[i])*sizeof(sorted_pair_node*);
554    else
555      size=(pn-a[i])*sizeof(sorted_pair_node*); //as indices begin with 0
556    memmove (p+a[i]+(1+i), p+a[i], size);
557    p[a[i]+i]=q[i];
558  }
559  omfree(a);
560  return p;
561}
562
563
564static BOOLEAN trivial_syzygie(int pos1,int pos2,poly bound,slimgb_alg* c){
565
566
567  poly p1=c->S->m[pos1];
568  poly p2=c->S->m[pos2];
569   
570
571  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
572    return FALSE;
573  int i = 1;
574  poly m=NULL;
575  poly gcd1=c->gcd_of_terms[pos1];
576  poly gcd2=c->gcd_of_terms[pos2];
577 
578  if((gcd1!=NULL) && (gcd2!=NULL)) 
579    {
580      gcd1->next=gcd2; //may ordered incorrect
581      m=gcd_of_terms(gcd1,c->r);
582      gcd1->next=NULL;
583     
584    } 
585
586  if (m==NULL) 
587  {
588     loop
589      {
590  if (pGetExp(p1, i)+ pGetExp(p2, i) > pGetExp(bound,i))   return FALSE;
591  if (i == pVariables){
592    //PrintS("trivial");
593    return TRUE;
594  }
595  i++;
596      }
597  }
598  else 
599  {
600    loop
601      {
602  if (pGetExp(p1, i)-pGetExp(m,i) + pGetExp(p2, i) > pGetExp(bound,i))   return FALSE;
603  if (i == pVariables){
604    pDelete(&m);
605    //PrintS("trivial");
606    return TRUE;
607  }
608  i++;
609      }
610  }
611
612 
613
614 
615}
616
617//! returns position sets w as weight
618int find_best(red_object* r,int l, int u, int &w, slimgb_alg* c){
619  int best=l;
620  int i;
621  w=r[l].guess_quality(c);
622  for(i=l+1;i<=u;i++){
623    int w2=r[i].guess_quality(c);
624    if(w2<w){
625      w=w2;
626      best=i;
627    }
628   
629  }
630 return best;
631}
632
633
634void red_object::canonicalize(){
635  kBucketCanonicalize(bucket);
636 
637 
638}
639BOOLEAN good_has_t_rep(int i, int j,slimgb_alg* c){
640  assume(i>=0);
641    assume(j>=0);
642  if (has_t_rep(i,j,c)) return TRUE;
643  //poly lm=pOne();
644  assume (c->tmp_lm!=NULL);
645  poly lm=c->tmp_lm;
646
647  pLcm(c->S->m[i], c->S->m[j], lm);
648  pSetm(lm);
649  assume(lm!=NULL);
650  //int deciding_deg= pTotaldegree(lm);
651  int* i_con =make_connections(i,j,lm,c);
652  //p_Delete(&lm,c->r);
653
654
655  for (int n=0;((n<c->n) && (i_con[n]>=0));n++){
656    if (i_con[n]==j){
657      now_t_rep(i,j,c);
658      omfree(i_con);
659
660      return TRUE;
661    }
662  }
663  omfree(i_con);
664
665  return FALSE;
666}
667BOOLEAN lenS_correct(kStrategy strat){
668  int i;
669  for(i=0;i<=strat->sl;i++){
670    if (strat->lenS[i]!=pLength(strat->S[i]))
671      return FALSE;
672  }
673  return TRUE;
674}
675
676
677static void cleanS(kStrategy strat, slimgb_alg* c){
678  int i=0;
679  LObject P;
680  while(i<=strat->sl){
681    P.p=strat->S[i];
682    P.sev=strat->sevS[i];
683    if(kFindDivisibleByInS(strat->S,strat->sevS,strat->sl,&P)!=i)
684    {
685      deleteInS(i,strat);
686      //remember destroying poly
687      BOOLEAN found=FALSE;
688      int j;
689      for(j=0;j<c->n;j++)
690  if(c->S->m[j]==P.p)
691  {
692    found=TRUE;
693    break;
694  }
695      if (!found)
696  pDelete(&P.p);
697      //remember additional reductors
698    }
699    else i++;
700  }
701}
702static int bucket_guess(kBucket* bucket){
703  int sum=0;
704  int i;
705  for (i=bucket->buckets_used;i>=0;i--){
706    if(bucket->buckets[i])
707       sum+=bucket->buckets_length[i];
708  }
709  return sum;
710}
711
712
713
714
715
716
717static int add_to_reductors(slimgb_alg* c, poly h, int len, BOOLEAN simplified){
718  //inDebug(h);
719  assume(lenS_correct(c->strat));
720  assume(len==pLength(h));
721  int i;
722//   if (c->is_char0)
723//        i=simple_posInS(c->strat,h,pSLength(h,len),c->is_char0);
724//   else
725//     i=simple_posInS(c->strat,h,len,c->is_char0);
726
727  LObject P; memset(&P,0,sizeof(P));
728  P.tailRing=c->r;
729  P.p=h; /*p_Copy(h,c->r);*/
730  P.FDeg=pFDeg(P.p,c->r);
731  if (!(simplified)){
732      if (!rField_is_Zp(c->r)){ 
733        pCleardenom(P.p);
734        pContent(P.p); //is a duplicate call, but belongs here
735       
736      }
737      else                     
738        pNorm(P.p);
739    pNormalize(P.p);
740  }
741  wlen_type pq=pQuality(h,c,len);
742  i=simple_posInS(c->strat,h,len,pq);
743  c->strat->enterS(P,i,c->strat,-1);
744 
745 
746
747  c->strat->lenS[i]=len;
748  assume(pLength(c->strat->S[i])==c->strat->lenS[i]);
749  if(c->strat->lenSw!=NULL)
750    c->strat->lenSw[i]=pq;
751 
752  return i;
753 
754}
755static void length_one_crit(slimgb_alg* c, int pos, int len)
756{
757  if (c->nc)
758    return;
759  if (len==1)
760  {
761    int i;
762    for ( i=0;i<pos;i++)
763    {
764      if (c->lengths[i]==1)
765        c->states[pos][i]=HASTREP;
766    }
767    for ( i=pos+1;i<c->n;i++){
768      if (c->lengths[i]==1)
769        c->states[i][pos]=HASTREP;
770    }
771    if (!c->nc)
772      shorten_tails(c,c->S->m[pos]);
773  }
774}
775
776
777static void move_forward_in_S(int old_pos, int new_pos,kStrategy strat)
778{
779  assume(old_pos>=new_pos);
780  poly p=strat->S[old_pos];
781  int ecart=strat->ecartS[old_pos];
782  long sev=strat->sevS[old_pos];
783  int s_2_r=strat->S_2_R[old_pos];
784  int length=strat->lenS[old_pos];
785  assume(length==pLength(strat->S[old_pos]));
786  wlen_type length_w;
787  if(strat->lenSw!=NULL)
788    length_w=strat->lenSw[old_pos];
789  int i;
790  for (i=old_pos; i>new_pos; i--)
791  {
792    strat->S[i] = strat->S[i-1];
793    strat->ecartS[i] = strat->ecartS[i-1];
794    strat->sevS[i] = strat->sevS[i-1];
795    strat->S_2_R[i] = strat->S_2_R[i-1];
796  }
797  if (strat->lenS!=NULL)
798    for (i=old_pos; i>new_pos; i--)
799      strat->lenS[i] = strat->lenS[i-1];
800  if (strat->lenSw!=NULL)
801    for (i=old_pos; i>new_pos; i--)
802      strat->lenSw[i] = strat->lenSw[i-1];
803
804  strat->S[new_pos]=p;
805  strat->ecartS[new_pos]=ecart;
806  strat->sevS[new_pos]=sev;
807  strat->S_2_R[new_pos]=s_2_r;
808  strat->lenS[new_pos]=length;
809  if(strat->lenSw!=NULL)
810    strat->lenSw[new_pos]=length_w;
811  //assume(lenS_correct(strat));
812}
813
814static int* make_connections(int from, int to, poly bound, slimgb_alg* c)
815{
816  ideal I=c->S;
817  int* cans=(int*) omalloc(c->n*sizeof(int));
818  int* connected=(int*) omalloc(c->n*sizeof(int));
819  cans[0]=to;
820  int cans_length=1;
821  connected[0]=from;
822  int last_cans_pos=-1;
823  int connected_length=1;
824  long neg_bounds_short= ~p_GetShortExpVector(bound,c->r);
825
826  int not_yet_found=cans_length;
827  int con_checked=0;
828  int pos;
829 
830  while(TRUE){
831    if ((con_checked<connected_length)&& (not_yet_found>0)){
832      pos=connected[con_checked];
833      for(int i=0;i<cans_length;i++){
834        if (cans[i]<0) continue;
835        //FIXME: triv. syz. does not hold on noncommutative, check it for modules
836        if ((has_t_rep(pos,cans[i],c)) ||((!rIsPluralRing(c->r))&&(trivial_syzygie(pos,cans[i],bound,c))))
837{
838
839          connected[connected_length]=cans[i];
840          connected_length++;
841          cans[i]=-1;
842          --not_yet_found;
843
844          if (connected[connected_length-1]==to){
845            if (connected_length<c->n){
846              connected[connected_length]=-1;
847            }
848            omfree(cans);
849            return connected;
850          }
851        }
852      }
853      con_checked++;
854    }
855    else
856    {
857      for(last_cans_pos++;last_cans_pos<=c->n;last_cans_pos++){
858        if (last_cans_pos==c->n){
859          if (connected_length<c->n){
860            connected[connected_length]=-1;
861          }
862          omfree(cans);
863          return connected;
864        }
865        if ((last_cans_pos==from)||(last_cans_pos==to))
866          continue;
867        if(p_LmShortDivisibleBy(I->m[last_cans_pos],c->short_Exps[last_cans_pos],bound,neg_bounds_short,c->r)){
868          cans[cans_length]=last_cans_pos;
869          cans_length++;
870          break;
871        }
872      }
873      not_yet_found++;
874      for (int i=0;i<con_checked;i++){
875        if (has_t_rep(connected[i],last_cans_pos,c)){
876
877          connected[connected_length]=last_cans_pos;
878          connected_length++;
879          cans[cans_length-1]=-1;
880
881          --not_yet_found;
882          if (connected[connected_length-1]==to){
883            if (connected_length<c->n){
884              connected[connected_length]=-1;
885            }
886
887            omfree(cans);
888            return connected;
889          }
890          break;
891        }
892      }
893    }
894  }
895  if (connected_length<c->n){
896    connected[connected_length]=-1;
897  }
898
899  omfree(cans);
900  return connected;
901}
902#ifdef HEAD_BIN
903static inline poly p_MoveHead(poly p, omBin b)
904{
905  poly np;
906  omTypeAllocBin(poly, np, b);
907  memmove(np, p, b->sizeW*sizeof(long));
908  omFreeBinAddr(p);
909  return np;
910}
911#endif
912
913
914
915static void add_later(poly p, char* prot, slimgb_alg* c){
916    int i=0;
917    //check, if it is already in the queue
918   
919   
920    while(c->add_later->m[i]!=NULL){
921        if (p_LmEqual(c->add_later->m[i],p,c->r))
922            return;
923        i++;
924    }
925    if (TEST_OPT_PROT)
926        PrintS(prot);
927    c->add_later->m[i]=p;
928}
929static int simple_posInS (kStrategy strat, poly p,int len, wlen_type wlen)
930{
931
932
933  if(strat->sl==-1) return 0;
934  if (strat->lenSw) return pos_helper(strat,p,(wlen_type) wlen,(wlen_set) strat->lenSw,strat->S);
935  return pos_helper(strat,p,len,strat->lenS,strat->S);
936
937}
938
939/*2
940 *if the leading term of p
941 *divides the leading term of some S[i] it will be canceled
942 */
943static inline void clearS (poly p, unsigned long p_sev,int l, int* at, int* k,
944                           kStrategy strat)
945{
946  assume(p_sev == pGetShortExpVector(p));
947  if (!pLmShortDivisibleBy(p,p_sev, strat->S[*at], ~ strat->sevS[*at])) return;
948  if (l>=strat->lenS[*at]) return;
949  if (TEST_OPT_PROT)
950    PrintS("!");
951  mflush();
952  //pDelete(&strat->S[*at]);
953  deleteInS((*at),strat);
954  (*at)--;
955  (*k)--;
956//  assume(lenS_correct(strat));
957}
958
959
960
961static int iq_crit(const void* ap,const void* bp){
962
963  sorted_pair_node* a=*((sorted_pair_node**)ap);
964  sorted_pair_node* b=*((sorted_pair_node**)bp);
965  assume(a->i>a->j);
966  assume(b->i>b->j);
967 
968 
969  if (a->deg<b->deg) return -1;
970  if (a->deg>b->deg) return 1;
971  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
972  if(comp!=0)
973    return comp;
974  if (a->expected_length<b->expected_length) return -1;
975  if (a->expected_length>b->expected_length) return 1;
976  if (a->j>b->j) return 1;
977  if (a->j<b->j) return -1;
978  return 0;
979}
980static wlen_type coeff_mult_size_estimate(int s1, int s2, ring r){
981    if (rField_is_Q(r)) return s1+s2;
982    else return s1*s2;
983}
984static wlen_type pair_weighted_length(int i, int j, slimgb_alg* c){
985    if ((c->is_char0) && (pLexOrder) &&(!c->is_homog))  {
986        int c1=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r);
987        int c2=slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
988        wlen_type el1=c->weighted_lengths[i]/c1;
989        assume(el1!=0);
990        assume(c->weighted_lengths[i] %c1==0);
991        wlen_type el2=c->weighted_lengths[j]/c2;
992        assume(el2!=0);
993        assume(c->weighted_lengths[j] %c2==0);
994        //should be * for function fields
995        //return (c1+c2) * (el1+el2-2);
996        wlen_type res=coeff_mult_size_estimate(c1,c2,c->r);
997        res*=el1+el2-2;
998        return res;
999       
1000    }
1001    if (c->is_char0) {
1002        //int cs=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r)+
1003        //    slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1004        wlen_type cs=
1005            coeff_mult_size_estimate(
1006                slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r),
1007                slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r),c->r);
1008        return (wlen_type)(c->lengths[i]+c->lengths[j]-2)*
1009            (wlen_type)cs;
1010    }
1011    if ((pLexOrder) &&(!c->is_homog)) {
1012
1013        return (c->weighted_lengths[i]+c->weighted_lengths[j]-2);
1014    }
1015    return c->lengths[i]+c->lengths[j]-2;
1016   
1017}
1018sorted_pair_node** add_to_basis_ideal_quotient(poly h, int i_pos, int j_pos,slimgb_alg* c, int* ip)
1019{
1020
1021  assume(h!=NULL);
1022  poly got=gcd_of_terms(h,c->r);
1023  if((got!=NULL) &&(TEST_V_UPTORADICAL)) {
1024    poly copy=p_Copy(got,c->r);
1025    //p_wrp(got,c->r);
1026    BOOLEAN changed=monomial_root(got,c->r);
1027    if (changed)
1028    {
1029         poly div_by=pDivide(copy, got);
1030         poly iter=h;
1031         while(iter){
1032            pExpVectorSub(iter,div_by);
1033            pIter(iter);
1034         }
1035         p_Delete(&div_by, c->r);
1036         PrintS("U");
1037    }
1038    p_Delete(&copy,c->r);
1039  }
1040
1041#define ENLARGE(pointer, type) pointer=(type*) omrealloc(pointer, c->array_lengths*sizeof(type))
1042//  BOOLEAN corr=lenS_correct(c->strat);
1043  BOOLEAN R_found=FALSE;
1044  void* hp;
1045
1046  ++(c->n);
1047  ++(c->S->ncols);
1048  int i,j;
1049  i=c->n-1;
1050  sorted_pair_node** nodes=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1051  int spc=0;
1052  if(c->n>c->array_lengths){
1053    c->array_lengths=c->array_lengths*2;
1054    assume(c->array_lengths>=c->n);
1055    ENLARGE(c->T_deg, int);
1056    ENLARGE(c->tmp_pair_lm,poly);
1057    ENLARGE(c->tmp_spn,sorted_pair_node*);
1058
1059    ENLARGE(c->short_Exps,long);
1060    ENLARGE(c->lengths,int);
1061    #ifndef HAVE_BOOST
1062    ENLARGE(c->states, char*);
1063    #endif
1064    ENLARGE(c->gcd_of_terms,poly);
1065    //if (c->weighted_lengths!=NULL) {
1066    ENLARGE(c->weighted_lengths,wlen_type);
1067    //}
1068    //ENLARGE(c->S->m,poly);
1069   
1070  }
1071  pEnlargeSet(&c->S->m,c->n-1,1);
1072  if (c->T_deg_full)
1073    ENLARGE(c->T_deg_full,int);
1074  c->T_deg[i]=pTotaldegree(h);
1075  if(c->T_deg_full){
1076    c->T_deg_full[i]=pTotaldegree_full(h);
1077  }
1078 
1079
1080  c->tmp_pair_lm[i]=pOne_Special(c->r);
1081
1082
1083  c->tmp_spn[i]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1084
1085
1086  c->lengths[i]=pLength(h);
1087 
1088  //necessary for correct weighted length
1089 
1090   if (!rField_is_Zp(c->r)){ 
1091    pCleardenom(h);
1092    pContent(h); //is a duplicate call, but belongs here
1093   
1094  }
1095  else                     
1096    pNorm(h);
1097  pNormalize(h);
1098 
1099  c->weighted_lengths[i]=pQuality(h, c, c->lengths[i]);
1100  c->gcd_of_terms[i]=got;
1101  #ifdef HAVE_BOOST
1102    c->states.push_back(dynamic_bitset<>(i));
1103 
1104  #else
1105  if (i>0)
1106    c->states[i]=(char*)  omalloc(i*sizeof(char));
1107  else
1108    c->states[i]=NULL;
1109  #endif
1110 
1111  c->S->m[i]=h;
1112  c->short_Exps[i]=p_GetShortExpVector(h,c->r);
1113 
1114#undef ENLARGE
1115  for (j=0;j<i;j++){
1116   
1117
1118    #ifndef HAVE_BOOST
1119    c->states[i][j]=UNCALCULATED;
1120    #endif
1121    assume(p_LmDivisibleBy(c->S->m[i],c->S->m[j],c->r)==
1122     p_LmShortDivisibleBy(c->S->m[i],c->short_Exps[i],c->S->m[j],~(c->short_Exps[j]),c->r));
1123
1124   
1125    if (_p_GetComp(c->S->m[i],c->r)!=_p_GetComp(c->S->m[j],c->r)){
1126      //c->states[i][j]=UNCALCULATED;
1127      //WARNUNG: be careful
1128      continue;
1129    } else
1130    if ((!c->nc) && (c->lengths[i]==1) && (c->lengths[j]==1)){
1131      c->states[i][j]=HASTREP;
1132     
1133      }
1134    else if ((!(c->nc)) &&  (pHasNotCF(c->S->m[i],c->S->m[j])))
1135    {
1136      c->easy_product_crit++;
1137      c->states[i][j]=HASTREP;
1138      continue;
1139    }
1140    else if(extended_product_criterion(c->S->m[i],c->gcd_of_terms[i],c->S->m[j],c->gcd_of_terms[j],c))
1141    {
1142      c->states[i][j]=HASTREP;
1143      c->extended_product_crit++;
1144     
1145      //PrintS("E");
1146    }
1147      //  if (c->states[i][j]==UNCALCULATED){
1148
1149    if ((TEST_V_FINDMONOM) &&(!c->nc)) {
1150        //PrintS("COMMU");
1151       //  if (c->lengths[i]==c->lengths[j]){
1152//             poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1153//             if (short_s==NULL){
1154//                 c->states[i][j]=HASTREP;
1155//             } else
1156//             {
1157//                 p_Delete(&short_s, currRing);
1158//             }
1159//         }
1160        if (c->lengths[i]+c->lengths[j]==3){
1161            poly                 short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1162            if (short_s==NULL){
1163                c->states[i][j]=HASTREP;
1164            } else
1165            {
1166                assume(pLength(short_s)==1);
1167                if (TEST_V_UPTORADICAL)
1168                   monomial_root(short_s,c->r);
1169                int iS=
1170                   kFindDivisibleByInS_easy(c->strat,short_s, p_GetShortExpVector(short_s,c->r));
1171                if (iS<0){
1172                    //PrintS("N");
1173                    c->states[i][j]=HASTREP;
1174                    add_later(short_s,"N",c);
1175                }
1176                else {
1177                    if (c->strat->lenS[iS]>1){
1178                        //PrintS("O");
1179                        c->states[i][j]=HASTREP;
1180                        add_later(short_s,"O",c);
1181                    }
1182                    else
1183                     p_Delete(&short_s, currRing);
1184                }
1185               
1186               
1187            }
1188        }
1189    }
1190      //    if (short_s)
1191      //    {
1192    assume(spc<=j);
1193    sorted_pair_node* s=c->tmp_spn[spc];//(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1194    s->i=si_max(i,j);
1195    s->j=si_min(i,j);
1196    assume(s->j==j);
1197    s->expected_length=pair_weighted_length(i,j,c);//c->lengths[i]+c->lengths[j]-2;
1198     
1199    poly lm=c->tmp_pair_lm[spc];//=pOne_Special();
1200     
1201    pLcm(c->S->m[i], c->S->m[j], lm);
1202    pSetm(lm);
1203    s->deg=pTotaldegree(lm);
1204
1205    if(c->T_deg_full)//Sugar
1206    {
1207      int t_i=c->T_deg_full[s->i]-c->T_deg[s->i];
1208      int t_j=c->T_deg_full[s->j]-c->T_deg[s->j];
1209      s->deg+=si_max(t_i,t_j);
1210      //Print("\n max: %d\n",max(t_i,t_j));
1211    }
1212    s->lcm_of_lm=lm;
1213    //          pDelete(&short_s);
1214    //assume(lm!=NULL);
1215    nodes[spc]=s;
1216    spc++;
1217 
1218  // }
1219  //else
1220  //{
1221        //c->states[i][j]=HASTREP;
1222  //}
1223  }
1224 
1225  assume(spc<=i);
1226  //now ideal quotient crit
1227  qsort(nodes,spc,sizeof(sorted_pair_node*),iq_crit);
1228 
1229    sorted_pair_node** nodes_final=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1230  int spc_final=0;
1231  j=0;
1232  while(j<spc)
1233  {
1234    int lower=j;
1235    int upper;
1236    BOOLEAN has=FALSE;
1237    for(upper=lower+1;upper<spc;upper++)
1238    {
1239     
1240      if(!pLmEqual(nodes[lower]->lcm_of_lm,nodes[upper]->lcm_of_lm))
1241      {
1242  break;
1243      }
1244      if (has_t_rep(nodes[upper]->i,nodes[upper]->j,c))
1245  has=TRUE;
1246
1247    }
1248    upper=upper-1;
1249    int z;
1250    assume(spc_final<=j);
1251    for(z=0;z<spc_final;z++)
1252    {
1253      if(p_LmDivisibleBy(nodes_final[z]->lcm_of_lm,nodes[lower]->lcm_of_lm,c->r))
1254      {
1255  has=TRUE;
1256  break;
1257      }
1258    }
1259   
1260    if(has)
1261    {
1262      for(;lower<=upper;lower++)
1263      {
1264  //free_sorted_pair_node(nodes[lower],c->r);
1265  //omfree(nodes[lower]);
1266  nodes[lower]=NULL;
1267      }
1268      j=upper+1;
1269      continue;
1270    }
1271    else
1272    {
1273      nodes[lower]->lcm_of_lm=pCopy(nodes[lower]->lcm_of_lm);
1274      assume(_p_GetComp(c->S->m[nodes[lower]->i],c->r)==_p_GetComp(c->S->m[nodes[lower]->j],c->r));
1275      nodes_final[spc_final]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1276     
1277      *(nodes_final[spc_final++])=*(nodes[lower]);
1278      //c->tmp_spn[nodes[lower]->j]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1279      nodes[lower]=NULL;
1280      for(lower=lower+1;lower<=upper;lower++)
1281      {
1282  //  free_sorted_pair_node(nodes[lower],c->r);
1283  //omfree(nodes[lower]);
1284  nodes[lower]=NULL;
1285      }
1286      j=upper+1;
1287      continue;
1288    }
1289  }
1290
1291  //  Print("i:%d,spc_final:%d",i,spc_final);
1292
1293
1294
1295
1296  assume(spc_final<=spc);
1297  omfree(nodes);
1298  nodes=NULL;
1299
1300  add_to_reductors(c, h, c->lengths[c->n-1], TRUE);
1301  //i=posInS(c->strat,c->strat->sl,h,0 ecart);
1302  if (!(c->nc)){
1303    if (c->lengths[c->n-1]==1)
1304      shorten_tails(c,c->S->m[c->n-1]);
1305  }
1306  //you should really update c->lengths, c->strat->lenS, and the oder of polys in strat if you sort after lengths
1307
1308  //for(i=c->strat->sl; i>0;i--)
1309  //  if(c->strat->lenS[i]<c->strat->lenS[i-1]) printf("fehler bei %d\n",i);
1310  if (c->Rcounter>50) {
1311    c->Rcounter=0;
1312    cleanS(c->strat,c);
1313  }
1314  if(!ip){
1315    qsort(nodes_final,spc_final,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1316 
1317   
1318    c->apairs=spn_merge(c->apairs,c->pair_top+1,nodes_final,spc_final,c);
1319    c->pair_top+=spc_final;
1320    clean_top_of_pair_list(c);
1321    omfree(nodes_final);
1322    return NULL;
1323  }
1324  {
1325    *ip=spc_final;
1326    return nodes_final;
1327  }
1328
1329 
1330
1331}
1332
1333
1334static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
1335{
1336  m=nInit(1);
1337  if (h==NULL) return NULL;
1338
1339  assume(len==pLength(h));
1340  kStrategy strat=c->strat;
1341  if (0 > strat->sl)
1342  {
1343    return h;
1344  }
1345  int j;
1346 
1347  LObject P(h);
1348  P.SetShortExpVector();
1349  P.bucket = kBucketCreate(currRing);
1350  // BOOLEAN corr=lenS_correct(strat);
1351  kBucketInit(P.bucket,P.p,len /*pLength(P.p)*/);
1352  //wlen_set lenSw=(wlen_set) c->strat->lenS;
1353  //FIXME: plainly wrong
1354  //strat->lenS;
1355  //if (strat->lenSw!=NULL)
1356  //  lenSw=strat->lenSw;
1357  //int max_pos=simple_posInS(strat,P.p);
1358  loop
1359    {
1360
1361      j=kFindDivisibleByInS(strat->S,strat->sevS,strat->sl,&P);
1362      if ((j>=0) && ((!n)||
1363        ((strat->lenS[j]<=n) &&
1364         ((strat->lenSw==NULL)||
1365         (strat->lenSw[j]<=n)))))
1366      {
1367
1368 
1369        nNormalize(pGetCoeff(P.p));
1370#ifdef KDEBUG
1371        if (TEST_OPT_DEBUG)
1372        {
1373          PrintS("red:");
1374          wrp(h);
1375          PrintS(" with ");
1376          wrp(strat->S[j]);
1377        }
1378#endif
1379       
1380        number coef=kBucketPolyRed(P.bucket,strat->S[j],
1381                                   strat->lenS[j]/*pLength(strat->S[j])*/,
1382                                   strat->kNoether);
1383  number m2=nMult(m,coef);
1384  nDelete(&m);
1385  m=m2;
1386        nDelete(&coef);
1387        h = kBucketGetLm(P.bucket);
1388 
1389  if (h==NULL) {
1390    len=0;
1391    kBucketDestroy(&P.bucket);
1392    return 
1393    NULL;}
1394        P.p=h;
1395        P.t_p=NULL;
1396        P.SetShortExpVector();
1397#ifdef KDEBUG
1398        if (TEST_OPT_DEBUG)
1399        {
1400          PrintS("\nto:");
1401          wrp(h);
1402          PrintLn();
1403        }
1404#endif
1405      }
1406      else
1407      {
1408        kBucketClear(P.bucket,&(P.p),&len);
1409        kBucketDestroy(&P.bucket);
1410        pNormalize(P.p);
1411  assume(len==(pLength(P.p)));
1412        return P.p;
1413      }
1414    }
1415}
1416
1417
1418
1419static poly redTailShort(poly h, kStrategy strat){
1420  if (h==NULL) return NULL;//n_Init(1,currRing);
1421  if (TEST_V_MODPSOLVSB){
1422    bit_reduce(pNext(h), strat->tailRing);
1423  }
1424  int sl=strat->sl;
1425  int i;
1426  int len=pLength(h);
1427  for(i=0;i<=strat->sl;i++){
1428    if((strat->lenS[i]>2) || ((strat->lenSw!=NULL) && (strat->lenSw[i]>2)))
1429      break;
1430  }
1431  return(redNFTail(h,i-1,strat, len));
1432}
1433
1434static void line_of_extended_prod(int fixpos,slimgb_alg* c){
1435    if (c->gcd_of_terms[fixpos]==NULL)
1436  {
1437    c->gcd_of_terms[fixpos]=gcd_of_terms(c->S->m[fixpos],c->r);
1438    if (c->gcd_of_terms[fixpos])
1439    {
1440      int i;
1441      for(i=0;i<fixpos;i++)
1442        if((c->states[fixpos][i]!=HASTREP)&& (extended_product_criterion(c->S->m[fixpos],c->gcd_of_terms[fixpos], c->S->m[i],c->gcd_of_terms[i],c)))
1443{
1444          c->states[fixpos][i]=HASTREP;
1445    c->extended_product_crit++;
1446}     
1447      for(i=fixpos+1;i<c->n;i++)
1448        if((c->states[i][fixpos]!=HASTREP)&& (extended_product_criterion(c->S->m[fixpos],c->gcd_of_terms[fixpos], c->S->m[i],c->gcd_of_terms[i],c)))
1449  {        c->states[i][fixpos]=HASTREP;
1450  c->extended_product_crit++;
1451  }
1452    }
1453  }
1454}
1455static void c_S_element_changed_hook(int pos, slimgb_alg* c){
1456  length_one_crit(c,pos, c->lengths[pos]);
1457  if (!c->nc)
1458    line_of_extended_prod(pos,c);
1459}
1460class poly_tree_node {
1461public:
1462  poly p;
1463  poly_tree_node* l;
1464  poly_tree_node* r;
1465  int n;
1466  poly_tree_node(int sn):l(NULL),r(NULL),n(sn){}
1467};
1468class exp_number_builder{
1469public:
1470  poly_tree_node* top_level;
1471  int n;
1472  int get_n(poly p);
1473  exp_number_builder():top_level(0),n(0){}
1474};
1475int exp_number_builder::get_n(poly p){
1476  poly_tree_node** node=&top_level;
1477  while(*node!=NULL){
1478    int c=pLmCmp(p,(*node)->p);
1479    if (c==0) return (*node)->n;
1480    if (c==-1) node=&((*node)->r);
1481    else
1482      node=&((*node)->l);
1483  }
1484  (*node)= new poly_tree_node(n);
1485  n++;
1486  (*node)->p=pLmInit(p);
1487  return (*node)->n;
1488}
1489
1490//mac_polys exp are smaller iff they are greater by monomial ordering
1491//corresponding to solving linear equations notation
1492
1493//! obsolete
1494struct int_poly_pair{
1495  poly p;
1496  int n;
1497};
1498
1499
1500//! obsolete
1501void t2ippa_rec(poly* ip,int* ia, poly_tree_node* k, int &offset){
1502    if(!k) return;
1503    t2ippa_rec(ip,ia,k->l,offset);
1504    ip[offset]=k->p;
1505    ia[k->n]=offset;
1506    ++offset;
1507
1508    t2ippa_rec(ip,ia,k->r,offset);
1509    delete k;
1510  }
1511
1512//! obsolete
1513void t2ippa(poly* ip,int* ia,exp_number_builder & e){
1514
1515  int o=0;
1516  t2ippa_rec(ip,ia,e.top_level,o);
1517}
1518int anti_poly_order(const void* a, const void* b){
1519  return -pLmCmp(((int_poly_pair*) a)->p,((int_poly_pair*) b)->p );
1520}
1521
1522BOOLEAN is_valid_ro(red_object & ro){
1523  red_object r2=ro;
1524  ro.validate();
1525  if ((r2.p!=ro.p)||(r2.sev!=ro.sev)) return FALSE;
1526  return TRUE;
1527}
1528
1529
1530
1531static void go_on (slimgb_alg* c){
1532  //set limit of 1000 for multireductions, at the moment for
1533  //programming reasons
1534  int i=0;
1535  c->average_length=0;
1536  for(i=0;i<c->n;i++){
1537    c->average_length+=c->lengths[i];
1538  }
1539  c->average_length=c->average_length/c->n;
1540  i=0;
1541  poly* p=(poly*) omalloc((PAR_N+1)*sizeof(poly));//nullterminated
1542
1543  int curr_deg=-1;
1544  while(i<PAR_N){
1545    sorted_pair_node* s=top_pair(c);//here is actually chain criterium done
1546    if (!s) break;
1547    if(curr_deg>=0){
1548      if (s->deg >curr_deg) break;
1549    }
1550
1551    else curr_deg=s->deg;
1552    quick_pop_pair(c);
1553    if(s->i>=0){
1554      //be careful replace_pair use createShortSpoly which is not noncommutative
1555      //replace_pair(s->i,s->j,c);
1556    if(s->i==s->j) {
1557      free_sorted_pair_node(s,c->r);
1558      continue;
1559  }
1560    }
1561    poly h;
1562    if(s->i>=0){
1563      if (!c->nc)
1564  h=ksOldCreateSpoly(c->S->m[s->i], c->S->m[s->j], NULL, c->r);
1565      else
1566  h= nc_CreateSpoly(c->S->m[s->i], c->S->m[s->j], NULL, c->r);
1567    } 
1568    else
1569      h=s->lcm_of_lm;
1570    if(s->i>=0)
1571      now_t_rep(s->j,s->i,c);
1572    number coef;
1573    int mlen=pLength(h);
1574    if (!c->nc){
1575      h=redNF2(h,c,mlen,coef,2);
1576      redTailShort(h,c->strat);
1577      nDelete(&coef);
1578    }
1579    free_sorted_pair_node(s,c->r);
1580    if(!h) continue;
1581    int len=pLength(h);
1582    p[i]=h;
1583   
1584    i++;
1585  }
1586  p[i]=NULL;
1587//  pre_comp(p,i,c);
1588  if(i==0){
1589    omfree(p);
1590    return;
1591  }
1592  #ifdef TGB_RESORT_PAIRS
1593  c->replaced=new bool[c->n];
1594  c->used_b=FALSE;
1595  #endif
1596  red_object* buf=(red_object*) omalloc(i*sizeof(red_object));
1597  c->normal_forms+=i;
1598  int j;
1599  for(j=0;j<i;j++){
1600    buf[j].p=p[j];
1601    buf[j].sev=pGetShortExpVector(p[j]);
1602    buf[j].bucket = kBucketCreate(currRing);
1603   
1604    int len=pLength(p[j]);
1605    kBucketInit(buf[j].bucket,buf[j].p,len);
1606    buf[j].initial_quality=buf[j].guess_quality(c);
1607    assume(buf[j].initial_quality>=0);
1608  }
1609  omfree(p);
1610  qsort(buf,i,sizeof(red_object),red_object_better_gen);
1611//    Print("\ncurr_deg:%i\n",curr_deg);
1612  if (TEST_OPT_PROT)
1613  {
1614    Print("%dM[%d,",curr_deg,i);
1615  }
1616#ifdef FIND_DETERMINISTIC
1617  c->modifiedS=(BOOLEAN*) omalloc((c->strat->sl+1)*sizeof(BOOLEAN));
1618  c->expandS=(poly*) omalloc((1)*sizeof(poly));
1619  c->expandS[0]=NULL;
1620  int z2;
1621  for(z2=0;z2<=c->strat->sl;z2++)
1622    c->modifiedS[z2]=FALSE;
1623#endif
1624  multi_reduction(buf, i, c);
1625  #ifdef TGB_RESORT_PAIRS
1626  if (c->used_b) {
1627    if (TEST_OPT_PROT)
1628        PrintS("B");
1629    int e;
1630    for(e=0;e<=c->pair_top;e++){
1631        if(c->apairs[e]->i<0) continue;
1632        assume(c->apairs[e]->j>=0);
1633        if ((c->replaced[c->apairs[e]->i])||(c->replaced[c->apairs[e]->j])) {
1634            sorted_pair_node* s=c->apairs[e];
1635            s->expected_length=pair_weighted_length(s->i,s->j,c);
1636        }
1637    }
1638    qsort(c->apairs,c->pair_top+1,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1639  }
1640  #endif
1641#ifdef TGB_DEBUG
1642 {
1643   int k;
1644   for(k=0;k<i;k++)
1645   {
1646     assume(kFindDivisibleByInS_easy(c->strat,buf[k])<0);
1647     int k2;
1648     for(k2=0;k2<i;k2++)
1649     {
1650       if(k==k2) continue;
1651       assume((!(p_LmDivisibleBy(buf[k].p,buf[k2].p,c->r)))||(wrp(buf[k].p),Print(" k %d k2 %d ",k,k2),wrp(buf[k2].p),FALSE));
1652     }
1653   }
1654 }
1655#endif
1656  //resort S
1657#ifdef FIND_DETERMINISTIC
1658  for(z2=0;z2<=c->strat->sl;z2++)
1659  {
1660    if (c->modifiedS[z2])
1661    {
1662      wlen_type qual;
1663      int new_pos;
1664      if (c->strat->lenSw!=NULL)
1665          new_pos=simple_posInS(c->strat,c->strat->S[z2],strat->lenS[z2],strat->Sw[z2]);
1666      else
1667          new_pos=simple_posInS(c->strat,c->strat->S[z2],strat->lenS[z2],lenS[z2]);
1668     
1669      if (new_pos<z2)
1670      { 
1671         move_forward_in_S(z2,new_pos,c->strat);
1672      }
1673     
1674      assume(new_pos<=z2);
1675    }
1676  }
1677  for(z2=0;c->expandS[z2]!=NULL;z2++)
1678  {
1679    add_to_reductors(c,c->expandS[z2],pLength(c->expandS[z2]));
1680    // PrintS("E");
1681  }
1682  omfree(c->modifiedS);
1683  c->modifiedS=NULL;
1684  omfree(c->expandS);
1685  c->expandS=NULL;
1686#endif
1687  if (TEST_OPT_PROT)
1688      Print("%i]",i); 
1689
1690  int* ibuf=(int*) omalloc(i*sizeof(int));
1691  sorted_pair_node*** sbuf=(sorted_pair_node***) omalloc(i*sizeof(sorted_pair_node**));
1692  for(j=0;j<i;j++)
1693  {
1694    int len;
1695    poly p;
1696    buf[j].flatten();
1697    kBucketClear(buf[j].bucket,&p, &len);
1698    kBucketDestroy(&buf[j].bucket);
1699
1700    if (!c->nc)
1701      p=redTailShort(p, c->strat);
1702    sbuf[j]=add_to_basis_ideal_quotient(p,-1,-1,c,ibuf+j);
1703    //sbuf[j]=add_to_basis(p,-1,-1,c,ibuf+j);
1704  }
1705  int sum=0;
1706  for(j=0;j<i;j++){
1707    sum+=ibuf[j];
1708  }
1709  sorted_pair_node** big_sbuf=(sorted_pair_node**) omalloc(sum*sizeof(sorted_pair_node*));
1710  int partsum=0;
1711  for(j=0;j<i;j++)
1712  {
1713    memmove(big_sbuf+partsum, sbuf[j],ibuf[j]*sizeof(sorted_pair_node*));
1714    omfree(sbuf[j]);
1715    partsum+=ibuf[j];
1716  }
1717
1718  qsort(big_sbuf,sum,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1719  c->apairs=spn_merge(c->apairs,c->pair_top+1,big_sbuf,sum,c);
1720  c->pair_top+=sum;
1721  clean_top_of_pair_list(c);
1722  omfree(big_sbuf);
1723  omfree(sbuf);
1724  omfree(ibuf);
1725  omfree(buf);
1726#ifdef TGB_DEBUG
1727  int z;
1728  for(z=1;z<=c->pair_top;z++)
1729  {
1730    assume(pair_better(c->apairs[z],c->apairs[z-1],c));
1731  }
1732#endif
1733  if (TEST_OPT_PROT)
1734      Print("(%d)",c->pair_top+1); 
1735  while(!(idIs0(c->add_later))){
1736    ideal add=c->add_later;
1737   
1738    ideal add2=kInterRed(add,NULL);
1739    id_Delete(&add,currRing);
1740    idSkipZeroes(add2);
1741    c->add_later=idInit(5000,c->S->rank);
1742    memset(c->add_later->m,0,5000*sizeof(poly));
1743    for(i=0;i<add2->idelems();i++){
1744      if (add2->m[i]!=NULL)
1745          add_to_basis_ideal_quotient(add2->m[i],-1,-1,c,NULL);
1746      add2->m[i]=NULL;
1747    }
1748    id_Delete(&add2, c->r);
1749  }
1750  #ifdef TGB_RESORT_PAIRS
1751  delete c->replaced;
1752  c->replaced=NULL;
1753  c->used_b=FALSE;
1754  #endif
1755  return;
1756}
1757
1758
1759
1760#ifdef REDTAIL_S
1761
1762static poly redNFTail (poly h,const int sl,kStrategy strat, int len)
1763{
1764  if (h==NULL) return NULL;
1765  pTest(h);
1766  if (0 > sl)
1767    return h;
1768  if (pNext(h)==NULL) return h;
1769
1770  int j;
1771  poly res=h;
1772  poly act=res;
1773  LObject P(pNext(h));
1774  pNext(res)=NULL;
1775  P.bucket = kBucketCreate(currRing);
1776  len--;
1777  h=P.p;
1778  if (len <=0) len=pLength(h);
1779  kBucketInit(P.bucket,h /*P.p*/,len /*pLength(P.p)*/);
1780  pTest(h);
1781  loop
1782  {
1783      P.p=h;
1784      P.t_p=NULL;
1785      P.SetShortExpVector();
1786      loop
1787      {
1788          j=kFindDivisibleByInS(strat->S,strat->sevS,sl,&P);
1789          if (j>=0)
1790          {
1791#ifdef REDTAIL_PROT
1792            PrintS("r");
1793#endif
1794            nNormalize(pGetCoeff(P.p));
1795#ifdef KDEBUG
1796            if (TEST_OPT_DEBUG)
1797            {
1798              PrintS("red tail:");
1799              wrp(h);
1800              PrintS(" with ");
1801              wrp(strat->S[j]);
1802            }
1803#endif
1804            number coef;
1805            pTest(strat->S[j]);
1806            coef=kBucketPolyRed(P.bucket,strat->S[j],
1807                                strat->lenS[j]/*pLength(strat->S[j])*/,strat->kNoether);
1808            pMult_nn(res,coef);
1809            nDelete(&coef);
1810            h = kBucketGetLm(P.bucket);
1811            pTest(h);
1812            if (h==NULL)
1813            {
1814#ifdef REDTAIL_PROT
1815              PrintS(" ");
1816#endif
1817        kBucketDestroy(&P.bucket);
1818              return res;
1819            }
1820            pTest(h);
1821            P.p=h;
1822            P.t_p=NULL;
1823            P.SetShortExpVector();
1824#ifdef KDEBUG
1825            if (TEST_OPT_DEBUG)
1826            {
1827              PrintS("\nto tail:");
1828              wrp(h);
1829              PrintLn();
1830            }
1831#endif
1832          }
1833          else
1834          {
1835#ifdef REDTAIL_PROT
1836            PrintS("n");
1837#endif
1838            break;
1839          }
1840      } /* end loop current mon */
1841      //   poly tmp=pHead(h /*kBucketGetLm(P.bucket)*/);
1842      //act->next=tmp;pIter(act);
1843      act->next=kBucketExtractLm(P.bucket);pIter(act);
1844      h = kBucketGetLm(P.bucket);
1845      if (h==NULL)
1846      {
1847#ifdef REDTAIL_PROT
1848        PrintS(" ");
1849#endif
1850  kBucketDestroy(&P.bucket);
1851        return res;
1852      }
1853      pTest(h);
1854  }
1855}
1856#endif
1857
1858
1859//try to fill, return FALSE iff queue is empty
1860
1861//transfers ownership of m to mat
1862void init_with_mac_poly(tgb_sparse_matrix* mat, int row, mac_poly m){
1863  assume(mat->mp[row]==NULL);
1864  mat->mp[row]=m;
1865#ifdef TGB_DEBUG
1866  mac_poly r=m;
1867  while(r){
1868    assume(r->exp<mat->columns);
1869    r=r->next;
1870  }
1871#endif
1872}
1873poly free_row_to_poly(tgb_sparse_matrix* mat, int row, poly* monoms, int monom_index){
1874  poly p=NULL;
1875  poly* set_this=&p;
1876  mac_poly r=mat->mp[row];
1877  mat->mp[row]=NULL;
1878  while(r)
1879  {
1880    (*set_this)=pLmInit(monoms[monom_index-1-r->exp]);
1881    pSetCoeff((*set_this),r->coef);
1882    set_this=&((*set_this)->next);
1883    mac_poly old=r;
1884    r=r->next;
1885    delete old;
1886   
1887  }
1888  return p;
1889
1890}
1891
1892
1893
1894
1895
1896
1897
1898
1899
1900
1901
1902static int poly_crit(const void* ap1, const void* ap2){
1903  poly p1,p2;
1904  p1=*((poly*) ap1);
1905  p2=*((poly*)ap2);
1906
1907  int c=pLmCmp(p1,p2);
1908  if (c !=0) return c;
1909  int l1=pLength(p1);
1910  int l2=pLength(p2);
1911  if (l1<l2) return -1;
1912  if (l1>l2) return 1;
1913  return 0;
1914}
1915slimgb_alg::slimgb_alg(ideal I, BOOLEAN F4){
1916 
1917 
1918  r=currRing;
1919  nc=rIsPluralRing(r);
1920 
1921  is_homog=TRUE;
1922  {
1923    int hz;
1924    for(hz=0;hz<IDELEMS(I);hz++){
1925      assume(I->m[hz]!=NULL);
1926      int d=pTotaldegree(I->m[hz]);
1927      poly t=I->m[hz]->next;
1928      while(t)
1929      {
1930        if (d!=pTotaldegree(t,r))
1931        {
1932          is_homog=FALSE;
1933          break;
1934        }
1935        t=t->next;
1936      }
1937      if(!(is_homog)) break;
1938    }
1939  }
1940  //  Print("is homog:%d",c->is_homog);
1941  void* h;
1942  poly hp;
1943  int i,j;
1944  to_destroy=NULL;
1945  easy_product_crit=0;
1946  extended_product_crit=0;
1947  if (rField_is_Zp(r))
1948    is_char0=FALSE;
1949  else
1950    is_char0=TRUE;
1951  //not fully correct
1952  //(rChar()==0);
1953  F4_mode=F4;
1954 
1955  if ((!F4_mode)&&(!is_homog) &&(pLexOrder)){
1956    this->doubleSugar=TRUE;
1957  }
1958  else this->doubleSugar=FALSE;
1959  reduction_steps=0;
1960  last_index=-1;
1961
1962  F=NULL;
1963  F_minus=NULL;
1964
1965  Rcounter=0;
1966
1967  soon_free=NULL;
1968
1969  tmp_lm=pOne();
1970
1971  normal_forms=0;
1972  current_degree=1;
1973 
1974  max_pairs=5*I->idelems();
1975 
1976  apairs=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*max_pairs);
1977  pair_top=-1;
1978
1979  int n=I->idelems();
1980  array_lengths=n;
1981
1982 
1983  i=0;
1984  this->n=0;
1985  T_deg=(int*) omalloc(n*sizeof(int));
1986  if((!(is_homog)) &&(pLexOrder))
1987    T_deg_full=(int*) omalloc(n*sizeof(int));
1988  else
1989    T_deg_full=NULL;
1990  tmp_pair_lm=(poly*) omalloc(n*sizeof(poly));
1991  tmp_spn=(sorted_pair_node**) omalloc(n*sizeof(sorted_pair_node*));
1992  lm_bin=omGetSpecBin(POLYSIZE + (r->ExpL_Size)*sizeof(long));
1993#ifdef HEAD_BIN
1994  HeadBin=omGetSpecBin(POLYSIZE + (currRing->ExpL_Size)*sizeof(long));
1995#endif
1996  /* omUnGetSpecBin(&(c->HeadBin)); */
1997  #ifndef HAVE_BOOST
1998  h=omalloc(n*sizeof(char*));
1999 
2000  states=(char**) h;
2001  #endif
2002  h=omalloc(n*sizeof(int));
2003  lengths=(int*) h;
2004  weighted_lengths=(wlen_type*)omalloc(n*sizeof(wlen_type));
2005  gcd_of_terms=(poly*) omalloc(n*sizeof(poly));
2006 
2007  short_Exps=(long*) omalloc(n*sizeof(long));
2008  if (F4_mode)
2009    S=idInit(n,I->rank);
2010  else
2011    S=idInit(1,I->rank);
2012  strat=new skStrategy;
2013  strat->syzComp = 0;
2014  initBuchMoraCrit(strat);
2015  initBuchMoraPos(strat);
2016  strat->initEcart = initEcartBBA;
2017  strat->tailRing=r;
2018  strat->enterS = enterSBba;
2019  strat->sl = -1;
2020  i=n;
2021  i=1;//some strange bug else
2022  /* initS(c->S,NULL,c->strat); */
2023  /* intS start: */
2024  // i=((i+IDELEMS(c->S)+15)/16)*16;
2025  strat->ecartS=(intset)omAlloc(i*sizeof(int)); /*initec(i);*/
2026  strat->sevS=(unsigned long*)omAlloc0(i*sizeof(unsigned long));
2027  /*initsevS(i);*/
2028  strat->S_2_R=(int*)omAlloc0(i*sizeof(int));/*initS_2_R(i);*/
2029  strat->fromQ=NULL;
2030  strat->Shdl=idInit(1,1);
2031  strat->S=strat->Shdl->m;
2032  strat->lenS=(int*)omAlloc0(i*sizeof(int));
2033  if((is_char0)||((pLexOrder) &&(!is_homog)))
2034    strat->lenSw=(wlen_type*)omAlloc0(i*sizeof(wlen_type));
2035  else
2036    strat->lenSw=NULL;
2037  sorted_pair_node* si;
2038  assume(n>0);
2039  add_to_basis_ideal_quotient(I->m[0],-1,-1,this,NULL);
2040
2041  assume(strat->sl==strat->Shdl->idelems()-1);
2042  if(!(F4_mode))
2043  {
2044    for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
2045    {
2046      //     add_to_basis(I->m[i],-1,-1,c);
2047      si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
2048      si->i=-1;
2049      si->j=-2;
2050      si->expected_length=pQuality(I->m[i],this,pLength(I->m[i]));
2051      si->deg=pTotaldegree(I->m[i]);
2052      if (!rField_is_Zp(r)){ 
2053        pCleardenom(I->m[i]);
2054      }
2055      si->lcm_of_lm=I->m[i];
2056     
2057      //      c->apairs[n-1-i]=si;
2058      apairs[n-i-1]=si;
2059      ++(pair_top);
2060    }
2061  }
2062  else
2063  {
2064    for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
2065      add_to_basis_ideal_quotient(I->m[i],-1,-1,this,NULL);
2066  }
2067  for(i=0;i<I->idelems();i++)
2068  {
2069    I->m[i]=NULL;
2070  }
2071  idDelete(&I);
2072  add_later=idInit(5000,S->rank);
2073  memset(add_later->m,0,5000*sizeof(poly));
2074}
2075slimgb_alg::~slimgb_alg(){
2076  id_Delete(&add_later,r);
2077  int i,j;
2078  slimgb_alg* c=this;
2079  while(c->to_destroy)
2080  {
2081    pDelete(&(c->to_destroy->p));
2082    poly_list_node* old=c->to_destroy;
2083    c->to_destroy=c->to_destroy->next;
2084    omfree(old);
2085  }
2086  while(c->F)
2087  {
2088    for(i=0;i<c->F->size;i++){
2089      pDelete(&(c->F->mp[i].m));
2090    }
2091    omfree(c->F->mp);
2092    c->F->mp=NULL;
2093    mp_array_list* old=c->F;
2094    c->F=c->F->next;
2095    omfree(old);
2096  }
2097  while(c->F_minus)
2098  {
2099    for(i=0;i<c->F_minus->size;i++){
2100      pDelete(&(c->F_minus->p[i]));
2101    }
2102    omfree(c->F_minus->p);
2103    c->F_minus->p=NULL;
2104    poly_array_list* old=c->F_minus;
2105    c->F_minus=c->F_minus->next;
2106    omfree(old);
2107  }
2108  #ifndef HAVE_BOOST
2109  for(int z=1 /* zero length at 0 */;z<c->n;z++){
2110    omfree(c->states[z]);
2111  }
2112  omfree(c->states);
2113  #endif
2114
2115  omfree(c->lengths);
2116  omfree(c->weighted_lengths);
2117  for(int z=0;z<c->n;z++)
2118  {
2119    pDelete(&c->tmp_pair_lm[z]);
2120    omfree(c->tmp_spn[z]);
2121  }
2122  omfree(c->tmp_pair_lm);
2123  omfree(c->tmp_spn);
2124 
2125  omfree(c->T_deg);
2126  if(c->T_deg_full)
2127    omfree(c->T_deg_full);
2128
2129  omFree(c->strat->ecartS);
2130  omFree(c->strat->sevS);
2131//   initsevS(i);
2132  omFree(c->strat->S_2_R);
2133   
2134
2135  omFree(c->strat->lenS);
2136
2137  if(c->strat->lenSw)  omFree(c->strat->lenSw);
2138
2139
2140
2141
2142  for(i=0;i<c->n;i++){
2143    if(c->gcd_of_terms[i])
2144      pDelete(&(c->gcd_of_terms[i]));
2145  }
2146  omfree(c->gcd_of_terms);
2147
2148  omfree(c->apairs);
2149  if (TEST_OPT_PROT)
2150  {
2151    Print("calculated %d NFs\n",c->normal_forms);
2152    Print("applied %i product crit, %i extended_product crit \n", c->easy_product_crit, c->extended_product_crit);
2153  }
2154  int deleted_form_c_s=0;
2155 
2156  for(i=0;i<=c->strat->sl;i++){
2157    if (!c->strat->S[i]) continue;
2158    BOOLEAN found=FALSE;
2159    for(j=0;j<c->n;j++){
2160      if (c->S->m[j]==c->strat->S[i]){
2161        found=TRUE;
2162        break;
2163      }
2164    }
2165    if(!found) pDelete(&c->strat->S[i]);
2166  }
2167//   for(i=0;i<c->n;i++){
2168//     if (c->rep[i]!=i){
2169// //       for(j=0;j<=c->strat->sl;j++){
2170// //   if(c->strat->S[j]==c->S->m[i]){
2171// //     c->strat->S[j]=NULL;
2172// //     break;
2173// //   }
2174// //       }
2175// //      PrintS("R_delete");
2176//       pDelete(&c->S->m[i]);
2177//     }
2178//   }
2179
2180 
2181  for(i=0;i<c->n;i++)
2182  {
2183    assume(c->S->m[i]!=NULL);
2184    for(j=0;j<c->n;j++)
2185    {
2186      if((c->S->m[j]==NULL)||(i==j)) 
2187        continue;
2188      assume(p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
2189             c->S->m[i],~c->short_Exps[i],
2190             c->r)==p_LmDivisibleBy(c->S->m[j],
2191             c->S->m[i],
2192             c->r));
2193      if (p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
2194          c->S->m[i],~c->short_Exps[i],
2195          c->r))
2196      {
2197        pDelete(&c->S->m[i]);
2198        break;
2199      }
2200    }
2201  }
2202  omfree(c->short_Exps);
2203 
2204
2205  ideal I=c->S;
2206 
2207  IDELEMS(I)=c->n;
2208
2209  idSkipZeroes(I);
2210  for(i=0;i<=c->strat->sl;i++)
2211    c->strat->S[i]=NULL;
2212  id_Delete(&c->strat->Shdl,c->r);
2213  pDelete(&c->tmp_lm);
2214  omUnGetSpecBin(&lm_bin);
2215  delete c->strat;
2216}
2217ideal t_rep_gb(ring r,ideal arg_I, BOOLEAN F4_mode){
2218 
2219  //  Print("QlogSize(0) %d, QlogSize(1) %d,QlogSize(-2) %d, QlogSize(5) %d\n", QlogSize(nlInit(0)),QlogSize(nlInit(1)),QlogSize(nlInit(-2)),QlogSize(nlInit(5)));
2220 
2221  if (TEST_OPT_PROT)
2222    if (F4_mode)
2223      PrintS("F4 Modus \n");
2224     
2225  //debug_Ideal=arg_debug_Ideal;
2226  //if (debug_Ideal) PrintS("DebugIdeal received\n");
2227  // Print("Idelems %i \n----------\n",IDELEMS(arg_I));
2228  ideal I=idCompactify(arg_I);
2229   int i;
2230  for(i=0;i<IDELEMS(I);i++){
2231    simplify_poly(I->m[i],currRing);
2232  }
2233  if (idIs0(I)) return I;
2234
2235  qsort(I->m,IDELEMS(I),sizeof(poly),poly_crit);
2236  //Print("Idelems %i \n----------\n",IDELEMS(I));
2237  //slimgb_alg* c=(slimgb_alg*) omalloc(sizeof(slimgb_alg));
2238  slimgb_alg* c=new slimgb_alg(I, F4_mode);
2239   
2240
2241  while(c->pair_top>=0)
2242  {
2243    if(F4_mode)
2244      go_on_F4(c);
2245    else
2246      go_on(c);
2247  }
2248  I=c->S;
2249  delete c;
2250  if (TEST_OPT_REDSB){
2251    ideal erg=kInterRed(I,NULL);
2252    assume(I!=erg);
2253    id_Delete(&I, currRing);
2254    return erg;
2255  }
2256  //qsort(I->m, IDELEMS(I),sizeof(poly),pLmCmp_func);
2257  assume(I->rank>=idRankFreeModule(I));
2258  return(I);
2259}
2260void now_t_rep(const int & arg_i, const int & arg_j, slimgb_alg* c){
2261  int i,j;
2262  if (arg_i==arg_j){
2263    return;
2264  }
2265  if (arg_i>arg_j){
2266    i=arg_j;
2267    j=arg_i;
2268  } else {
2269    i=arg_i;
2270    j=arg_j;
2271  }
2272  c->states[j][i]=HASTREP;
2273}
2274
2275static BOOLEAN has_t_rep(const int & arg_i, const  int & arg_j, slimgb_alg* state){
2276  assume(0<=arg_i);
2277  assume(0<=arg_j);
2278  assume(arg_i<state->n);
2279  assume(arg_j<state->n);
2280  if (arg_i==arg_j)
2281  {
2282    return (TRUE);
2283  }
2284  if (arg_i>arg_j)
2285  {
2286    return (state->states[arg_i][arg_j]==HASTREP);
2287  } else
2288  {
2289    return (state->states[arg_j][arg_i]==HASTREP);
2290  }
2291}
2292static int pLcmDeg(poly a, poly b)
2293{
2294  int i;
2295  int n=0;
2296  for (i=pVariables; i; i--)
2297  {
2298    n+=si_max( pGetExp(a,i), pGetExp(b,i));
2299  }
2300  return n;
2301
2302}
2303
2304
2305
2306static void shorten_tails(slimgb_alg* c, poly monom)
2307{
2308  return;
2309// BOOLEAN corr=lenS_correct(c->strat);
2310  for(int i=0;i<c->n;i++)
2311  {
2312    //enter tail
2313   
2314    if (c->S->m[i]==NULL) continue;
2315    poly tail=c->S->m[i]->next;
2316    poly prev=c->S->m[i];
2317    BOOLEAN did_something=FALSE;
2318    while((tail!=NULL)&& (pLmCmp(tail, monom)>=0))
2319    {
2320      if (p_LmDivisibleBy(monom,tail,c->r))
2321      {
2322        did_something=TRUE;
2323        prev->next=tail->next;
2324        tail->next=NULL;
2325        p_Delete(& tail,c->r);
2326        tail=prev;
2327        //PrintS("Shortened");
2328        c->lengths[i]--;
2329      }
2330      prev=tail;
2331      tail=tail->next;
2332    }
2333    if (did_something)
2334    {
2335      int new_pos;
2336      wlen_type q;
2337      q=pQuality(c->S->m[i],c,c->lengths[i]);
2338      new_pos=simple_posInS(c->strat,c->S->m[i],c->lengths[i],q);
2339
2340      int old_pos=-1;
2341      //assume new_pos<old_pos
2342      for (int z=0;z<=c->strat->sl;z++)
2343      {
2344        if (c->strat->S[z]==c->S->m[i])
2345        {
2346          old_pos=z;
2347          break;
2348        }
2349      }
2350      if (old_pos== -1)
2351        for (int z=new_pos-1;z>=0;z--)
2352        {
2353          if (c->strat->S[z]==c->S->m[i])
2354          {
2355            old_pos=z;
2356            break;
2357          }
2358        }
2359      assume(old_pos>=0);
2360      assume(new_pos<=old_pos);
2361      assume(pLength(c->strat->S[old_pos])==c->lengths[i]);
2362      c->strat->lenS[old_pos]=c->lengths[i];
2363      if (c->strat->lenSw)
2364        c->strat->lenSw[old_pos]=q;
2365
2366      if (new_pos<old_pos)
2367        move_forward_in_S(old_pos,new_pos,c->strat);
2368
2369      length_one_crit(c,i,c->lengths[i]);
2370    }
2371  }
2372}
2373static sorted_pair_node* pop_pair(slimgb_alg* c){
2374  clean_top_of_pair_list(c);
2375
2376  if(c->pair_top<0) return NULL;
2377  else return (c->apairs[c->pair_top--]);
2378}
2379sorted_pair_node* top_pair(slimgb_alg* c){
2380  super_clean_top_of_pair_list(c);//yeah, I know, it's odd that I use a different proc here
2381
2382  if(c->pair_top<0) return NULL;
2383  else return (c->apairs[c->pair_top]);
2384}
2385sorted_pair_node* quick_pop_pair(slimgb_alg* c){
2386  if(c->pair_top<0) return NULL;
2387  else return (c->apairs[c->pair_top--]);
2388}
2389
2390
2391
2392static void super_clean_top_of_pair_list(slimgb_alg* c){
2393  while((c->pair_top>=0)
2394  && (c->apairs[c->pair_top]->i>=0)
2395  && (good_has_t_rep(c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,c)))
2396  {
2397
2398    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
2399    c->pair_top--;
2400
2401  }
2402}
2403void clean_top_of_pair_list(slimgb_alg* c){
2404  while((c->pair_top>=0) && (c->apairs[c->pair_top]->i>=0) && (!state_is(UNCALCULATED,c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,c))){
2405
2406    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
2407    c->pair_top--;
2408
2409  }
2410}
2411static BOOLEAN state_is(calc_state state, const int & arg_i, const  int & arg_j, slimgb_alg* c){
2412  assume(0<=arg_i);
2413  assume(0<=arg_j);
2414  assume(arg_i<c->n);
2415  assume(arg_j<c->n);
2416  if (arg_i==arg_j)
2417  {
2418    return (TRUE);
2419  }
2420  if (arg_i>arg_j)
2421  {
2422    return (c->states[arg_i][arg_j]==state);
2423  }
2424  else return(c->states[arg_j][arg_i]==state);
2425}
2426
2427
2428void free_sorted_pair_node(sorted_pair_node* s, ring r){
2429  if (s->i>=0)
2430    p_Delete(&s->lcm_of_lm,r);
2431  omfree(s);
2432}
2433static BOOLEAN pair_better(sorted_pair_node* a,sorted_pair_node* b, slimgb_alg* c){
2434  if (a->deg<b->deg) return TRUE;
2435  if (a->deg>b->deg) return FALSE;
2436
2437//  if (a->expected_length<b->expected_length) return TRUE;
2438  //  if (a->expected_length>b->expected_length) return FALSE;
2439  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
2440  if (comp==1) return FALSE;
2441  if (-1==comp) return TRUE;
2442  if (a->i+a->j<b->i+b->j) return TRUE;
2443   if (a->i+a->j>b->i+b->j) return FALSE;
2444  if (a->i<b->i) return TRUE;
2445  if (a->i>b->i) return FALSE;
2446  return TRUE;
2447}
2448
2449static int tgb_pair_better_gen(const void* ap,const void* bp){
2450
2451  sorted_pair_node* a=*((sorted_pair_node**)ap);
2452  sorted_pair_node* b=*((sorted_pair_node**)bp);
2453  assume(a->i>a->j);
2454  assume(b->i>b->j);
2455  if (a->deg<b->deg) return -1;
2456  if (a->deg>b->deg) return 1;
2457
2458
2459//  if (a->expected_length<b->expected_length) return -1;
2460  // if (a->expected_length>b->expected_length) return 1;
2461 int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
2462 
2463  if (comp==1) return 1;
2464  if (-1==comp) return -1;
2465  if (a->i+a->j<b->i+b->j) return -1;
2466   if (a->i+a->j>b->i+b->j) return 1;
2467  if (a->i<b->i) return -1;
2468   if (a->i>b->i) return 1;
2469  return 0;
2470}
2471
2472
2473static poly gcd_of_terms(poly p, ring r){
2474  int max_g_0=0;
2475  assume(p!=NULL);
2476  int i;
2477  poly m=pOne();
2478  poly t;
2479  for (i=pVariables; i; i--)
2480  {
2481      pSetExp(m,i, pGetExp(p,i));
2482      if (max_g_0==0)
2483  if (pGetExp(m,i)>0)
2484    max_g_0=i;
2485  }
2486 
2487  t=p->next;
2488  while (t!=NULL){
2489   
2490    if (max_g_0==0) break;
2491    for (i=max_g_0; i; i--)
2492    {
2493      pSetExp(m,i, si_min(pGetExp(t,i),pGetExp(m,i)));
2494      if (max_g_0==i)
2495  if (pGetExp(m,i)==0)
2496    max_g_0=0;
2497      if ((max_g_0==0) && (pGetExp(m,i)>0)){
2498  max_g_0=i;
2499      }
2500    }
2501    t=t->next;
2502  }
2503
2504  if (max_g_0>0)
2505    return m;
2506  pDelete(&m);
2507  return NULL;
2508}
2509static inline BOOLEAN pHasNotCFExtended(poly p1, poly p2, poly m)
2510{
2511
2512  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
2513    return FALSE;
2514  int i = 1;
2515  loop
2516  {
2517    if ((pGetExp(p1, i)-pGetExp(m,i) >0) && (pGetExp(p2, i) -pGetExp(m,i)> 0))   return FALSE;
2518    if (i == pVariables)                                return TRUE;
2519    i++;
2520  }
2521}
2522
2523
2524//for impl reasons may return false if the the normal product criterion matches
2525static inline BOOLEAN extended_product_criterion(poly p1, poly gcd1, poly p2, poly gcd2, slimgb_alg* c){
2526  if (c->nc)
2527    return FALSE;
2528  if(gcd1==NULL) return FALSE;
2529        if(gcd2==NULL) return FALSE;
2530        gcd1->next=gcd2; //may ordered incorrect
2531        poly m=gcd_of_terms(gcd1,c->r);
2532        gcd1->next=NULL;
2533        if (m==NULL) return FALSE;
2534
2535        BOOLEAN erg=pHasNotCFExtended(p1,p2,m);
2536        pDelete(&m);
2537        return erg;
2538}
2539static poly kBucketGcd(kBucket* b, ring r)
2540{
2541  int s=0;
2542  int i;
2543  poly m, n;
2544  BOOLEAN initialized=FALSE;
2545  for (i=MAX_BUCKET-1;i>=0;i--)
2546  { 
2547    if (b->buckets[i]!=NULL){
2548      if (!initialized){
2549  m=gcd_of_terms(b->buckets[i],r);
2550  initialized=TRUE;
2551  if (m==NULL) return NULL;
2552      }
2553      else
2554  {
2555    n=gcd_of_terms(b->buckets[i],r);
2556    if (n==NULL) {
2557      pDelete(&m);
2558      return NULL;   
2559    }
2560    n->next=m;
2561    poly t=gcd_of_terms(n,r);
2562    n->next=NULL;
2563    pDelete(&m);
2564    pDelete(&n);
2565    m=t;
2566    if (m==NULL) return NULL;
2567   
2568  }
2569    }
2570  }
2571  return m;
2572}
2573
2574
2575
2576
2577static inline int quality_of_pos_in_strat_S(int pos, slimgb_alg* c){
2578  if (c->strat->lenSw!=NULL) return c->strat->lenSw[pos];
2579  return c->strat->lenS[pos];
2580}
2581static inline int quality_of_pos_in_strat_S_mult_high(int pos, poly high, slimgb_alg* c)
2582  //meant only for nc
2583{
2584  poly m=pOne();
2585  pExpVectorDiff(m,high ,c->strat->S[pos]);
2586  poly product=nc_mm_Mult_p(m, pCopy(c->strat->S[pos]), c->r);
2587  int erg=pQuality(product,c);
2588  pDelete(&m);
2589  pDelete(&product);
2590  return erg;
2591}
2592
2593static void multi_reduction_lls_trick(red_object* los, int losl,slimgb_alg* c,find_erg & erg){
2594  erg.expand=NULL;
2595  BOOLEAN swap_roles; //from reduce_by, to_reduce_u if fromS
2596  if(erg.fromS){
2597    if(pLmEqual(c->strat->S[erg.reduce_by],los[erg.to_reduce_u].p))
2598    {
2599      int i;
2600      int quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
2601      int best=erg.to_reduce_u+1;
2602/*
2603      for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
2604  int qc=los[i].guess_quality(c);
2605  if (qc<quality_a){
2606    best=i;
2607    quality_a=qc;
2608  }
2609      }
2610      if(best!=erg.to_reduce_u+1){*/
2611      int qc;
2612      best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
2613      if(qc<quality_a){
2614  los[best].flatten();
2615  int b_pos=kBucketCanonicalize(los[best].bucket);
2616  los[best].p=los[best].bucket->buckets[b_pos];
2617  qc=pQuality(los[best].bucket->buckets[b_pos],c);
2618  if(qc<quality_a){
2619    red_object h=los[erg.to_reduce_u];
2620    los[erg.to_reduce_u]=los[best];
2621    los[best]=h;
2622    swap_roles=TRUE;
2623  }
2624  else
2625    swap_roles=FALSE;
2626      }
2627      else{
2628 
2629  swap_roles=FALSE;
2630      }
2631 
2632    }
2633      else
2634    {
2635      if (erg.to_reduce_u>erg.to_reduce_l){
2636
2637  int i;
2638  int quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
2639  if (c->nc)
2640    quality_a=quality_of_pos_in_strat_S_mult_high(erg.reduce_by, los[erg.to_reduce_u].p, c);
2641  int best=erg.to_reduce_u+1;
2642  int qc;
2643  best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
2644  assume(qc==los[best].guess_quality(c));
2645  if(qc<quality_a){
2646    los[best].flatten();
2647    int b_pos=kBucketCanonicalize(los[best].bucket);
2648    los[best].p=los[best].bucket->buckets[b_pos];
2649    qc==pQuality(los[best].bucket->buckets[b_pos],c);
2650    //(best!=erg.to_reduce_u+1)
2651    if(qc<quality_a){
2652    red_object h=los[erg.to_reduce_u];
2653    los[erg.to_reduce_u]=los[best];
2654    los[best]=h;
2655    erg.reduce_by=erg.to_reduce_u;
2656    erg.fromS=FALSE;
2657    erg.to_reduce_u--;
2658    }
2659  }
2660      }
2661      else 
2662      {
2663  assume(erg.to_reduce_u==erg.to_reduce_l);
2664  wlen_type quality_a=
2665        quality_of_pos_in_strat_S(erg.reduce_by,c);
2666  wlen_type qc=los[erg.to_reduce_u].guess_quality(c);
2667  if (qc<0) PrintS("Wrong wlen_type");
2668  if(qc<quality_a){
2669    int best=erg.to_reduce_u;
2670    los[best].flatten();
2671    int b_pos=kBucketCanonicalize(los[best].bucket);
2672    los[best].p=los[best].bucket->buckets[b_pos];
2673    qc=pQuality(los[best].bucket->buckets[b_pos],c);
2674    assume(qc>=0);
2675    if(qc<quality_a){
2676      BOOLEAN exp=FALSE;
2677      if(qc<=2){
2678         //Print("\n qc is %lld \n",qc);
2679         exp=TRUE;
2680      }
2681       
2682      else {
2683         if (qc<quality_a/2)
2684          exp=TRUE;
2685         else
2686       if(erg.reduce_by<c->n/4)
2687          exp=TRUE;
2688      }
2689      if (exp){
2690        poly clear_into;
2691        los[erg.to_reduce_u].flatten();
2692        kBucketClear(los[erg.to_reduce_u].bucket,&clear_into,&erg.expand_length);
2693        erg.expand=pCopy(clear_into);
2694        kBucketInit(los[erg.to_reduce_u].bucket,clear_into,erg.expand_length);
2695        if (TEST_OPT_PROT)
2696    PrintS("e");
2697       
2698      }
2699    }
2700  }
2701
2702 
2703      }
2704     
2705      swap_roles=FALSE;
2706      return;
2707      }
2708   
2709  }
2710  else{
2711    if(erg.reduce_by>erg.to_reduce_u){
2712      //then lm(rb)>= lm(tru) so =
2713      assume(erg.reduce_by==erg.to_reduce_u+1);
2714      int best=erg.reduce_by;
2715      wlen_type quality_a=los[erg.reduce_by].guess_quality(c);
2716      int qc;
2717      best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
2718     
2719      int i;
2720      if(qc<quality_a){
2721    red_object h=los[erg.reduce_by];
2722    los[erg.reduce_by]=los[best];
2723    los[best]=h;
2724  }
2725  swap_roles=FALSE;
2726  return;
2727 
2728   
2729    }
2730    else
2731    {
2732      assume(!pLmEqual(los[erg.reduce_by].p,los[erg.to_reduce_l].p));
2733      assume(erg.to_reduce_u==erg.to_reduce_l);
2734      //further assume, that reduce_by is the above all other polys
2735      //with same leading term
2736      int il=erg.reduce_by;
2737      int quality_a =los[erg.reduce_by].guess_quality(c);
2738      int qc;
2739      while((il>0) && pLmEqual(los[il-1].p,los[il].p)){
2740  il--;
2741  qc=los[il].guess_quality(c);
2742  if (qc<quality_a){
2743    quality_a=qc;
2744    erg.reduce_by=il;
2745  }
2746      }
2747      swap_roles=FALSE;
2748    }
2749 
2750  }
2751  if(swap_roles)
2752  {
2753    if (TEST_OPT_PROT)
2754      PrintS("b");
2755    poly clear_into;
2756    int dummy_len;
2757    int new_length;
2758    int bp=erg.to_reduce_u;//bucket_positon
2759    //kBucketClear(los[bp].bucket,&clear_into,&new_length);
2760    new_length=los[bp].clear_to_poly();
2761    clear_into=los[bp].p;
2762    poly p=c->strat->S[erg.reduce_by];
2763    int j=erg.reduce_by;
2764    int old_length=c->strat->lenS[j];// in view of S
2765    los[bp].p=p;
2766    kBucketInit(los[bp].bucket,p,old_length);
2767    wlen_type qal=pQuality(clear_into,c,new_length);
2768    int pos_in_c=-1;   
2769    int z;
2770    int new_pos;
2771    new_pos=simple_posInS(c->strat,clear_into,new_length, qal);
2772    assume(new_pos<=j);
2773    for (z=c->n;z;z--)
2774    {
2775      if(p==c->S->m[z-1])
2776      {
2777  pos_in_c=z-1;
2778  break;
2779      }
2780    }
2781    if(pos_in_c>=0)
2782    {
2783      #ifdef TGB_RESORT_PAIRS
2784      c->used_b=TRUE;
2785      c->replaced[pos_in_c]=TRUE;
2786      #endif
2787      c->S->m[pos_in_c]=clear_into;
2788      c->lengths[pos_in_c]=new_length;
2789      c->weighted_lengths[pos_in_c]=qal;
2790      if (c->gcd_of_terms[pos_in_c]==NULL)
2791        c->gcd_of_terms[pos_in_c]=gcd_of_terms(clear_into,c->r);
2792      if (c->T_deg_full)
2793  c->T_deg_full[pos_in_c]=pTotaldegree_full(clear_into);
2794      c_S_element_changed_hook(pos_in_c,c);
2795    }
2796    c->strat->S[j]=clear_into;
2797    c->strat->lenS[j]=new_length;
2798   
2799    assume(pLength(clear_into)==new_length);
2800    if(c->strat->lenSw!=NULL)
2801      c->strat->lenSw[j]=qal;
2802    if (!rField_is_Zp(c->r))
2803    {
2804      pContent(clear_into);
2805      pCleardenom(clear_into);//should be unnecessary
2806    }
2807    else                     
2808      pNorm(clear_into);
2809#ifdef FIND_DETERMINISTIC
2810    erg.reduce_by=j;
2811    //resort later see diploma thesis, find_in_S must be deterministic
2812    //during multireduction if spolys are only in the span of the
2813    //input polys
2814#else
2815   
2816    if (new_pos<j)
2817    { 
2818      move_forward_in_S(j,new_pos,c->strat);
2819      erg.reduce_by=new_pos;
2820    }
2821#endif
2822  }
2823}
2824static int fwbw(red_object* los, int i){
2825   int i2=i;
2826   int step=1;
2827   
2828   BOOLEAN bw=FALSE;
2829   BOOLEAN incr=TRUE;
2830   
2831   while(1)
2832   {
2833     if(!bw)
2834     {
2835       step=si_min(i2,step);
2836       if (step==0) break;
2837       i2-=step;
2838   
2839       if(!pLmEqual(los[i].p,los[i2].p))
2840       {
2841   bw=TRUE;
2842   incr=FALSE;
2843       }
2844       else
2845       {
2846   if ((!incr) &&(step==1)) break;
2847       }
2848       
2849       
2850     }
2851     else
2852     {
2853       
2854       step=si_min(i-i2,step);
2855       if (step==0) break;
2856       i2+=step;
2857       if(pLmEqual(los[i].p,los[i2].p)){
2858   if(step==1) break;
2859   else
2860   {
2861     bw=FALSE;
2862   }
2863       }
2864       
2865     }
2866     if (incr)
2867       step*=2;
2868     else
2869     {
2870       if (step%2==1)
2871   step=(step+1)/2;
2872       else
2873   step/=2;
2874       
2875     }
2876   }
2877   return i2;
2878}
2879static void canonicalize_region(red_object* los, int l, int u,slimgb_alg* c){
2880    assume(l<=u+1);
2881    int i;
2882    for(i=l;i<=u;i++){
2883        kBucketCanonicalize(los[i].bucket);
2884    }
2885
2886}
2887static void multi_reduction_find(red_object* los, int losl,slimgb_alg* c,int startf,find_erg & erg){
2888  kStrategy strat=c->strat;
2889
2890  assume(startf<=losl);
2891  assume((startf==losl-1)||(pLmCmp(los[startf].p,los[startf+1].p)==-1));
2892  int i=startf;
2893 
2894  int j;
2895  while(i>=0){
2896    assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)<=0));
2897    assume(is_valid_ro(los[i]));
2898    j=kFindDivisibleByInS_easy(strat,los[i]);
2899    if(j>=0){
2900     
2901      erg.to_reduce_u=i;
2902      erg.reduce_by=j;
2903      erg.fromS=TRUE;
2904      int i2=fwbw(los,i);
2905      assume(pLmEqual(los[i].p,los[i2].p));
2906      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
2907      assume(i>=i2);
2908
2909
2910      erg.to_reduce_l=i2;
2911      assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2912      canonicalize_region(los,erg.to_reduce_u+1,startf,c);
2913      return;
2914    }
2915    if (j<0){
2916     
2917      //not reduceable, try to use this for reducing higher terms
2918      int i2=fwbw(los,i);
2919      assume(pLmEqual(los[i].p,los[i2].p));
2920      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
2921      assume(i>=i2);
2922      if(i2!=i){
2923 
2924 
2925  erg.to_reduce_u=i-1;
2926  erg.to_reduce_l=i2;
2927  erg.reduce_by=i;
2928  erg.fromS=FALSE;
2929  assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2930  canonicalize_region(los,erg.to_reduce_u+1,startf,c);
2931  return;
2932      }
2933      if((!(c->is_homog)) &&(!(c->doubleSugar)))
2934      {
2935
2936  for (i2=i+1;i2<losl;i2++){
2937    if (p_LmShortDivisibleBy(los[i].p,los[i].sev,los[i2].p,~los[i2].sev,
2938           c->r)){
2939      int i3=i2;
2940      while((i3+1<losl) && (pLmEqual(los[i2].p, los[i3+1].p)))
2941        i3++;
2942      erg.to_reduce_u=i3;
2943      erg.to_reduce_l=i2;
2944      erg.reduce_by=i;
2945      erg.fromS=FALSE;
2946      assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2947      canonicalize_region(los,erg.to_reduce_u+1,startf,c);
2948      return;
2949    }
2950    //  else {assume(!p_LmDivisibleBy(los[i].p, los[i2].p,c->r));}
2951  }
2952      }
2953      i--;
2954    }
2955  }
2956  erg.reduce_by=-1;//error code
2957  return;
2958}
2959
2960 //  nicht reduzierbare eintraege in ergebnisliste schreiben
2961//   nullen loeschen
2962//   while(finde_groessten leitterm reduzierbar(c,erg)){
2963 
2964static int multi_reduction_clear_zeroes(red_object* los, int  losl, int l, int u)
2965{
2966
2967
2968  int deleted=0;
2969  int  i=l;
2970  int last=-1;
2971  while(i<=u)
2972  {
2973   
2974    if(los[i].p==NULL){
2975      kBucketDestroy(&los[i].bucket);
2976//      delete los[i];//here we assume los are constructed with new
2977      //destroy resources, must be added here   
2978     if (last>=0)
2979     {
2980       memmove(los+(int)(last+1-deleted),los+(last+1),sizeof(red_object)*(i-1-last));
2981     }
2982     last=i;
2983     deleted++;
2984    }
2985    i++;
2986  }
2987  if((last>=0)&&(last!=losl-1))
2988      memmove(los+(int)(last+1-deleted),los+last+1,sizeof(red_object)*(losl-1-last));
2989  return deleted;
2990 
2991}
2992
2993static void sort_region_down(red_object* los, int l, int u, slimgb_alg* c)
2994{
2995  qsort(los+l,u-l+1,sizeof(red_object),red_object_better_gen);
2996  int i;
2997
2998  for(i=l;i<=u;i++)
2999  {
3000    BOOLEAN moved=FALSE;
3001    int j;
3002    for(j=i;j;j--)
3003    {
3004      if(pLmCmp(los[j].p,los[j-1].p)==-1){
3005  red_object h=los[j];
3006  los[j]=los[j-1];
3007  los[j-1]=h;
3008  moved=TRUE;
3009      }
3010      else break;
3011    }
3012    if(!moved) return;
3013  }
3014}
3015
3016//assume that los is ordered ascending by leading term, all non zero
3017static void multi_reduction(red_object* los, int & losl, slimgb_alg* c)
3018{
3019  poly* delay=(poly*) omalloc(losl*sizeof(poly));
3020  int delay_s=0;
3021  //initialize;
3022  assume(c->strat->sl>=0);
3023  assume(losl>0);
3024  int i;
3025  wlen_type max_initial_quality=0;
3026 
3027  for(i=0;i<losl;i++){
3028    los[i].sev=pGetShortExpVector(los[i].p);
3029//SetShortExpVector();
3030    los[i].p=kBucketGetLm(los[i].bucket);
3031    if (los[i].initial_quality>max_initial_quality)
3032        max_initial_quality=los[i].initial_quality;
3033    // else
3034//         Print("init2_qal=%lld;", los[i].initial_quality);
3035//     Print("initial_quality=%lld;",max_initial_quality);
3036  }
3037
3038  kStrategy strat=c->strat;
3039  int curr_pos=losl-1;
3040
3041
3042//  nicht reduzierbare einträge in ergebnisliste schreiben
3043  // nullen loeschen
3044  while(curr_pos>=0){
3045   
3046    find_erg erg;
3047    multi_reduction_find(los, losl,c,curr_pos,erg);//last argument should be curr_pos
3048    if(erg.reduce_by<0) break;
3049
3050
3051
3052    erg.expand=NULL;
3053    int d=erg.to_reduce_u-erg.to_reduce_l+1;
3054   
3055   
3056    multi_reduction_lls_trick(los,losl,c,erg);
3057   
3058   
3059    int i;
3060    int len;
3061    //    wrp(los[erg.to_reduce_u].p);
3062    //Print("\n");
3063    multi_reduce_step(erg,los,c);
3064   
3065
3066    if(!K_TEST_OPT_REDTHROUGH){
3067  for(i=erg.to_reduce_l;i<=erg.to_reduce_u;i++){
3068     if  (los[i].p!=NULL)  //the check (los[i].p!=NULL) might be invalid
3069     {
3070         //
3071         assume(los[i].initial_quality>0);
3072         
3073               if(los[i].guess_quality(c)
3074                  >1.5*delay_factor*max_initial_quality){
3075                       if (TEST_OPT_PROT)
3076                           PrintS("v");
3077                       los[i].canonicalize();
3078                       if(los[i].guess_quality(c)
3079                           >delay_factor*max_initial_quality){
3080                               if (TEST_OPT_PROT)
3081                                   PrintS(".");
3082                               los[i].clear_to_poly();
3083                               //delay.push_back(los[i].p);
3084                               delay[delay_s]=los[i].p;
3085                               delay_s++;
3086
3087                               los[i].p=NULL;
3088                     
3089                      }
3090                  }
3091           
3092            }
3093     }
3094  }
3095    int deleted=multi_reduction_clear_zeroes(los, losl, erg.to_reduce_l, erg.to_reduce_u);
3096    if(erg.fromS==FALSE)
3097      curr_pos=si_max(erg.to_reduce_u,erg.reduce_by);
3098    else
3099      curr_pos=erg.to_reduce_u;
3100    losl -= deleted;
3101    curr_pos -= deleted;
3102
3103    //Print("deleted %i \n",deleted);
3104    if ((TEST_V_UPTORADICAL) &&(!(erg.fromS)))
3105        sort_region_down(los,si_min(erg.to_reduce_l,erg.reduce_by),(si_max(erg.to_reduce_u,erg.reduce_by))-deleted,c);
3106    else   
3107    sort_region_down(los, erg.to_reduce_l, erg.to_reduce_u-deleted, c);
3108
3109
3110    if(erg.expand)
3111    {
3112#ifdef FIND_DETERMINISTIC
3113      int i;
3114      for(i=0;c->expandS[i];i++);
3115      c->expandS=(poly*) omrealloc(c->expandS,(i+2)*sizeof(poly));
3116      c->expandS[i]=erg.expand;
3117      c->expandS[i+1]=NULL;
3118#else
3119      add_to_reductors(c,erg.expand,erg.expand_length);
3120#endif
3121    }
3122     
3123  }
3124 
3125
3126  sorted_pair_node** pairs=(sorted_pair_node**)
3127    omalloc(delay_s*sizeof(sorted_pair_node*)); 
3128  for(i=0;i<delay_s;i++){
3129       
3130      poly p=delay[i];
3131      //if (rPar(c->r)==0)
3132      simplify_poly(p,c->r);
3133      sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
3134      si->i=-1;
3135      si->j=-1;
3136       if (!rField_is_Zp(c->r)){
3137        if (!c->nc)
3138            p=redTailShort(p, c->strat);
3139        pCleardenom(p);
3140        pContent(p);
3141      }
3142      si->expected_length=pQuality(p,c,pLength(p));
3143      si->deg=pTotaldegree(p);
3144     
3145      si->lcm_of_lm=p;
3146      pairs[i]=si;
3147  }
3148  qsort(pairs,delay_s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
3149  c->apairs=spn_merge(c->apairs,c->pair_top+1,pairs,delay_s,c);
3150  c->pair_top+=delay_s;
3151  omfree(delay);
3152  return;
3153}
3154void red_object::flatten(){
3155  assume(p==kBucketGetLm(bucket));
3156}
3157void red_object::validate(){
3158  p=kBucketGetLm(bucket);
3159  if(p)
3160    sev=pGetShortExpVector(p);
3161}
3162int red_object::clear_to_poly(){
3163  flatten();
3164  int l;
3165  kBucketClear(bucket,&p,&l);
3166  return l;
3167}
3168
3169 
3170
3171
3172
3173void reduction_step::reduce(red_object* r, int l, int u){}
3174void simple_reducer::do_reduce(red_object & ro){
3175  number coef;
3176  if (!c->nc)
3177    coef=kBucketPolyRed(ro.bucket,p,
3178       p_len,
3179       c->strat->kNoether);
3180  else
3181    nc_kBucketPolyRed_Z(ro.bucket, p, &coef);
3182  nDelete(&coef);
3183}
3184
3185
3186void simple_reducer::reduce(red_object* r, int l, int u){
3187  this->pre_reduce(r,l,u);
3188  int i;
3189//debug start
3190  int im;
3191
3192
3193  for(i=l;i<=u;i++){
3194 
3195
3196
3197    this->do_reduce(r[i]);
3198 
3199  }
3200  for(i=l;i<=u;i++){
3201 
3202    kBucketSimpleContent(r[i].bucket);
3203    r[i].validate();
3204    #ifdef TGB_DEBUG
3205    #endif
3206  }
3207}
3208reduction_step::~reduction_step(){}
3209simple_reducer::~simple_reducer(){
3210  if(fill_back!=NULL)
3211  {
3212    kBucketInit(fill_back,p,p_len);
3213  }
3214  fill_back=NULL;
3215   
3216}
3217 
3218void multi_reduce_step(find_erg & erg, red_object* r, slimgb_alg* c){
3219  static int id=0;
3220  id++;
3221  unsigned long sev;
3222    BOOLEAN lt_changed=FALSE;
3223  int rn=erg.reduce_by;
3224  poly red;
3225  int red_len;
3226  simple_reducer* pointer;
3227  BOOLEAN work_on_copy=FALSE;
3228  if(erg.fromS){
3229    red=c->strat->S[rn];
3230    red_len=c->strat->lenS[rn];
3231    assume(red_len==pLength(red));
3232  }
3233  else
3234  {
3235    r[rn].flatten();
3236    kBucketClear(r[rn].bucket,&red,&red_len);
3237   
3238    if (!rField_is_Zp(c->r))
3239    {
3240      pContent(red);
3241      pCleardenom(red);//should be unnecessary
3242     
3243    }
3244    pNormalize(red);
3245   
3246
3247    if ((!(erg.fromS))&&(TEST_V_UPTORADICAL)){
3248         
3249         if (polynomial_root(red,c->r))
3250            lt_changed=TRUE;
3251            sev=p_GetShortExpVector(red,c->r);}
3252    red_len=pLength(red);
3253  }
3254  if (((TEST_V_MODPSOLVSB)&&(red_len>1))||((c->nc)||(erg.to_reduce_u-erg.to_reduce_l>5))){
3255    work_on_copy=TRUE;
3256    // poly m=pOne();
3257    poly m=c->tmp_lm;
3258    pSetCoeff(m,nInit(1));
3259    for(int i=1;i<=pVariables;i++)
3260      pSetExp(m,i,(pGetExp(r[erg.to_reduce_l].p, i)-pGetExp(red,i)));
3261    pSetm(m);
3262    poly red_cp;
3263    if (!c->nc)
3264      red_cp=ppMult_mm(red,m);
3265    else
3266      red_cp=nc_mm_Mult_p(m, pCopy(red), c->r);
3267    if(!erg.fromS){
3268      kBucketInit(r[rn].bucket,red,red_len);
3269    }
3270    //now reduce the copy
3271    //static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
3272
3273    if (!c->nc)
3274      redTailShort(red_cp,c->strat);
3275    //number mul;
3276    // red_len--;
3277//     red_cp->next=redNF2(red_cp->next,c,red_len,mul,c->average_length);
3278//     pSetCoeff(red_cp,nMult(red_cp->coef,mul));
3279//     nDelete(&mul);
3280//     red_len++;
3281    red=red_cp;
3282    red_len=pLength(red);
3283    // pDelete(&m);
3284   
3285  }
3286  int i;
3287
3288
3289
3290  assume(red_len==pLength(red));
3291 
3292  pointer=new simple_reducer(red,red_len,c);
3293
3294  if ((!work_on_copy) && (!erg.fromS))
3295    pointer->fill_back=r[rn].bucket;
3296  else
3297    pointer->fill_back=NULL;
3298  pointer->reduction_id=id;
3299  pointer->c=c;
3300
3301  pointer->reduce(r,erg.to_reduce_l, erg.to_reduce_u);
3302  if(work_on_copy) pDelete(&pointer->p);
3303  delete pointer;
3304  if (lt_changed){
3305    assume(!erg.fromS);
3306    r[erg.reduce_by].sev=sev;
3307  }
3308 
3309};
3310
3311
3312
3313 
3314void simple_reducer:: pre_reduce(red_object* r, int l, int u){}
3315
Note: See TracBrowser for help on using the repository browser.