source: git/kernel/tgb.cc @ 690e21e

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