source: git/kernel/tgb.cc @ 4a40cb0

spielwiese
Last change on this file since 4a40cb0 was 4a40cb0, checked in by Hans Schoenemann <hannes@…>, 12 years ago
typo git-svn-id: file:///usr/local/Singular/svn/trunk@14143 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 104.7 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$ */
8/*
9* ABSTRACT: slimgb and F4 implementation
10*/
11//#include <vector>
12//using namespace std;
13
14
15///@TODO: delay nur auf Sugarvergr?erung
16///@TODO: grade aus ecartS, setze dazu strat->honey; und nutze p.ecart
17///@TODO: no tail reductions in syz comp
18#include <kernel/mod2.h>
19#include <kernel/options.h>
20#include <kernel/tgb.h>
21#include <kernel/tgb_internal.h>
22#include <kernel/tgbgauss.h>
23
24#include <kernel/digitech.h>
25#include <kernel/gring.h>
26#include <kernel/sca.h>
27#include <kernel/prCopy.h>
28#include <kernel/longrat.h>
29#include <kernel/modulop.h>
30#include <stdlib.h>
31#include <stdio.h>
32#include <queue>
33#define BUCKETS_FOR_NORO_RED 1
34#define SR_HDL(A) ((long)(A))
35static const int bundle_size=100;
36static const int bundle_size_noro=10000;
37static const int delay_factor=3;
38int QlogSize(number n);
39#define ADD_LATER_SIZE 500
40#if 1
41static omBin lm_bin=NULL;
42
43static void simplify_poly(poly p, ring r)
44{
45     assume(r==currRing);
46     if (!rField_is_Zp(r))
47     {
48        p_Cleardenom(p,r);
49        //p_Content(p,r); //is a duplicate call, but belongs here
50     }
51     else
52       pNorm(p);
53}
54//static const BOOLEAN up_to_radical=TRUE;
55
56int slim_nsize(number n, ring r)
57{
58  if (rField_is_Zp(r))
59  {
60    return 1;
61  }
62  if (rField_is_Q(r))
63  {
64    return QlogSize(n);
65  }
66  else
67  {
68    return n_Size(n,r);
69  }
70}
71static BOOLEAN monomial_root(poly m, ring r){
72    BOOLEAN changed=FALSE;
73    int i;
74    for(i=1;i<=rVar(r);i++){
75        int e=p_GetExp(m,i,r);
76        if (e>1){
77            p_SetExp(m,i,1,r);
78            changed=TRUE;
79        }
80    }
81    if (changed) {
82        p_Setm(m,r);
83    }
84    return changed;
85}
86static BOOLEAN polynomial_root(poly h, ring r){
87  poly got=gcd_of_terms(h,r);
88  BOOLEAN changed=FALSE;
89  if((got!=NULL) &&(TEST_V_UPTORADICAL)) {
90    poly copy=p_Copy(got,r);
91    //p_wrp(got,c->r);
92    changed=monomial_root(got,r);
93    if (changed)
94    {
95         poly div_by=pDivide(copy, got);
96         poly iter=h;
97         while(iter){
98            pExpVectorSub(iter,div_by);
99            pIter(iter);
100         }
101         p_Delete(&div_by, r);
102         if (TEST_OPT_PROT)
103             PrintS("U");
104    }
105    p_Delete(&copy,r);
106  }
107  p_Delete(&got,r);
108  return changed;
109}
110static inline poly p_Init_Special(const ring r)
111{
112  return p_Init(r,lm_bin);
113}
114static inline poly pOne_Special(const ring r=currRing)
115{
116  poly rc = p_Init_Special(r);
117  pSetCoeff0(rc,n_Init(1,r));
118  return rc;
119}
120// zum Initialiseren: in t_rep_gb plazieren:
121
122#endif
123#define LEN_VAR3
124#define degbound(p) assume(pTotaldegree(p)<10)
125//#define inDebug(p) assume((debug_Ideal==NULL)||(kNF(debug_Ideal,NULL,p,0,0)==0))
126
127//die meisten Varianten stossen sich an coef_buckets
128
129
130
131#ifdef LEN_VAR1
132// erste Variante: Laenge: Anzahl der Monome
133static inline int pSLength(poly p, int l) { return l; }
134static inline int kSBucketLength(kBucket* bucket, poly lm)
135  {return bucket_guess(bucket);}
136#endif
137
138#ifdef LEN_VAR2
139// 2. Variante: Laenge: Platz fuer die Koeff.
140int pSLength(poly p,int l)
141{
142  int s=0;
143  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
144  return s;
145}
146int kSBucketLength(kBucket* b, poly lm)
147{
148  int s=0;
149  int i;
150  for (i=MAX_BUCKET;i>=0;i--)
151  {
152    s+=pSLength(b->buckets[i],0);
153  }
154  return s;
155}
156#endif
157
158
159
160
161
162
163int QlogSize(number n){
164
165    if (SR_HDL(n) & SR_INT){
166       long i=SR_TO_INT(n);
167       if (i==0) return 0;
168
169       unsigned long v;
170       v=(i>=0)?i:-i;
171       int r = 0;
172
173       while (v >>= 1)
174       {
175        r++;
176       }
177       return r+1;
178    }
179    //assume denominator is 0
180    return mpz_sizeinbase(n->z,2);
181}
182
183
184#ifdef LEN_VAR3
185static inline wlen_type pSLength(poly p,int l)
186{
187  wlen_type c;
188  number coef=pGetCoeff(p);
189  if (rField_is_Q(currRing)){
190    c=QlogSize(coef);
191  }
192  else
193    c=nSize(coef);
194  if (!(TEST_V_COEFSTRAT))
195      return (wlen_type)c*(wlen_type)l /*pLength(p)*/;
196  else {
197    wlen_type res=l;
198    res*=c;
199    res*=c;
200    return res;
201  }
202}
203//! TODO CoefBuckets bercksichtigen
204wlen_type kSBucketLength(kBucket* b, poly lm=NULL)
205{
206  int s=0;
207  wlen_type c;
208  number coef;
209  if(lm==NULL)
210    coef=pGetCoeff(kBucketGetLm(b));
211    //c=nSize(pGetCoeff(kBucketGetLm(b)));
212  else
213    coef=pGetCoeff(lm);
214    //c=nSize(pGetCoeff(lm));
215  if (rField_is_Q(currRing)){
216    c=QlogSize(coef);
217  }
218  else
219    c=nSize(coef);
220
221  int i;
222  for (i=b->buckets_used;i>=0;i--)
223  {
224    assume((b->buckets_length[i]==0)||(b->buckets[i]!=NULL));
225    s+=b->buckets_length[i] /*pLength(b->buckets[i])*/;
226  }
227  #ifdef HAVE_COEF_BUCKETS
228  assume(b->buckets[0]==kBucketGetLm(b));
229  if (b->coef[0]!=NULL){
230
231    if (rField_is_Q(currRing)){
232      int modifier=QlogSize(pGetCoeff(b->coef[0]));
233      c+=modifier;
234  }
235    else{
236      int modifier=nSize(pGetCoeff(b->coef[0]));
237      c*=modifier;}
238    }
239  #endif
240  if (!(TEST_V_COEFSTRAT)){
241  return s*c;
242  } else
243  {
244    wlen_type res=s;
245    res*=c;
246    res*=c;
247    return res;
248  }
249}
250#endif
251#ifdef LEN_VAR5
252static inline wlen_type pSLength(poly p,int l)
253{
254  int c;
255  number coef=pGetCoeff(p);
256  if (rField_is_Q(currRing)){
257    c=QlogSize(coef);
258  }
259  else
260    c=nSize(coef);
261  wlen_type erg=l;
262  erg*=c;
263  erg*=c;
264  //PrintS("lenvar 5");
265  assume(erg>=0);
266  return erg; /*pLength(p)*/;
267}
268//! TODO CoefBuckets bercksichtigen
269wlen_type kSBucketLength(kBucket* b, poly lm=NULL)
270{
271  wlen_type s=0;
272  int c;
273  number coef;
274  if(lm==NULL)
275    coef=pGetCoeff(kBucketGetLm(b));
276    //c=nSize(pGetCoeff(kBucketGetLm(b)));
277  else
278    coef=pGetCoeff(lm);
279    //c=nSize(pGetCoeff(lm));
280  if (rField_is_Q(currRing)){
281    c=QlogSize(coef);
282  }
283  else
284    c=nSize(coef);
285
286  int i;
287  for (i=b->buckets_used;i>=0;i--)
288  {
289    assume((b->buckets_length[i]==0)||(b->buckets[i]!=NULL));
290    s+=b->buckets_length[i] /*pLength(b->buckets[i])*/;
291  }
292  #ifdef HAVE_COEF_BUCKETS
293  assume(b->buckets[0]==kBucketGetLm(b));
294  if (b->coef[0]!=NULL){
295
296    if (rField_is_Q(currRing)){
297      int modifier=QlogSize(pGetCoeff(b->coef[0]));
298      c+=modifier;
299  }
300    else{
301      int modifier=nSize(pGetCoeff(b->coef[0]));
302      c*=modifier;}
303    }
304  #endif
305  wlen_type erg=s;
306  erg*=c;
307  erg*=c;
308  return erg;
309}
310#endif
311
312#ifdef LEN_VAR4
313// 4.Variante: Laenge: Platz fuer Leitk * (1+Platz fuer andere Koeff.)
314int pSLength(poly p, int l)
315{
316  int s=1;
317  int c=nSize(pGetCoeff(p));
318  pIter(p);
319  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
320  return s*c;
321}
322int kSBucketLength(kBucket* b)
323{
324  int s=1;
325  int c=nSize(pGetCoeff(kBucketGetLm(b)));
326  int i;
327  for (i=MAX_BUCKET;i>0;i--)
328  {
329    if(b->buckets[i]==NULL) continue;
330    s+=pSLength(b->buckets[i],0);
331  }
332  return s*c;
333}
334#endif
335//BUG/TODO this stuff will fail on internal Schreyer orderings
336static BOOLEAN elength_is_normal_length(poly p, slimgb_alg* c){
337    ring r=c->r;
338    if (p_GetComp(p,r)!=0) return FALSE;
339    if (c->lastDpBlockStart<=pVariables){
340        int i;
341        for(i=1;i<c->lastDpBlockStart;i++){
342            if (p_GetExp(p,i,r)!=0){
343                break;
344            }
345        }
346        if (i>=c->lastDpBlockStart) {
347        //wrp(p);
348        //PrintS("\n");
349        return TRUE;
350        }
351        else return FALSE;
352    }else
353    return FALSE;
354}
355
356static BOOLEAN lies_in_last_dp_block(poly p, slimgb_alg* c){
357    ring r=c->r;
358    if (p_GetComp(p,r)!=0) return FALSE;
359    if (c->lastDpBlockStart<=pVariables){
360        int i;
361        for(i=1;i<c->lastDpBlockStart;i++){
362            if (p_GetExp(p,i,r)!=0){
363                break;
364            }
365        }
366        if (i>=c->lastDpBlockStart) {
367        //wrp(p);
368        //PrintS("\n");
369        return TRUE;
370        }
371        else return FALSE;
372    }else
373    return FALSE;
374}
375
376static int get_last_dp_block_start(ring r){
377    //ring r=c->r;
378    int last_block;
379
380    if (rRing_has_CompLastBlock(r)){
381        last_block=rBlocks(r) - 3;
382    }
383    else {last_block=rBlocks(r)-2;}
384    assume(last_block>=0);
385    if (r->order[last_block]==ringorder_dp)
386        return r->block0[last_block];
387    return pVariables+1;
388
389}
390
391static wlen_type do_pELength(poly p, slimgb_alg* c, int dlm=-1){
392
393  if(p==NULL) return 0;
394  wlen_type s=0;
395  poly pi=p;
396  if(dlm<0){
397    dlm=c->pTotaldegree(p);
398    s=1;
399    pi=p->next;
400  }
401
402  while(pi){
403    int d=c->pTotaldegree(pi);
404    if(d>dlm)
405      s+=1+d-dlm;
406    else
407      ++s;
408    pi=pi->next;
409  }
410  return s;
411}
412
413wlen_type pELength(poly p, slimgb_alg* c, ring r){
414  if(p==NULL) return 0;
415  wlen_type s=0;
416  poly pi=p;
417  int dlm;
418    dlm=c->pTotaldegree(p);
419    s=1;
420    pi=p->next;
421
422
423  while(pi){
424    int d=c->pTotaldegree(pi);
425    if(d>dlm)
426      s+=1+d-dlm;
427    else
428      ++s;
429    pi=pi->next;
430  }
431  return s;
432}
433
434wlen_type kEBucketLength(kBucket* b, poly lm,int sugar,slimgb_alg* ca)
435{
436  wlen_type s=0;
437  if(lm==NULL){
438    lm=kBucketGetLm(b);
439  }
440  if(lm==NULL) return 0;
441  if(elength_is_normal_length(lm,ca)) {
442    return bucket_guess(b);
443  }
444  int d=ca->pTotaldegree(lm);
445  #if 0
446  assume(sugar>=d);
447  s=1+(bucket_guess(b)-1)*(sugar-d+1);
448  return s;
449  #else
450
451
452  //int d=pTotaldegree(lm,ca->r);
453  int i;
454  for (i=b->buckets_used;i>=0;i--)
455  {
456    if(b->buckets[i]==NULL) continue;
457
458    if ((ca->pTotaldegree(b->buckets[i])<=d) &&(elength_is_normal_length(b->buckets[i],ca))){
459        s+=b->buckets_length[i];
460    } else
461    {
462    s+=do_pELength(b->buckets[i],ca,d);
463    }
464  }
465  return s;
466  #endif
467}
468
469static inline int pELength(poly p, slimgb_alg* c,int l){
470  if (p==NULL) return 0;
471  if ((l>0) &&(elength_is_normal_length(p,c)))
472    return l;
473  return do_pELength(p,c);
474}
475
476
477
478
479static inline wlen_type pQuality(poly p, slimgb_alg* c, int l=-1){
480
481  if(l<0)
482    l=pLength(p);
483  if(c->isDifficultField) {
484    if(c->eliminationProblem){
485      wlen_type cs;
486      number coef=pGetCoeff(p);
487      if (rField_is_Q(currRing)){
488         cs=QlogSize(coef);
489      }
490      else
491        cs=nSize(coef);
492     wlen_type erg=cs;
493     if(TEST_V_COEFSTRAT)
494        erg*=cs;
495     //erg*=cs;//for quadratic
496     erg*=pELength(p,c,l);
497    //FIXME: not quadratic coeff size
498      //return cs*pELength(p,c,l);
499      return erg;
500    }
501    //PrintS("I am here");
502    wlen_type r=pSLength(p,l);
503    assume(r>=0);
504    return r;
505  }
506  if(c->eliminationProblem) return pELength(p,c,l);
507  return l;
508}
509
510static inline int pTotaldegree_full(poly p){
511  int r=0;
512  while(p){
513    int d=pTotaldegree(p);
514    r=si_max(r,d);
515    pIter(p);
516  }
517  return r;
518}
519
520wlen_type red_object::guess_quality(slimgb_alg* c){
521    //works at the moment only for lenvar 1, because in different
522    //case, you have to look on coefs
523    wlen_type s=0;
524    if (c->isDifficultField){
525      //s=kSBucketLength(bucket,this->p);
526      if(c->eliminationProblem){
527    wlen_type cs;
528    number coef;
529
530    coef=pGetCoeff(kBucketGetLm(bucket));
531    //c=nSize(pGetCoeff(kBucketGetLm(b)));
532
533    //c=nSize(pGetCoeff(lm));
534    if (rField_is_Q(currRing)){
535      cs=QlogSize(coef);
536    }
537    else
538      cs=nSize(coef);
539    #ifdef HAVE_COEF_BUCKETS
540    if (bucket->coef[0]!=NULL){
541      if (rField_is_Q(currRing)){
542        int modifier=QlogSize(pGetCoeff(bucket->coef[0]));
543        cs+=modifier;
544      }
545      else{
546        int modifier=nSize(pGetCoeff(bucket->coef[0]));
547        cs*=modifier;}
548    }
549    #endif
550    //FIXME:not quadratic
551    wlen_type erg=kEBucketLength(this->bucket,this->p,this->sugar,c);
552    //erg*=cs;//for quadratic
553    erg*=cs;
554    if (TEST_V_COEFSTRAT)
555        erg*=cs;
556    //return cs*kEBucketLength(this->bucket,this->p,c);
557    return erg;
558      }
559      s=kSBucketLength(bucket,NULL);
560    }
561    else
562    {
563      if(c->eliminationProblem)
564  //if (false)
565  s=kEBucketLength(this->bucket,this->p,this->sugar,c);
566      else s=bucket_guess(bucket);
567    }
568
569    return s;
570}
571
572
573
574static void finalize_reduction_step(reduction_step* r){
575  delete r;
576}
577static int LObject_better_gen(const void* ap, const void* bp)
578{
579  LObject* a=*(LObject**)ap;
580  LObject* b=*(LObject**)bp;
581  return(pLmCmp(a->p,b->p));
582}
583static int red_object_better_gen(const void* ap, const void* bp)
584{
585
586
587  return(pLmCmp(((red_object*) ap)->p,((red_object*) bp)->p));
588}
589
590
591static int pLmCmp_func_inverted(const void* ap1, const void* ap2){
592    poly p1,p2;
593  p1=*((poly*) ap1);
594  p2=*((poly*)ap2);
595
596  return -pLmCmp(p1,p2);
597}
598
599int tgb_pair_better_gen2(const void* ap,const void* bp){
600  return(-tgb_pair_better_gen(ap,bp));
601}
602int kFindDivisibleByInS_easy(kStrategy strat,const red_object & obj){
603  int i;
604  long not_sev=~obj.sev;
605  poly p=obj.p;
606  for(i=0;i<=strat->sl;i++){
607    if (pLmShortDivisibleBy(strat->S[i],strat->sevS[i],p,not_sev))
608      return i;
609  }
610  return -1;
611}
612int kFindDivisibleByInS_easy(kStrategy strat,poly p, long sev){
613  int i;
614  long not_sev=~sev;
615  for(i=0;i<=strat->sl;i++){
616    if (pLmShortDivisibleBy(strat->S[i],strat->sevS[i],p,not_sev))
617      return i;
618  }
619  return -1;
620}
621static int posInPairs (sorted_pair_node**  p, int pn, sorted_pair_node* qe,slimgb_alg* c,int an=0)
622{
623  if(pn==0) return 0;
624
625  int length=pn-1;
626  int i;
627  //int an = 0;
628  int en= length;
629
630  if (pair_better(qe,p[en],c))
631    return length+1;
632
633  while(1)
634    {
635      //if (an >= en-1)
636      if(en-1<=an)
637      {
638        if (pair_better(p[an],qe,c)) return an;
639        return en;
640      }
641      i=(an+en) / 2;
642        if (pair_better(p[i],qe,c))
643          en=i;
644      else an=i;
645    }
646}
647
648static BOOLEAN  ascending(int* i,int top){
649  if(top<1) return TRUE;
650  if(i[top]<i[top-1]) return FALSE;
651  return ascending(i,top-1);
652}
653
654sorted_pair_node**  spn_merge(sorted_pair_node** p, int pn,sorted_pair_node **q, int qn,slimgb_alg* c){
655  int i;
656  int* a= (int*) omalloc(qn*sizeof(int));
657//   int mc;
658//   PrintS("Debug\n");
659//   for(mc=0;mc<qn;mc++)
660// {
661
662//     wrp(q[mc]->lcm_of_lm);
663//     PrintS("\n");
664// }
665//    PrintS("Debug they are in\n");
666//   for(mc=0;mc<pn;mc++)
667// {
668
669//     wrp(p[mc]->lcm_of_lm);
670//     PrintS("\n");
671// }
672  int lastpos=0;
673  for(i=0;i<qn;i++){
674    lastpos=posInPairs(p,pn,q[i],c, si_max(lastpos-1,0));
675    //   cout<<lastpos<<"\n";
676    a[i]=lastpos;
677
678  }
679  if((pn+qn)>c->max_pairs){
680    p=(sorted_pair_node**) omrealloc(p,2*(pn+qn)*sizeof(sorted_pair_node*));
681    c->max_pairs=2*(pn+qn);
682  }
683  for(i=qn-1;i>=0;i--){
684    size_t size;
685    if(qn-1>i)
686      size=(a[i+1]-a[i])*sizeof(sorted_pair_node*);
687    else
688      size=(pn-a[i])*sizeof(sorted_pair_node*); //as indices begin with 0
689    memmove (p+a[i]+(1+i), p+a[i], size);
690    p[a[i]+i]=q[i];
691  }
692  omfree(a);
693  return p;
694}
695
696
697static BOOLEAN trivial_syzygie(int pos1,int pos2,poly bound,slimgb_alg* c){
698
699
700  poly p1=c->S->m[pos1];
701  poly p2=c->S->m[pos2];
702
703
704  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
705    return FALSE;
706  int i = 1;
707  poly m=NULL;
708  poly gcd1=c->gcd_of_terms[pos1];
709  poly gcd2=c->gcd_of_terms[pos2];
710
711  if((gcd1!=NULL) && (gcd2!=NULL))
712    {
713      gcd1->next=gcd2; //may ordered incorrect
714      m=gcd_of_terms(gcd1,c->r);
715      gcd1->next=NULL;
716
717    }
718
719  if (m==NULL)
720  {
721     loop
722      {
723  if (pGetExp(p1, i)+ pGetExp(p2, i) > pGetExp(bound,i))   return FALSE;
724  if (i == pVariables){
725    //PrintS("trivial");
726    return TRUE;
727  }
728  i++;
729      }
730  }
731  else
732  {
733    loop
734      {
735  if (pGetExp(p1, i)-pGetExp(m,i) + pGetExp(p2, i) > pGetExp(bound,i))  {
736    pDelete(&m);
737    return FALSE;}
738  if (i == pVariables){
739    pDelete(&m);
740    //PrintS("trivial");
741    return TRUE;
742  }
743  i++;
744      }
745  }
746
747
748
749
750}
751
752//! returns position sets w as weight
753int find_best(red_object* r,int l, int u, wlen_type &w, slimgb_alg* c){
754  int best=l;
755  int i;
756  w=r[l].guess_quality(c);
757  for(i=l+1;i<=u;i++){
758    wlen_type w2=r[i].guess_quality(c);
759    if(w2<w){
760      w=w2;
761      best=i;
762    }
763
764  }
765 return best;
766}
767
768
769void red_object::canonicalize(){
770  kBucketCanonicalize(bucket);
771
772
773}
774BOOLEAN good_has_t_rep(int i, int j,slimgb_alg* c){
775  assume(i>=0);
776    assume(j>=0);
777  if (has_t_rep(i,j,c)) return TRUE;
778  //poly lm=pOne();
779  assume (c->tmp_lm!=NULL);
780  poly lm=c->tmp_lm;
781
782  pLcm(c->S->m[i], c->S->m[j], lm);
783  pSetm(lm);
784  assume(lm!=NULL);
785  //int deciding_deg= pTotaldegree(lm);
786  int* i_con =make_connections(i,j,lm,c);
787  //p_Delete(&lm,c->r);
788
789
790  for (int n=0;((n<c->n) && (i_con[n]>=0));n++){
791    if (i_con[n]==j){
792      now_t_rep(i,j,c);
793      omfree(i_con);
794
795      return TRUE;
796    }
797  }
798  omfree(i_con);
799
800  return FALSE;
801}
802BOOLEAN lenS_correct(kStrategy strat){
803  int i;
804  for(i=0;i<=strat->sl;i++){
805    if (strat->lenS[i]!=pLength(strat->S[i]))
806      return FALSE;
807  }
808  return TRUE;
809}
810
811
812static void cleanS(kStrategy strat, slimgb_alg* c)
813{
814  int i=0;
815  LObject P;
816  while(i<=strat->sl)
817  {
818    P.p=strat->S[i];
819    P.sev=strat->sevS[i];
820    int dummy=strat->sl;
821    //if(kFindDivisibleByInS(strat,&dummy,&P)!=i)
822    if (kFindDivisibleByInS_easy(strat,P.p,P.sev)!=i)
823    {
824      deleteInS(i,strat);
825      //remember destroying poly
826      BOOLEAN found=FALSE;
827      int j;
828      for(j=0;j<c->n;j++)
829      {
830        if(c->S->m[j]==P.p)
831        {
832          found=TRUE;
833          break;
834        }
835      }
836      if (!found)
837        pDelete(&P.p);
838      //remember additional reductors
839    }
840    else i++;
841  }
842}
843static int bucket_guess(kBucket* bucket){
844  int sum=0;
845  int i;
846  for (i=bucket->buckets_used;i>=0;i--){
847    if(bucket->buckets[i])
848       sum+=bucket->buckets_length[i];
849  }
850  return sum;
851}
852
853
854
855
856
857
858static int add_to_reductors(slimgb_alg* c, poly h, int len, int ecart, BOOLEAN simplified)
859{
860  //inDebug(h);
861  assume(lenS_correct(c->strat));
862  assume(len==pLength(h));
863  int i;
864//   if (c->isDifficultField)
865//        i=simple_posInS(c->strat,h,pSLength(h,len),c->isDifficultField);
866//   else
867//     i=simple_posInS(c->strat,h,len,c->isDifficultField);
868
869  LObject P; memset(&P,0,sizeof(P));
870  P.tailRing=c->r;
871  P.p=h; /*p_Copy(h,c->r);*/
872  P.ecart=ecart;
873  P.FDeg=pFDeg(P.p,c->r);
874  if (!(simplified))
875  {
876      if (!rField_is_Zp(c->r))
877      {
878        p_Cleardenom(P.p,c->r);
879        //p_Content(P.p,c->r ); //is a duplicate call, but belongs here
880
881      }
882      else
883        pNorm(P.p);
884    pNormalize(P.p);
885  }
886  wlen_type pq=pQuality(h,c,len);
887  i=simple_posInS(c->strat,h,len,pq);
888  c->strat->enterS(P,i,c->strat,-1);
889
890
891
892  c->strat->lenS[i]=len;
893  assume(pLength(c->strat->S[i])==c->strat->lenS[i]);
894  if(c->strat->lenSw!=NULL)
895    c->strat->lenSw[i]=pq;
896
897  return i;
898
899}
900static void length_one_crit(slimgb_alg* c, int pos, int len)
901{
902  if (c->nc)
903    return;
904  if (len==1)
905  {
906    int i;
907    for ( i=0;i<pos;i++)
908    {
909      if (c->lengths[i]==1)
910        c->states[pos][i]=HASTREP;
911    }
912    for ( i=pos+1;i<c->n;i++){
913      if (c->lengths[i]==1)
914        c->states[i][pos]=HASTREP;
915    }
916    if (!c->nc)
917      shorten_tails(c,c->S->m[pos]);
918  }
919}
920
921
922static void move_forward_in_S(int old_pos, int new_pos,kStrategy strat)
923{
924  assume(old_pos>=new_pos);
925  poly p=strat->S[old_pos];
926  int ecart=strat->ecartS[old_pos];
927  long sev=strat->sevS[old_pos];
928  int s_2_r=strat->S_2_R[old_pos];
929  int length=strat->lenS[old_pos];
930  assume(length==pLength(strat->S[old_pos]));
931  wlen_type length_w;
932  if(strat->lenSw!=NULL)
933    length_w=strat->lenSw[old_pos];
934  int i;
935  for (i=old_pos; i>new_pos; i--)
936  {
937    strat->S[i] = strat->S[i-1];
938    strat->ecartS[i] = strat->ecartS[i-1];
939    strat->sevS[i] = strat->sevS[i-1];
940    strat->S_2_R[i] = strat->S_2_R[i-1];
941  }
942  if (strat->lenS!=NULL)
943    for (i=old_pos; i>new_pos; i--)
944      strat->lenS[i] = strat->lenS[i-1];
945  if (strat->lenSw!=NULL)
946    for (i=old_pos; i>new_pos; i--)
947      strat->lenSw[i] = strat->lenSw[i-1];
948
949  strat->S[new_pos]=p;
950  strat->ecartS[new_pos]=ecart;
951  strat->sevS[new_pos]=sev;
952  strat->S_2_R[new_pos]=s_2_r;
953  strat->lenS[new_pos]=length;
954  if(strat->lenSw!=NULL)
955    strat->lenSw[new_pos]=length_w;
956  //assume(lenS_correct(strat));
957}
958
959static void move_backward_in_S(int old_pos, int new_pos,kStrategy strat)
960{
961  assume(old_pos<=new_pos);
962  poly p=strat->S[old_pos];
963  int ecart=strat->ecartS[old_pos];
964  long sev=strat->sevS[old_pos];
965  int s_2_r=strat->S_2_R[old_pos];
966  int length=strat->lenS[old_pos];
967  assume(length==pLength(strat->S[old_pos]));
968  wlen_type length_w;
969  if(strat->lenSw!=NULL)
970    length_w=strat->lenSw[old_pos];
971  int i;
972  for (i=old_pos; i<new_pos; i++)
973  {
974    strat->S[i] = strat->S[i+1];
975    strat->ecartS[i] = strat->ecartS[i+1];
976    strat->sevS[i] = strat->sevS[i+1];
977    strat->S_2_R[i] = strat->S_2_R[i+1];
978  }
979  if (strat->lenS!=NULL)
980    for (i=old_pos; i<new_pos; i++)
981      strat->lenS[i] = strat->lenS[i+1];
982  if (strat->lenSw!=NULL)
983    for (i=old_pos; i<new_pos; i++)
984      strat->lenSw[i] = strat->lenSw[i+1];
985
986  strat->S[new_pos]=p;
987  strat->ecartS[new_pos]=ecart;
988  strat->sevS[new_pos]=sev;
989  strat->S_2_R[new_pos]=s_2_r;
990  strat->lenS[new_pos]=length;
991  if(strat->lenSw!=NULL)
992    strat->lenSw[new_pos]=length_w;
993  //assume(lenS_correct(strat));
994}
995
996static int* make_connections(int from, int to, poly bound, slimgb_alg* c)
997{
998  ideal I=c->S;
999  int* cans=(int*) omalloc(c->n*sizeof(int));
1000  int* connected=(int*) omalloc(c->n*sizeof(int));
1001  cans[0]=to;
1002  int cans_length=1;
1003  connected[0]=from;
1004  int last_cans_pos=-1;
1005  int connected_length=1;
1006  long neg_bounds_short= ~p_GetShortExpVector(bound,c->r);
1007
1008  int not_yet_found=cans_length;
1009  int con_checked=0;
1010  int pos;
1011
1012  while(TRUE){
1013    if ((con_checked<connected_length)&& (not_yet_found>0)){
1014      pos=connected[con_checked];
1015      for(int i=0;i<cans_length;i++){
1016        if (cans[i]<0) continue;
1017        //FIXME: triv. syz. does not hold on noncommutative, check it for modules
1018        if ((has_t_rep(pos,cans[i],c)) ||((!rIsPluralRing(c->r))&&(trivial_syzygie(pos,cans[i],bound,c))))
1019{
1020
1021          connected[connected_length]=cans[i];
1022          connected_length++;
1023          cans[i]=-1;
1024          --not_yet_found;
1025
1026          if (connected[connected_length-1]==to){
1027            if (connected_length<c->n){
1028              connected[connected_length]=-1;
1029            }
1030            omfree(cans);
1031            return connected;
1032          }
1033        }
1034      }
1035      con_checked++;
1036    }
1037    else
1038    {
1039      for(last_cans_pos++;last_cans_pos<=c->n;last_cans_pos++){
1040        if (last_cans_pos==c->n){
1041          if (connected_length<c->n){
1042            connected[connected_length]=-1;
1043          }
1044          omfree(cans);
1045          return connected;
1046        }
1047        if ((last_cans_pos==from)||(last_cans_pos==to))
1048          continue;
1049        if(p_LmShortDivisibleBy(I->m[last_cans_pos],c->short_Exps[last_cans_pos],bound,neg_bounds_short,c->r)){
1050          cans[cans_length]=last_cans_pos;
1051          cans_length++;
1052          break;
1053        }
1054      }
1055      not_yet_found++;
1056      for (int i=0;i<con_checked;i++){
1057        if (has_t_rep(connected[i],last_cans_pos,c)){
1058
1059          connected[connected_length]=last_cans_pos;
1060          connected_length++;
1061          cans[cans_length-1]=-1;
1062
1063          --not_yet_found;
1064          if (connected[connected_length-1]==to){
1065            if (connected_length<c->n){
1066              connected[connected_length]=-1;
1067            }
1068
1069            omfree(cans);
1070            return connected;
1071          }
1072          break;
1073        }
1074      }
1075    }
1076  }
1077  if (connected_length<c->n){
1078    connected[connected_length]=-1;
1079  }
1080
1081  omfree(cans);
1082  return connected;
1083}
1084#ifdef HEAD_BIN
1085static inline poly p_MoveHead(poly p, omBin b)
1086{
1087  poly np;
1088  omTypeAllocBin(poly, np, b);
1089  memmove(np, p, b->sizeW*sizeof(long));
1090  omFreeBinAddr(p);
1091  return np;
1092}
1093#endif
1094
1095static void replace_pair(int & i, int & j,slimgb_alg* c)
1096{
1097  if (i<0) return;
1098  c->soon_free=NULL;
1099  int syz_deg;
1100  poly lm=pOne();
1101
1102  pLcm(c->S->m[i], c->S->m[j], lm);
1103  pSetm(lm);
1104
1105  int* i_con =make_connections(i,j,lm,c);
1106
1107  for (int n=0;((n<c->n) && (i_con[n]>=0));n++)
1108  {
1109    if (i_con[n]==j)
1110    {
1111      now_t_rep(i,j,c);
1112      omfree(i_con);
1113      p_Delete(&lm,c->r);
1114      return;
1115    }
1116  }
1117
1118  int* j_con =make_connections(j,i,lm,c);
1119
1120//   if(c->n>1){
1121//     if (i_con[1]>=0)
1122//       i=i_con[1];
1123//     else {
1124//       if (j_con[1]>=0)
1125//         j=j_con[1];
1126//     }
1127 // }
1128
1129  int sugar=c->pTotaldegree(lm);
1130  p_Delete(&lm, c->r);
1131  if(c->T_deg_full)//Sugar
1132  {
1133    int t_i=c->T_deg_full[i]-c->T_deg[i];
1134    int t_j=c->T_deg_full[j]-c->T_deg[j];
1135    sugar+=si_max(t_i,t_j);
1136    //Print("\n max: %d\n",max(t_i,t_j));
1137  }
1138
1139  for (int m=0;((m<c->n) && (i_con[m]>=0));m++)
1140  {
1141    if(c->T_deg_full!=NULL)
1142    {
1143        int s1=c->T_deg_full[i_con[m]]+syz_deg-c->T_deg[i_con[m]];
1144        if (s1>sugar) continue;
1145    }
1146    if (c->weighted_lengths[i_con[m]]<c->weighted_lengths[i])
1147        i=i_con[m];
1148  }
1149  for (int m=0;((m<c->n) && (j_con[m]>=0));m++)
1150  {
1151    if (c->T_deg_full!=NULL)
1152    {
1153      int s1=c->T_deg_full[j_con[m]]+syz_deg-c->T_deg[j_con[m]];
1154      if (s1>sugar) continue;
1155    }
1156    if (c->weighted_lengths[j_con[m]]<c->weighted_lengths[j])
1157      j=j_con[m];
1158  }
1159
1160  //can also try dependend search
1161  omfree(i_con);
1162  omfree(j_con);
1163}
1164
1165
1166static void add_later(poly p, const char* prot, slimgb_alg* c)
1167{
1168    int i=0;
1169    //check, if it is already in the queue
1170
1171
1172    while(c->add_later->m[i]!=NULL)
1173    {
1174        if (p_LmEqual(c->add_later->m[i],p,c->r))
1175            return;
1176        i++;
1177    }
1178    if (TEST_OPT_PROT)
1179        PrintS(prot);
1180    c->add_later->m[i]=p;
1181}
1182static int simple_posInS (kStrategy strat, poly p,int len, wlen_type wlen)
1183{
1184
1185
1186  if(strat->sl==-1) return 0;
1187  if (strat->lenSw) return pos_helper(strat,p,(wlen_type) wlen,(wlen_set) strat->lenSw,strat->S);
1188  return pos_helper(strat,p,len,strat->lenS,strat->S);
1189
1190}
1191
1192/*2
1193 *if the leading term of p
1194 *divides the leading term of some S[i] it will be canceled
1195 */
1196static inline void clearS (poly p, unsigned long p_sev,int l, int* at, int* k,
1197                           kStrategy strat)
1198{
1199  assume(p_sev == pGetShortExpVector(p));
1200  if (!pLmShortDivisibleBy(p,p_sev, strat->S[*at], ~ strat->sevS[*at])) return;
1201  if (l>=strat->lenS[*at]) return;
1202  if (TEST_OPT_PROT)
1203    PrintS("!");
1204  mflush();
1205  //pDelete(&strat->S[*at]);
1206  deleteInS((*at),strat);
1207  (*at)--;
1208  (*k)--;
1209//  assume(lenS_correct(strat));
1210}
1211
1212
1213
1214static int iq_crit(const void* ap,const void* bp)
1215{
1216  sorted_pair_node* a=*((sorted_pair_node**)ap);
1217  sorted_pair_node* b=*((sorted_pair_node**)bp);
1218  assume(a->i>a->j);
1219  assume(b->i>b->j);
1220
1221
1222  if (a->deg<b->deg) return -1;
1223  if (a->deg>b->deg) return 1;
1224  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
1225  if(comp!=0)
1226    return comp;
1227  if (a->expected_length<b->expected_length) return -1;
1228  if (a->expected_length>b->expected_length) return 1;
1229  if (a->j>b->j) return 1;
1230  if (a->j<b->j) return -1;
1231  return 0;
1232}
1233static wlen_type coeff_mult_size_estimate(int s1, int s2, ring r){
1234    if (rField_is_Q(r)) return s1+s2;
1235    else return s1*s2;
1236}
1237static wlen_type pair_weighted_length(int i, int j, slimgb_alg* c){
1238    if ((c->isDifficultField) && (c->eliminationProblem))  {
1239        int c1=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r);
1240        int c2=slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1241        wlen_type el1=c->weighted_lengths[i]/c1;
1242        assume(el1!=0);
1243        assume(c->weighted_lengths[i] %c1==0);
1244        wlen_type el2=c->weighted_lengths[j]/c2;
1245        assume(el2!=0);
1246        assume(c->weighted_lengths[j] %c2==0);
1247        //should be * for function fields
1248        //return (c1+c2) * (el1+el2-2);
1249        wlen_type res=coeff_mult_size_estimate(c1,c2,c->r);
1250        res*=el1+el2-2;
1251        return res;
1252
1253    }
1254    if (c->isDifficultField) {
1255        //int cs=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r)+
1256        //    slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1257        if(!(TEST_V_COEFSTRAT)){
1258        wlen_type cs=
1259            coeff_mult_size_estimate(
1260                slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r),
1261                slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r),c->r);
1262        return (wlen_type)(c->lengths[i]+c->lengths[j]-2)*
1263            (wlen_type)cs;}
1264            else {
1265
1266            wlen_type cs=
1267            coeff_mult_size_estimate(
1268                slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r),
1269                slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r),c->r);
1270            cs*=cs;
1271        return (wlen_type)(c->lengths[i]+c->lengths[j]-2)*
1272            (wlen_type)cs;
1273            }
1274    }
1275    if (c->eliminationProblem) {
1276
1277        return (c->weighted_lengths[i]+c->weighted_lengths[j]-2);
1278    }
1279    return c->lengths[i]+c->lengths[j]-2;
1280
1281}
1282sorted_pair_node** add_to_basis_ideal_quotient(poly h, slimgb_alg* c, int* ip)
1283{
1284  p_Test(h,c->r);
1285  assume(h!=NULL);
1286  poly got=gcd_of_terms(h,c->r);
1287  if((got!=NULL) &&(TEST_V_UPTORADICAL)) {
1288    poly copy=p_Copy(got,c->r);
1289    //p_wrp(got,c->r);
1290    BOOLEAN changed=monomial_root(got,c->r);
1291    if (changed)
1292    {
1293         poly div_by=pDivide(copy, got);
1294         poly iter=h;
1295         while(iter){
1296            pExpVectorSub(iter,div_by);
1297            pIter(iter);
1298         }
1299         p_Delete(&div_by, c->r);
1300         PrintS("U");
1301    }
1302    p_Delete(&copy,c->r);
1303  }
1304
1305#define ENLARGE(pointer, type) pointer=(type*) omrealloc(pointer, c->array_lengths*sizeof(type))
1306//  BOOLEAN corr=lenS_correct(c->strat);
1307  BOOLEAN R_found=FALSE;
1308  void* hp;
1309  int sugar;
1310  int ecart=0;
1311  ++(c->n);
1312  ++(c->S->ncols);
1313  int i,j;
1314  i=c->n-1;
1315  sorted_pair_node** nodes=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1316  int spc=0;
1317  if(c->n>c->array_lengths)
1318  {
1319    c->array_lengths=c->array_lengths*2;
1320    assume(c->array_lengths>=c->n);
1321    ENLARGE(c->T_deg, int);
1322    ENLARGE(c->tmp_pair_lm,poly);
1323    ENLARGE(c->tmp_spn,sorted_pair_node*);
1324
1325    ENLARGE(c->short_Exps,long);
1326    ENLARGE(c->lengths,int);
1327    #ifndef HAVE_BOOST
1328    #ifndef USE_STDVECBOOL
1329
1330    ENLARGE(c->states, char*);
1331    #endif
1332    #endif
1333    ENLARGE(c->gcd_of_terms,poly);
1334    //if (c->weighted_lengths!=NULL) {
1335    ENLARGE(c->weighted_lengths,wlen_type);
1336    //}
1337    //ENLARGE(c->S->m,poly);
1338
1339  }
1340  pEnlargeSet(&c->S->m,c->n-1,1);
1341  if (c->T_deg_full)
1342    ENLARGE(c->T_deg_full,int);
1343  sugar=c->T_deg[i]=c->pTotaldegree(h);
1344  if(c->T_deg_full)
1345  {
1346    sugar=c->T_deg_full[i]=c->pTotaldegree_full(h);
1347    ecart=sugar-c->T_deg[i];
1348    assume(ecart>=0);
1349  }
1350
1351
1352  c->tmp_pair_lm[i]=pOne_Special(c->r);
1353
1354
1355  c->tmp_spn[i]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1356
1357
1358  c->lengths[i]=pLength(h);
1359
1360  //necessary for correct weighted length
1361
1362  if (!rField_is_Zp(c->r))
1363  {
1364    p_Cleardenom(h, c->r);
1365    //p_Content(h,c->r); //is a duplicate call, but belongs here
1366
1367  }
1368  else
1369    pNorm(h);
1370  pNormalize(h);
1371
1372  c->weighted_lengths[i]=pQuality(h, c, c->lengths[i]);
1373  c->gcd_of_terms[i]=got;
1374  #ifdef HAVE_BOOST
1375    c->states.push_back(dynamic_bitset<>(i));
1376
1377  #else
1378  #ifdef USE_STDVECBOOL
1379
1380    c->states.push_back(vector<bool>(i));
1381
1382
1383  #else
1384  if (i>0)
1385    c->states[i]=(char*)  omalloc(i*sizeof(char));
1386  else
1387    c->states[i]=NULL;
1388  #endif
1389  #endif
1390
1391  c->S->m[i]=h;
1392  c->short_Exps[i]=p_GetShortExpVector(h,c->r);
1393
1394#undef ENLARGE
1395  if (p_GetComp(h,currRing)<=c->syz_comp){
1396  for (j=0;j<i;j++){
1397
1398
1399    #ifndef HAVE_BOOST
1400    c->states[i][j]=UNCALCULATED;
1401    #endif
1402    assume(p_LmDivisibleBy(c->S->m[i],c->S->m[j],c->r)==
1403     p_LmShortDivisibleBy(c->S->m[i],c->short_Exps[i],c->S->m[j],~(c->short_Exps[j]),c->r));
1404
1405
1406    if (_p_GetComp(c->S->m[i],c->r)!=_p_GetComp(c->S->m[j],c->r)){
1407      //c->states[i][j]=UNCALCULATED;
1408      //WARNUNG: be careful
1409      continue;
1410    } else
1411    if ((!c->nc) && (c->lengths[i]==1) && (c->lengths[j]==1)){
1412      c->states[i][j]=HASTREP;
1413
1414      }
1415    else if (( (!c->nc) || (c->is_homog && rIsSCA(c->r) ) ) &&  (pHasNotCF(c->S->m[i],c->S->m[j])))
1416//     else if ((!(c->nc)) &&  (pHasNotCF(c->S->m[i],c->S->m[j])))
1417    {
1418      c->easy_product_crit++;
1419      c->states[i][j]=HASTREP;
1420      continue;
1421    }
1422    else if(extended_product_criterion(c->S->m[i],c->gcd_of_terms[i],c->S->m[j],c->gcd_of_terms[j],c))
1423    {
1424      c->states[i][j]=HASTREP;
1425      c->extended_product_crit++;
1426
1427      //PrintS("E");
1428    }
1429      //  if (c->states[i][j]==UNCALCULATED){
1430
1431    if ((TEST_V_FINDMONOM) &&(!c->nc)) {
1432        //PrintS("COMMU");
1433       //  if (c->lengths[i]==c->lengths[j]){
1434//             poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1435//             if (short_s==NULL){
1436//                 c->states[i][j]=HASTREP;
1437//             } else
1438//             {
1439//                 p_Delete(&short_s, currRing);
1440//             }
1441//         }
1442        if (c->lengths[i]+c->lengths[j]==3){
1443
1444
1445             poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1446            if (short_s==NULL)
1447            {
1448                c->states[i][j]=HASTREP;
1449            }
1450            else
1451            {
1452                assume(pLength(short_s)==1);
1453                if (TEST_V_UPTORADICAL)
1454                   monomial_root(short_s,c->r);
1455                int iS=
1456                   kFindDivisibleByInS_easy(c->strat,short_s, p_GetShortExpVector(short_s,c->r));
1457                if (iS<0)
1458                {
1459                    //PrintS("N");
1460                    if (TRUE)
1461                    {
1462                      c->states[i][j]=HASTREP;
1463                      add_later(short_s,"N",c);
1464                    }
1465                    else
1466                      p_Delete(&short_s,currRing);
1467                }
1468                else
1469                {
1470                    if (c->strat->lenS[iS]>1){
1471                        //PrintS("O");
1472                        if (TRUE) {
1473                        c->states[i][j]=HASTREP;
1474                        add_later(short_s,"O",c);
1475                        } else p_Delete(&short_s,currRing);
1476                    }
1477                    else
1478                     p_Delete(&short_s, currRing);
1479                     c->states[i][j]=HASTREP;
1480                }
1481
1482
1483            }
1484        }
1485    }
1486      //    if (short_s)
1487      //    {
1488    assume(spc<=j);
1489    sorted_pair_node* s=c->tmp_spn[spc];//(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1490    s->i=si_max(i,j);
1491    s->j=si_min(i,j);
1492    assume(s->j==j);
1493    s->expected_length=pair_weighted_length(i,j,c);//c->lengths[i]+c->lengths[j]-2;
1494
1495    poly lm=c->tmp_pair_lm[spc];//=pOne_Special();
1496
1497    pLcm(c->S->m[i], c->S->m[j], lm);
1498    pSetm(lm);
1499    p_Test(lm,c->r);
1500    s->deg=c->pTotaldegree(lm);
1501
1502    if(c->T_deg_full)//Sugar
1503    {
1504      int t_i=c->T_deg_full[s->i]-c->T_deg[s->i];
1505      int t_j=c->T_deg_full[s->j]-c->T_deg[s->j];
1506      s->deg+=si_max(t_i,t_j);
1507      //Print("\n max: %d\n",max(t_i,t_j));
1508    }
1509    p_Test(lm,c->r);
1510    s->lcm_of_lm=lm;
1511    //          pDelete(&short_s);
1512    //assume(lm!=NULL);
1513    nodes[spc]=s;
1514    spc++;
1515
1516  // }
1517  //else
1518  //{
1519        //c->states[i][j]=HASTREP;
1520  //}
1521  }
1522  }//if syz_comp end
1523
1524
1525
1526
1527  assume(spc<=i);
1528  //now ideal quotient crit
1529  qsort(nodes,spc,sizeof(sorted_pair_node*),iq_crit);
1530
1531    sorted_pair_node** nodes_final=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1532  int spc_final=0;
1533  j=0;
1534  while(j<spc)
1535  {
1536    int lower=j;
1537    int upper;
1538    BOOLEAN has=FALSE;
1539    for(upper=lower+1;upper<spc;upper++)
1540    {
1541      if(!pLmEqual(nodes[lower]->lcm_of_lm,nodes[upper]->lcm_of_lm))
1542      {
1543        break;
1544      }
1545      if (has_t_rep(nodes[upper]->i,nodes[upper]->j,c))
1546        has=TRUE;
1547    }
1548    upper=upper-1;
1549    int z;
1550    assume(spc_final<=j);
1551    for(z=0;z<spc_final;z++)
1552    {
1553      if(p_LmDivisibleBy(nodes_final[z]->lcm_of_lm,nodes[lower]->lcm_of_lm,c->r))
1554      {
1555  has=TRUE;
1556  break;
1557      }
1558    }
1559
1560    if(has)
1561    {
1562      for(;lower<=upper;lower++)
1563      {
1564  //free_sorted_pair_node(nodes[lower],c->r);
1565  //omfree(nodes[lower]);
1566  nodes[lower]=NULL;
1567      }
1568      j=upper+1;
1569      continue;
1570    }
1571    else
1572    {
1573      p_Test(nodes[lower]->lcm_of_lm,c->r);
1574      nodes[lower]->lcm_of_lm=pCopy(nodes[lower]->lcm_of_lm);
1575      assume(_p_GetComp(c->S->m[nodes[lower]->i],c->r)==_p_GetComp(c->S->m[nodes[lower]->j],c->r));
1576      nodes_final[spc_final]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1577
1578      *(nodes_final[spc_final++])=*(nodes[lower]);
1579      //c->tmp_spn[nodes[lower]->j]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1580      nodes[lower]=NULL;
1581      for(lower=lower+1;lower<=upper;lower++)
1582      {
1583  //  free_sorted_pair_node(nodes[lower],c->r);
1584  //omfree(nodes[lower]);
1585  nodes[lower]=NULL;
1586      }
1587      j=upper+1;
1588      continue;
1589    }
1590  }
1591
1592  //  Print("i:%d,spc_final:%d",i,spc_final);
1593
1594
1595
1596
1597  assume(spc_final<=spc);
1598  omfree(nodes);
1599  nodes=NULL;
1600
1601  add_to_reductors(c, h, c->lengths[c->n-1], ecart,TRUE);
1602  //i=posInS(c->strat,c->strat->sl,h,0 ecart);
1603  if (!(c->nc)){
1604    if (c->lengths[c->n-1]==1)
1605      shorten_tails(c,c->S->m[c->n-1]);
1606  }
1607  //you should really update c->lengths, c->strat->lenS, and the oder of polys in strat if you sort after lengths
1608
1609  //for(i=c->strat->sl; i>0;i--)
1610  //  if(c->strat->lenS[i]<c->strat->lenS[i-1]) printf("fehler bei %d\n",i);
1611  if (c->Rcounter>50) {
1612    c->Rcounter=0;
1613    cleanS(c->strat,c);
1614  }
1615
1616#ifdef HAVE_PLURAL
1617  // for SCA:
1618  // here write at the end of nodes_final[spc_final,...,spc_final+lmdeg-1]
1619  if(rIsSCA(c->r))
1620  {
1621    const poly pNext = pNext(h);
1622
1623    if(pNext != NULL)
1624    {
1625      // for additional polynomials
1626      const unsigned int m_iFirstAltVar = scaFirstAltVar(c->r);
1627      const unsigned int m_iLastAltVar  = scaLastAltVar(c->r);
1628
1629      int N = // c->r->N;
1630              m_iLastAltVar - m_iFirstAltVar + 1; // should be enough
1631      // TODO: but we may also use got = gcd({m}_{m\in f}))!
1632
1633       poly* array_arg=(poly*)omalloc(N*sizeof(poly)); // !
1634       int j = 0;
1635
1636
1637      for( unsigned short v = m_iFirstAltVar; v <= m_iLastAltVar; v++ )
1638      // for all x_v | Ann(lm(h))
1639      if( p_GetExp(h, v, c->r) ) // TODO: use 'got' here!
1640      {
1641        assume(p_GetExp(h, v, c->r)==1);
1642
1643        poly p = sca_pp_Mult_xi_pp(v, pNext, c->r); // x_v * h;
1644
1645        if(p != NULL) // if (x_v * h != 0)
1646          array_arg[j++] = p;
1647      } // for all x_v | Ann(lm(h))
1648
1649      c->introduceDelayedPairs(array_arg, j);
1650
1651      omfree(array_arg); // !!!
1652    }
1653//     PrintS("Saturation - done!!!\n");
1654  }
1655#endif // if SCAlgebra
1656
1657
1658  if(!ip){
1659    qsort(nodes_final,spc_final,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1660
1661
1662    c->apairs=spn_merge(c->apairs,c->pair_top+1,nodes_final,spc_final,c);
1663    c->pair_top+=spc_final;
1664    clean_top_of_pair_list(c);
1665    omfree(nodes_final);
1666    return NULL;
1667  }
1668  {
1669    *ip=spc_final;
1670    return nodes_final;
1671  }
1672
1673
1674
1675}
1676
1677
1678static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
1679{
1680  m=nInit(1);
1681  if (h==NULL) return NULL;
1682
1683  assume(len==pLength(h));
1684  kStrategy strat=c->strat;
1685  if (0 > strat->sl)
1686  {
1687    return h;
1688  }
1689  int j;
1690
1691  LObject P(h);
1692  P.SetShortExpVector();
1693  P.bucket = kBucketCreate(currRing);
1694  // BOOLEAN corr=lenS_correct(strat);
1695  kBucketInit(P.bucket,P.p,len /*pLength(P.p)*/);
1696  //wlen_set lenSw=(wlen_set) c->strat->lenS;
1697  //FIXME: plainly wrong
1698  //strat->lenS;
1699  //if (strat->lenSw!=NULL)
1700  //  lenSw=strat->lenSw;
1701  //int max_pos=simple_posInS(strat,P.p);
1702  loop
1703  {
1704      int dummy=strat->sl;
1705      j=kFindDivisibleByInS_easy(strat,P.p,P.sev);
1706      //j=kFindDivisibleByInS(strat,&dummy,&P);
1707      if ((j>=0) && ((!n)||
1708        ((strat->lenS[j]<=n) &&
1709         ((strat->lenSw==NULL)||
1710         (strat->lenSw[j]<=n)))))
1711      {
1712        nNormalize(pGetCoeff(P.p));
1713#ifdef KDEBUG
1714        if (TEST_OPT_DEBUG)
1715        {
1716          PrintS("red:");
1717          wrp(h);
1718          PrintS(" with ");
1719          wrp(strat->S[j]);
1720        }
1721#endif
1722
1723        number coef=kBucketPolyRed(P.bucket,strat->S[j],
1724                                   strat->lenS[j]/*pLength(strat->S[j])*/,
1725                                   strat->kNoether);
1726  number m2=nMult(m,coef);
1727  nDelete(&m);
1728  m=m2;
1729        nDelete(&coef);
1730        h = kBucketGetLm(P.bucket);
1731
1732  if (h==NULL) {
1733    len=0;
1734    kBucketDestroy(&P.bucket);
1735    return
1736    NULL;}
1737        P.p=h;
1738        P.t_p=NULL;
1739        P.SetShortExpVector();
1740#ifdef KDEBUG
1741        if (TEST_OPT_DEBUG)
1742        {
1743          PrintS("\nto:");
1744          wrp(h);
1745          PrintLn();
1746        }
1747#endif
1748      }
1749      else
1750      {
1751        kBucketClear(P.bucket,&(P.p),&len);
1752        kBucketDestroy(&P.bucket);
1753        pNormalize(P.p);
1754  assume(len==(pLength(P.p)));
1755        return P.p;
1756      }
1757    }
1758}
1759
1760
1761
1762static poly redTailShort(poly h, kStrategy strat){
1763  if (h==NULL) return NULL;//n_Init(1,currRing);
1764  if (TEST_V_MODPSOLVSB){
1765    bit_reduce(pNext(h), strat->tailRing);
1766  }
1767  int sl=strat->sl;
1768  int i;
1769  int len=pLength(h);
1770  for(i=0;i<=strat->sl;i++){
1771    if((strat->lenS[i]>2) || ((strat->lenSw!=NULL) && (strat->lenSw[i]>2)))
1772      break;
1773  }
1774  return(redNFTail(h,i-1,strat, len));
1775}
1776
1777static void line_of_extended_prod(int fixpos,slimgb_alg* c){
1778    if (c->gcd_of_terms[fixpos]==NULL)
1779  {
1780    c->gcd_of_terms[fixpos]=gcd_of_terms(c->S->m[fixpos],c->r);
1781    if (c->gcd_of_terms[fixpos])
1782    {
1783      int i;
1784      for(i=0;i<fixpos;i++)
1785        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)))
1786{
1787          c->states[fixpos][i]=HASTREP;
1788    c->extended_product_crit++;
1789}
1790      for(i=fixpos+1;i<c->n;i++)
1791        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)))
1792  {        c->states[i][fixpos]=HASTREP;
1793  c->extended_product_crit++;
1794  }
1795    }
1796  }
1797}
1798static void c_S_element_changed_hook(int pos, slimgb_alg* c){
1799  length_one_crit(c,pos, c->lengths[pos]);
1800  if (!c->nc)
1801    line_of_extended_prod(pos,c);
1802}
1803class poly_tree_node {
1804public:
1805  poly p;
1806  poly_tree_node* l;
1807  poly_tree_node* r;
1808  int n;
1809  poly_tree_node(int sn):l(NULL),r(NULL),n(sn){}
1810};
1811class exp_number_builder{
1812public:
1813  poly_tree_node* top_level;
1814  int n;
1815  int get_n(poly p);
1816  exp_number_builder():top_level(0),n(0){}
1817};
1818int exp_number_builder::get_n(poly p){
1819  poly_tree_node** node=&top_level;
1820  while(*node!=NULL){
1821    int c=pLmCmp(p,(*node)->p);
1822    if (c==0) return (*node)->n;
1823    if (c==-1) node=&((*node)->r);
1824    else
1825      node=&((*node)->l);
1826  }
1827  (*node)= new poly_tree_node(n);
1828  n++;
1829  (*node)->p=pLmInit(p);
1830  return (*node)->n;
1831}
1832
1833//mac_polys exp are smaller iff they are greater by monomial ordering
1834//corresponding to solving linear equations notation
1835
1836//! obsolete
1837struct int_poly_pair{
1838  poly p;
1839  int n;
1840};
1841
1842
1843//! obsolete
1844void t2ippa_rec(poly* ip,int* ia, poly_tree_node* k, int &offset){
1845    if(!k) return;
1846    t2ippa_rec(ip,ia,k->l,offset);
1847    ip[offset]=k->p;
1848    ia[k->n]=offset;
1849    ++offset;
1850
1851    t2ippa_rec(ip,ia,k->r,offset);
1852    delete k;
1853  }
1854
1855//! obsolete
1856void t2ippa(poly* ip,int* ia,exp_number_builder & e){
1857
1858  int o=0;
1859  t2ippa_rec(ip,ia,e.top_level,o);
1860}
1861int anti_poly_order(const void* a, const void* b){
1862  return -pLmCmp(((int_poly_pair*) a)->p,((int_poly_pair*) b)->p );
1863}
1864
1865BOOLEAN is_valid_ro(red_object & ro){
1866  red_object r2=ro;
1867  ro.validate();
1868  if ((r2.p!=ro.p)||(r2.sev!=ro.sev)) return FALSE;
1869  return TRUE;
1870}
1871int terms_sort_crit(const void* a, const void* b){
1872  return -pLmCmp(*((poly*) a),*((poly*) b));
1873}
1874static void unify_terms(poly* terms,int & sum){
1875  if (sum==0) return;
1876  int last=0;
1877  int curr=1;
1878  while(curr<sum){
1879    if (!(pLmEqual(terms[curr],terms[last]))){
1880      terms[++last]=terms[curr];
1881    }
1882    ++curr;
1883  }
1884  sum=last+1;
1885}
1886static void export_mat(number* number_array,int pn, int tn,const char* format_str, int mat_nr){
1887  char matname[20];
1888  sprintf(matname,format_str,mat_nr);
1889  FILE* out=fopen(matname,"w");
1890  int i,j;
1891  fprintf(out,"mat=[\n");
1892  for(i=0;i<pn;i++){
1893    fprintf(out,"[\n");
1894    for(j=0;j<tn;j++)
1895    {
1896      if (j>0)
1897      {
1898        fprintf(out,", ");
1899      }
1900      fprintf(out,"%i",npInt(number_array[i*tn+j],currRing));
1901    }
1902    if (i<pn-1)
1903      fprintf(out,"],\n");
1904    else
1905      fprintf(out,"],\n");
1906  }
1907  fprintf(out,"]\n");
1908  fclose(out);
1909}
1910//typedef unsigned short number_type;
1911
1912
1913#ifdef USE_NORO
1914#ifndef NORO_CACHE
1915static void linalg_step_modp(poly *p, poly* p_out, int& pn, poly* terms,int tn, slimgb_alg* c){
1916  static int export_n=0;
1917  assume(terms[tn-1]!=NULL);
1918  assume(rField_is_Zp(c->r));
1919  //I don't do deletes, copies of number_types ...
1920  const number_type zero=0;//npInit(0);
1921  int array_size=pn*tn;
1922  number_type* number_array=(number_type*) omalloc(pn*tn*sizeof(number_type));
1923  int i;
1924  for(i=0;i<array_size;i++){
1925    number_array[i]=zero;
1926  }
1927  for(i=0;i<pn;i++){
1928    poly h=p[i];
1929    //int base=tn*i;
1930    write_poly_to_row(number_array+tn*i,h,terms,tn,c->r);
1931
1932  }
1933#if 0
1934  //export matrix
1935  export_mat(number_array,pn,tn,"mat%i.py",++export_n);
1936#endif
1937  int rank=pn;
1938  simplest_gauss_modp(number_array,rank,tn);
1939  int act_row=0;
1940  int p_pos=0;
1941  for(i=0;i<pn;i++){
1942    poly h=NULL;
1943    int j;
1944    int base=tn*i;
1945    number* row=number_array+base;
1946    h=row_to_poly(row,terms,tn,c->r);
1947
1948   if (h!=NULL)
1949   {
1950     p_out[p_pos++]=h;
1951   }
1952  }
1953  pn=p_pos;
1954  //assert(p_pos==rank)
1955  while(p_pos<pn){
1956    p_out[p_pos++]=NULL;
1957  }
1958#if 0
1959  export_mat(number_array,pn,tn,"mat%i.py",++export_n);
1960#endif
1961}
1962#endif
1963#endif
1964static void mass_add(poly* p, int pn,slimgb_alg* c){
1965    int j;
1966    int* ibuf=(int*) omalloc(pn*sizeof(int));
1967    sorted_pair_node*** sbuf=(sorted_pair_node***) omalloc(pn*sizeof(sorted_pair_node**));
1968    for(j=0;j<pn;j++){
1969      p_Test(p[j],c->r);
1970      sbuf[j]=add_to_basis_ideal_quotient(p[j],c,ibuf+j);
1971    }
1972    int sum=0;
1973    for(j=0;j<pn;j++){
1974      sum+=ibuf[j];
1975    }
1976    sorted_pair_node** big_sbuf=(sorted_pair_node**) omalloc(sum*sizeof(sorted_pair_node*));
1977    int partsum=0;
1978    for(j=0;j<pn;j++)
1979    {
1980      memmove(big_sbuf+partsum, sbuf[j],ibuf[j]*sizeof(sorted_pair_node*));
1981      omfree(sbuf[j]);
1982      partsum+=ibuf[j];
1983    }
1984
1985    qsort(big_sbuf,sum,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1986    c->apairs=spn_merge(c->apairs,c->pair_top+1,big_sbuf,sum,c);
1987    c->pair_top+=sum;
1988    clean_top_of_pair_list(c);
1989    omfree(big_sbuf);
1990    omfree(sbuf);
1991    omfree(ibuf);
1992    //omfree(buf);
1993  #ifdef TGB_DEBUG
1994    int z;
1995    for(z=1;z<=c->pair_top;z++)
1996    {
1997      assume(pair_better(c->apairs[z],c->apairs[z-1],c));
1998    }
1999  #endif
2000
2001}
2002
2003#ifdef NORO_CACHE
2004#ifndef NORO_NON_POLY
2005void NoroCache::evaluateRows(){
2006  //after that can evaluate placeholders
2007  int i;
2008  buffer=(number*) omalloc(nIrreducibleMonomials*sizeof(number));
2009  for(i=0;i<root.branches_len;i++){
2010    evaluateRows(1,root.branches[i]);
2011  }
2012  omfree(buffer);
2013  buffer=NULL;
2014}
2015void NoroCache::evaluateRows(int level, NoroCacheNode* node){
2016  assume(level>=0);
2017  if (node==NULL) return;
2018  if (level<pVariables){
2019    int i,sum;
2020    for(i=0;i<node->branches_len;i++){
2021      evaluateRows(level+1,node->branches[i]);
2022    }
2023  } else {
2024    DataNoroCacheNode* dn=(DataNoroCacheNode*) node;
2025    if (dn->value_len!=backLinkCode){
2026      poly p=dn->value_poly;
2027      #ifndef NORO_SPARSE_ROWS_PRE
2028      dn->row=new DenseRow();
2029      DenseRow* row=dn->row;
2030      memset(buffer,0,sizeof(number)*nIrreducibleMonomials);
2031
2032      if (p==NULL) {row->array=NULL;row->begin=0;row->end=0; return;}
2033      int i=0;
2034      int idx;
2035      number* a=buffer;
2036      while(p){
2037        DataNoroCacheNode* ref=getCacheReference(p);
2038
2039        idx=ref->term_index;
2040        assume(idx>=0);
2041        a[idx]=p_GetCoeff(p,currRing);
2042        if (i==0) row->begin=idx;
2043        i++;
2044        pIter(p);
2045      }
2046      row->end=idx+1;
2047      assume(row->end>row->begin);
2048      int len=row->end-row->begin;
2049      row->array=(number*) omalloc((len)*sizeof(number));
2050      memcpy(row->array,a+row->begin,len*sizeof(number));
2051      #else
2052      assume(dn->value_len==pLength(dn->value_poly));
2053      dn->row=new SparseRow(dn->value_len);
2054      SparseRow* row=dn->row;
2055      int i=0;
2056      while(p){
2057        DataNoroCacheNode* ref=getCacheReference(p);
2058
2059        int idx=ref->term_index;
2060        assume(idx>=0);
2061        row->idx_array[i]=idx;
2062        row->coef_array[i]=p_GetCoeff(p,currRing);
2063
2064        i++;
2065        pIter(p);
2066      }
2067      if (i!=dn->value_len){
2068        PrintS("F4 calc wrong, as poly len was wrong\n");
2069      }
2070      assume(i==dn->value_len);
2071      #endif
2072
2073    }
2074  }
2075}
2076
2077void NoroCache::evaluatePlaceHolder(number* row,std::vector<NoroPlaceHolder>& place_holders){
2078  int i;
2079  int s=place_holders.size();
2080  for(i=0;i<s;i++){
2081    DataNoroCacheNode* ref=place_holders[i].ref;
2082    number coef=place_holders[i].coef;
2083    if (ref->value_len==backLinkCode){
2084      row[ref->term_index]=npAddM(row[ref->term_index],coef);
2085    } else {
2086      #ifndef NORO_SPARSE_ROWS_PRE
2087      DenseRow* ref_row=ref->row;
2088      if (ref_row==NULL) continue;
2089      number* ref_begin=ref_row->array;
2090      number* ref_end=ref_row->array+(ref_row->end-ref_row->begin);
2091      number* my_pos=row+ref_row->begin;
2092      //TODO npisOne distinction
2093      if (!(npIsOne(coef))){
2094        while(ref_begin!=ref_end){
2095
2096          *my_pos=npAddM(*my_pos,npMult(coef,*ref_begin));
2097          ++ref_begin;
2098          ++my_pos;
2099        }
2100      }
2101      else{
2102        while(ref_begin!=ref_end){
2103
2104          *my_pos=npAddM(*my_pos,*ref_begin);
2105          ++ref_begin;
2106          ++my_pos;
2107        }
2108      }
2109
2110
2111
2112    #else
2113    SparseRow* ref_row=ref->row;
2114    if (ref_row==NULL) continue;
2115    int n=ref_row->len;
2116    int j;
2117    int* idx_array=ref_row->idx_array;
2118    number* coef_array=ref_row->coef_array;
2119    for(j=0;j<n;j++){
2120      int idx=idx_array[j];
2121      number ref_coef=coef_array[j];
2122      row[idx]=npAddM(row[idx],npMult(coef,ref_coef));
2123    }
2124    #endif
2125  }
2126  }
2127
2128}
2129#endif
2130
2131
2132
2133//poly noro_red_non_unique(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2134
2135
2136#ifndef NORO_NON_POLY
2137MonRedRes noro_red_mon(poly t, BOOLEAN force_unique, NoroCache* cache,slimgb_alg* c){
2138  MonRedRes res_holder;
2139
2140  //wrp(t);
2141  res_holder.changed=TRUE;
2142  if (force_unique){
2143    DataNoroCacheNode* ref=cache->getCacheReference(t);
2144    if (ref!=NULL){
2145      res_holder.len=ref->value_len;
2146      if (res_holder.len==NoroCache::backLinkCode){
2147        res_holder.len=1;
2148
2149
2150      }
2151      res_holder.coef=p_GetCoeff(t,c->r);
2152      res_holder.p=ref->value_poly;
2153      res_holder.ref=ref;
2154      res_holder.onlyBorrowed=TRUE;
2155      res_holder.changed=TRUE;
2156      p_Delete(&t,c->r);
2157      return res_holder;
2158    }
2159  } else{
2160    BOOLEAN succ;
2161    poly cache_lookup=cache->lookup(t,succ, res_holder.len);//don't own this yet
2162    if (succ){
2163      if (cache_lookup==t){
2164      //know they are equal
2165      //res_holder.len=1;
2166
2167        res_holder.changed=FALSE;
2168        res_holder.p=t;
2169        res_holder.coef=npInit(1);
2170
2171        res_holder.onlyBorrowed=FALSE;
2172        return res_holder;
2173      }
2174
2175      res_holder.coef=p_GetCoeff(t,c->r);
2176      p_Delete(&t,c->r);
2177
2178      res_holder.p=cache_lookup;
2179
2180      res_holder.onlyBorrowed=TRUE;
2181      return res_holder;
2182
2183    }
2184  }
2185
2186  unsigned long sev=p_GetShortExpVector(t,currRing);
2187  int i=kFindDivisibleByInS_easy(c->strat,t,sev);
2188  if (i>=0)
2189  {
2190    number coef_bak=p_GetCoeff(t,c->r);
2191
2192    p_SetCoeff(t,npInit(1),c->r);
2193    assume(npIsOne(p_GetCoeff(c->strat->S[i],c->r)));
2194    number coefstrat=p_GetCoeff(c->strat->S[i],c->r);
2195
2196      //poly t_copy_mon=p_Copy(t,c->r);
2197    poly exp_diff=cache->temp_term;
2198    p_ExpVectorDiff(exp_diff,t,c->strat->S[i],c->r);
2199    p_SetCoeff(exp_diff,npNeg(nInvers(coefstrat)),c->r);
2200      // nInvers may be npInvers or nvInvers
2201    p_Setm(exp_diff,c->r);
2202    assume(c->strat->S[i]!=NULL);
2203      //poly t_to_del=t;
2204    poly res;
2205    res=pp_Mult_mm(pNext(c->strat->S[i]),exp_diff,c->r);
2206
2207    res_holder.len=c->strat->lenS[i]-1;
2208    res=noro_red_non_unique(res,res_holder.len,cache,c);
2209
2210    DataNoroCacheNode* ref=cache->insert(t,res,res_holder.len);
2211    p_Delete(&t,c->r);
2212      //p_Delete(&t_copy_mon,c->r);
2213      //res=pMult_nn(res,coef_bak);
2214
2215    res_holder.changed=TRUE;
2216    res_holder.p=res;
2217    res_holder.coef=coef_bak;
2218    res_holder.onlyBorrowed=TRUE;
2219    res_holder.ref=ref;
2220    return res_holder;
2221
2222  } else {
2223    number coef_bak=p_GetCoeff(t,c->r);
2224    number one=npInit(1);
2225    p_SetCoeff(t,one,c->r);
2226    res_holder.len=1;
2227    if (!(force_unique)){
2228      res_holder.ref=cache->insert(t,t,res_holder.len);
2229      p_SetCoeff(t,coef_bak,c->r);
2230      //return t;
2231
2232      //we need distinction
2233      res_holder.changed=FALSE;
2234      res_holder.p=t;
2235
2236      res_holder.coef=npInit(1);
2237      res_holder.onlyBorrowed=FALSE;
2238      return res_holder;
2239    } else {
2240      res_holder.ref=cache->insertAndTransferOwnerShip(t,c->r);
2241      res_holder.coef=coef_bak;
2242      res_holder.onlyBorrowed=TRUE;
2243      res_holder.changed=FALSE;
2244      res_holder.p=t;
2245      return res_holder;
2246    }
2247  }
2248
2249}
2250#endif
2251//SparseRow* noro_red_to_non_poly(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2252#ifndef NORO_NON_POLY
2253//len input and out: Idea: reverse addition
2254poly noro_red_non_unique(poly p, int &len, NoroCache* cache,slimgb_alg* c){
2255  assume(len==pLength(p));
2256  poly orig_p=p;
2257  if (p==NULL) {
2258    len=0;
2259    return NULL;
2260  }
2261  kBucket_pt bucket=kBucketCreate(currRing);
2262  kBucketInit(bucket,NULL,0);
2263  poly unchanged_head=NULL;
2264  poly unchanged_tail=NULL;
2265  int unchanged_size=0;
2266
2267  while(p){
2268    poly t=p;
2269    pIter(p);
2270    pNext(t)=NULL;
2271#ifndef NDEBUG
2272    number coef_debug=p_GetCoeff(t,currRing);
2273#endif
2274    MonRedRes red=noro_red_mon(t,FALSE,cache,c);
2275    if ((!(red.changed))&&(!(red.onlyBorrowed))){
2276      unchanged_size++;
2277      assume(npIsOne(red.coef));
2278      assume(p_GetCoeff(red.p,currRing)==coef_debug);
2279      if (unchanged_head){
2280        pNext(unchanged_tail)=red.p;
2281        pIter(unchanged_tail);
2282      } else{
2283        unchanged_tail=red.p;
2284        unchanged_head=red.p;
2285      }
2286    } else{
2287      assume(red.len==pLength(red.p));
2288      if (red.onlyBorrowed){
2289        if (npIsOne(red.coef)){
2290          t=p_Copy(red.p,currRing);
2291        }else
2292        t=pp_Mult_nn(red.p,red.coef,currRing);
2293      } else {
2294        if (npIsOne(red.coef))
2295        t=red.p;
2296        else
2297         t=p_Mult_nn(red.p,red.coef,currRing);
2298      }
2299      kBucket_Add_q(bucket,t,&red.len);
2300    }
2301  }
2302  poly res=NULL;
2303  len=0;
2304  kBucket_Add_q(bucket,unchanged_head,&unchanged_size);
2305  kBucketClear(bucket,&res,&len);
2306  kBucketDestroy(&bucket);
2307
2308
2309
2310
2311  return res;
2312}
2313#endif
2314#ifdef NORO_SPARSE_ROWS_PRE
2315//len input and out: Idea: reverse addition
2316
2317/*template <class number_type> SparseRow<number_type>* noro_red_to_non_poly(poly p, int &len, NoroCache<number_type>* cache,slimgb_alg* c){
2318  if (npPrimeM<255){
2319    return noro_red_to_non_poly_t<tgb_uint8>(p,len,cache,c);
2320  } else {
2321    if (npPrimeM<65000){
2322      return noro_red_to_non_poly_t<tgb_uint16>(p,len,cache,c);
2323    } else{
2324      return noro_red_to_non_poly_t<tgb_uint32>(p,len,cache,c);
2325    }
2326  }
2327}*/
2328#endif
2329//len input and out: Idea: reverse addition
2330#ifndef NORO_NON_POLY
2331std::vector<NoroPlaceHolder> noro_red(poly p, int &len, NoroCache* cache,slimgb_alg* c){
2332  std::vector<NoroPlaceHolder> res;
2333   while(p){
2334      poly t=p;
2335      pIter(p);
2336      pNext(t)=NULL;
2337
2338      MonRedRes red=noro_red_mon(t,TRUE,cache,c);
2339      assume(red.onlyBorrowed);
2340      assume(red.coef);
2341      assume(red.ref);
2342      NoroPlaceHolder h;
2343      h.ref=red.ref;
2344      h.coef=red.coef;
2345      assume(!((h.ref->value_poly==NULL) &&(h.ref->value_len!=0)));
2346      if (h.ref->value_poly)
2347        res.push_back(h);
2348    }
2349    return res;
2350}
2351#endif
2352
2353
2354
2355
2356
2357#endif
2358#ifdef USE_NORO
2359#ifndef NORO_CACHE
2360void noro_step(poly*p,int &pn,slimgb_alg* c){
2361  poly* reduced=(poly*) omalloc(pn*sizeof(poly));
2362  int j;
2363  int* reduced_len=(int*) omalloc(pn*sizeof(int));
2364  int reduced_c=0;
2365  //if (TEST_OPT_PROT)
2366  //  PrintS("reduced system:\n");
2367#ifdef NORO_CACHE
2368  NoroCache cache;
2369#endif
2370  for(j=0;j<pn;j++){
2371
2372    poly h=p[j];
2373    int h_len=pLength(h);
2374
2375    number coef;
2376#ifndef NORO_CACHE
2377    h=redNF2(p_Copy(h,c->r),c,h_len,coef,0);
2378#else
2379    h=noro_red(p_Copy(h,c->r),h_len,&cache,c);
2380    assume(pLength(h)==h_len);
2381#endif
2382    if (h!=NULL){
2383#ifndef NORO_CACHE
2384
2385      h=redNFTail(h,c->strat->sl,c->strat,h_len);
2386      h_len=pLength(h);
2387#endif
2388
2389
2390      reduced[reduced_c]=h;
2391      reduced_len[reduced_c]=h_len;
2392      reduced_c++;
2393      if (TEST_OPT_PROT)
2394        Print("%d ",h_len);
2395    }
2396  }
2397  int reduced_sum=0;
2398  for(j=0;j<reduced_c;j++){
2399    reduced_sum+=reduced_len[j];
2400  }
2401  poly* terms=(poly*) omalloc(reduced_sum*sizeof(poly));
2402  int tc=0;
2403  for(j=0;j<reduced_c;j++){
2404    poly h=reduced[j];
2405
2406    while(h!=NULL){
2407      terms[tc++]=h;
2408      pIter(h);
2409      assume(tc<=reduced_sum);
2410    }
2411  }
2412  assume(tc==reduced_sum);
2413  qsort(terms,reduced_sum,sizeof(poly),terms_sort_crit);
2414  int nterms=reduced_sum;
2415  //if (TEST_OPT_PROT)
2416  //Print("orig estimation:%i\n",reduced_sum);
2417  unify_terms(terms,nterms);
2418  //if (TEST_OPT_PROT)
2419  //    Print("actual number of columns:%i\n",nterms);
2420  int rank=reduced_c;
2421  linalg_step_modp(reduced, p,rank,terms,nterms,c);
2422  omfree(terms);
2423
2424  pn=rank;
2425  omfree(reduced);
2426
2427  if (TEST_OPT_PROT)
2428    PrintS("\n");
2429}
2430#else
2431
2432
2433
2434#endif
2435#endif
2436static void go_on (slimgb_alg* c){
2437  //set limit of 1000 for multireductions, at the moment for
2438  //programming reasons
2439  #ifdef USE_NORO
2440  //Print("module rank%d\n",c->S->rank);
2441    const BOOLEAN use_noro=c->use_noro;
2442  #else
2443  const BOOLEAN use_noro=FALSE;
2444  #endif
2445  int i=0;
2446  c->average_length=0;
2447  for(i=0;i<c->n;i++){
2448    c->average_length+=c->lengths[i];
2449  }
2450  c->average_length=c->average_length/c->n;
2451  i=0;
2452  int max_pairs=bundle_size;
2453
2454  #ifdef USE_NORO
2455  if ((use_noro)||(c->use_noro_last_block))
2456    max_pairs=bundle_size_noro;
2457  #endif
2458  poly* p=(poly*) omalloc((max_pairs+1)*sizeof(poly));//nullterminated
2459
2460  int curr_deg=-1;
2461  while(i<max_pairs){
2462    sorted_pair_node* s=top_pair(c);//here is actually chain criterium done
2463
2464
2465    if (!s) break;
2466
2467    if(curr_deg>=0){
2468      if (s->deg >curr_deg) break;
2469    }
2470
2471    else curr_deg=s->deg;
2472    quick_pop_pair(c);
2473    if(s->i>=0){
2474      //be careful replace_pair use createShortSpoly which is not noncommutative
2475      now_t_rep(s->i,s->j,c);
2476    replace_pair(s->i,s->j,c);
2477
2478    if(s->i==s->j) {
2479      free_sorted_pair_node(s,c->r);
2480      continue;
2481    }
2482    now_t_rep(s->i,s->j,c);
2483    }
2484    poly h;
2485    if(s->i>=0)
2486    {
2487#ifdef HAVE_PLURAL
2488      if (c->nc)
2489      {
2490        h= nc_CreateSpoly(c->S->m[s->i], c->S->m[s->j]/*, NULL*/, c->r);
2491
2492        if (h!=NULL)
2493          p_Cleardenom(h, c->r);
2494      }
2495      else
2496#endif
2497        h=ksOldCreateSpoly(c->S->m[s->i], c->S->m[s->j], NULL, c->r);
2498    p_Test(h,c->r);
2499    }
2500    else{
2501      h=s->lcm_of_lm;
2502      p_Test(h,c->r);
2503  }
2504    // if(s->i>=0)
2505//       now_t_rep(s->j,s->i,c);
2506    number coef;
2507    int mlen=pLength(h);
2508    p_Test(h,c->r);
2509    if ((!c->nc)&(!(use_noro))){
2510      h=redNF2(h,c,mlen,coef,2);
2511      redTailShort(h,c->strat);
2512      nDelete(&coef);
2513    }
2514    p_Test(h,c->r);
2515    free_sorted_pair_node(s,c->r);
2516    if(!h) continue;
2517    int len=pLength(h);
2518    p[i]=h;
2519
2520    i++;
2521  }
2522  p[i]=NULL;
2523//  pre_comp(p,i,c);
2524  if(i==0){
2525    omfree(p);
2526    return;
2527  }
2528  #ifdef TGB_RESORT_PAIRS
2529  c->replaced=new bool[c->n];
2530  c->used_b=FALSE;
2531  #endif
2532
2533  c->normal_forms+=i;
2534  int j;
2535#ifdef USE_NORO
2536  //if ((!(c->nc))&&(rField_is_Zp(c->r))){
2537  if (use_noro){
2538    int pn=i;
2539    if (pn==0) {omfree(p);return;}
2540
2541    {
2542
2543      if (npPrimeM<255){
2544        noro_step<tgb_uint8>(p,pn,c);
2545      } else {
2546        if (npPrimeM<65000){
2547          noro_step<tgb_uint16>(p,pn,c);
2548        } else{
2549          noro_step<tgb_uint32>(p,pn,c);
2550        }
2551      }
2552
2553
2554
2555    }
2556
2557    //if (TEST_OPT_PROT){
2558    //  Print("reported rank:%i\n",pn);
2559    //}
2560    mass_add(p,pn,c);
2561    omfree(p);
2562    return;
2563    /*if (TEST_OPT_PROT)
2564      for(j=0;j<pn;j++){
2565        p_wrp(p[j],c->r);
2566      }*/
2567  }
2568#endif
2569  red_object* buf=(red_object*) omalloc(i*sizeof(red_object));
2570  for(j=0;j<i;j++){
2571    p_Test(p[j],c->r);
2572    buf[j].p=p[j];
2573    buf[j].sev=pGetShortExpVector(p[j]);
2574    buf[j].bucket = kBucketCreate(currRing);
2575    p_Test(p[j],c->r);
2576    if (c->eliminationProblem){
2577        buf[j].sugar=c->pTotaldegree_full(p[j]);
2578    }
2579    int len=pLength(p[j]);
2580    kBucketInit(buf[j].bucket,buf[j].p,len);
2581    buf[j].initial_quality=buf[j].guess_quality(c);
2582    assume(buf[j].initial_quality>=0);
2583  }
2584  omfree(p);
2585  qsort(buf,i,sizeof(red_object),red_object_better_gen);
2586//    Print("\ncurr_deg:%i\n",curr_deg);
2587  if (TEST_OPT_PROT)
2588  {
2589    Print("%dM[%d,",curr_deg,i);
2590  }
2591
2592  multi_reduction(buf, i, c);
2593  #ifdef TGB_RESORT_PAIRS
2594  if (c->used_b) {
2595    if (TEST_OPT_PROT)
2596        PrintS("B");
2597    int e;
2598    for(e=0;e<=c->pair_top;e++){
2599        if(c->apairs[e]->i<0) continue;
2600        assume(c->apairs[e]->j>=0);
2601        if ((c->replaced[c->apairs[e]->i])||(c->replaced[c->apairs[e]->j])) {
2602            sorted_pair_node* s=c->apairs[e];
2603            s->expected_length=pair_weighted_length(s->i,s->j,c);
2604        }
2605    }
2606    qsort(c->apairs,c->pair_top+1,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
2607  }
2608  #endif
2609#ifdef TGB_DEBUG
2610 {
2611   int k;
2612   for(k=0;k<i;k++)
2613   {
2614     assume(kFindDivisibleByInS_easy(c->strat,buf[k])<0);
2615     int k2;
2616     for(k2=0;k2<i;k2++)
2617     {
2618       if(k==k2) continue;
2619       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));
2620     }
2621   }
2622 }
2623#endif
2624  //resort S
2625
2626  if (TEST_OPT_PROT)
2627      Print("%i]",i);
2628
2629  poly* add_those=(poly*) omalloc(i*sizeof(poly));
2630  for(j=0;j<i;j++)
2631  {
2632    int len;
2633    poly p;
2634    buf[j].flatten();
2635    kBucketClear(buf[j].bucket,&p, &len);
2636    kBucketDestroy(&buf[j].bucket);
2637    p_Test(p,c->r);
2638    //if (!c->nc) {
2639      if ((c->tailReductions) ||(lies_in_last_dp_block(p,c))){
2640      p=redNFTail(p,c->strat->sl,c->strat, 0);
2641      } else {
2642      p=redTailShort(p, c->strat);
2643      }
2644      //}
2645      p_Test(p,c->r);
2646      add_those[j]=p;
2647
2648    //sbuf[j]=add_to_basis(p,-1,-1,c,ibuf+j);
2649  }
2650  mass_add(add_those,i,c);
2651  omfree(add_those);
2652  omfree(buf);
2653
2654
2655
2656
2657
2658
2659
2660  if (TEST_OPT_PROT)
2661      Print("(%d)",c->pair_top+1);
2662  //TODO: implement that while(!(idIs0(c->add_later)))
2663  #ifdef TGB_RESORT_PAIRS
2664  delete c->replaced;
2665  c->replaced=NULL;
2666  c->used_b=FALSE;
2667  #endif
2668  return;
2669}
2670
2671
2672
2673#ifdef REDTAIL_S
2674
2675static poly redNFTail (poly h,const int sl,kStrategy strat, int len)
2676{
2677  BOOLEAN nc=rIsPluralRing(currRing);
2678  if (h==NULL) return NULL;
2679  pTest(h);
2680  if (0 > sl)
2681    return h;
2682  if (pNext(h)==NULL) return h;
2683
2684  int j;
2685  poly res=h;
2686  poly act=res;
2687  LObject P(pNext(h));
2688  pNext(res)=NULL;
2689  P.bucket = kBucketCreate(currRing);
2690  len--;
2691  h=P.p;
2692  if (len <=0) len=pLength(h);
2693  kBucketInit(P.bucket,h /*P.p*/,len /*pLength(P.p)*/);
2694  pTest(h);
2695  loop
2696  {
2697      P.p=h;
2698      P.t_p=NULL;
2699      P.SetShortExpVector();
2700      loop
2701      {
2702          int dummy=strat->sl;
2703          j=kFindDivisibleByInS_easy(strat,P.p,P.sev);//kFindDivisibleByInS(strat,&dummy,&P);
2704          if (j>=0)
2705          {
2706#ifdef REDTAIL_PROT
2707            PrintS("r");
2708#endif
2709            nNormalize(pGetCoeff(P.p));
2710#ifdef KDEBUG
2711            if (TEST_OPT_DEBUG)
2712            {
2713              PrintS("red tail:");
2714              wrp(h);
2715              PrintS(" with ");
2716              wrp(strat->S[j]);
2717            }
2718#endif
2719            number coef;
2720            pTest(strat->S[j]);
2721#ifdef HAVE_PLURAL
2722            if (nc){
2723              nc_BucketPolyRed_Z(P.bucket, strat->S[j], &coef);
2724            } else
2725#endif
2726              coef=kBucketPolyRed(P.bucket,strat->S[j],
2727                                strat->lenS[j]/*pLength(strat->S[j])*/,strat->kNoether);
2728            pMult_nn(res,coef);
2729            nDelete(&coef);
2730            h = kBucketGetLm(P.bucket);
2731            pTest(h);
2732            if (h==NULL)
2733            {
2734#ifdef REDTAIL_PROT
2735              PrintS(" ");
2736#endif
2737        kBucketDestroy(&P.bucket);
2738              return res;
2739            }
2740            pTest(h);
2741            P.p=h;
2742            P.t_p=NULL;
2743            P.SetShortExpVector();
2744#ifdef KDEBUG
2745            if (TEST_OPT_DEBUG)
2746            {
2747              PrintS("\nto tail:");
2748              wrp(h);
2749              PrintLn();
2750            }
2751#endif
2752          }
2753          else
2754          {
2755#ifdef REDTAIL_PROT
2756            PrintS("n");
2757#endif
2758            break;
2759          }
2760      } /* end loop current mon */
2761      //   poly tmp=pHead(h /*kBucketGetLm(P.bucket)*/);
2762      //act->next=tmp;pIter(act);
2763      act->next=kBucketExtractLm(P.bucket);pIter(act);
2764      h = kBucketGetLm(P.bucket);
2765      if (h==NULL)
2766      {
2767#ifdef REDTAIL_PROT
2768        PrintS(" ");
2769#endif
2770  kBucketDestroy(&P.bucket);
2771        return res;
2772      }
2773      pTest(h);
2774  }
2775}
2776#endif
2777
2778
2779//try to fill, return FALSE iff queue is empty
2780
2781//transfers ownership of m to mat
2782void init_with_mac_poly(tgb_sparse_matrix* mat, int row, mac_poly m){
2783  assume(mat->mp[row]==NULL);
2784  mat->mp[row]=m;
2785#ifdef TGB_DEBUG
2786  mac_poly r=m;
2787  while(r){
2788    assume(r->exp<mat->columns);
2789    r=r->next;
2790  }
2791#endif
2792}
2793poly free_row_to_poly(tgb_sparse_matrix* mat, int row, poly* monoms, int monom_index){
2794  poly p=NULL;
2795  poly* set_this=&p;
2796  mac_poly r=mat->mp[row];
2797  mat->mp[row]=NULL;
2798  while(r)
2799  {
2800    (*set_this)=pLmInit(monoms[monom_index-1-r->exp]);
2801    pSetCoeff((*set_this),r->coef);
2802    set_this=&((*set_this)->next);
2803    mac_poly old=r;
2804    r=r->next;
2805    delete old;
2806
2807  }
2808  return p;
2809
2810}
2811
2812
2813
2814
2815
2816
2817
2818
2819
2820
2821
2822static int poly_crit(const void* ap1, const void* ap2){
2823  poly p1,p2;
2824  p1=*((poly*) ap1);
2825  p2=*((poly*)ap2);
2826
2827  int c=pLmCmp(p1,p2);
2828  if (c !=0) return c;
2829  int l1=pLength(p1);
2830  int l2=pLength(p2);
2831  if (l1<l2) return -1;
2832  if (l1>l2) return 1;
2833  return 0;
2834}
2835void slimgb_alg::introduceDelayedPairs(poly* pa,int s){
2836    if (s==0) return;
2837    sorted_pair_node** si_array=(sorted_pair_node**) omalloc(s* sizeof(sorted_pair_node*));
2838
2839    for( unsigned int i = 0; i < s; i++ )
2840    {
2841        sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
2842        si->i=-1;
2843        si->j=-2;
2844        poly p=pa[i];
2845        simplify_poly(p,r);
2846        si->expected_length=pQuality(p,this,pLength(p));
2847        p_Test(p,r);
2848        si->deg=this->pTotaldegree_full(p);
2849        /*if (!rField_is_Zp(r))
2850       {
2851          p_Content(p,r);
2852          p_Cleardenom(p,r);
2853        }*/
2854
2855        si->lcm_of_lm=p;
2856
2857        //      c->apairs[n-1-i]=si;
2858        si_array[i]=si;
2859
2860  }
2861
2862  qsort(si_array,s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
2863    apairs=spn_merge(apairs,pair_top+1,si_array,s,this);
2864  pair_top+=s;
2865  omfree(si_array);
2866}
2867slimgb_alg::slimgb_alg(ideal I, int syz_comp,BOOLEAN F4,int deg_pos)
2868{
2869  this->deg_pos=deg_pos;
2870  lastCleanedDeg=-1;
2871  completed=FALSE;
2872  this->syz_comp=syz_comp;
2873  r=currRing;
2874  nc=rIsPluralRing(r);
2875  this->lastDpBlockStart=get_last_dp_block_start(r);
2876  //Print("last dp Block start: %i\n", this->lastDpBlockStart);
2877  is_homog=TRUE;
2878  {
2879    int hzz;
2880    for(hzz=0;hzz<IDELEMS(I);hzz++)
2881    {
2882      assume(I->m[hzz]!=NULL);
2883      int d=this->pTotaldegree(I->m[hzz]);
2884      poly t=I->m[hzz]->next;
2885      while(t)
2886      {
2887        if (d!=this->pTotaldegree(t))
2888        {
2889          is_homog=FALSE;
2890          break;
2891        }
2892        t=t->next;
2893      }
2894      if(!(is_homog)) break;
2895    }
2896  }
2897  eliminationProblem=((!(is_homog))&&((pLexOrder)||(I->rank>1)));
2898  tailReductions=((is_homog)||((TEST_OPT_REDTAIL)&&(!(I->rank>1))));
2899  //  Print("is homog:%d",c->is_homog);
2900  void* h;
2901  poly hp;
2902  int i,j;
2903  to_destroy=NULL;
2904  easy_product_crit=0;
2905  extended_product_crit=0;
2906  if (rField_is_Zp(r))
2907    isDifficultField=FALSE;
2908  else
2909    isDifficultField=TRUE;
2910  //not fully correct
2911  //(rChar()==0);
2912  F4_mode=F4;
2913
2914  reduction_steps=0;
2915  last_index=-1;
2916
2917  F=NULL;
2918  F_minus=NULL;
2919
2920  Rcounter=0;
2921
2922  soon_free=NULL;
2923
2924  tmp_lm=pOne();
2925
2926  normal_forms=0;
2927  current_degree=1;
2928
2929  max_pairs=5*IDELEMS(I);
2930
2931  apairs=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*max_pairs);
2932  pair_top=-1;
2933
2934  int nn=IDELEMS(I);
2935  array_lengths=nn;
2936
2937  i=0;
2938  this->n=0;
2939  T_deg=(int*) omalloc(nn*sizeof(int));
2940  if(eliminationProblem)
2941    T_deg_full=(int*) omalloc(nn*sizeof(int));
2942  else
2943    T_deg_full=NULL;
2944  tmp_pair_lm=(poly*) omalloc(nn*sizeof(poly));
2945  tmp_spn=(sorted_pair_node**) omalloc(nn*sizeof(sorted_pair_node*));
2946  lm_bin=omGetSpecBin(POLYSIZE + (r->ExpL_Size)*sizeof(long));
2947#ifdef HEAD_BIN
2948  HeadBin=omGetSpecBin(POLYSIZE + (currRing->ExpL_Size)*sizeof(long));
2949#endif
2950  /* omUnGetSpecBin(&(c->HeadBin)); */
2951  #ifndef HAVE_BOOST
2952  #ifdef USE_STDVECBOOL
2953  #else
2954  h=omalloc(nn*sizeof(char*));
2955
2956  states=(char**) h;
2957  #endif
2958  #endif
2959  h=omalloc(n*sizeof(int));
2960  lengths=(int*) h;
2961  weighted_lengths=(wlen_type*)omalloc(nn*sizeof(wlen_type));
2962  gcd_of_terms=(poly*) omalloc(nn*sizeof(poly));
2963
2964  short_Exps=(long*) omalloc(nn*sizeof(long));
2965  if (F4_mode)
2966    S=idInit(nn,I->rank);
2967  else
2968    S=idInit(1,I->rank);
2969  strat=new skStrategy;
2970  if (eliminationProblem)
2971    strat->honey=TRUE;
2972  strat->syzComp = 0;
2973  initBuchMoraCrit(strat);
2974  initBuchMoraPos(strat);
2975  strat->initEcart = initEcartBBA;
2976  strat->tailRing=r;
2977  strat->enterS = enterSBba;
2978  strat->sl = -1;
2979  i=nn;
2980  i=1;//some strange bug else
2981  /* initS(c->S,NULL,c->strat); */
2982  /* intS start: */
2983  // i=((i+IDELEMS(c->S)+15)/16)*16;
2984  strat->ecartS=(intset)omAlloc(i*sizeof(int)); /*initec(i);*/
2985  strat->sevS=(unsigned long*)omAlloc0(i*sizeof(unsigned long));
2986  /*initsevS(i);*/
2987  strat->S_2_R=(int*)omAlloc0(i*sizeof(int));/*initS_2_R(i);*/
2988  strat->fromQ=NULL;
2989  strat->Shdl=idInit(1,1);
2990  strat->S=strat->Shdl->m;
2991  strat->lenS=(int*)omAlloc0(i*sizeof(int));
2992  if((isDifficultField)||(eliminationProblem))
2993    strat->lenSw=(wlen_type*)omAlloc0(i*sizeof(wlen_type));
2994  else
2995    strat->lenSw=NULL;
2996  sorted_pair_node* si;
2997  assume(nn>0);
2998  add_to_basis_ideal_quotient(I->m[0],this,NULL);
2999
3000  assume(strat->sl==IDELEMS(strat->Shdl)-1);
3001  if(!(F4_mode))
3002  {
3003    poly* array_arg=I->m;
3004    array_arg++;
3005    introduceDelayedPairs(array_arg,n-1);
3006    /*
3007    for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
3008    {
3009      //     add_to_basis(I->m[i],-1,-1,c);
3010      si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
3011      si->i=-1;
3012      si->j=-2;
3013      si->expected_length=pQuality(I->m[i],this,pLength(I->m[i]));
3014      si->deg=pTotaldegree(I->m[i]);
3015      if (!rField_is_Zp(r))
3016      {
3017        p_Cleardenom(I->m[i], r);
3018      }
3019      si->lcm_of_lm=I->m[i];
3020
3021      //      c->apairs[nn-1-i]=si;
3022      apairs[nn-i-1]=si;
3023      ++(pair_top);
3024    }
3025    */
3026  }
3027  else
3028  {
3029    for (i=1;i<nn;i++)//the 1 is wanted, because first element is added to basis
3030      add_to_basis_ideal_quotient(I->m[i],this,NULL);
3031  }
3032  for(i=0;i<IDELEMS(I);i++)
3033  {
3034    I->m[i]=NULL;
3035  }
3036  idDelete(&I);
3037  add_later=idInit(ADD_LATER_SIZE,S->rank);
3038  #ifdef USE_NORO
3039  use_noro=((!(nc))&&(S->rank<=1)&&(rField_is_Zp(r)) &&(!(eliminationProblem))&&(npPrimeM<=32003));
3040  use_noro_last_block=false;
3041  if ((!(use_noro))&&(lastDpBlockStart<=pVariables))
3042  {
3043      use_noro_last_block=((!(nc))&&(S->rank<=1)&&(rField_is_Zp(r)) &&(npPrimeM<=32003));
3044  }
3045  #else
3046  use_noro=false;
3047  use_noro_last_block=false;
3048  #endif
3049  //Print("NORO last block %i",use_noro_last_block);
3050  memset(add_later->m,0,ADD_LATER_SIZE*sizeof(poly));
3051}
3052slimgb_alg::~slimgb_alg()
3053{
3054  if (!(completed))
3055  {
3056    poly* add=(poly*) omalloc((pair_top+2)*sizeof(poly));
3057    int piter;
3058    int pos=0;
3059    for(piter=0;piter<=pair_top;piter++)
3060    {
3061      sorted_pair_node* s=apairs[piter];
3062      if (s->i<0)
3063      {
3064        //delayed element
3065        if (s->lcm_of_lm!=NULL)
3066        {
3067          add[pos]=s->lcm_of_lm;
3068          pos++;
3069        }
3070      }
3071      free_sorted_pair_node(s,r);
3072      apairs[piter]=NULL;
3073    }
3074    pair_top=-1;
3075    add[pos]=NULL;
3076    pos=0;
3077    while(add[pos]!=NULL)
3078    {
3079      add_to_basis_ideal_quotient(add[pos],this,NULL);
3080      pos++;
3081    }
3082    for(piter=0;piter<=pair_top;piter++)
3083    {
3084      sorted_pair_node* s=apairs[piter];
3085      assume(s->i>=0);
3086      free_sorted_pair_node(s,r);
3087      apairs[piter]=NULL;
3088    }
3089    pair_top=-1;
3090  }
3091  id_Delete(&add_later,r);
3092  int i,j;
3093  slimgb_alg* c=this;
3094  while(c->to_destroy)
3095  {
3096    pDelete(&(c->to_destroy->p));
3097    poly_list_node* old=c->to_destroy;
3098    c->to_destroy=c->to_destroy->next;
3099    omfree(old);
3100  }
3101  while(c->F)
3102  {
3103    for(i=0;i<c->F->size;i++)
3104    {
3105      pDelete(&(c->F->mp[i].m));
3106    }
3107    omfree(c->F->mp);
3108    c->F->mp=NULL;
3109    mp_array_list* old=c->F;
3110    c->F=c->F->next;
3111    omfree(old);
3112  }
3113  while(c->F_minus)
3114  {
3115    for(i=0;i<c->F_minus->size;i++)
3116    {
3117      pDelete(&(c->F_minus->p[i]));
3118    }
3119    omfree(c->F_minus->p);
3120    c->F_minus->p=NULL;
3121    poly_array_list* old=c->F_minus;
3122    c->F_minus=c->F_minus->next;
3123    omfree(old);
3124  }
3125  #ifndef HAVE_BOOST
3126  #ifndef USE_STDVECBOOL
3127  for(int z=1 /* zero length at 0 */;z<c->n;z++)
3128  {
3129    omfree(c->states[z]);
3130  }
3131  omfree(c->states);
3132  #endif
3133  #endif
3134
3135  omfree(c->lengths);
3136  omfree(c->weighted_lengths);
3137  for(int z=0;z<c->n;z++)
3138  {
3139    pDelete(&c->tmp_pair_lm[z]);
3140    omfree(c->tmp_spn[z]);
3141  }
3142  omfree(c->tmp_pair_lm);
3143  omfree(c->tmp_spn);
3144
3145  omfree(c->T_deg);
3146  if(c->T_deg_full)
3147    omfree(c->T_deg_full);
3148
3149  omFree(c->strat->ecartS);
3150  omFree(c->strat->sevS);
3151//   initsevS(i);
3152  omFree(c->strat->S_2_R);
3153  omFree(c->strat->lenS);
3154
3155  if(c->strat->lenSw)  omFree(c->strat->lenSw);
3156
3157  for(i=0;i<c->n;i++)
3158  {
3159    if(c->gcd_of_terms[i])
3160      pDelete(&(c->gcd_of_terms[i]));
3161  }
3162  omfree(c->gcd_of_terms);
3163
3164  omfree(c->apairs);
3165  if (TEST_OPT_PROT)
3166  {
3167    //Print("calculated %d NFs\n",c->normal_forms);
3168      Print("\nNF:%i product criterion:%i, ext_product criterion:%i \n", c->normal_forms, c->easy_product_crit, c->extended_product_crit);
3169  }
3170  int deleted_form_c_s=0;
3171
3172  for(i=0;i<=c->strat->sl;i++)
3173  {
3174    if (!c->strat->S[i]) continue;
3175    BOOLEAN found=FALSE;
3176    for(j=0;j<c->n;j++)
3177    {
3178      if (c->S->m[j]==c->strat->S[i])
3179      {
3180        found=TRUE;
3181        break;
3182      }
3183    }
3184    if(!found) pDelete(&c->strat->S[i]);
3185  }
3186//   for(i=0;i<c->n;i++){
3187//     if (c->rep[i]!=i){
3188// //       for(j=0;j<=c->strat->sl;j++){
3189// //   if(c->strat->S[j]==c->S->m[i]){
3190// //     c->strat->S[j]=NULL;
3191// //     break;
3192// //   }
3193// //       }
3194// //      PrintS("R_delete");
3195//       pDelete(&c->S->m[i]);
3196//     }
3197//   }
3198
3199  if (completed)
3200  {
3201    for(i=0;i<c->n;i++)
3202    {
3203      assume(c->S->m[i]!=NULL);
3204      if (p_GetComp(c->S->m[i],currRing)>this->syz_comp) continue;
3205      for(j=0;j<c->n;j++)
3206      {
3207        if((c->S->m[j]==NULL)||(i==j))
3208          continue;
3209        assume(p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
3210             c->S->m[i],~c->short_Exps[i],
3211             c->r)==p_LmDivisibleBy(c->S->m[j],
3212             c->S->m[i],
3213             c->r));
3214        if (p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
3215          c->S->m[i],~c->short_Exps[i],
3216          c->r))
3217        {
3218          pDelete(&c->S->m[i]);
3219          break;
3220        }
3221      }
3222    }
3223  }
3224  omfree(c->short_Exps);
3225
3226  ideal I=c->S;
3227
3228  IDELEMS(I)=c->n;
3229
3230  idSkipZeroes(I);
3231  for(i=0;i<=c->strat->sl;i++)
3232    c->strat->S[i]=NULL;
3233  id_Delete(&c->strat->Shdl,c->r);
3234  pDelete(&c->tmp_lm);
3235  omUnGetSpecBin(&lm_bin);
3236  delete c->strat;
3237}
3238ideal t_rep_gb(ring r,ideal arg_I, int syz_comp, BOOLEAN F4_mode)
3239{
3240  assume(r==currRing);
3241  ring orig_ring=r;
3242  int pos;
3243  ring new_ring=rAssure_TDeg(orig_ring,1,rVar(orig_ring),pos);
3244
3245  ideal s_h;
3246  if (orig_ring != new_ring)
3247  {
3248    rChangeCurrRing(new_ring);
3249    s_h=idrCopyR_NoSort(arg_I,orig_ring);
3250    idTest(s_h);
3251    /*int i;
3252
3253    for(i=0;i<IDELEMS(s_h);i++)
3254    {
3255      poly p=s_h->m[i];
3256      while(p)
3257      {
3258                p_Setm(p,new_ring);
3259                pIter(p);
3260      }
3261    }
3262    */
3263  }
3264  else
3265  {
3266    s_h = id_Copy(arg_I,orig_ring);
3267  }
3268
3269  ideal s_result=do_t_rep_gb(new_ring,s_h,syz_comp,F4_mode,pos);
3270  ideal result;
3271  if(orig_ring != new_ring)
3272  {
3273    idTest(s_result);
3274    rChangeCurrRing(orig_ring);
3275    result = idrMoveR_NoSort(s_result, new_ring);
3276
3277    idTest(result);
3278    //rChangeCurrRing(new_ring);
3279    rKill(new_ring);
3280    //rChangeCurrRing(orig_ring);
3281  }
3282  else
3283    result=s_result;
3284  idTest(result);
3285  return result;
3286}
3287
3288ideal do_t_rep_gb(ring r,ideal arg_I, int syz_comp, BOOLEAN F4_mode,int deg_pos){
3289  //  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)));
3290
3291  if (TEST_OPT_PROT)
3292    if (F4_mode)
3293      PrintS("F4 Modus \n");
3294
3295  //debug_Ideal=arg_debug_Ideal;
3296  //if (debug_Ideal) PrintS("DebugIdeal received\n");
3297  // Print("Idelems %i \n----------\n",IDELEMS(arg_I));
3298  ideal I=arg_I;
3299  idCompactify(I);
3300  if (idIs0(I)) return I;
3301  int i;
3302  for(i=0;i<IDELEMS(I);i++)
3303  {
3304    assume(I->m[i]!=NULL);
3305    simplify_poly(I->m[i],currRing);
3306  }
3307
3308  qsort(I->m,IDELEMS(I),sizeof(poly),poly_crit);
3309  //Print("Idelems %i \n----------\n",IDELEMS(I));
3310  //slimgb_alg* c=(slimgb_alg*) omalloc(sizeof(slimgb_alg));
3311  //int syz_comp=arg_I->rank;
3312  slimgb_alg* c=new slimgb_alg(I, syz_comp,F4_mode,deg_pos);
3313
3314  while ((c->pair_top>=0) && ((!(TEST_OPT_DEGBOUND)) || (c->apairs[c->pair_top]->deg<=Kstd1_deg)))
3315  {
3316    #ifdef HAVE_F4
3317    if(F4_mode)
3318      go_on_F4(c);
3319    else
3320    #endif
3321      go_on(c);
3322  }
3323  if (c->pair_top<0)
3324    c->completed=TRUE;
3325  I=c->S;
3326  delete c;
3327  if (TEST_OPT_REDSB)
3328  {
3329    ideal erg=kInterRed(I,NULL);
3330    assume(I!=erg);
3331    id_Delete(&I, currRing);
3332    return erg;
3333  }
3334  //qsort(I->m, IDELEMS(I),sizeof(poly),pLmCmp_func);
3335  assume(I->rank>=idRankFreeModule(I));
3336  return(I);
3337}
3338void now_t_rep(const int & arg_i, const int & arg_j, slimgb_alg* c)
3339{
3340  int i,j;
3341  if (arg_i==arg_j)
3342  {
3343    return;
3344  }
3345  if (arg_i>arg_j)
3346  {
3347    i=arg_j;
3348    j=arg_i;
3349  }
3350  else
3351  {
3352    i=arg_i;
3353    j=arg_j;
3354  }
3355  c->states[j][i]=HASTREP;
3356}
3357
3358static BOOLEAN has_t_rep(const int & arg_i, const  int & arg_j, slimgb_alg* state)
3359{
3360  assume(0<=arg_i);
3361  assume(0<=arg_j);
3362  assume(arg_i<state->n);
3363  assume(arg_j<state->n);
3364  if (arg_i==arg_j)
3365  {
3366    return (TRUE);
3367  }
3368  if (arg_i>arg_j)
3369  {
3370    return (state->states[arg_i][arg_j]==HASTREP);
3371  }
3372  else
3373  {
3374    return (state->states[arg_j][arg_i]==HASTREP);
3375  }
3376}
3377static int pLcmDeg(poly a, poly b)
3378{
3379  int i;
3380  int n=0;
3381  for (i=pVariables; i; i--)
3382  {
3383    n+=si_max( pGetExp(a,i), pGetExp(b,i));
3384  }
3385  return n;
3386}
3387
3388static void shorten_tails(slimgb_alg* c, poly monom)
3389{
3390  return;
3391// BOOLEAN corr=lenS_correct(c->strat);
3392  for(int i=0;i<c->n;i++)
3393  {
3394    //enter tail
3395    if (c->S->m[i]==NULL) continue;
3396    poly tail=c->S->m[i]->next;
3397    poly prev=c->S->m[i];
3398    BOOLEAN did_something=FALSE;
3399    while((tail!=NULL)&& (pLmCmp(tail, monom)>=0))
3400    {
3401      if (p_LmDivisibleBy(monom,tail,c->r))
3402      {
3403        did_something=TRUE;
3404        prev->next=tail->next;
3405        tail->next=NULL;
3406        p_Delete(& tail,c->r);
3407        tail=prev;
3408        //PrintS("Shortened");
3409        c->lengths[i]--;
3410      }
3411      prev=tail;
3412      tail=tail->next;
3413    }
3414    if (did_something)
3415    {
3416      int new_pos;
3417      wlen_type q;
3418      q=pQuality(c->S->m[i],c,c->lengths[i]);
3419      new_pos=simple_posInS(c->strat,c->S->m[i],c->lengths[i],q);
3420
3421      int old_pos=-1;
3422      //assume new_pos<old_pos
3423      for (int z=0;z<=c->strat->sl;z++)
3424      {
3425        if (c->strat->S[z]==c->S->m[i])
3426        {
3427          old_pos=z;
3428          break;
3429        }
3430      }
3431      if (old_pos== -1)
3432        for (int z=new_pos-1;z>=0;z--)
3433        {
3434          if (c->strat->S[z]==c->S->m[i])
3435          {
3436            old_pos=z;
3437            break;
3438          }
3439        }
3440      assume(old_pos>=0);
3441      assume(new_pos<=old_pos);
3442      assume(pLength(c->strat->S[old_pos])==c->lengths[i]);
3443      c->strat->lenS[old_pos]=c->lengths[i];
3444      if (c->strat->lenSw)
3445        c->strat->lenSw[old_pos]=q;
3446
3447      if (new_pos<old_pos)
3448        move_forward_in_S(old_pos,new_pos,c->strat);
3449
3450      length_one_crit(c,i,c->lengths[i]);
3451    }
3452  }
3453}
3454static sorted_pair_node* pop_pair(slimgb_alg* c)
3455{
3456  clean_top_of_pair_list(c);
3457
3458  if(c->pair_top<0) return NULL;
3459  else return (c->apairs[c->pair_top--]);
3460}
3461void slimgb_alg::cleanDegs(int lower, int upper)
3462{
3463  assume(is_homog);
3464  int deg;
3465  if (TEST_OPT_PROT)
3466  {
3467    PrintS("C");
3468  }
3469  for(deg=lower;deg<=upper;deg++)
3470  {
3471    int i;
3472    for(i=0;i<n;i++)
3473    {
3474      if (T_deg[i]==deg)
3475      {
3476        poly h;
3477        h=S->m[i];
3478        h=redNFTail(h,strat->sl,strat,lengths[i]);
3479        if (!rField_is_Zp(r))
3480        {
3481          p_Cleardenom(h,r);
3482          //p_Content(h,r);
3483        }
3484        else pNorm(h);
3485          //TODO:GCD of TERMS
3486        poly got=::gcd_of_terms(h,r);
3487        p_Delete(&gcd_of_terms[i],r);
3488        gcd_of_terms[i]=got;
3489        int len=pLength(h);
3490        wlen_type wlen=pQuality(h,this,len);
3491        if (weighted_lengths)
3492          weighted_lengths[i]=wlen;
3493        lengths[i]=len;
3494        assume(h==S->m[i]);
3495        int j;
3496        for(j=0;j<=strat->sl;j++)
3497        {
3498          if (h==strat->S[j])
3499          {
3500            int new_pos=simple_posInS(strat, h,len, wlen);
3501            if (strat->lenS)
3502            {
3503              strat->lenS[j]=len;
3504            }
3505            if (strat->lenSw)
3506            {
3507              strat->lenSw[j]=wlen;
3508            }
3509            if (new_pos<j)
3510            {
3511              move_forward_in_S(j,new_pos,strat);
3512            }
3513            else
3514            {
3515              if (new_pos>j)
3516                new_pos=new_pos-1;//is identical with one element
3517              if (new_pos>j)
3518                move_backward_in_S(j,new_pos,strat);
3519            }
3520            break;
3521          }
3522        }
3523      }
3524    }
3525  }
3526  {
3527    int i,j;
3528    for(i=0;i<this->n;i++)
3529    {
3530      for(j=0;j<i;j++)
3531      {
3532        if (T_deg[i]+T_deg[j]<=upper)
3533        {
3534          now_t_rep(i,j,this);
3535        }
3536      }
3537    }
3538  }
3539  //TODO resort and update strat->S,strat->lenSw
3540  //TODO mark pairs
3541}
3542sorted_pair_node* top_pair(slimgb_alg* c)
3543{
3544  while(c->pair_top>=0)
3545  {
3546    super_clean_top_of_pair_list(c);//yeah, I know, it's odd that I use a different proc here
3547    if ((c->is_homog)&&(c->pair_top>=0)&&(c->apairs[c->pair_top]->deg>=c->lastCleanedDeg+2))
3548    {
3549      int upper=c->apairs[c->pair_top]->deg-1;
3550      c->cleanDegs(c->lastCleanedDeg+1,upper);
3551      c->lastCleanedDeg=upper;
3552    }
3553    else
3554    {
3555      break;
3556    }
3557  }
3558  if(c->pair_top<0) return NULL;
3559  else return (c->apairs[c->pair_top]);
3560}
3561sorted_pair_node* quick_pop_pair(slimgb_alg* c)
3562{
3563  if(c->pair_top<0) return NULL;
3564  else return (c->apairs[c->pair_top--]);
3565}
3566
3567static void super_clean_top_of_pair_list(slimgb_alg* c)
3568{
3569  while((c->pair_top>=0)
3570  && (c->apairs[c->pair_top]->i>=0)
3571  && (good_has_t_rep(c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,c)))
3572  {
3573    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
3574    c->pair_top--;
3575  }
3576}
3577void clean_top_of_pair_list(slimgb_alg* c)
3578{
3579  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)))
3580  {
3581    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
3582    c->pair_top--;
3583  }
3584}
3585static BOOLEAN state_is(calc_state state, const int & arg_i, const  int & arg_j, slimgb_alg* c)
3586{
3587  assume(0<=arg_i);
3588  assume(0<=arg_j);
3589  assume(arg_i<c->n);
3590  assume(arg_j<c->n);
3591  if (arg_i==arg_j)
3592  {
3593    return (TRUE);
3594  }
3595  if (arg_i>arg_j)
3596  {
3597    return (c->states[arg_i][arg_j]==state);
3598  }
3599  else return(c->states[arg_j][arg_i]==state);
3600}
3601
3602void free_sorted_pair_node(sorted_pair_node* s, ring r)
3603{
3604  if (s->i>=0)
3605    p_Delete(&s->lcm_of_lm,r);
3606  omfree(s);
3607}
3608static BOOLEAN pair_better(sorted_pair_node* a,sorted_pair_node* b, slimgb_alg* c)
3609{
3610  if (a->deg<b->deg) return TRUE;
3611  if (a->deg>b->deg) return FALSE;
3612
3613  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
3614  if (comp==1) return FALSE;
3615  if (-1==comp) return TRUE;
3616  if (a->expected_length<b->expected_length) return TRUE;
3617  if (a->expected_length>b->expected_length) return FALSE;
3618  if (a->i+a->j<b->i+b->j) return TRUE;
3619   if (a->i+a->j>b->i+b->j) return FALSE;
3620  if (a->i<b->i) return TRUE;
3621  if (a->i>b->i) return FALSE;
3622  return TRUE;
3623}
3624
3625static int tgb_pair_better_gen(const void* ap,const void* bp)
3626{
3627  sorted_pair_node* a=*((sorted_pair_node**)ap);
3628  sorted_pair_node* b=*((sorted_pair_node**)bp);
3629  assume((a->i>a->j) || (a->i < 0));
3630  assume((b->i>b->j) || (b->i < 0));
3631  if (a->deg<b->deg) return -1;
3632  if (a->deg>b->deg) return 1;
3633
3634  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
3635
3636  if (comp==1) return 1;
3637  if (-1==comp) return -1;
3638  if (a->expected_length<b->expected_length) return -1;
3639  if (a->expected_length>b->expected_length) return 1;
3640  if (a->i+a->j<b->i+b->j) return -1;
3641  if (a->i+a->j>b->i+b->j) return 1;
3642  if (a->i<b->i) return -1;
3643  if (a->i>b->i) return 1;
3644  return 0;
3645}
3646
3647static poly gcd_of_terms(poly p, ring r)
3648{
3649  int max_g_0=0;
3650  assume(p!=NULL);
3651  int i;
3652  poly m=pOne();
3653  poly t;
3654  for (i=pVariables; i; i--)
3655  {
3656    pSetExp(m,i, pGetExp(p,i));
3657    if (max_g_0==0)
3658      if (pGetExp(m,i)>0)
3659        max_g_0=i;
3660  }
3661  t=p->next;
3662  while (t!=NULL)
3663  {
3664    if (max_g_0==0) break;
3665    for (i=max_g_0; i; i--)
3666    {
3667      pSetExp(m,i, si_min(pGetExp(t,i),pGetExp(m,i)));
3668      if (max_g_0==i)
3669        if (pGetExp(m,i)==0)
3670          max_g_0=0;
3671      if ((max_g_0==0) && (pGetExp(m,i)>0))
3672      {
3673        max_g_0=i;
3674      }
3675    }
3676    t=t->next;
3677  }
3678  p_Setm(m,r);
3679  if (max_g_0>0)
3680    return m;
3681  pDelete(&m);
3682  return NULL;
3683}
3684static inline BOOLEAN pHasNotCFExtended(poly p1, poly p2, poly m)
3685{
3686  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
3687    return FALSE;
3688  int i = 1;
3689  loop
3690  {
3691    if ((pGetExp(p1, i)-pGetExp(m,i) >0) && (pGetExp(p2, i) -pGetExp(m,i)> 0))   return FALSE;
3692    if (i == pVariables)                                return TRUE;
3693    i++;
3694  }
3695}
3696
3697
3698//for impl reasons may return false if the the normal product criterion matches
3699static inline BOOLEAN extended_product_criterion(poly p1, poly gcd1, poly p2, poly gcd2, slimgb_alg* c)
3700{
3701  if (c->nc)
3702    return FALSE;
3703  if(gcd1==NULL) return FALSE;
3704  if(gcd2==NULL) return FALSE;
3705  gcd1->next=gcd2; //may ordered incorrect
3706  poly m=gcd_of_terms(gcd1,c->r);
3707  gcd1->next=NULL;
3708  if (m==NULL) return FALSE;
3709
3710  BOOLEAN erg=pHasNotCFExtended(p1,p2,m);
3711  pDelete(&m);
3712  return erg;
3713}
3714static poly kBucketGcd(kBucket* b, ring r)
3715{
3716  int s=0;
3717  int i;
3718  poly m, n;
3719  BOOLEAN initialized=FALSE;
3720  for (i=MAX_BUCKET-1;i>=0;i--)
3721  {
3722    if (b->buckets[i]!=NULL)
3723    {
3724      if (!initialized)
3725      {
3726        m=gcd_of_terms(b->buckets[i],r);
3727        initialized=TRUE;
3728        if (m==NULL) return NULL;
3729      }
3730      else
3731      {
3732        n=gcd_of_terms(b->buckets[i],r);
3733        if (n==NULL)
3734        {
3735          pDelete(&m);
3736          return NULL;
3737        }
3738        n->next=m;
3739        poly t=gcd_of_terms(n,r);
3740        n->next=NULL;
3741        pDelete(&m);
3742        pDelete(&n);
3743        m=t;
3744        if (m==NULL) return NULL;
3745      }
3746    }
3747  }
3748  return m;
3749}
3750
3751static inline wlen_type quality_of_pos_in_strat_S(int pos, slimgb_alg* c)
3752{
3753  if (c->strat->lenSw!=NULL) return c->strat->lenSw[pos];
3754  return c->strat->lenS[pos];
3755}
3756#ifdef HAVE_PLURAL
3757static inline wlen_type quality_of_pos_in_strat_S_mult_high(int pos, poly high, slimgb_alg* c)
3758  //meant only for nc
3759{
3760  poly m=pOne();
3761  pExpVectorDiff(m,high ,c->strat->S[pos]);
3762  poly product = nc_mm_Mult_pp(m, c->strat->S[pos], c->r);
3763  wlen_type erg=pQuality(product,c);
3764  pDelete(&m);
3765  pDelete(&product);
3766  return erg;
3767}
3768#endif
3769
3770static void multi_reduction_lls_trick(red_object* los, int losl,slimgb_alg* c,find_erg & erg)
3771{
3772  erg.expand=NULL;
3773  BOOLEAN swap_roles; //from reduce_by, to_reduce_u if fromS
3774  if(erg.fromS)
3775  {
3776    if(pLmEqual(c->strat->S[erg.reduce_by],los[erg.to_reduce_u].p))
3777    {
3778      int i;
3779      wlen_type quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
3780      int best=erg.to_reduce_u+1;
3781/*
3782      for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
3783  int qc=los[i].guess_quality(c);
3784  if (qc<quality_a){
3785    best=i;
3786    quality_a=qc;
3787  }
3788      }
3789      if(best!=erg.to_reduce_u+1){
3790*/
3791      wlen_type qc;
3792      best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3793      if(qc<quality_a)
3794      {
3795        los[best].flatten();
3796        int b_pos=kBucketCanonicalize(los[best].bucket);
3797        los[best].p=los[best].bucket->buckets[b_pos];
3798        qc=pQuality(los[best].bucket->buckets[b_pos],c);
3799        if(qc<quality_a)
3800        {
3801          red_object h=los[erg.to_reduce_u];
3802          los[erg.to_reduce_u]=los[best];
3803          los[best]=h;
3804          swap_roles=TRUE;
3805        }
3806        else
3807          swap_roles=FALSE;
3808      }
3809      else
3810      {
3811        swap_roles=FALSE;
3812      }
3813    }
3814    else
3815    {
3816      if (erg.to_reduce_u>erg.to_reduce_l)
3817      {
3818        int i;
3819        wlen_type quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
3820#ifdef HAVE_PLURAL
3821        if ((c->nc) && (!(rIsSCA(c->r))))
3822         quality_a=quality_of_pos_in_strat_S_mult_high(erg.reduce_by, los[erg.to_reduce_u].p, c);
3823#endif
3824        int best=erg.to_reduce_u+1;
3825        wlen_type qc;
3826        best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3827        assume(qc==los[best].guess_quality(c));
3828        if(qc<quality_a)
3829        {
3830          los[best].flatten();
3831          int b_pos=kBucketCanonicalize(los[best].bucket);
3832          los[best].p=los[best].bucket->buckets[b_pos];
3833          qc==pQuality(los[best].bucket->buckets[b_pos],c);
3834          //(best!=erg.to_reduce_u+1)
3835          if(qc<quality_a)
3836          {
3837            red_object h=los[erg.to_reduce_u];
3838            los[erg.to_reduce_u]=los[best];
3839            los[best]=h;
3840            erg.reduce_by=erg.to_reduce_u;
3841            erg.fromS=FALSE;
3842            erg.to_reduce_u--;
3843          }
3844        }
3845      }
3846      else
3847      {
3848        assume(erg.to_reduce_u==erg.to_reduce_l);
3849        wlen_type quality_a=
3850          quality_of_pos_in_strat_S(erg.reduce_by,c);
3851         wlen_type qc=los[erg.to_reduce_u].guess_quality(c);
3852         if (qc<0) PrintS("Wrong wlen_type");
3853         if(qc<quality_a)
3854         {
3855           int best=erg.to_reduce_u;
3856           los[best].flatten();
3857           int b_pos=kBucketCanonicalize(los[best].bucket);
3858           los[best].p=los[best].bucket->buckets[b_pos];
3859           qc=pQuality(los[best].bucket->buckets[b_pos],c);
3860           assume(qc>=0);
3861           if(qc<quality_a)
3862           {
3863             BOOLEAN exp=FALSE;
3864             if(qc<=2)
3865             {
3866               //Print("\n qc is %lld \n",qc);
3867               exp=TRUE;
3868             }
3869             else
3870             {
3871               if (qc<quality_a/2)
3872                 exp=TRUE;
3873               else
3874                 if(erg.reduce_by<c->n/4)
3875                   exp=TRUE;
3876             }
3877             if (exp)
3878             {
3879               poly clear_into;
3880               los[erg.to_reduce_u].flatten();
3881               kBucketClear(los[erg.to_reduce_u].bucket,&clear_into,&erg.expand_length);
3882               erg.expand=pCopy(clear_into);
3883               kBucketInit(los[erg.to_reduce_u].bucket,clear_into,erg.expand_length);
3884               if (TEST_OPT_PROT) PrintS("e");
3885             }
3886           }
3887         }
3888       }
3889       swap_roles=FALSE;
3890       return;
3891     }
3892   }
3893   else
3894   {
3895     if(erg.reduce_by>erg.to_reduce_u)
3896     {
3897       //then lm(rb)>= lm(tru) so =
3898       assume(erg.reduce_by==erg.to_reduce_u+1);
3899       int best=erg.reduce_by;
3900       wlen_type quality_a=los[erg.reduce_by].guess_quality(c);
3901       wlen_type qc;
3902       best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3903
3904       int i;
3905       if(qc<quality_a)
3906       {
3907         red_object h=los[erg.reduce_by];
3908         los[erg.reduce_by]=los[best];
3909         los[best]=h;
3910       }
3911       swap_roles=FALSE;
3912       return;
3913     }
3914     else
3915     {
3916       assume(!pLmEqual(los[erg.reduce_by].p,los[erg.to_reduce_l].p));
3917       assume(erg.to_reduce_u==erg.to_reduce_l);
3918       //further assume, that reduce_by is the above all other polys
3919       //with same leading term
3920       int il=erg.reduce_by;
3921       wlen_type quality_a =los[erg.reduce_by].guess_quality(c);
3922       wlen_type qc;
3923       while((il>0) && pLmEqual(los[il-1].p,los[il].p))
3924       {
3925         il--;
3926         qc=los[il].guess_quality(c);
3927         if (qc<quality_a)
3928         {
3929           quality_a=qc;
3930           erg.reduce_by=il;
3931         }
3932       }
3933       swap_roles=FALSE;
3934     }
3935   }
3936   if(swap_roles)
3937   {
3938     if (TEST_OPT_PROT) PrintS("b");
3939     poly clear_into;
3940     int dummy_len;
3941     int new_length;
3942     int bp=erg.to_reduce_u;//bucket_positon
3943     //kBucketClear(los[bp].bucket,&clear_into,&new_length);
3944     new_length=los[bp].clear_to_poly();
3945     clear_into=los[bp].p;
3946     poly p=c->strat->S[erg.reduce_by];
3947     int j=erg.reduce_by;
3948     int old_length=c->strat->lenS[j];// in view of S
3949     los[bp].p=p;
3950     if (c->eliminationProblem)
3951     {
3952        los[bp].sugar=c->pTotaldegree_full(p);
3953     }
3954     kBucketInit(los[bp].bucket,p,old_length);
3955     wlen_type qal=pQuality(clear_into,c,new_length);
3956     int pos_in_c=-1;
3957     int z;
3958     int new_pos;
3959     new_pos=simple_posInS(c->strat,clear_into,new_length, qal);
3960     assume(new_pos<=j);
3961     for (z=c->n;z;z--)
3962     {
3963       if(p==c->S->m[z-1])
3964       {
3965         pos_in_c=z-1;
3966         break;
3967       }
3968     }
3969
3970     int tdeg_full=-1;
3971     int tdeg=-1;
3972     if(pos_in_c>=0)
3973     {
3974       #ifdef TGB_RESORT_PAIRS
3975       c->used_b=TRUE;
3976       c->replaced[pos_in_c]=TRUE;
3977       #endif
3978       tdeg=c->T_deg[pos_in_c];
3979       c->S->m[pos_in_c]=clear_into;
3980       c->lengths[pos_in_c]=new_length;
3981       c->weighted_lengths[pos_in_c]=qal;
3982       if (c->gcd_of_terms[pos_in_c]==NULL)
3983         c->gcd_of_terms[pos_in_c]=gcd_of_terms(clear_into,c->r);
3984       if (c->T_deg_full)
3985         tdeg_full=c->T_deg_full[pos_in_c]=c->pTotaldegree_full(clear_into);
3986       else tdeg_full=tdeg;
3987       c_S_element_changed_hook(pos_in_c,c);
3988     }
3989     else
3990     {
3991       if (c->eliminationProblem)
3992       {
3993         tdeg_full=c->pTotaldegree_full(clear_into);
3994         tdeg=c->pTotaldegree(clear_into);
3995       }
3996     }
3997     c->strat->S[j]=clear_into;
3998     c->strat->lenS[j]=new_length;
3999
4000     assume(pLength(clear_into)==new_length);
4001     if(c->strat->lenSw!=NULL)
4002       c->strat->lenSw[j]=qal;
4003     if (!rField_is_Zp(c->r))
4004     {
4005       p_Cleardenom(clear_into,c->r);//should be unnecessary
4006       //p_Content(clear_into, c->r);
4007     }
4008     else
4009       pNorm(clear_into);
4010#ifdef FIND_DETERMINISTIC
4011     erg.reduce_by=j;
4012     //resort later see diploma thesis, find_in_S must be deterministic
4013     //during multireduction if spolys are only in the span of the
4014     //input polys
4015#else
4016     if (new_pos<j)
4017     {
4018       if (c->strat->honey) c->strat->ecartS[j]=tdeg_full-tdeg;
4019       move_forward_in_S(j,new_pos,c->strat);
4020       erg.reduce_by=new_pos;
4021     }
4022#endif
4023   }
4024}
4025static int fwbw(red_object* los, int i){
4026   int i2=i;
4027   int step=1;
4028
4029   BOOLEAN bw=FALSE;
4030   BOOLEAN incr=TRUE;
4031
4032   while(1)
4033   {
4034     if(!bw)
4035     {
4036       step=si_min(i2,step);
4037       if (step==0) break;
4038       i2-=step;
4039
4040       if(!pLmEqual(los[i].p,los[i2].p))
4041       {
4042   bw=TRUE;
4043   incr=FALSE;
4044       }
4045       else
4046       {
4047   if ((!incr) &&(step==1)) break;
4048       }
4049
4050
4051     }
4052     else
4053     {
4054
4055       step=si_min(i-i2,step);
4056       if (step==0) break;
4057       i2+=step;
4058       if(pLmEqual(los[i].p,los[i2].p)){
4059   if(step==1) break;
4060   else
4061   {
4062     bw=FALSE;
4063   }
4064       }
4065
4066     }
4067     if (incr)
4068       step*=2;
4069     else
4070     {
4071       if (step%2==1)
4072   step=(step+1)/2;
4073       else
4074   step/=2;
4075
4076     }
4077   }
4078   return i2;
4079}
4080static void canonicalize_region(red_object* los, int l, int u,slimgb_alg* c){
4081    assume(l<=u+1);
4082    int i;
4083    for(i=l;i<=u;i++){
4084        kBucketCanonicalize(los[i].bucket);
4085    }
4086
4087}
4088static void multi_reduction_find(red_object* los, int losl,slimgb_alg* c,int startf,find_erg & erg){
4089  kStrategy strat=c->strat;
4090
4091  assume(startf<=losl);
4092  assume((startf==losl-1)||(pLmCmp(los[startf].p,los[startf+1].p)==-1));
4093  int i=startf;
4094
4095  int j;
4096  while(i>=0){
4097    assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)<=0));
4098    assume(is_valid_ro(los[i]));
4099    assume((!(c->eliminationProblem))||(los[i].sugar>=c->pTotaldegree(los[i].p)));
4100    j=kFindDivisibleByInS_easy(strat,los[i]);
4101    if(j>=0){
4102
4103      erg.to_reduce_u=i;
4104      erg.reduce_by=j;
4105      erg.fromS=TRUE;
4106      int i2=fwbw(los,i);
4107      assume(pLmEqual(los[i].p,los[i2].p));
4108      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
4109      assume(i>=i2);
4110
4111
4112      erg.to_reduce_l=i2;
4113      assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
4114      canonicalize_region(los,erg.to_reduce_u+1,startf,c);
4115      return;
4116    }
4117    if (j<0){
4118
4119      //not reduceable, try to use this for reducing higher terms
4120      int i2=fwbw(los,i);
4121      assume(pLmEqual(los[i].p,los[i2].p));
4122      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
4123      assume(i>=i2);
4124      if(i2!=i){
4125
4126
4127  erg.to_reduce_u=i-1;
4128  erg.to_reduce_l=i2;
4129  erg.reduce_by=i;
4130  erg.fromS=FALSE;
4131  assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
4132  canonicalize_region(los,erg.to_reduce_u+1,startf,c);
4133  return;
4134      }
4135
4136      i--;
4137    }
4138  }
4139  erg.reduce_by=-1;//error code
4140  return;
4141}
4142
4143 //  nicht reduzierbare eintraege in ergebnisliste schreiben
4144//   nullen loeschen
4145//   while(finde_groessten leitterm reduzierbar(c,erg)){
4146
4147static int multi_reduction_clear_zeroes(red_object* los, int  losl, int l, int u)
4148{
4149
4150
4151  int deleted=0;
4152  int  i=l;
4153  int last=-1;
4154  while(i<=u)
4155  {
4156
4157    if(los[i].p==NULL){
4158      kBucketDestroy(&los[i].bucket);
4159//      delete los[i];//here we assume los are constructed with new
4160      //destroy resources, must be added here
4161     if (last>=0)
4162     {
4163       memmove(los+(int)(last+1-deleted),los+(last+1),sizeof(red_object)*(i-1-last));
4164     }
4165     last=i;
4166     deleted++;
4167    }
4168    i++;
4169  }
4170  if((last>=0)&&(last!=losl-1))
4171      memmove(los+(int)(last+1-deleted),los+last+1,sizeof(red_object)*(losl-1-last));
4172  return deleted;
4173
4174}
4175int search_red_object_pos(red_object* a, int top, red_object* key ){
4176
4177    int an = 0;
4178    int en= top;
4179    if (top==-1) return 0;
4180    if (pLmCmp(key->p,a[top].p)==1)
4181      return top+1;
4182    int i;
4183    loop
4184    {
4185      if (an >= en-1)
4186      {
4187        if (pLmCmp(key->p,a[an].p)==-1)
4188           return an;
4189        return en;
4190      }
4191      i=(an+en) / 2;
4192      if (pLmCmp(key->p,a[i].p)==-1)
4193        en=i;
4194      else
4195        an=i;
4196    }
4197
4198}
4199static void sort_region_down(red_object* los, int l, int u, slimgb_alg* c)
4200{
4201  int r_size=u-l+1;
4202  qsort(los+l,r_size,sizeof(red_object),red_object_better_gen);
4203  int i;
4204  int * new_indices=(int*) omalloc((r_size)*sizeof(int));
4205  int bound=0;
4206  BOOLEAN at_end=FALSE;
4207  for(i=l;i<=u;i++){
4208    if (!(at_end)){
4209      bound=new_indices[i-l]=bound+search_red_object_pos(los+bound,l-bound-1,los+i);
4210      if (bound==l) at_end=TRUE;
4211    }
4212    else{
4213      new_indices[i-l]=l;
4214    }
4215  }
4216  red_object* los_region=(red_object*) omalloc(sizeof(red_object)*(u-l+1));
4217  for (int i=0;i<r_size;i++){
4218    new_indices[i]+=i;
4219    los_region[i]=los[l+i];
4220    assume((i==0)||(new_indices[i]>new_indices[i-1]));
4221
4222  }
4223
4224  i=r_size-1;
4225  int j=u;
4226  int j2=l-1;
4227  while(i>=0){
4228    if (new_indices[i]==j){
4229      los[j]=los_region[i];
4230      i--;
4231      j--;
4232    } else{
4233      assume(new_indices[i]<j);
4234      los[j]=los[j2];
4235      assume(j2>=0);
4236      j2--;
4237      j--;
4238    }
4239  }
4240  omfree(los_region);
4241
4242  omfree(new_indices);
4243
4244}
4245
4246//assume that los is ordered ascending by leading term, all non zero
4247static void multi_reduction(red_object* los, int & losl, slimgb_alg* c)
4248{
4249  poly* delay=(poly*) omalloc(losl*sizeof(poly));
4250  int delay_s=0;
4251  //initialize;
4252  assume(c->strat->sl>=0);
4253  assume(losl>0);
4254  int i;
4255  wlen_type max_initial_quality=0;
4256
4257  for(i=0;i<losl;i++){
4258    los[i].sev=pGetShortExpVector(los[i].p);
4259//SetShortExpVector();
4260    los[i].p=kBucketGetLm(los[i].bucket);
4261    if (los[i].initial_quality>max_initial_quality)
4262        max_initial_quality=los[i].initial_quality;
4263    // else
4264//         Print("init2_qal=%lld;", los[i].initial_quality);
4265//     Print("initial_quality=%lld;",max_initial_quality);
4266  }
4267
4268  kStrategy strat=c->strat;
4269  int curr_pos=losl-1;
4270
4271
4272//  nicht reduzierbare eintrᅵe in ergebnisliste schreiben
4273  // nullen loeschen
4274  while(curr_pos>=0){
4275    if ((c->use_noro_last_block)&&(lies_in_last_dp_block(los[curr_pos].p,c))){
4276        int pn_noro=curr_pos+1;
4277        poly* p_noro=(poly*) omalloc(pn_noro*sizeof(poly));
4278        for(i=0;i<pn_noro;i++){
4279            int dummy_len;
4280            poly p;
4281            los[i].p=NULL;
4282            kBucketClear(los[i].bucket,&p,&dummy_len);
4283            p_noro[i]=p;
4284        }
4285
4286
4287        if (npPrimeM<255){
4288          noro_step<tgb_uint8>(p_noro,pn_noro,c);
4289        } else {
4290          if (npPrimeM<65000){
4291            noro_step<tgb_uint16>(p_noro,pn_noro,c);
4292          } else{
4293            noro_step<tgb_uint32>(p_noro,pn_noro,c);
4294          }
4295        }
4296        for(i=0;i<pn_noro;i++){
4297            los[i].p=p_noro[i];
4298            los[i].sev=pGetShortExpVector(los[i].p);
4299            //ignore quality
4300            kBucketInit(los[i].bucket,los[i].p,pLength(los[i].p));
4301
4302        }
4303        qsort(los,pn_noro,sizeof(red_object),red_object_better_gen);
4304        int deleted=multi_reduction_clear_zeroes(los, losl, pn_noro, curr_pos);
4305        losl -= deleted;
4306        curr_pos -= deleted;
4307        break;
4308    }
4309    find_erg erg;
4310
4311    multi_reduction_find(los, losl,c,curr_pos,erg);//last argument should be curr_pos
4312    if(erg.reduce_by<0) break;
4313
4314
4315
4316    erg.expand=NULL;
4317    int d=erg.to_reduce_u-erg.to_reduce_l+1;
4318
4319
4320    multi_reduction_lls_trick(los,losl,c,erg);
4321
4322
4323    int i;
4324    int len;
4325    //    wrp(los[erg.to_reduce_u].p);
4326    //PrintLn();
4327    multi_reduce_step(erg,los,c);
4328
4329
4330    if(!TEST_OPT_REDTHROUGH)
4331    {
4332     for(i=erg.to_reduce_l;i<=erg.to_reduce_u;i++)
4333     {
4334     if  (los[i].p!=NULL)  //the check (los[i].p!=NULL) might be invalid
4335     {
4336         //
4337         assume(los[i].initial_quality>0);
4338
4339               if(los[i].guess_quality(c)
4340                  >1.5*delay_factor*max_initial_quality){
4341                       if (TEST_OPT_PROT)
4342                           PrintS("v");
4343                       los[i].canonicalize();
4344                       if(los[i].guess_quality(c)
4345                           >delay_factor*max_initial_quality){
4346                               if (TEST_OPT_PROT)
4347                                   PrintS(".");
4348                               los[i].clear_to_poly();
4349                               //delay.push_back(los[i].p);
4350                               delay[delay_s]=los[i].p;
4351                               delay_s++;
4352
4353                               los[i].p=NULL;
4354
4355                      }
4356                  }
4357
4358            }
4359     }
4360  }
4361    int deleted=multi_reduction_clear_zeroes(los, losl, erg.to_reduce_l, erg.to_reduce_u);
4362    if(erg.fromS==FALSE)
4363      curr_pos=si_max(erg.to_reduce_u,erg.reduce_by);
4364    else
4365      curr_pos=erg.to_reduce_u;
4366    losl -= deleted;
4367    curr_pos -= deleted;
4368
4369    //Print("deleted %i \n",deleted);
4370    if ((TEST_V_UPTORADICAL) &&(!(erg.fromS)))
4371        sort_region_down(los,si_min(erg.to_reduce_l,erg.reduce_by),(si_max(erg.to_reduce_u,erg.reduce_by))-deleted,c);
4372    else
4373    sort_region_down(los, erg.to_reduce_l, erg.to_reduce_u-deleted, c);
4374
4375
4376    if(erg.expand)
4377    {
4378#ifdef FIND_DETERMINISTIC
4379      int i;
4380      for(i=0;c->expandS[i];i++);
4381      c->expandS=(poly*) omrealloc(c->expandS,(i+2)*sizeof(poly));
4382      c->expandS[i]=erg.expand;
4383      c->expandS[i+1]=NULL;
4384#else
4385      int ecart=0;
4386      if (c->eliminationProblem){
4387        ecart=c->pTotaldegree_full(erg.expand)-c->pTotaldegree(erg.expand);
4388      }
4389      add_to_reductors(c,erg.expand,erg.expand_length,ecart);
4390#endif
4391    }
4392
4393  }
4394
4395
4396  //sorted_pair_node** pairs=(sorted_pair_node**)
4397  //  omalloc(delay_s*sizeof(sorted_pair_node*));
4398  c->introduceDelayedPairs(delay,delay_s);
4399  /*
4400  for(i=0;i<delay_s;i++){
4401
4402      poly p=delay[i];
4403      //if (rPar(c->r)==0)
4404      simplify_poly(p,c->r);
4405      sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
4406      si->i=-1;
4407      si->j=-1;
4408       if (!rField_is_Zp(c->r)){
4409        if (!c->nc)
4410            p=redTailShort(p, c->strat);
4411        p_Cleardenom(p, c->r);
4412        p_Content(p, c->r);
4413      }
4414      si->expected_length=pQuality(p,c,pLength(p));
4415      si->deg=pTotaldegree(p);
4416
4417      si->lcm_of_lm=p;
4418      pairs[i]=si;
4419  }
4420  qsort(pairs,delay_s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
4421  c->apairs=spn_merge(c->apairs,c->pair_top+1,pairs,delay_s,c);
4422  c->pair_top+=delay_s;*/
4423  omfree(delay);
4424  //omfree(pairs);
4425  return;
4426}
4427void red_object::flatten(){
4428  assume(p==kBucketGetLm(bucket));
4429}
4430void red_object::validate(){
4431  p=kBucketGetLm(bucket);
4432  if(p)
4433    sev=pGetShortExpVector(p);
4434}
4435int red_object::clear_to_poly(){
4436  flatten();
4437  int l;
4438  kBucketClear(bucket,&p,&l);
4439  return l;
4440}
4441
4442
4443
4444
4445
4446void reduction_step::reduce(red_object* r, int l, int u){}
4447void simple_reducer::do_reduce(red_object & ro)
4448{
4449  number coef;
4450#ifdef HAVE_PLURAL
4451  if (c->nc)
4452    nc_BucketPolyRed_Z(ro.bucket, p, &coef);
4453  else
4454#endif
4455    coef=kBucketPolyRed(ro.bucket,p,
4456       p_len,
4457       c->strat->kNoether);
4458  nDelete(&coef);
4459}
4460
4461
4462void simple_reducer::reduce(red_object* r, int l, int u){
4463  this->pre_reduce(r,l,u);
4464  int i;
4465//debug start
4466  int im;
4467
4468
4469  if(c->eliminationProblem){
4470    assume(p_LmEqual(r[l].p,r[u].p,c->r));
4471    /*int lm_deg=pTotaldegree(r[l].p);
4472    reducer_deg=lm_deg+pTotaldegree_full(p)-pTotaldegree(p);*/
4473  }
4474
4475  for(i=l;i<=u;i++){
4476
4477
4478
4479    this->do_reduce(r[i]);
4480    if (c->eliminationProblem){
4481        r[i].sugar=si_max(r[i].sugar,reducer_deg);
4482    }
4483  }
4484  for(i=l;i<=u;i++){
4485
4486    kBucketSimpleContent(r[i].bucket);
4487    r[i].validate();
4488    #ifdef TGB_DEBUG
4489    #endif
4490  }
4491}
4492reduction_step::~reduction_step(){}
4493simple_reducer::~simple_reducer(){
4494  if(fill_back!=NULL)
4495  {
4496    kBucketInit(fill_back,p,p_len);
4497  }
4498  fill_back=NULL;
4499
4500}
4501
4502void multi_reduce_step(find_erg & erg, red_object* r, slimgb_alg* c){
4503  static int id=0;
4504  id++;
4505  unsigned long sev;
4506    BOOLEAN lt_changed=FALSE;
4507  int rn=erg.reduce_by;
4508  poly red;
4509  int red_len;
4510  simple_reducer* pointer;
4511  BOOLEAN work_on_copy=FALSE;
4512  if(erg.fromS){
4513    red=c->strat->S[rn];
4514    red_len=c->strat->lenS[rn];
4515    assume(red_len==pLength(red));
4516  }
4517  else
4518  {
4519    r[rn].flatten();
4520    kBucketClear(r[rn].bucket,&red,&red_len);
4521
4522    if (!rField_is_Zp(c->r))
4523    {
4524      p_Cleardenom(red, c->r);//should be unnecessary
4525      //p_Content(red, c->r);
4526    }
4527    pNormalize(red);
4528    if (c->eliminationProblem)
4529    {
4530        r[rn].sugar=c->pTotaldegree_full(red);
4531    }
4532
4533    if ((!(erg.fromS))&&(TEST_V_UPTORADICAL))
4534    {
4535         if (polynomial_root(red,c->r))
4536            lt_changed=TRUE;
4537            sev=p_GetShortExpVector(red,c->r);
4538    }
4539    red_len=pLength(red);
4540  }
4541  if (((TEST_V_MODPSOLVSB)&&(red_len>1))||((c->nc)||(erg.to_reduce_u-erg.to_reduce_l>5))){
4542    work_on_copy=TRUE;
4543    // poly m=pOne();
4544    poly m=c->tmp_lm;
4545    pSetCoeff(m,nInit(1));
4546    pSetComp(m,0);
4547    for(int i=1;i<=pVariables;i++)
4548      pSetExp(m,i,(pGetExp(r[erg.to_reduce_l].p, i)-pGetExp(red,i)));
4549    pSetm(m);
4550    poly red_cp;
4551    #ifdef HAVE_PLURAL
4552    if (c->nc)
4553      red_cp = nc_mm_Mult_pp(m, red, c->r);
4554    else
4555    #endif
4556      red_cp=ppMult_mm(red,m);
4557    if(!erg.fromS){
4558      kBucketInit(r[rn].bucket,red,red_len);
4559    }
4560    //now reduce the copy
4561    //static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
4562
4563    if (!c->nc)
4564      redTailShort(red_cp,c->strat);
4565    //number mul;
4566    // red_len--;
4567//     red_cp->next=redNF2(red_cp->next,c,red_len,mul,c->average_length);
4568//     pSetCoeff(red_cp,nMult(red_cp->coef,mul));
4569//     nDelete(&mul);
4570//     red_len++;
4571    red=red_cp;
4572    red_len=pLength(red);
4573    // pDelete(&m);
4574
4575  }
4576  int i;
4577
4578
4579
4580  assume(red_len==pLength(red));
4581
4582  int reducer_deg=0;
4583  if (c->eliminationProblem){
4584     int lm_deg=c->pTotaldegree(r[erg.to_reduce_l].p);
4585     int ecart;
4586     if (erg.fromS){
4587       ecart=c->strat->ecartS[erg.reduce_by];
4588     } else {
4589       ecart=c->pTotaldegree_full(red)-lm_deg;
4590     }
4591     reducer_deg=lm_deg+ecart;
4592  }
4593  pointer=new simple_reducer(red,red_len,reducer_deg,c);
4594
4595  if ((!work_on_copy) && (!erg.fromS))
4596    pointer->fill_back=r[rn].bucket;
4597  else
4598    pointer->fill_back=NULL;
4599  pointer->reduction_id=id;
4600  pointer->c=c;
4601
4602  pointer->reduce(r,erg.to_reduce_l, erg.to_reduce_u);
4603  if(work_on_copy) pDelete(&pointer->p);
4604  delete pointer;
4605  if (lt_changed){
4606    assume(!erg.fromS);
4607    r[erg.reduce_by].sev=sev;
4608  }
4609
4610};
4611
4612
4613
4614
4615void simple_reducer:: pre_reduce(red_object* r, int l, int u){}
4616
Note: See TracBrowser for help on using the repository browser.