source: git/kernel/tgb.cc @ bdabc8

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