source: git/kernel/tgb.cc @ 338842d

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