source: git/kernel/tgb.cc @ 3a0e1a

spielwiese
Last change on this file since 3a0e1a was 6a2e9c, checked in by Hans Schönemann <hannes@…>, 15 years ago
*hannes: gcc 4.3 git-svn-id: file:///usr/local/Singular/svn/trunk@11871 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 103.6 KB
Line 
1//! \file tgb.cc
2//       multiple rings
3//       shorten_tails und dessen Aufrufe pruefen wlength!!!
4/****************************************
5*  Computer Algebra System SINGULAR     *
6****************************************/
7/* $Id: tgb.cc,v 1.163 2009-06-04 08:55:36 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, const char* prot, slimgb_alg* c)
1159{
1160    int i=0;
1161    //check, if it is already in the queue
1162
1163
1164    while(c->add_later->m[i]!=NULL)
1165    {
1166        if (p_LmEqual(c->add_later->m[i],p,c->r))
1167            return;
1168        i++;
1169    }
1170    if (TEST_OPT_PROT)
1171        PrintS(prot);
1172    c->add_later->m[i]=p;
1173}
1174static int simple_posInS (kStrategy strat, poly p,int len, wlen_type wlen)
1175{
1176
1177
1178  if(strat->sl==-1) return 0;
1179  if (strat->lenSw) return pos_helper(strat,p,(wlen_type) wlen,(wlen_set) strat->lenSw,strat->S);
1180  return pos_helper(strat,p,len,strat->lenS,strat->S);
1181
1182}
1183
1184/*2
1185 *if the leading term of p
1186 *divides the leading term of some S[i] it will be canceled
1187 */
1188static inline void clearS (poly p, unsigned long p_sev,int l, int* at, int* k,
1189                           kStrategy strat)
1190{
1191  assume(p_sev == pGetShortExpVector(p));
1192  if (!pLmShortDivisibleBy(p,p_sev, strat->S[*at], ~ strat->sevS[*at])) return;
1193  if (l>=strat->lenS[*at]) return;
1194  if (TEST_OPT_PROT)
1195    PrintS("!");
1196  mflush();
1197  //pDelete(&strat->S[*at]);
1198  deleteInS((*at),strat);
1199  (*at)--;
1200  (*k)--;
1201//  assume(lenS_correct(strat));
1202}
1203
1204
1205
1206static int iq_crit(const void* ap,const void* bp){
1207
1208  sorted_pair_node* a=*((sorted_pair_node**)ap);
1209  sorted_pair_node* b=*((sorted_pair_node**)bp);
1210  assume(a->i>a->j);
1211  assume(b->i>b->j);
1212
1213
1214  if (a->deg<b->deg) return -1;
1215  if (a->deg>b->deg) return 1;
1216  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
1217  if(comp!=0)
1218    return comp;
1219  if (a->expected_length<b->expected_length) return -1;
1220  if (a->expected_length>b->expected_length) return 1;
1221  if (a->j>b->j) return 1;
1222  if (a->j<b->j) return -1;
1223  return 0;
1224}
1225static wlen_type coeff_mult_size_estimate(int s1, int s2, ring r){
1226    if (rField_is_Q(r)) return s1+s2;
1227    else return s1*s2;
1228}
1229static wlen_type pair_weighted_length(int i, int j, slimgb_alg* c){
1230    if ((c->isDifficultField) && (c->eliminationProblem))  {
1231        int c1=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r);
1232        int c2=slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1233        wlen_type el1=c->weighted_lengths[i]/c1;
1234        assume(el1!=0);
1235        assume(c->weighted_lengths[i] %c1==0);
1236        wlen_type el2=c->weighted_lengths[j]/c2;
1237        assume(el2!=0);
1238        assume(c->weighted_lengths[j] %c2==0);
1239        //should be * for function fields
1240        //return (c1+c2) * (el1+el2-2);
1241        wlen_type res=coeff_mult_size_estimate(c1,c2,c->r);
1242        res*=el1+el2-2;
1243        return res;
1244
1245    }
1246    if (c->isDifficultField) {
1247        //int cs=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r)+
1248        //    slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1249        if(!(TEST_V_COEFSTRAT)){
1250        wlen_type cs=
1251            coeff_mult_size_estimate(
1252                slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r),
1253                slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r),c->r);
1254        return (wlen_type)(c->lengths[i]+c->lengths[j]-2)*
1255            (wlen_type)cs;}
1256            else {
1257
1258            wlen_type cs=
1259            coeff_mult_size_estimate(
1260                slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r),
1261                slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r),c->r);
1262            cs*=cs;
1263        return (wlen_type)(c->lengths[i]+c->lengths[j]-2)*
1264            (wlen_type)cs;
1265            }
1266    }
1267    if (c->eliminationProblem) {
1268
1269        return (c->weighted_lengths[i]+c->weighted_lengths[j]-2);
1270    }
1271    return c->lengths[i]+c->lengths[j]-2;
1272
1273}
1274sorted_pair_node** add_to_basis_ideal_quotient(poly h, slimgb_alg* c, int* ip)
1275{
1276  p_Test(h,c->r);
1277  assume(h!=NULL);
1278  poly got=gcd_of_terms(h,c->r);
1279  if((got!=NULL) &&(TEST_V_UPTORADICAL)) {
1280    poly copy=p_Copy(got,c->r);
1281    //p_wrp(got,c->r);
1282    BOOLEAN changed=monomial_root(got,c->r);
1283    if (changed)
1284    {
1285         poly div_by=pDivide(copy, got);
1286         poly iter=h;
1287         while(iter){
1288            pExpVectorSub(iter,div_by);
1289            pIter(iter);
1290         }
1291         p_Delete(&div_by, c->r);
1292         PrintS("U");
1293    }
1294    p_Delete(&copy,c->r);
1295  }
1296
1297#define ENLARGE(pointer, type) pointer=(type*) omrealloc(pointer, c->array_lengths*sizeof(type))
1298//  BOOLEAN corr=lenS_correct(c->strat);
1299  BOOLEAN R_found=FALSE;
1300  void* hp;
1301  int sugar;
1302  int ecart=0;
1303  ++(c->n);
1304  ++(c->S->ncols);
1305  int i,j;
1306  i=c->n-1;
1307  sorted_pair_node** nodes=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1308  int spc=0;
1309  if(c->n>c->array_lengths){
1310    c->array_lengths=c->array_lengths*2;
1311    assume(c->array_lengths>=c->n);
1312    ENLARGE(c->T_deg, int);
1313    ENLARGE(c->tmp_pair_lm,poly);
1314    ENLARGE(c->tmp_spn,sorted_pair_node*);
1315
1316    ENLARGE(c->short_Exps,long);
1317    ENLARGE(c->lengths,int);
1318    #ifndef HAVE_BOOST
1319    #ifndef USE_STDVECBOOL
1320
1321    ENLARGE(c->states, char*);
1322    #endif
1323    #endif
1324    ENLARGE(c->gcd_of_terms,poly);
1325    //if (c->weighted_lengths!=NULL) {
1326    ENLARGE(c->weighted_lengths,wlen_type);
1327    //}
1328    //ENLARGE(c->S->m,poly);
1329
1330  }
1331  pEnlargeSet(&c->S->m,c->n-1,1);
1332  if (c->T_deg_full)
1333    ENLARGE(c->T_deg_full,int);
1334  sugar=c->T_deg[i]=c->pTotaldegree(h);
1335  if(c->T_deg_full){
1336    sugar=c->T_deg_full[i]=c->pTotaldegree_full(h);
1337    ecart=sugar-c->T_deg[i];
1338    assume(ecart>=0);
1339  }
1340
1341
1342  c->tmp_pair_lm[i]=pOne_Special(c->r);
1343
1344
1345  c->tmp_spn[i]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1346
1347
1348  c->lengths[i]=pLength(h);
1349
1350  //necessary for correct weighted length
1351
1352   if (!rField_is_Zp(c->r)){
1353    pCleardenom(h);
1354    pContent(h); //is a duplicate call, but belongs here
1355
1356  }
1357  else
1358    pNorm(h);
1359  pNormalize(h);
1360
1361  c->weighted_lengths[i]=pQuality(h, c, c->lengths[i]);
1362  c->gcd_of_terms[i]=got;
1363  #ifdef HAVE_BOOST
1364    c->states.push_back(dynamic_bitset<>(i));
1365
1366  #else
1367  #ifdef USE_STDVECBOOL
1368
1369    c->states.push_back(vector<bool>(i));
1370
1371
1372  #else
1373  if (i>0)
1374    c->states[i]=(char*)  omalloc(i*sizeof(char));
1375  else
1376    c->states[i]=NULL;
1377  #endif
1378  #endif
1379
1380  c->S->m[i]=h;
1381  c->short_Exps[i]=p_GetShortExpVector(h,c->r);
1382
1383#undef ENLARGE
1384  if (p_GetComp(h,currRing)<=c->syz_comp){
1385  for (j=0;j<i;j++){
1386
1387
1388    #ifndef HAVE_BOOST
1389    c->states[i][j]=UNCALCULATED;
1390    #endif
1391    assume(p_LmDivisibleBy(c->S->m[i],c->S->m[j],c->r)==
1392     p_LmShortDivisibleBy(c->S->m[i],c->short_Exps[i],c->S->m[j],~(c->short_Exps[j]),c->r));
1393
1394
1395    if (_p_GetComp(c->S->m[i],c->r)!=_p_GetComp(c->S->m[j],c->r)){
1396      //c->states[i][j]=UNCALCULATED;
1397      //WARNUNG: be careful
1398      continue;
1399    } else
1400    if ((!c->nc) && (c->lengths[i]==1) && (c->lengths[j]==1)){
1401      c->states[i][j]=HASTREP;
1402
1403      }
1404    else if (( (!c->nc) || (c->is_homog && rIsSCA(c->r) ) ) &&  (pHasNotCF(c->S->m[i],c->S->m[j])))
1405//     else if ((!(c->nc)) &&  (pHasNotCF(c->S->m[i],c->S->m[j])))
1406    {
1407      c->easy_product_crit++;
1408      c->states[i][j]=HASTREP;
1409      continue;
1410    }
1411    else if(extended_product_criterion(c->S->m[i],c->gcd_of_terms[i],c->S->m[j],c->gcd_of_terms[j],c))
1412    {
1413      c->states[i][j]=HASTREP;
1414      c->extended_product_crit++;
1415
1416      //PrintS("E");
1417    }
1418      //  if (c->states[i][j]==UNCALCULATED){
1419
1420    if ((TEST_V_FINDMONOM) &&(!c->nc)) {
1421        //PrintS("COMMU");
1422       //  if (c->lengths[i]==c->lengths[j]){
1423//             poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1424//             if (short_s==NULL){
1425//                 c->states[i][j]=HASTREP;
1426//             } else
1427//             {
1428//                 p_Delete(&short_s, currRing);
1429//             }
1430//         }
1431        if (c->lengths[i]+c->lengths[j]==3){
1432
1433
1434             poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1435            if (short_s==NULL)
1436            {
1437                c->states[i][j]=HASTREP;
1438            }
1439            else
1440            {
1441                assume(pLength(short_s)==1);
1442                if (TEST_V_UPTORADICAL)
1443                   monomial_root(short_s,c->r);
1444                int iS=
1445                   kFindDivisibleByInS_easy(c->strat,short_s, p_GetShortExpVector(short_s,c->r));
1446                if (iS<0)
1447                {
1448                    //PrintS("N");
1449                    if (TRUE)
1450                    {
1451                      c->states[i][j]=HASTREP;
1452                      add_later(short_s,"N",c);
1453                    }
1454                    else
1455                      p_Delete(&short_s,currRing);
1456                }
1457                else
1458                {
1459                    if (c->strat->lenS[iS]>1){
1460                        //PrintS("O");
1461                        if (TRUE) {
1462                        c->states[i][j]=HASTREP;
1463                        add_later(short_s,"O",c);
1464                        } else p_Delete(&short_s,currRing);
1465                    }
1466                    else
1467                     p_Delete(&short_s, currRing);
1468                     c->states[i][j]=HASTREP;
1469                }
1470
1471
1472            }
1473        }
1474    }
1475      //    if (short_s)
1476      //    {
1477    assume(spc<=j);
1478    sorted_pair_node* s=c->tmp_spn[spc];//(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1479    s->i=si_max(i,j);
1480    s->j=si_min(i,j);
1481    assume(s->j==j);
1482    s->expected_length=pair_weighted_length(i,j,c);//c->lengths[i]+c->lengths[j]-2;
1483
1484    poly lm=c->tmp_pair_lm[spc];//=pOne_Special();
1485
1486    pLcm(c->S->m[i], c->S->m[j], lm);
1487    pSetm(lm);
1488    p_Test(lm,c->r);
1489    s->deg=c->pTotaldegree(lm);
1490
1491    if(c->T_deg_full)//Sugar
1492    {
1493      int t_i=c->T_deg_full[s->i]-c->T_deg[s->i];
1494      int t_j=c->T_deg_full[s->j]-c->T_deg[s->j];
1495      s->deg+=si_max(t_i,t_j);
1496      //Print("\n max: %d\n",max(t_i,t_j));
1497    }
1498    p_Test(lm,c->r);
1499    s->lcm_of_lm=lm;
1500    //          pDelete(&short_s);
1501    //assume(lm!=NULL);
1502    nodes[spc]=s;
1503    spc++;
1504
1505  // }
1506  //else
1507  //{
1508        //c->states[i][j]=HASTREP;
1509  //}
1510  }
1511  }//if syz_comp end
1512
1513
1514
1515
1516  assume(spc<=i);
1517  //now ideal quotient crit
1518  qsort(nodes,spc,sizeof(sorted_pair_node*),iq_crit);
1519
1520    sorted_pair_node** nodes_final=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
1521  int spc_final=0;
1522  j=0;
1523  while(j<spc)
1524  {
1525    int lower=j;
1526    int upper;
1527    BOOLEAN has=FALSE;
1528    for(upper=lower+1;upper<spc;upper++)
1529    {
1530
1531      if(!pLmEqual(nodes[lower]->lcm_of_lm,nodes[upper]->lcm_of_lm))
1532      {
1533  break;
1534      }
1535      if (has_t_rep(nodes[upper]->i,nodes[upper]->j,c))
1536  has=TRUE;
1537
1538    }
1539    upper=upper-1;
1540    int z;
1541    assume(spc_final<=j);
1542    for(z=0;z<spc_final;z++)
1543    {
1544      if(p_LmDivisibleBy(nodes_final[z]->lcm_of_lm,nodes[lower]->lcm_of_lm,c->r))
1545      {
1546  has=TRUE;
1547  break;
1548      }
1549    }
1550
1551    if(has)
1552    {
1553      for(;lower<=upper;lower++)
1554      {
1555  //free_sorted_pair_node(nodes[lower],c->r);
1556  //omfree(nodes[lower]);
1557  nodes[lower]=NULL;
1558      }
1559      j=upper+1;
1560      continue;
1561    }
1562    else
1563    {
1564      p_Test(nodes[lower]->lcm_of_lm,c->r);
1565      nodes[lower]->lcm_of_lm=pCopy(nodes[lower]->lcm_of_lm);
1566      assume(_p_GetComp(c->S->m[nodes[lower]->i],c->r)==_p_GetComp(c->S->m[nodes[lower]->j],c->r));
1567      nodes_final[spc_final]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1568
1569      *(nodes_final[spc_final++])=*(nodes[lower]);
1570      //c->tmp_spn[nodes[lower]->j]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1571      nodes[lower]=NULL;
1572      for(lower=lower+1;lower<=upper;lower++)
1573      {
1574  //  free_sorted_pair_node(nodes[lower],c->r);
1575  //omfree(nodes[lower]);
1576  nodes[lower]=NULL;
1577      }
1578      j=upper+1;
1579      continue;
1580    }
1581  }
1582
1583  //  Print("i:%d,spc_final:%d",i,spc_final);
1584
1585
1586
1587
1588  assume(spc_final<=spc);
1589  omfree(nodes);
1590  nodes=NULL;
1591
1592  add_to_reductors(c, h, c->lengths[c->n-1], ecart,TRUE);
1593  //i=posInS(c->strat,c->strat->sl,h,0 ecart);
1594  if (!(c->nc)){
1595    if (c->lengths[c->n-1]==1)
1596      shorten_tails(c,c->S->m[c->n-1]);
1597  }
1598  //you should really update c->lengths, c->strat->lenS, and the oder of polys in strat if you sort after lengths
1599
1600  //for(i=c->strat->sl; i>0;i--)
1601  //  if(c->strat->lenS[i]<c->strat->lenS[i-1]) printf("fehler bei %d\n",i);
1602  if (c->Rcounter>50) {
1603    c->Rcounter=0;
1604    cleanS(c->strat,c);
1605  }
1606
1607#ifdef HAVE_PLURAL
1608  // for SCA:
1609  // here write at the end of nodes_final[spc_final,...,spc_final+lmdeg-1]
1610  if(rIsSCA(c->r))
1611  {
1612    const poly pNext = pNext(h);
1613
1614    if(pNext != NULL)
1615    {
1616      // for additional polynomials
1617      const unsigned int m_iFirstAltVar = scaFirstAltVar(c->r);
1618      const unsigned int m_iLastAltVar  = scaLastAltVar(c->r);
1619
1620      int N = // c->r->N;
1621              m_iLastAltVar - m_iFirstAltVar + 1; // should be enough
1622      // TODO: but we may also use got = gcd({m}_{m\in f}))!
1623
1624       poly* array_arg=(poly*)omalloc(N*sizeof(poly)); // !
1625       int j = 0;
1626
1627
1628      for( unsigned short v = m_iFirstAltVar; v <= m_iLastAltVar; v++ )
1629      // for all x_v | Ann(lm(h))
1630      if( p_GetExp(h, v, c->r) ) // TODO: use 'got' here!
1631      {
1632        assume(p_GetExp(h, v, c->r)==1);
1633
1634        poly p = sca_pp_Mult_xi_pp(v, pNext, c->r); // x_v * h;
1635
1636        if(p != NULL) // if (x_v * h != 0)
1637          array_arg[j++] = p;
1638      } // for all x_v | Ann(lm(h))
1639
1640      c->introduceDelayedPairs(array_arg, j);
1641
1642      omfree(array_arg); // !!!
1643    }
1644//     PrintS("Saturation - done!!!\n");
1645  }
1646#endif // if SCAlgebra
1647
1648
1649  if(!ip){
1650    qsort(nodes_final,spc_final,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1651
1652
1653    c->apairs=spn_merge(c->apairs,c->pair_top+1,nodes_final,spc_final,c);
1654    c->pair_top+=spc_final;
1655    clean_top_of_pair_list(c);
1656    omfree(nodes_final);
1657    return NULL;
1658  }
1659  {
1660    *ip=spc_final;
1661    return nodes_final;
1662  }
1663
1664
1665
1666}
1667
1668
1669static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
1670{
1671  m=nInit(1);
1672  if (h==NULL) return NULL;
1673
1674  assume(len==pLength(h));
1675  kStrategy strat=c->strat;
1676  if (0 > strat->sl)
1677  {
1678    return h;
1679  }
1680  int j;
1681
1682  LObject P(h);
1683  P.SetShortExpVector();
1684  P.bucket = kBucketCreate(currRing);
1685  // BOOLEAN corr=lenS_correct(strat);
1686  kBucketInit(P.bucket,P.p,len /*pLength(P.p)*/);
1687  //wlen_set lenSw=(wlen_set) c->strat->lenS;
1688  //FIXME: plainly wrong
1689  //strat->lenS;
1690  //if (strat->lenSw!=NULL)
1691  //  lenSw=strat->lenSw;
1692  //int max_pos=simple_posInS(strat,P.p);
1693  loop
1694  {
1695      int dummy=strat->sl;
1696      j=kFindDivisibleByInS_easy(strat,P.p,P.sev);
1697      //j=kFindDivisibleByInS(strat,&dummy,&P);
1698      if ((j>=0) && ((!n)||
1699        ((strat->lenS[j]<=n) &&
1700         ((strat->lenSw==NULL)||
1701         (strat->lenSw[j]<=n)))))
1702      {
1703        nNormalize(pGetCoeff(P.p));
1704#ifdef KDEBUG
1705        if (TEST_OPT_DEBUG)
1706        {
1707          PrintS("red:");
1708          wrp(h);
1709          PrintS(" with ");
1710          wrp(strat->S[j]);
1711        }
1712#endif
1713
1714        number coef=kBucketPolyRed(P.bucket,strat->S[j],
1715                                   strat->lenS[j]/*pLength(strat->S[j])*/,
1716                                   strat->kNoether);
1717  number m2=nMult(m,coef);
1718  nDelete(&m);
1719  m=m2;
1720        nDelete(&coef);
1721        h = kBucketGetLm(P.bucket);
1722
1723  if (h==NULL) {
1724    len=0;
1725    kBucketDestroy(&P.bucket);
1726    return
1727    NULL;}
1728        P.p=h;
1729        P.t_p=NULL;
1730        P.SetShortExpVector();
1731#ifdef KDEBUG
1732        if (TEST_OPT_DEBUG)
1733        {
1734          PrintS("\nto:");
1735          wrp(h);
1736          PrintLn();
1737        }
1738#endif
1739      }
1740      else
1741      {
1742        kBucketClear(P.bucket,&(P.p),&len);
1743        kBucketDestroy(&P.bucket);
1744        pNormalize(P.p);
1745  assume(len==(pLength(P.p)));
1746        return P.p;
1747      }
1748    }
1749}
1750
1751
1752
1753static poly redTailShort(poly h, kStrategy strat){
1754  if (h==NULL) return NULL;//n_Init(1,currRing);
1755  if (TEST_V_MODPSOLVSB){
1756    bit_reduce(pNext(h), strat->tailRing);
1757  }
1758  int sl=strat->sl;
1759  int i;
1760  int len=pLength(h);
1761  for(i=0;i<=strat->sl;i++){
1762    if((strat->lenS[i]>2) || ((strat->lenSw!=NULL) && (strat->lenSw[i]>2)))
1763      break;
1764  }
1765  return(redNFTail(h,i-1,strat, len));
1766}
1767
1768static void line_of_extended_prod(int fixpos,slimgb_alg* c){
1769    if (c->gcd_of_terms[fixpos]==NULL)
1770  {
1771    c->gcd_of_terms[fixpos]=gcd_of_terms(c->S->m[fixpos],c->r);
1772    if (c->gcd_of_terms[fixpos])
1773    {
1774      int i;
1775      for(i=0;i<fixpos;i++)
1776        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)))
1777{
1778          c->states[fixpos][i]=HASTREP;
1779    c->extended_product_crit++;
1780}
1781      for(i=fixpos+1;i<c->n;i++)
1782        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)))
1783  {        c->states[i][fixpos]=HASTREP;
1784  c->extended_product_crit++;
1785  }
1786    }
1787  }
1788}
1789static void c_S_element_changed_hook(int pos, slimgb_alg* c){
1790  length_one_crit(c,pos, c->lengths[pos]);
1791  if (!c->nc)
1792    line_of_extended_prod(pos,c);
1793}
1794class poly_tree_node {
1795public:
1796  poly p;
1797  poly_tree_node* l;
1798  poly_tree_node* r;
1799  int n;
1800  poly_tree_node(int sn):l(NULL),r(NULL),n(sn){}
1801};
1802class exp_number_builder{
1803public:
1804  poly_tree_node* top_level;
1805  int n;
1806  int get_n(poly p);
1807  exp_number_builder():top_level(0),n(0){}
1808};
1809int exp_number_builder::get_n(poly p){
1810  poly_tree_node** node=&top_level;
1811  while(*node!=NULL){
1812    int c=pLmCmp(p,(*node)->p);
1813    if (c==0) return (*node)->n;
1814    if (c==-1) node=&((*node)->r);
1815    else
1816      node=&((*node)->l);
1817  }
1818  (*node)= new poly_tree_node(n);
1819  n++;
1820  (*node)->p=pLmInit(p);
1821  return (*node)->n;
1822}
1823
1824//mac_polys exp are smaller iff they are greater by monomial ordering
1825//corresponding to solving linear equations notation
1826
1827//! obsolete
1828struct int_poly_pair{
1829  poly p;
1830  int n;
1831};
1832
1833
1834//! obsolete
1835void t2ippa_rec(poly* ip,int* ia, poly_tree_node* k, int &offset){
1836    if(!k) return;
1837    t2ippa_rec(ip,ia,k->l,offset);
1838    ip[offset]=k->p;
1839    ia[k->n]=offset;
1840    ++offset;
1841
1842    t2ippa_rec(ip,ia,k->r,offset);
1843    delete k;
1844  }
1845
1846//! obsolete
1847void t2ippa(poly* ip,int* ia,exp_number_builder & e){
1848
1849  int o=0;
1850  t2ippa_rec(ip,ia,e.top_level,o);
1851}
1852int anti_poly_order(const void* a, const void* b){
1853  return -pLmCmp(((int_poly_pair*) a)->p,((int_poly_pair*) b)->p );
1854}
1855
1856BOOLEAN is_valid_ro(red_object & ro){
1857  red_object r2=ro;
1858  ro.validate();
1859  if ((r2.p!=ro.p)||(r2.sev!=ro.sev)) return FALSE;
1860  return TRUE;
1861}
1862int terms_sort_crit(const void* a, const void* b){
1863  return -pLmCmp(*((poly*) a),*((poly*) b));
1864}
1865static void unify_terms(poly* terms,int & sum){
1866  if (sum==0) return;
1867  int last=0;
1868  int curr=1;
1869  while(curr<sum){
1870    if (!(pLmEqual(terms[curr],terms[last]))){
1871      terms[++last]=terms[curr];
1872    }
1873    ++curr;
1874  }
1875  sum=last+1;
1876}
1877static void export_mat(number* number_array,int pn, int tn,const char* format_str, int mat_nr){
1878  char matname[20];
1879  sprintf(matname,format_str,mat_nr);
1880  FILE* out=fopen(matname,"w");
1881  int i,j;
1882  fprintf(out,"mat=[\n");
1883  for(i=0;i<pn;i++){
1884    fprintf(out,"[\n");
1885    for(j=0;j<tn;j++)
1886    {
1887      if (j>0){
1888        fprintf(out,", ");
1889      }
1890      fprintf(out,"%i",npInt(number_array[i*tn+j]));
1891    }
1892    if (i<pn-1)
1893      fprintf(out,"],\n");
1894    else
1895      fprintf(out,"],\n");
1896  }
1897  fprintf(out,"]\n");
1898  fclose(out);
1899}
1900//typedef unsigned short number_type;
1901
1902
1903#ifdef USE_NORO
1904#ifndef NORO_CACHE
1905static void linalg_step_modp(poly *p, poly* p_out, int& pn, poly* terms,int tn, slimgb_alg* c){
1906  static int export_n=0;
1907  assume(terms[tn-1]!=NULL);
1908  assume(rField_is_Zp(c->r));
1909  //I don't do deletes, copies of number_types ...
1910  const number_type zero=0;//npInit(0);
1911  int array_size=pn*tn;
1912  number_type* number_array=(number_type*) omalloc(pn*tn*sizeof(number_type));
1913  int i;
1914  for(i=0;i<array_size;i++){
1915    number_array[i]=zero;
1916  }
1917  for(i=0;i<pn;i++){
1918    poly h=p[i];
1919    //int base=tn*i;
1920    write_poly_to_row(number_array+tn*i,h,terms,tn,c->r);
1921
1922  }
1923#if 0
1924  //export matrix
1925  export_mat(number_array,pn,tn,"mat%i.py",++export_n);
1926#endif
1927  int rank=pn;
1928  simplest_gauss_modp(number_array,rank,tn);
1929  int act_row=0;
1930  int p_pos=0;
1931  for(i=0;i<pn;i++){
1932    poly h=NULL;
1933    int j;
1934    int base=tn*i;
1935    number* row=number_array+base;
1936    h=row_to_poly(row,terms,tn,c->r);
1937
1938   if (h!=NULL)
1939   {
1940     p_out[p_pos++]=h;
1941   }
1942  }
1943  pn=p_pos;
1944  //assert(p_pos==rank)
1945  while(p_pos<pn){
1946    p_out[p_pos++]=NULL;
1947  }
1948#if 0
1949  export_mat(number_array,pn,tn,"mat%i.py",++export_n);
1950#endif
1951}
1952#endif
1953#endif
1954static void mass_add(poly* p, int pn,slimgb_alg* c){
1955    int j;
1956    int* ibuf=(int*) omalloc(pn*sizeof(int));
1957    sorted_pair_node*** sbuf=(sorted_pair_node***) omalloc(pn*sizeof(sorted_pair_node**));
1958    for(j=0;j<pn;j++){
1959      p_Test(p[j],c->r);
1960      sbuf[j]=add_to_basis_ideal_quotient(p[j],c,ibuf+j);
1961    }
1962    int sum=0;
1963    for(j=0;j<pn;j++){
1964      sum+=ibuf[j];
1965    }
1966    sorted_pair_node** big_sbuf=(sorted_pair_node**) omalloc(sum*sizeof(sorted_pair_node*));
1967    int partsum=0;
1968    for(j=0;j<pn;j++)
1969    {
1970      memmove(big_sbuf+partsum, sbuf[j],ibuf[j]*sizeof(sorted_pair_node*));
1971      omfree(sbuf[j]);
1972      partsum+=ibuf[j];
1973    }
1974
1975    qsort(big_sbuf,sum,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
1976    c->apairs=spn_merge(c->apairs,c->pair_top+1,big_sbuf,sum,c);
1977    c->pair_top+=sum;
1978    clean_top_of_pair_list(c);
1979    omfree(big_sbuf);
1980    omfree(sbuf);
1981    omfree(ibuf);
1982    //omfree(buf);
1983  #ifdef TGB_DEBUG
1984    int z;
1985    for(z=1;z<=c->pair_top;z++)
1986    {
1987      assume(pair_better(c->apairs[z],c->apairs[z-1],c));
1988    }
1989  #endif
1990
1991}
1992
1993#ifdef NORO_CACHE
1994#ifndef NORO_NON_POLY
1995void NoroCache::evaluateRows(){
1996  //after that can evaluate placeholders
1997  int i;
1998  buffer=(number*) omalloc(nIrreducibleMonomials*sizeof(number));
1999  for(i=0;i<root.branches_len;i++){
2000    evaluateRows(1,root.branches[i]);
2001  }
2002  omfree(buffer);
2003  buffer=NULL;
2004}
2005void NoroCache::evaluateRows(int level, NoroCacheNode* node){
2006  assume(level>=0);
2007  if (node==NULL) return;
2008  if (level<pVariables){
2009    int i,sum;
2010    for(i=0;i<node->branches_len;i++){
2011      evaluateRows(level+1,node->branches[i]);
2012    }
2013  } else {
2014    DataNoroCacheNode* dn=(DataNoroCacheNode*) node;
2015    if (dn->value_len!=backLinkCode){
2016      poly p=dn->value_poly;
2017      #ifndef NORO_SPARSE_ROWS_PRE
2018      dn->row=new DenseRow();
2019      DenseRow* row=dn->row;
2020      memset(buffer,0,sizeof(number)*nIrreducibleMonomials);
2021
2022      if (p==NULL) {row->array=NULL;row->begin=0;row->end=0; return;}
2023      int i=0;
2024      int idx;
2025      number* a=buffer;
2026      while(p){
2027        DataNoroCacheNode* ref=getCacheReference(p);
2028
2029        idx=ref->term_index;
2030        assume(idx>=0);
2031        a[idx]=p_GetCoeff(p,currRing);
2032        if (i==0) row->begin=idx;
2033        i++;
2034        pIter(p);
2035      }
2036      row->end=idx+1;
2037      assume(row->end>row->begin);
2038      int len=row->end-row->begin;
2039      row->array=(number*) omalloc((len)*sizeof(number));
2040      memcpy(row->array,a+row->begin,len*sizeof(number));
2041      #else
2042      assume(dn->value_len==pLength(dn->value_poly));
2043      dn->row=new SparseRow(dn->value_len);
2044      SparseRow* row=dn->row;
2045      int i=0;
2046      while(p){
2047        DataNoroCacheNode* ref=getCacheReference(p);
2048
2049        int idx=ref->term_index;
2050        assume(idx>=0);
2051        row->idx_array[i]=idx;
2052        row->coef_array[i]=p_GetCoeff(p,currRing);
2053
2054        i++;
2055        pIter(p);
2056      }
2057      if (i!=dn->value_len){
2058        PrintS("F4 calc wrong, as poly len was wrong\n");
2059      }
2060      assume(i==dn->value_len);
2061      #endif
2062
2063    }
2064  }
2065}
2066
2067void NoroCache::evaluatePlaceHolder(number* row,std::vector<NoroPlaceHolder>& place_holders){
2068  int i;
2069  int s=place_holders.size();
2070  for(i=0;i<s;i++){
2071    DataNoroCacheNode* ref=place_holders[i].ref;
2072    number coef=place_holders[i].coef;
2073    if (ref->value_len==backLinkCode){
2074      row[ref->term_index]=npAddM(row[ref->term_index],coef);
2075    } else {
2076      #ifndef NORO_SPARSE_ROWS_PRE
2077      DenseRow* ref_row=ref->row;
2078      if (ref_row==NULL) continue;
2079      number* ref_begin=ref_row->array;
2080      number* ref_end=ref_row->array+(ref_row->end-ref_row->begin);
2081      number* my_pos=row+ref_row->begin;
2082      //TODO npisOne distinction
2083      if (!(npIsOne(coef))){
2084        while(ref_begin!=ref_end){
2085
2086          *my_pos=npAddM(*my_pos,npMult(coef,*ref_begin));
2087          ++ref_begin;
2088          ++my_pos;
2089        }
2090      }
2091      else{
2092        while(ref_begin!=ref_end){
2093
2094          *my_pos=npAddM(*my_pos,*ref_begin);
2095          ++ref_begin;
2096          ++my_pos;
2097        }
2098      }
2099
2100
2101
2102    #else
2103    SparseRow* ref_row=ref->row;
2104    if (ref_row==NULL) continue;
2105    int n=ref_row->len;
2106    int j;
2107    int* idx_array=ref_row->idx_array;
2108    number* coef_array=ref_row->coef_array;
2109    for(j=0;j<n;j++){
2110      int idx=idx_array[j];
2111      number ref_coef=coef_array[j];
2112      row[idx]=npAddM(row[idx],npMult(coef,ref_coef));
2113    }
2114    #endif
2115  }
2116  }
2117
2118}
2119#endif
2120
2121
2122
2123//poly noro_red_non_unique(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2124
2125
2126#ifndef NORO_NON_POLY
2127MonRedRes noro_red_mon(poly t, BOOLEAN force_unique, NoroCache* cache,slimgb_alg* c){
2128  MonRedRes res_holder;
2129
2130  //wrp(t);
2131  res_holder.changed=TRUE;
2132  if (force_unique){
2133    DataNoroCacheNode* ref=cache->getCacheReference(t);
2134    if (ref!=NULL){
2135      res_holder.len=ref->value_len;
2136      if (res_holder.len==NoroCache::backLinkCode){
2137        res_holder.len=1;
2138
2139
2140      }
2141      res_holder.coef=p_GetCoeff(t,c->r);
2142      res_holder.p=ref->value_poly;
2143      res_holder.ref=ref;
2144      res_holder.onlyBorrowed=TRUE;
2145      res_holder.changed=TRUE;
2146      p_Delete(&t,c->r);
2147      return res_holder;
2148    }
2149  } else{
2150    BOOLEAN succ;
2151    poly cache_lookup=cache->lookup(t,succ, res_holder.len);//don't own this yet
2152    if (succ){
2153      if (cache_lookup==t){
2154      //know they are equal
2155      //res_holder.len=1;
2156
2157        res_holder.changed=FALSE;
2158        res_holder.p=t;
2159        res_holder.coef=npInit(1);
2160
2161        res_holder.onlyBorrowed=FALSE;
2162        return res_holder;
2163      }
2164
2165      res_holder.coef=p_GetCoeff(t,c->r);
2166      p_Delete(&t,c->r);
2167
2168      res_holder.p=cache_lookup;
2169
2170      res_holder.onlyBorrowed=TRUE;
2171      return res_holder;
2172
2173    }
2174  }
2175
2176  unsigned long sev=p_GetShortExpVector(t,currRing);
2177  int i=kFindDivisibleByInS_easy(c->strat,t,sev);
2178  if (i>=0){
2179    number coef_bak=p_GetCoeff(t,c->r);
2180
2181    p_SetCoeff(t,npInit(1),c->r);
2182    assume(npIsOne(p_GetCoeff(c->strat->S[i],c->r)));
2183    number coefstrat=p_GetCoeff(c->strat->S[i],c->r);
2184
2185      //poly t_copy_mon=p_Copy(t,c->r);
2186    poly exp_diff=cache->temp_term;
2187    p_ExpVectorDiff(exp_diff,t,c->strat->S[i],c->r);
2188    p_SetCoeff(exp_diff,npNeg(npInvers(coefstrat)),c->r);
2189    p_Setm(exp_diff,c->r);
2190    assume(c->strat->S[i]!=NULL);
2191      //poly t_to_del=t;
2192    poly res;
2193    res=pp_Mult_mm(pNext(c->strat->S[i]),exp_diff,c->r);
2194
2195    res_holder.len=c->strat->lenS[i]-1;
2196    res=noro_red_non_unique(res,res_holder.len,cache,c);
2197
2198    DataNoroCacheNode* ref=cache->insert(t,res,res_holder.len);
2199    p_Delete(&t,c->r);
2200      //p_Delete(&t_copy_mon,c->r);
2201      //res=pMult_nn(res,coef_bak);
2202
2203    res_holder.changed=TRUE;
2204    res_holder.p=res;
2205    res_holder.coef=coef_bak;
2206    res_holder.onlyBorrowed=TRUE;
2207    res_holder.ref=ref;
2208    return res_holder;
2209
2210  } else {
2211    number coef_bak=p_GetCoeff(t,c->r);
2212    number one=npInit(1);
2213    p_SetCoeff(t,one,c->r);
2214    res_holder.len=1;
2215    if (!(force_unique)){
2216      res_holder.ref=cache->insert(t,t,res_holder.len);
2217      p_SetCoeff(t,coef_bak,c->r);
2218      //return t;
2219
2220      //we need distinction
2221      res_holder.changed=FALSE;
2222      res_holder.p=t;
2223
2224      res_holder.coef=npInit(1);
2225      res_holder.onlyBorrowed=FALSE;
2226      return res_holder;
2227    } else {
2228      res_holder.ref=cache->insertAndTransferOwnerShip(t,c->r);
2229      res_holder.coef=coef_bak;
2230      res_holder.onlyBorrowed=TRUE;
2231      res_holder.changed=FALSE;
2232      res_holder.p=t;
2233      return res_holder;
2234    }
2235  }
2236
2237}
2238#endif
2239//SparseRow* noro_red_to_non_poly(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2240#ifndef NORO_NON_POLY
2241//len input and out: Idea: reverse addition
2242poly noro_red_non_unique(poly p, int &len, NoroCache* cache,slimgb_alg* c){
2243  assume(len==pLength(p));
2244  poly orig_p=p;
2245  if (p==NULL) {
2246    len=0;
2247    return NULL;
2248  }
2249  kBucket_pt bucket=kBucketCreate(currRing);
2250  kBucketInit(bucket,NULL,0);
2251  poly unchanged_head=NULL;
2252  poly unchanged_tail=NULL;
2253  int unchanged_size=0;
2254
2255  while(p){
2256    poly t=p;
2257    pIter(p);
2258    pNext(t)=NULL;
2259#ifndef NDEBUG
2260    number coef_debug=p_GetCoeff(t,currRing);
2261#endif
2262    MonRedRes red=noro_red_mon(t,FALSE,cache,c);
2263    if ((!(red.changed))&&(!(red.onlyBorrowed))){
2264      unchanged_size++;
2265      assume(npIsOne(red.coef));
2266      assume(p_GetCoeff(red.p,currRing)==coef_debug);
2267      if (unchanged_head){
2268        pNext(unchanged_tail)=red.p;
2269        pIter(unchanged_tail);
2270      } else{
2271        unchanged_tail=red.p;
2272        unchanged_head=red.p;
2273      }
2274    } else{
2275      assume(red.len==pLength(red.p));
2276      if (red.onlyBorrowed){
2277        if (npIsOne(red.coef)){
2278          t=p_Copy(red.p,currRing);
2279        }else
2280        t=pp_Mult_nn(red.p,red.coef,currRing);
2281      } else {
2282        if (npIsOne(red.coef))
2283        t=red.p;
2284        else
2285         t=p_Mult_nn(red.p,red.coef,currRing);
2286      }
2287      kBucket_Add_q(bucket,t,&red.len);
2288    }
2289  }
2290  poly res=NULL;
2291  len=0;
2292  kBucket_Add_q(bucket,unchanged_head,&unchanged_size);
2293  kBucketClear(bucket,&res,&len);
2294  kBucketDestroy(&bucket);
2295
2296
2297
2298
2299  return res;
2300}
2301#endif
2302#ifdef NORO_SPARSE_ROWS_PRE
2303//len input and out: Idea: reverse addition
2304
2305/*template <class number_type> SparseRow<number_type>* noro_red_to_non_poly(poly p, int &len, NoroCache<number_type>* cache,slimgb_alg* c){
2306  if (npPrimeM<255){
2307    return noro_red_to_non_poly_t<tgb_uint8>(p,len,cache,c);
2308  } else {
2309    if (npPrimeM<65000){
2310      return noro_red_to_non_poly_t<tgb_uint16>(p,len,cache,c);
2311    } else{
2312      return noro_red_to_non_poly_t<tgb_uint32>(p,len,cache,c);
2313    }
2314  }
2315}*/
2316#endif
2317//len input and out: Idea: reverse addition
2318#ifndef NORO_NON_POLY
2319std::vector<NoroPlaceHolder> noro_red(poly p, int &len, NoroCache* cache,slimgb_alg* c){
2320  std::vector<NoroPlaceHolder> res;
2321   while(p){
2322      poly t=p;
2323      pIter(p);
2324      pNext(t)=NULL;
2325
2326      MonRedRes red=noro_red_mon(t,TRUE,cache,c);
2327      assume(red.onlyBorrowed);
2328      assume(red.coef);
2329      assume(red.ref);
2330      NoroPlaceHolder h;
2331      h.ref=red.ref;
2332      h.coef=red.coef;
2333      assume(!((h.ref->value_poly==NULL) &&(h.ref->value_len!=0)));
2334      if (h.ref->value_poly)
2335        res.push_back(h);
2336    }
2337    return res;
2338}
2339#endif
2340
2341
2342
2343
2344
2345#endif
2346#ifdef USE_NORO
2347#ifndef NORO_CACHE
2348void noro_step(poly*p,int &pn,slimgb_alg* c){
2349  poly* reduced=(poly*) omalloc(pn*sizeof(poly));
2350  int j;
2351  int* reduced_len=(int*) omalloc(pn*sizeof(int));
2352  int reduced_c=0;
2353  //if (TEST_OPT_PROT)
2354  //  PrintS("reduced system:\n");
2355#ifdef NORO_CACHE
2356  NoroCache cache;
2357#endif
2358  for(j=0;j<pn;j++){
2359
2360    poly h=p[j];
2361    int h_len=pLength(h);
2362
2363    number coef;
2364#ifndef NORO_CACHE
2365    h=redNF2(p_Copy(h,c->r),c,h_len,coef,0);
2366#else
2367    h=noro_red(p_Copy(h,c->r),h_len,&cache,c);
2368    assume(pLength(h)==h_len);
2369#endif
2370    if (h!=NULL){
2371#ifndef NORO_CACHE
2372
2373      h=redNFTail(h,c->strat->sl,c->strat,h_len);
2374      h_len=pLength(h);
2375#endif
2376
2377
2378      reduced[reduced_c]=h;
2379      reduced_len[reduced_c]=h_len;
2380      reduced_c++;
2381      if (TEST_OPT_PROT)
2382        Print("%d ",h_len);
2383    }
2384  }
2385  int reduced_sum=0;
2386  for(j=0;j<reduced_c;j++){
2387    reduced_sum+=reduced_len[j];
2388  }
2389  poly* terms=(poly*) omalloc(reduced_sum*sizeof(poly));
2390  int tc=0;
2391  for(j=0;j<reduced_c;j++){
2392    poly h=reduced[j];
2393
2394    while(h!=NULL){
2395      terms[tc++]=h;
2396      pIter(h);
2397      assume(tc<=reduced_sum);
2398    }
2399  }
2400  assume(tc==reduced_sum);
2401  qsort(terms,reduced_sum,sizeof(poly),terms_sort_crit);
2402  int nterms=reduced_sum;
2403  //if (TEST_OPT_PROT)
2404  //Print("orig estimation:%i\n",reduced_sum);
2405  unify_terms(terms,nterms);
2406  //if (TEST_OPT_PROT)
2407  //    Print("actual number of columns:%i\n",nterms);
2408  int rank=reduced_c;
2409  linalg_step_modp(reduced, p,rank,terms,nterms,c);
2410  omfree(terms);
2411
2412  pn=rank;
2413  omfree(reduced);
2414
2415  if (TEST_OPT_PROT)
2416    PrintS("\n");
2417}
2418#else
2419
2420
2421
2422#endif
2423#endif
2424static void go_on (slimgb_alg* c){
2425  //set limit of 1000 for multireductions, at the moment for
2426  //programming reasons
2427  #ifdef USE_NORO
2428  //Print("module rank%d\n",c->S->rank);
2429    const BOOLEAN use_noro=c->use_noro;
2430  #else
2431  const BOOLEAN use_noro=FALSE;
2432  #endif
2433  int i=0;
2434  c->average_length=0;
2435  for(i=0;i<c->n;i++){
2436    c->average_length+=c->lengths[i];
2437  }
2438  c->average_length=c->average_length/c->n;
2439  i=0;
2440  int max_pairs=bundle_size;
2441
2442  #ifdef USE_NORO
2443  if ((use_noro)||(c->use_noro_last_block))
2444    max_pairs=bundle_size_noro;
2445  #endif
2446  poly* p=(poly*) omalloc((max_pairs+1)*sizeof(poly));//nullterminated
2447
2448  int curr_deg=-1;
2449  while(i<max_pairs){
2450    sorted_pair_node* s=top_pair(c);//here is actually chain criterium done
2451
2452
2453    if (!s) break;
2454
2455    if(curr_deg>=0){
2456      if (s->deg >curr_deg) break;
2457    }
2458
2459    else curr_deg=s->deg;
2460    quick_pop_pair(c);
2461    if(s->i>=0){
2462      //be careful replace_pair use createShortSpoly which is not noncommutative
2463      now_t_rep(s->i,s->j,c);
2464    replace_pair(s->i,s->j,c);
2465
2466    if(s->i==s->j) {
2467      free_sorted_pair_node(s,c->r);
2468      continue;
2469    }
2470    now_t_rep(s->i,s->j,c);
2471    }
2472    poly h;
2473    if(s->i>=0)
2474    {
2475#ifdef HAVE_PLURAL
2476      if (c->nc){
2477        h= nc_CreateSpoly(c->S->m[s->i], c->S->m[s->j]/*, NULL*/, c->r);
2478
2479        if (h!=NULL)
2480          pCleardenom(h);
2481      }
2482      else
2483#endif
2484        h=ksOldCreateSpoly(c->S->m[s->i], c->S->m[s->j], NULL, c->r);
2485    p_Test(h,c->r);
2486    }
2487    else{
2488      h=s->lcm_of_lm;
2489      p_Test(h,c->r);
2490  }
2491    // if(s->i>=0)
2492//       now_t_rep(s->j,s->i,c);
2493    number coef;
2494    int mlen=pLength(h);
2495    p_Test(h,c->r);
2496    if ((!c->nc)&(!(use_noro))){
2497      h=redNF2(h,c,mlen,coef,2);
2498      redTailShort(h,c->strat);
2499      nDelete(&coef);
2500    }
2501    p_Test(h,c->r);
2502    free_sorted_pair_node(s,c->r);
2503    if(!h) continue;
2504    int len=pLength(h);
2505    p[i]=h;
2506
2507    i++;
2508  }
2509  p[i]=NULL;
2510//  pre_comp(p,i,c);
2511  if(i==0){
2512    omfree(p);
2513    return;
2514  }
2515  #ifdef TGB_RESORT_PAIRS
2516  c->replaced=new bool[c->n];
2517  c->used_b=FALSE;
2518  #endif
2519
2520  c->normal_forms+=i;
2521  int j;
2522#ifdef USE_NORO
2523  //if ((!(c->nc))&&(rField_is_Zp(c->r))){
2524  if (use_noro){
2525    int pn=i;
2526    if (pn==0) {omfree(p);return;}
2527
2528    {
2529
2530      if (npPrimeM<255){
2531        noro_step<tgb_uint8>(p,pn,c);
2532      } else {
2533        if (npPrimeM<65000){
2534          noro_step<tgb_uint16>(p,pn,c);
2535        } else{
2536          noro_step<tgb_uint32>(p,pn,c);
2537        }
2538      }
2539
2540
2541
2542    }
2543
2544    //if (TEST_OPT_PROT){
2545    //  Print("reported rank:%i\n",pn);
2546    //}
2547    mass_add(p,pn,c);
2548    omfree(p);
2549    return;
2550    /*if (TEST_OPT_PROT)
2551      for(j=0;j<pn;j++){
2552        p_wrp(p[j],c->r);
2553      }*/
2554  }
2555#endif
2556  red_object* buf=(red_object*) omalloc(i*sizeof(red_object));
2557  for(j=0;j<i;j++){
2558    p_Test(p[j],c->r);
2559    buf[j].p=p[j];
2560    buf[j].sev=pGetShortExpVector(p[j]);
2561    buf[j].bucket = kBucketCreate(currRing);
2562    p_Test(p[j],c->r);
2563    if (c->eliminationProblem){
2564        buf[j].sugar=c->pTotaldegree_full(p[j]);
2565    }
2566    int len=pLength(p[j]);
2567    kBucketInit(buf[j].bucket,buf[j].p,len);
2568    buf[j].initial_quality=buf[j].guess_quality(c);
2569    assume(buf[j].initial_quality>=0);
2570  }
2571  omfree(p);
2572  qsort(buf,i,sizeof(red_object),red_object_better_gen);
2573//    Print("\ncurr_deg:%i\n",curr_deg);
2574  if (TEST_OPT_PROT)
2575  {
2576    Print("%dM[%d,",curr_deg,i);
2577  }
2578
2579  multi_reduction(buf, i, c);
2580  #ifdef TGB_RESORT_PAIRS
2581  if (c->used_b) {
2582    if (TEST_OPT_PROT)
2583        PrintS("B");
2584    int e;
2585    for(e=0;e<=c->pair_top;e++){
2586        if(c->apairs[e]->i<0) continue;
2587        assume(c->apairs[e]->j>=0);
2588        if ((c->replaced[c->apairs[e]->i])||(c->replaced[c->apairs[e]->j])) {
2589            sorted_pair_node* s=c->apairs[e];
2590            s->expected_length=pair_weighted_length(s->i,s->j,c);
2591        }
2592    }
2593    qsort(c->apairs,c->pair_top+1,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
2594  }
2595  #endif
2596#ifdef TGB_DEBUG
2597 {
2598   int k;
2599   for(k=0;k<i;k++)
2600   {
2601     assume(kFindDivisibleByInS_easy(c->strat,buf[k])<0);
2602     int k2;
2603     for(k2=0;k2<i;k2++)
2604     {
2605       if(k==k2) continue;
2606       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));
2607     }
2608   }
2609 }
2610#endif
2611  //resort S
2612
2613  if (TEST_OPT_PROT)
2614      Print("%i]",i);
2615
2616  poly* add_those=(poly*) omalloc(i*sizeof(poly));
2617  for(j=0;j<i;j++)
2618  {
2619    int len;
2620    poly p;
2621    buf[j].flatten();
2622    kBucketClear(buf[j].bucket,&p, &len);
2623    kBucketDestroy(&buf[j].bucket);
2624    p_Test(p,c->r);
2625    //if (!c->nc) {
2626      if ((c->tailReductions) ||(lies_in_last_dp_block(p,c))){
2627      p=redNFTail(p,c->strat->sl,c->strat, 0);
2628      } else {
2629      p=redTailShort(p, c->strat);
2630      }
2631      //}
2632      p_Test(p,c->r);
2633      add_those[j]=p;
2634
2635    //sbuf[j]=add_to_basis(p,-1,-1,c,ibuf+j);
2636  }
2637  mass_add(add_those,i,c);
2638  omfree(add_those);
2639  omfree(buf);
2640
2641
2642
2643
2644
2645
2646
2647  if (TEST_OPT_PROT)
2648      Print("(%d)",c->pair_top+1);
2649  //TODO: implement that while(!(idIs0(c->add_later)))
2650  #ifdef TGB_RESORT_PAIRS
2651  delete c->replaced;
2652  c->replaced=NULL;
2653  c->used_b=FALSE;
2654  #endif
2655  return;
2656}
2657
2658
2659
2660#ifdef REDTAIL_S
2661
2662static poly redNFTail (poly h,const int sl,kStrategy strat, int len)
2663{
2664  BOOLEAN nc=rIsPluralRing(currRing);
2665  if (h==NULL) return NULL;
2666  pTest(h);
2667  if (0 > sl)
2668    return h;
2669  if (pNext(h)==NULL) return h;
2670
2671  int j;
2672  poly res=h;
2673  poly act=res;
2674  LObject P(pNext(h));
2675  pNext(res)=NULL;
2676  P.bucket = kBucketCreate(currRing);
2677  len--;
2678  h=P.p;
2679  if (len <=0) len=pLength(h);
2680  kBucketInit(P.bucket,h /*P.p*/,len /*pLength(P.p)*/);
2681  pTest(h);
2682  loop
2683  {
2684      P.p=h;
2685      P.t_p=NULL;
2686      P.SetShortExpVector();
2687      loop
2688      {
2689          int dummy=strat->sl;
2690          j=kFindDivisibleByInS_easy(strat,P.p,P.sev);//kFindDivisibleByInS(strat,&dummy,&P);
2691          if (j>=0)
2692          {
2693#ifdef REDTAIL_PROT
2694            PrintS("r");
2695#endif
2696            nNormalize(pGetCoeff(P.p));
2697#ifdef KDEBUG
2698            if (TEST_OPT_DEBUG)
2699            {
2700              PrintS("red tail:");
2701              wrp(h);
2702              PrintS(" with ");
2703              wrp(strat->S[j]);
2704            }
2705#endif
2706            number coef;
2707            pTest(strat->S[j]);
2708#ifdef HAVE_PLURAL
2709            if (nc){
2710              nc_BucketPolyRed_Z(P.bucket, strat->S[j], &coef);
2711            } else
2712#endif
2713              coef=kBucketPolyRed(P.bucket,strat->S[j],
2714                                strat->lenS[j]/*pLength(strat->S[j])*/,strat->kNoether);
2715            pMult_nn(res,coef);
2716            nDelete(&coef);
2717            h = kBucketGetLm(P.bucket);
2718            pTest(h);
2719            if (h==NULL)
2720            {
2721#ifdef REDTAIL_PROT
2722              PrintS(" ");
2723#endif
2724        kBucketDestroy(&P.bucket);
2725              return res;
2726            }
2727            pTest(h);
2728            P.p=h;
2729            P.t_p=NULL;
2730            P.SetShortExpVector();
2731#ifdef KDEBUG
2732            if (TEST_OPT_DEBUG)
2733            {
2734              PrintS("\nto tail:");
2735              wrp(h);
2736              PrintLn();
2737            }
2738#endif
2739          }
2740          else
2741          {
2742#ifdef REDTAIL_PROT
2743            PrintS("n");
2744#endif
2745            break;
2746          }
2747      } /* end loop current mon */
2748      //   poly tmp=pHead(h /*kBucketGetLm(P.bucket)*/);
2749      //act->next=tmp;pIter(act);
2750      act->next=kBucketExtractLm(P.bucket);pIter(act);
2751      h = kBucketGetLm(P.bucket);
2752      if (h==NULL)
2753      {
2754#ifdef REDTAIL_PROT
2755        PrintS(" ");
2756#endif
2757  kBucketDestroy(&P.bucket);
2758        return res;
2759      }
2760      pTest(h);
2761  }
2762}
2763#endif
2764
2765
2766//try to fill, return FALSE iff queue is empty
2767
2768//transfers ownership of m to mat
2769void init_with_mac_poly(tgb_sparse_matrix* mat, int row, mac_poly m){
2770  assume(mat->mp[row]==NULL);
2771  mat->mp[row]=m;
2772#ifdef TGB_DEBUG
2773  mac_poly r=m;
2774  while(r){
2775    assume(r->exp<mat->columns);
2776    r=r->next;
2777  }
2778#endif
2779}
2780poly free_row_to_poly(tgb_sparse_matrix* mat, int row, poly* monoms, int monom_index){
2781  poly p=NULL;
2782  poly* set_this=&p;
2783  mac_poly r=mat->mp[row];
2784  mat->mp[row]=NULL;
2785  while(r)
2786  {
2787    (*set_this)=pLmInit(monoms[monom_index-1-r->exp]);
2788    pSetCoeff((*set_this),r->coef);
2789    set_this=&((*set_this)->next);
2790    mac_poly old=r;
2791    r=r->next;
2792    delete old;
2793
2794  }
2795  return p;
2796
2797}
2798
2799
2800
2801
2802
2803
2804
2805
2806
2807
2808
2809static int poly_crit(const void* ap1, const void* ap2){
2810  poly p1,p2;
2811  p1=*((poly*) ap1);
2812  p2=*((poly*)ap2);
2813
2814  int c=pLmCmp(p1,p2);
2815  if (c !=0) return c;
2816  int l1=pLength(p1);
2817  int l2=pLength(p2);
2818  if (l1<l2) return -1;
2819  if (l1>l2) return 1;
2820  return 0;
2821}
2822void slimgb_alg::introduceDelayedPairs(poly* pa,int s){
2823    if (s==0) return;
2824    sorted_pair_node** si_array=(sorted_pair_node**) omalloc(s* sizeof(sorted_pair_node*));
2825
2826    for( unsigned int i = 0; i < s; i++ )
2827    {
2828        sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
2829        si->i=-1;
2830        si->j=-2;
2831        poly p=pa[i];
2832        simplify_poly(p,r);
2833        si->expected_length=pQuality(p,this,pLength(p));
2834        p_Test(p,r);
2835        si->deg=this->pTotaldegree_full(p);
2836        /*if (!rField_is_Zp(r)){
2837          pContent(p);
2838          pCleardenom(p);
2839        }*/
2840
2841        si->lcm_of_lm=p;
2842
2843        //      c->apairs[n-1-i]=si;
2844        si_array[i]=si;
2845
2846  }
2847
2848  qsort(si_array,s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
2849    apairs=spn_merge(apairs,pair_top+1,si_array,s,this);
2850  pair_top+=s;
2851  omfree(si_array);
2852}
2853slimgb_alg::slimgb_alg(ideal I, int syz_comp,BOOLEAN F4,int deg_pos){
2854  this->deg_pos=deg_pos;
2855  lastCleanedDeg=-1;
2856  completed=FALSE;
2857  this->syz_comp=syz_comp;
2858  r=currRing;
2859  nc=rIsPluralRing(r);
2860  this->lastDpBlockStart=get_last_dp_block_start(r);
2861  //Print("last dp Block start: %i\n", this->lastDpBlockStart);
2862  is_homog=TRUE;
2863  {
2864    int hzz;
2865    for(hzz=0;hzz<IDELEMS(I);hzz++){
2866      assume(I->m[hzz]!=NULL);
2867      int d=this->pTotaldegree(I->m[hzz]);
2868      poly t=I->m[hzz]->next;
2869      while(t)
2870      {
2871        if (d!=this->pTotaldegree(t))
2872        {
2873          is_homog=FALSE;
2874          break;
2875        }
2876        t=t->next;
2877      }
2878      if(!(is_homog)) break;
2879    }
2880  }
2881  eliminationProblem=((!(is_homog))&&((pLexOrder)||(I->rank>1)));
2882  tailReductions=((is_homog)||((TEST_OPT_REDTAIL)&&(!(I->rank>1))));
2883  //  Print("is homog:%d",c->is_homog);
2884  void* h;
2885  poly hp;
2886  int i,j;
2887  to_destroy=NULL;
2888  easy_product_crit=0;
2889  extended_product_crit=0;
2890  if (rField_is_Zp(r))
2891    isDifficultField=FALSE;
2892  else
2893    isDifficultField=TRUE;
2894  //not fully correct
2895  //(rChar()==0);
2896  F4_mode=F4;
2897
2898  reduction_steps=0;
2899  last_index=-1;
2900
2901  F=NULL;
2902  F_minus=NULL;
2903
2904  Rcounter=0;
2905
2906  soon_free=NULL;
2907
2908  tmp_lm=pOne();
2909
2910  normal_forms=0;
2911  current_degree=1;
2912
2913  max_pairs=5*IDELEMS(I);
2914
2915  apairs=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*max_pairs);
2916  pair_top=-1;
2917
2918  int n=IDELEMS(I);
2919  array_lengths=n;
2920
2921
2922  i=0;
2923  this->n=0;
2924  T_deg=(int*) omalloc(n*sizeof(int));
2925  if(eliminationProblem)
2926    T_deg_full=(int*) omalloc(n*sizeof(int));
2927  else
2928    T_deg_full=NULL;
2929  tmp_pair_lm=(poly*) omalloc(n*sizeof(poly));
2930  tmp_spn=(sorted_pair_node**) omalloc(n*sizeof(sorted_pair_node*));
2931  lm_bin=omGetSpecBin(POLYSIZE + (r->ExpL_Size)*sizeof(long));
2932#ifdef HEAD_BIN
2933  HeadBin=omGetSpecBin(POLYSIZE + (currRing->ExpL_Size)*sizeof(long));
2934#endif
2935  /* omUnGetSpecBin(&(c->HeadBin)); */
2936  #ifndef HAVE_BOOST
2937  #ifdef USE_STDVECBOOL
2938  #else
2939  h=omalloc(n*sizeof(char*));
2940
2941  states=(char**) h;
2942  #endif
2943  #endif
2944  h=omalloc(n*sizeof(int));
2945  lengths=(int*) h;
2946  weighted_lengths=(wlen_type*)omalloc(n*sizeof(wlen_type));
2947  gcd_of_terms=(poly*) omalloc(n*sizeof(poly));
2948
2949  short_Exps=(long*) omalloc(n*sizeof(long));
2950  if (F4_mode)
2951    S=idInit(n,I->rank);
2952  else
2953    S=idInit(1,I->rank);
2954  strat=new skStrategy;
2955  if (eliminationProblem)
2956    strat->honey=TRUE;
2957  strat->syzComp = 0;
2958  initBuchMoraCrit(strat);
2959  initBuchMoraPos(strat);
2960  strat->initEcart = initEcartBBA;
2961  strat->tailRing=r;
2962  strat->enterS = enterSBba;
2963  strat->sl = -1;
2964  i=n;
2965  i=1;//some strange bug else
2966  /* initS(c->S,NULL,c->strat); */
2967  /* intS start: */
2968  // i=((i+IDELEMS(c->S)+15)/16)*16;
2969  strat->ecartS=(intset)omAlloc(i*sizeof(int)); /*initec(i);*/
2970  strat->sevS=(unsigned long*)omAlloc0(i*sizeof(unsigned long));
2971  /*initsevS(i);*/
2972  strat->S_2_R=(int*)omAlloc0(i*sizeof(int));/*initS_2_R(i);*/
2973  strat->fromQ=NULL;
2974  strat->Shdl=idInit(1,1);
2975  strat->S=strat->Shdl->m;
2976  strat->lenS=(int*)omAlloc0(i*sizeof(int));
2977  if((isDifficultField)||(eliminationProblem))
2978    strat->lenSw=(wlen_type*)omAlloc0(i*sizeof(wlen_type));
2979  else
2980    strat->lenSw=NULL;
2981  sorted_pair_node* si;
2982  assume(n>0);
2983  add_to_basis_ideal_quotient(I->m[0],this,NULL);
2984
2985  assume(strat->sl==IDELEMS(strat->Shdl)-1);
2986  if(!(F4_mode))
2987  {
2988        poly* array_arg=I->m;
2989        array_arg++;
2990        introduceDelayedPairs(array_arg,n-1);
2991        /*
2992    for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
2993    {
2994      //     add_to_basis(I->m[i],-1,-1,c);
2995      si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
2996      si->i=-1;
2997      si->j=-2;
2998      si->expected_length=pQuality(I->m[i],this,pLength(I->m[i]));
2999      si->deg=pTotaldegree(I->m[i]);
3000      if (!rField_is_Zp(r)){
3001        pCleardenom(I->m[i]);
3002      }
3003      si->lcm_of_lm=I->m[i];
3004
3005      //      c->apairs[n-1-i]=si;
3006      apairs[n-i-1]=si;
3007      ++(pair_top);
3008    }*/
3009  }
3010  else
3011  {
3012    for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
3013      add_to_basis_ideal_quotient(I->m[i],this,NULL);
3014  }
3015  for(i=0;i<IDELEMS(I);i++)
3016  {
3017    I->m[i]=NULL;
3018  }
3019  idDelete(&I);
3020  add_later=idInit(ADD_LATER_SIZE,S->rank);
3021  #ifdef USE_NORO
3022  use_noro=((!(nc))&&(S->rank<=1)&&(rField_is_Zp(r)) &&(!(eliminationProblem))&&(npPrimeM<=32003));
3023  use_noro_last_block=false;
3024  if ((!(use_noro))&&(lastDpBlockStart<=pVariables)){
3025      use_noro_last_block=((!(nc))&&(S->rank<=1)&&(rField_is_Zp(r)) &&(npPrimeM<=32003));
3026  }
3027  #else
3028  use_noro=false;
3029  use_noro_last_block=false;
3030  #endif
3031  //Print("NORO last block %i",use_noro_last_block);
3032  memset(add_later->m,0,ADD_LATER_SIZE*sizeof(poly));
3033}
3034slimgb_alg::~slimgb_alg(){
3035
3036  if (!(completed)){
3037      poly* add=(poly*) omalloc((pair_top+2)*sizeof(poly));
3038      int piter;
3039      int pos=0;
3040      for(piter=0;piter<=pair_top;piter++){
3041        sorted_pair_node* s=apairs[piter];
3042        if (s->i<0){
3043            //delayed element
3044            if (s->lcm_of_lm!=NULL){
3045                add[pos]=s->lcm_of_lm;
3046                pos++;
3047
3048            }
3049        }
3050
3051        free_sorted_pair_node(s,r);
3052        apairs[piter]=NULL;
3053      }
3054      pair_top=-1;
3055      add[pos]=NULL;
3056      pos=0;
3057      while(add[pos]!=NULL){
3058        add_to_basis_ideal_quotient(add[pos],this,NULL);
3059        pos++;
3060      }
3061      for(piter=0;piter<=pair_top;piter++){
3062        sorted_pair_node* s=apairs[piter];
3063        assume(s->i>=0);
3064        free_sorted_pair_node(s,r);
3065        apairs[piter]=NULL;
3066      }
3067      pair_top=-1;
3068  }
3069  id_Delete(&add_later,r);
3070  int i,j;
3071  slimgb_alg* c=this;
3072  while(c->to_destroy)
3073  {
3074    pDelete(&(c->to_destroy->p));
3075    poly_list_node* old=c->to_destroy;
3076    c->to_destroy=c->to_destroy->next;
3077    omfree(old);
3078  }
3079  while(c->F)
3080  {
3081    for(i=0;i<c->F->size;i++){
3082      pDelete(&(c->F->mp[i].m));
3083    }
3084    omfree(c->F->mp);
3085    c->F->mp=NULL;
3086    mp_array_list* old=c->F;
3087    c->F=c->F->next;
3088    omfree(old);
3089  }
3090  while(c->F_minus)
3091  {
3092    for(i=0;i<c->F_minus->size;i++){
3093      pDelete(&(c->F_minus->p[i]));
3094    }
3095    omfree(c->F_minus->p);
3096    c->F_minus->p=NULL;
3097    poly_array_list* old=c->F_minus;
3098    c->F_minus=c->F_minus->next;
3099    omfree(old);
3100  }
3101  #ifndef HAVE_BOOST
3102  #ifndef USE_STDVECBOOL
3103  for(int z=1 /* zero length at 0 */;z<c->n;z++){
3104    omfree(c->states[z]);
3105  }
3106  omfree(c->states);
3107  #endif
3108  #endif
3109
3110  omfree(c->lengths);
3111  omfree(c->weighted_lengths);
3112  for(int z=0;z<c->n;z++)
3113  {
3114    pDelete(&c->tmp_pair_lm[z]);
3115    omfree(c->tmp_spn[z]);
3116  }
3117  omfree(c->tmp_pair_lm);
3118  omfree(c->tmp_spn);
3119
3120  omfree(c->T_deg);
3121  if(c->T_deg_full)
3122    omfree(c->T_deg_full);
3123
3124  omFree(c->strat->ecartS);
3125  omFree(c->strat->sevS);
3126//   initsevS(i);
3127  omFree(c->strat->S_2_R);
3128
3129
3130  omFree(c->strat->lenS);
3131
3132  if(c->strat->lenSw)  omFree(c->strat->lenSw);
3133
3134
3135
3136
3137  for(i=0;i<c->n;i++){
3138    if(c->gcd_of_terms[i])
3139      pDelete(&(c->gcd_of_terms[i]));
3140  }
3141  omfree(c->gcd_of_terms);
3142
3143  omfree(c->apairs);
3144  if (TEST_OPT_PROT)
3145  {
3146    //Print("calculated %d NFs\n",c->normal_forms);
3147      Print("\nNF:%i product criterion:%i, ext_product criterion:%i \n", c->normal_forms, c->easy_product_crit, c->extended_product_crit);
3148
3149
3150  }
3151  int deleted_form_c_s=0;
3152
3153  for(i=0;i<=c->strat->sl;i++){
3154    if (!c->strat->S[i]) continue;
3155    BOOLEAN found=FALSE;
3156    for(j=0;j<c->n;j++){
3157      if (c->S->m[j]==c->strat->S[i]){
3158        found=TRUE;
3159        break;
3160      }
3161    }
3162    if(!found) pDelete(&c->strat->S[i]);
3163  }
3164//   for(i=0;i<c->n;i++){
3165//     if (c->rep[i]!=i){
3166// //       for(j=0;j<=c->strat->sl;j++){
3167// //   if(c->strat->S[j]==c->S->m[i]){
3168// //     c->strat->S[j]=NULL;
3169// //     break;
3170// //   }
3171// //       }
3172// //      PrintS("R_delete");
3173//       pDelete(&c->S->m[i]);
3174//     }
3175//   }
3176
3177  if (completed){
3178  for(i=0;i<c->n;i++)
3179  {
3180    assume(c->S->m[i]!=NULL);
3181    if (p_GetComp(c->S->m[i],currRing)>this->syz_comp) continue;
3182    for(j=0;j<c->n;j++)
3183    {
3184      if((c->S->m[j]==NULL)||(i==j))
3185        continue;
3186      assume(p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
3187             c->S->m[i],~c->short_Exps[i],
3188             c->r)==p_LmDivisibleBy(c->S->m[j],
3189             c->S->m[i],
3190             c->r));
3191      if (p_LmShortDivisibleBy(c->S->m[j],c->short_Exps[j],
3192          c->S->m[i],~c->short_Exps[i],
3193          c->r))
3194      {
3195        pDelete(&c->S->m[i]);
3196        break;
3197      }
3198    }
3199  }
3200  }
3201  omfree(c->short_Exps);
3202
3203
3204  ideal I=c->S;
3205
3206  IDELEMS(I)=c->n;
3207
3208  idSkipZeroes(I);
3209  for(i=0;i<=c->strat->sl;i++)
3210    c->strat->S[i]=NULL;
3211  id_Delete(&c->strat->Shdl,c->r);
3212  pDelete(&c->tmp_lm);
3213  omUnGetSpecBin(&lm_bin);
3214  delete c->strat;
3215}
3216ideal t_rep_gb(ring r,ideal arg_I, int syz_comp, BOOLEAN F4_mode){
3217
3218    assume(r==currRing);
3219    ring orig_ring=r;
3220    int pos;
3221    ring new_ring=rAssure_TDeg(orig_ring,1,rVar(orig_ring),pos);
3222
3223
3224    ideal s_h;
3225    if (orig_ring != new_ring)
3226    {
3227        rChangeCurrRing(new_ring);
3228        s_h=idrCopyR_NoSort(arg_I,orig_ring);
3229        idTest(s_h);
3230        /*int i;
3231
3232        for(i=0;i<IDELEMS(s_h);i++){
3233            poly p=s_h->m[i];
3234            while(p){
3235                p_Setm(p,new_ring);
3236                pIter(p);
3237            }
3238        }*/
3239    }
3240    else
3241    {
3242        s_h = id_Copy(arg_I,orig_ring);
3243    }
3244
3245    ideal s_result=do_t_rep_gb(new_ring,s_h,syz_comp,F4_mode,pos);
3246    ideal result;
3247    if(orig_ring != new_ring)
3248    {
3249
3250        idTest(s_result);
3251        rChangeCurrRing(orig_ring);
3252        result = idrMoveR_NoSort(s_result, new_ring);
3253
3254        idTest(result);
3255        //rChangeCurrRing(new_ring);
3256        rKill(new_ring);
3257        //rChangeCurrRing(orig_ring);
3258    }
3259    else
3260        result=s_result;
3261    idTest(result);
3262    return result;
3263}
3264
3265
3266
3267ideal do_t_rep_gb(ring r,ideal arg_I, int syz_comp, BOOLEAN F4_mode,int deg_pos){
3268
3269  //  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)));
3270
3271  if (TEST_OPT_PROT)
3272    if (F4_mode)
3273      PrintS("F4 Modus \n");
3274
3275  //debug_Ideal=arg_debug_Ideal;
3276  //if (debug_Ideal) PrintS("DebugIdeal received\n");
3277  // Print("Idelems %i \n----------\n",IDELEMS(arg_I));
3278  ideal I=arg_I;
3279  idCompactify(I);
3280  if (idIs0(I)) return I;
3281  int i;
3282  for(i=0;i<IDELEMS(I);i++)
3283  {
3284    assume(I->m[i]!=NULL);
3285    simplify_poly(I->m[i],currRing);
3286  }
3287
3288
3289  qsort(I->m,IDELEMS(I),sizeof(poly),poly_crit);
3290  //Print("Idelems %i \n----------\n",IDELEMS(I));
3291  //slimgb_alg* c=(slimgb_alg*) omalloc(sizeof(slimgb_alg));
3292  //int syz_comp=arg_I->rank;
3293  slimgb_alg* c=new slimgb_alg(I, syz_comp,F4_mode,deg_pos);
3294
3295
3296  while ((c->pair_top>=0) && ((!(TEST_OPT_DEGBOUND)) || (c->apairs[c->pair_top]->deg<=Kstd1_deg)))
3297  {
3298    #ifdef HAVE_F4
3299    if(F4_mode)
3300      go_on_F4(c);
3301
3302    else
3303    #endif
3304      go_on(c);
3305  }
3306  if (c->pair_top<0)
3307    c->completed=TRUE;
3308  I=c->S;
3309  delete c;
3310  if (TEST_OPT_REDSB)
3311  {
3312    ideal erg=kInterRed(I,NULL);
3313    assume(I!=erg);
3314    id_Delete(&I, currRing);
3315    return erg;
3316  }
3317  //qsort(I->m, IDELEMS(I),sizeof(poly),pLmCmp_func);
3318  assume(I->rank>=idRankFreeModule(I));
3319  return(I);
3320}
3321void now_t_rep(const int & arg_i, const int & arg_j, slimgb_alg* c){
3322  int i,j;
3323  if (arg_i==arg_j){
3324    return;
3325  }
3326  if (arg_i>arg_j){
3327    i=arg_j;
3328    j=arg_i;
3329  } else {
3330    i=arg_i;
3331    j=arg_j;
3332  }
3333  c->states[j][i]=HASTREP;
3334}
3335
3336static BOOLEAN has_t_rep(const int & arg_i, const  int & arg_j, slimgb_alg* state){
3337  assume(0<=arg_i);
3338  assume(0<=arg_j);
3339  assume(arg_i<state->n);
3340  assume(arg_j<state->n);
3341  if (arg_i==arg_j)
3342  {
3343    return (TRUE);
3344  }
3345  if (arg_i>arg_j)
3346  {
3347    return (state->states[arg_i][arg_j]==HASTREP);
3348  } else
3349  {
3350    return (state->states[arg_j][arg_i]==HASTREP);
3351  }
3352}
3353static int pLcmDeg(poly a, poly b)
3354{
3355  int i;
3356  int n=0;
3357  for (i=pVariables; i; i--)
3358  {
3359    n+=si_max( pGetExp(a,i), pGetExp(b,i));
3360  }
3361  return n;
3362
3363}
3364
3365
3366
3367static void shorten_tails(slimgb_alg* c, poly monom)
3368{
3369  return;
3370// BOOLEAN corr=lenS_correct(c->strat);
3371  for(int i=0;i<c->n;i++)
3372  {
3373    //enter tail
3374
3375    if (c->S->m[i]==NULL) continue;
3376    poly tail=c->S->m[i]->next;
3377    poly prev=c->S->m[i];
3378    BOOLEAN did_something=FALSE;
3379    while((tail!=NULL)&& (pLmCmp(tail, monom)>=0))
3380    {
3381      if (p_LmDivisibleBy(monom,tail,c->r))
3382      {
3383        did_something=TRUE;
3384        prev->next=tail->next;
3385        tail->next=NULL;
3386        p_Delete(& tail,c->r);
3387        tail=prev;
3388        //PrintS("Shortened");
3389        c->lengths[i]--;
3390      }
3391      prev=tail;
3392      tail=tail->next;
3393    }
3394    if (did_something)
3395    {
3396      int new_pos;
3397      wlen_type q;
3398      q=pQuality(c->S->m[i],c,c->lengths[i]);
3399      new_pos=simple_posInS(c->strat,c->S->m[i],c->lengths[i],q);
3400
3401      int old_pos=-1;
3402      //assume new_pos<old_pos
3403      for (int z=0;z<=c->strat->sl;z++)
3404      {
3405        if (c->strat->S[z]==c->S->m[i])
3406        {
3407          old_pos=z;
3408          break;
3409        }
3410      }
3411      if (old_pos== -1)
3412        for (int z=new_pos-1;z>=0;z--)
3413        {
3414          if (c->strat->S[z]==c->S->m[i])
3415          {
3416            old_pos=z;
3417            break;
3418          }
3419        }
3420      assume(old_pos>=0);
3421      assume(new_pos<=old_pos);
3422      assume(pLength(c->strat->S[old_pos])==c->lengths[i]);
3423      c->strat->lenS[old_pos]=c->lengths[i];
3424      if (c->strat->lenSw)
3425        c->strat->lenSw[old_pos]=q;
3426
3427      if (new_pos<old_pos)
3428        move_forward_in_S(old_pos,new_pos,c->strat);
3429
3430      length_one_crit(c,i,c->lengths[i]);
3431    }
3432  }
3433}
3434static sorted_pair_node* pop_pair(slimgb_alg* c){
3435  clean_top_of_pair_list(c);
3436
3437  if(c->pair_top<0) return NULL;
3438  else return (c->apairs[c->pair_top--]);
3439}
3440void slimgb_alg::cleanDegs(int lower, int upper){
3441  assume(is_homog);
3442  int deg;
3443  if (TEST_OPT_PROT){
3444    PrintS("C");
3445  }
3446  for(deg=lower;deg<=upper;deg++){
3447    int i;
3448    for(i=0;i<n;i++){
3449      if (T_deg[i]==deg){
3450          poly h;
3451          h=S->m[i];
3452          h=redNFTail(h,strat->sl,strat,lengths[i]);
3453          if (!rField_is_Zp(r))
3454          {
3455            pCleardenom(h);
3456            pContent(h);
3457
3458          } else pNorm(h);
3459          //TODO:GCD of TERMS
3460          poly got=::gcd_of_terms(h,r);
3461          p_Delete(&gcd_of_terms[i],r);
3462          gcd_of_terms[i]=got;
3463          int len=pLength(h);
3464          wlen_type wlen=pQuality(h,this,len);
3465          if (weighted_lengths)
3466            weighted_lengths[i]=wlen;
3467          lengths[i]=len;
3468          assume(h==S->m[i]);
3469          int j;
3470          for(j=0;j<=strat->sl;j++){
3471            if (h==strat->S[j]){
3472              int new_pos=simple_posInS(strat, h,len, wlen);
3473              if (strat->lenS){
3474                strat->lenS[j]=len;
3475              }
3476              if (strat->lenSw){
3477                strat->lenSw[j]=wlen;
3478              }
3479              if (new_pos<j){
3480                move_forward_in_S(j,new_pos,strat);
3481              }else{
3482                if (new_pos>j)
3483                 new_pos=new_pos-1;//is identical with one element
3484                if (new_pos>j)
3485                  move_backward_in_S(j,new_pos,strat);
3486              }
3487              break;
3488            }
3489          }
3490        }
3491
3492    }
3493
3494  }
3495  {
3496    int i,j;
3497    for(i=0;i<this->n;i++){
3498      for(j=0;j<i;j++){
3499        if (T_deg[i]+T_deg[j]<=upper){
3500          now_t_rep(i,j,this);
3501        }
3502      }
3503    }
3504  }
3505  //TODO resort and update strat->S,strat->lenSw
3506  //TODO mark pairs
3507}
3508sorted_pair_node* top_pair(slimgb_alg* c){
3509  while(c->pair_top>=0){
3510    super_clean_top_of_pair_list(c);//yeah, I know, it's odd that I use a different proc here
3511    if ((c->is_homog)&&(c->pair_top>=0)&&(c->apairs[c->pair_top]->deg>=c->lastCleanedDeg+2)){
3512      int upper=c->apairs[c->pair_top]->deg-1;
3513      c->cleanDegs(c->lastCleanedDeg+1,upper);
3514      c->lastCleanedDeg=upper;
3515    } else{
3516      break;
3517    }
3518
3519  }
3520
3521
3522  if(c->pair_top<0) return NULL;
3523  else return (c->apairs[c->pair_top]);
3524}
3525sorted_pair_node* quick_pop_pair(slimgb_alg* c){
3526  if(c->pair_top<0) return NULL;
3527  else return (c->apairs[c->pair_top--]);
3528}
3529
3530
3531
3532static void super_clean_top_of_pair_list(slimgb_alg* c){
3533  while((c->pair_top>=0)
3534  && (c->apairs[c->pair_top]->i>=0)
3535  && (good_has_t_rep(c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,c)))
3536  {
3537
3538    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
3539    c->pair_top--;
3540
3541  }
3542}
3543void clean_top_of_pair_list(slimgb_alg* c){
3544  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))){
3545
3546    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
3547    c->pair_top--;
3548
3549  }
3550}
3551static BOOLEAN state_is(calc_state state, const int & arg_i, const  int & arg_j, slimgb_alg* c){
3552  assume(0<=arg_i);
3553  assume(0<=arg_j);
3554  assume(arg_i<c->n);
3555  assume(arg_j<c->n);
3556  if (arg_i==arg_j)
3557  {
3558    return (TRUE);
3559  }
3560  if (arg_i>arg_j)
3561  {
3562    return (c->states[arg_i][arg_j]==state);
3563  }
3564  else return(c->states[arg_j][arg_i]==state);
3565}
3566
3567
3568void free_sorted_pair_node(sorted_pair_node* s, ring r){
3569  if (s->i>=0)
3570    p_Delete(&s->lcm_of_lm,r);
3571  omfree(s);
3572}
3573static BOOLEAN pair_better(sorted_pair_node* a,sorted_pair_node* b, slimgb_alg* c){
3574  if (a->deg<b->deg) return TRUE;
3575  if (a->deg>b->deg) return FALSE;
3576
3577
3578  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
3579  if (comp==1) return FALSE;
3580  if (-1==comp) return TRUE;
3581  if (a->expected_length<b->expected_length) return TRUE;
3582  if (a->expected_length>b->expected_length) return FALSE;
3583  if (a->i+a->j<b->i+b->j) return TRUE;
3584   if (a->i+a->j>b->i+b->j) return FALSE;
3585  if (a->i<b->i) return TRUE;
3586  if (a->i>b->i) return FALSE;
3587  return TRUE;
3588}
3589
3590static int tgb_pair_better_gen(const void* ap,const void* bp){
3591
3592  sorted_pair_node* a=*((sorted_pair_node**)ap);
3593  sorted_pair_node* b=*((sorted_pair_node**)bp);
3594  assume((a->i>a->j) || (a->i < 0));
3595  assume((b->i>b->j) || (b->i < 0));
3596  if (a->deg<b->deg) return -1;
3597  if (a->deg>b->deg) return 1;
3598
3599
3600
3601 int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
3602
3603  if (comp==1) return 1;
3604  if (-1==comp) return -1;
3605   if (a->expected_length<b->expected_length) return -1;
3606  if (a->expected_length>b->expected_length) return 1;
3607  if (a->i+a->j<b->i+b->j) return -1;
3608   if (a->i+a->j>b->i+b->j) return 1;
3609  if (a->i<b->i) return -1;
3610   if (a->i>b->i) return 1;
3611  return 0;
3612}
3613
3614
3615static poly gcd_of_terms(poly p, ring r){
3616  int max_g_0=0;
3617  assume(p!=NULL);
3618  int i;
3619  poly m=pOne();
3620  poly t;
3621  for (i=pVariables; i; i--)
3622  {
3623      pSetExp(m,i, pGetExp(p,i));
3624      if (max_g_0==0)
3625  if (pGetExp(m,i)>0)
3626    max_g_0=i;
3627  }
3628
3629  t=p->next;
3630  while (t!=NULL){
3631
3632    if (max_g_0==0) break;
3633    for (i=max_g_0; i; i--)
3634    {
3635      pSetExp(m,i, si_min(pGetExp(t,i),pGetExp(m,i)));
3636      if (max_g_0==i)
3637  if (pGetExp(m,i)==0)
3638    max_g_0=0;
3639      if ((max_g_0==0) && (pGetExp(m,i)>0)){
3640  max_g_0=i;
3641      }
3642    }
3643    t=t->next;
3644  }
3645  p_Setm(m,r);
3646  if (max_g_0>0)
3647    return m;
3648  pDelete(&m);
3649  return NULL;
3650}
3651static inline BOOLEAN pHasNotCFExtended(poly p1, poly p2, poly m)
3652{
3653
3654  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
3655    return FALSE;
3656  int i = 1;
3657  loop
3658  {
3659    if ((pGetExp(p1, i)-pGetExp(m,i) >0) && (pGetExp(p2, i) -pGetExp(m,i)> 0))   return FALSE;
3660    if (i == pVariables)                                return TRUE;
3661    i++;
3662  }
3663}
3664
3665
3666//for impl reasons may return false if the the normal product criterion matches
3667static inline BOOLEAN extended_product_criterion(poly p1, poly gcd1, poly p2, poly gcd2, slimgb_alg* c){
3668  if (c->nc)
3669    return FALSE;
3670  if(gcd1==NULL) return FALSE;
3671        if(gcd2==NULL) return FALSE;
3672        gcd1->next=gcd2; //may ordered incorrect
3673        poly m=gcd_of_terms(gcd1,c->r);
3674        gcd1->next=NULL;
3675        if (m==NULL) return FALSE;
3676
3677        BOOLEAN erg=pHasNotCFExtended(p1,p2,m);
3678        pDelete(&m);
3679        return erg;
3680}
3681static poly kBucketGcd(kBucket* b, ring r)
3682{
3683  int s=0;
3684  int i;
3685  poly m, n;
3686  BOOLEAN initialized=FALSE;
3687  for (i=MAX_BUCKET-1;i>=0;i--)
3688  {
3689    if (b->buckets[i]!=NULL){
3690      if (!initialized){
3691  m=gcd_of_terms(b->buckets[i],r);
3692  initialized=TRUE;
3693  if (m==NULL) return NULL;
3694      }
3695      else
3696  {
3697    n=gcd_of_terms(b->buckets[i],r);
3698    if (n==NULL) {
3699      pDelete(&m);
3700      return NULL;
3701    }
3702    n->next=m;
3703    poly t=gcd_of_terms(n,r);
3704    n->next=NULL;
3705    pDelete(&m);
3706    pDelete(&n);
3707    m=t;
3708    if (m==NULL) return NULL;
3709
3710  }
3711    }
3712  }
3713  return m;
3714}
3715
3716
3717
3718
3719static inline wlen_type quality_of_pos_in_strat_S(int pos, slimgb_alg* c){
3720  if (c->strat->lenSw!=NULL) return c->strat->lenSw[pos];
3721  return c->strat->lenS[pos];
3722}
3723#ifdef HAVE_PLURAL
3724static inline wlen_type quality_of_pos_in_strat_S_mult_high(int pos, poly high, slimgb_alg* c)
3725  //meant only for nc
3726{
3727  poly m=pOne();
3728  pExpVectorDiff(m,high ,c->strat->S[pos]);
3729  poly product = nc_mm_Mult_pp(m, c->strat->S[pos], c->r);
3730  wlen_type erg=pQuality(product,c);
3731  pDelete(&m);
3732  pDelete(&product);
3733  return erg;
3734}
3735#endif
3736
3737static void multi_reduction_lls_trick(red_object* los, int losl,slimgb_alg* c,find_erg & erg){
3738  erg.expand=NULL;
3739  BOOLEAN swap_roles; //from reduce_by, to_reduce_u if fromS
3740  if(erg.fromS){
3741    if(pLmEqual(c->strat->S[erg.reduce_by],los[erg.to_reduce_u].p))
3742    {
3743      int i;
3744      wlen_type quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
3745      int best=erg.to_reduce_u+1;
3746/*
3747      for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
3748  int qc=los[i].guess_quality(c);
3749  if (qc<quality_a){
3750    best=i;
3751    quality_a=qc;
3752  }
3753      }
3754      if(best!=erg.to_reduce_u+1){*/
3755      wlen_type qc;
3756      best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3757      if(qc<quality_a){
3758  los[best].flatten();
3759  int b_pos=kBucketCanonicalize(los[best].bucket);
3760  los[best].p=los[best].bucket->buckets[b_pos];
3761  qc=pQuality(los[best].bucket->buckets[b_pos],c);
3762  if(qc<quality_a){
3763    red_object h=los[erg.to_reduce_u];
3764    los[erg.to_reduce_u]=los[best];
3765    los[best]=h;
3766    swap_roles=TRUE;
3767  }
3768  else
3769    swap_roles=FALSE;
3770      }
3771      else{
3772
3773  swap_roles=FALSE;
3774      }
3775
3776    }
3777      else
3778    {
3779      if (erg.to_reduce_u>erg.to_reduce_l){
3780
3781  int i;
3782  wlen_type quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
3783  #ifdef HAVE_PLURAL
3784  if ((c->nc) && (!(rIsSCA(c->r))))
3785    quality_a=quality_of_pos_in_strat_S_mult_high(erg.reduce_by, los[erg.to_reduce_u].p, c);
3786  #endif
3787  int best=erg.to_reduce_u+1;
3788  wlen_type qc;
3789  best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3790  assume(qc==los[best].guess_quality(c));
3791  if(qc<quality_a){
3792    los[best].flatten();
3793    int b_pos=kBucketCanonicalize(los[best].bucket);
3794    los[best].p=los[best].bucket->buckets[b_pos];
3795    qc==pQuality(los[best].bucket->buckets[b_pos],c);
3796    //(best!=erg.to_reduce_u+1)
3797    if(qc<quality_a){
3798    red_object h=los[erg.to_reduce_u];
3799    los[erg.to_reduce_u]=los[best];
3800    los[best]=h;
3801    erg.reduce_by=erg.to_reduce_u;
3802    erg.fromS=FALSE;
3803    erg.to_reduce_u--;
3804    }
3805  }
3806      }
3807      else
3808      {
3809  assume(erg.to_reduce_u==erg.to_reduce_l);
3810  wlen_type quality_a=
3811        quality_of_pos_in_strat_S(erg.reduce_by,c);
3812  wlen_type qc=los[erg.to_reduce_u].guess_quality(c);
3813  if (qc<0) PrintS("Wrong wlen_type");
3814  if(qc<quality_a){
3815    int best=erg.to_reduce_u;
3816    los[best].flatten();
3817    int b_pos=kBucketCanonicalize(los[best].bucket);
3818    los[best].p=los[best].bucket->buckets[b_pos];
3819    qc=pQuality(los[best].bucket->buckets[b_pos],c);
3820    assume(qc>=0);
3821    if(qc<quality_a){
3822      BOOLEAN exp=FALSE;
3823      if(qc<=2){
3824         //Print("\n qc is %lld \n",qc);
3825         exp=TRUE;
3826      }
3827
3828      else {
3829         if (qc<quality_a/2)
3830          exp=TRUE;
3831         else
3832       if(erg.reduce_by<c->n/4)
3833          exp=TRUE;
3834      }
3835      if (exp){
3836        poly clear_into;
3837        los[erg.to_reduce_u].flatten();
3838        kBucketClear(los[erg.to_reduce_u].bucket,&clear_into,&erg.expand_length);
3839        erg.expand=pCopy(clear_into);
3840        kBucketInit(los[erg.to_reduce_u].bucket,clear_into,erg.expand_length);
3841        if (TEST_OPT_PROT)
3842    PrintS("e");
3843
3844      }
3845    }
3846  }
3847
3848
3849      }
3850
3851      swap_roles=FALSE;
3852      return;
3853      }
3854
3855  }
3856  else{
3857    if(erg.reduce_by>erg.to_reduce_u){
3858      //then lm(rb)>= lm(tru) so =
3859      assume(erg.reduce_by==erg.to_reduce_u+1);
3860      int best=erg.reduce_by;
3861      wlen_type quality_a=los[erg.reduce_by].guess_quality(c);
3862      wlen_type qc;
3863      best=find_best(los,erg.to_reduce_l,erg.to_reduce_u,qc,c);
3864
3865      int i;
3866      if(qc<quality_a){
3867    red_object h=los[erg.reduce_by];
3868    los[erg.reduce_by]=los[best];
3869    los[best]=h;
3870  }
3871  swap_roles=FALSE;
3872  return;
3873
3874
3875    }
3876    else
3877    {
3878      assume(!pLmEqual(los[erg.reduce_by].p,los[erg.to_reduce_l].p));
3879      assume(erg.to_reduce_u==erg.to_reduce_l);
3880      //further assume, that reduce_by is the above all other polys
3881      //with same leading term
3882      int il=erg.reduce_by;
3883      wlen_type quality_a =los[erg.reduce_by].guess_quality(c);
3884      wlen_type qc;
3885      while((il>0) && pLmEqual(los[il-1].p,los[il].p)){
3886  il--;
3887  qc=los[il].guess_quality(c);
3888  if (qc<quality_a){
3889    quality_a=qc;
3890    erg.reduce_by=il;
3891  }
3892      }
3893      swap_roles=FALSE;
3894    }
3895
3896  }
3897  if(swap_roles)
3898  {
3899    if (TEST_OPT_PROT)
3900      PrintS("b");
3901    poly clear_into;
3902    int dummy_len;
3903    int new_length;
3904    int bp=erg.to_reduce_u;//bucket_positon
3905    //kBucketClear(los[bp].bucket,&clear_into,&new_length);
3906    new_length=los[bp].clear_to_poly();
3907    clear_into=los[bp].p;
3908    poly p=c->strat->S[erg.reduce_by];
3909    int j=erg.reduce_by;
3910    int old_length=c->strat->lenS[j];// in view of S
3911    los[bp].p=p;
3912    if (c->eliminationProblem){
3913        los[bp].sugar=c->pTotaldegree_full(p);
3914    }
3915    kBucketInit(los[bp].bucket,p,old_length);
3916    wlen_type qal=pQuality(clear_into,c,new_length);
3917    int pos_in_c=-1;
3918    int z;
3919    int new_pos;
3920    new_pos=simple_posInS(c->strat,clear_into,new_length, qal);
3921    assume(new_pos<=j);
3922    for (z=c->n;z;z--)
3923    {
3924      if(p==c->S->m[z-1])
3925      {
3926  pos_in_c=z-1;
3927  break;
3928      }
3929    }
3930
3931    int tdeg_full=-1;
3932    int tdeg=-1;
3933    if(pos_in_c>=0)
3934    {
3935      #ifdef TGB_RESORT_PAIRS
3936      c->used_b=TRUE;
3937      c->replaced[pos_in_c]=TRUE;
3938      #endif
3939      tdeg=c->T_deg[pos_in_c];
3940      c->S->m[pos_in_c]=clear_into;
3941      c->lengths[pos_in_c]=new_length;
3942      c->weighted_lengths[pos_in_c]=qal;
3943      if (c->gcd_of_terms[pos_in_c]==NULL)
3944        c->gcd_of_terms[pos_in_c]=gcd_of_terms(clear_into,c->r);
3945      if (c->T_deg_full)
3946        tdeg_full=c->T_deg_full[pos_in_c]=c->pTotaldegree_full(clear_into);
3947      else tdeg_full=tdeg;
3948      c_S_element_changed_hook(pos_in_c,c);
3949    } else {
3950      if (c->eliminationProblem){
3951        tdeg_full=c->pTotaldegree_full(clear_into);
3952        tdeg=c->pTotaldegree(clear_into);
3953      }
3954    }
3955    c->strat->S[j]=clear_into;
3956    c->strat->lenS[j]=new_length;
3957
3958    assume(pLength(clear_into)==new_length);
3959    if(c->strat->lenSw!=NULL)
3960      c->strat->lenSw[j]=qal;
3961    if (!rField_is_Zp(c->r))
3962    {
3963      pCleardenom(clear_into);//should be unnecessary
3964      pContent(clear_into);
3965
3966    }
3967    else
3968      pNorm(clear_into);
3969#ifdef FIND_DETERMINISTIC
3970    erg.reduce_by=j;
3971    //resort later see diploma thesis, find_in_S must be deterministic
3972    //during multireduction if spolys are only in the span of the
3973    //input polys
3974#else
3975
3976    if (new_pos<j)
3977    {
3978      if (c->strat->honey) c->strat->ecartS[j]=tdeg_full-tdeg;
3979      move_forward_in_S(j,new_pos,c->strat);
3980      erg.reduce_by=new_pos;
3981    }
3982#endif
3983  }
3984}
3985static int fwbw(red_object* los, int i){
3986   int i2=i;
3987   int step=1;
3988
3989   BOOLEAN bw=FALSE;
3990   BOOLEAN incr=TRUE;
3991
3992   while(1)
3993   {
3994     if(!bw)
3995     {
3996       step=si_min(i2,step);
3997       if (step==0) break;
3998       i2-=step;
3999
4000       if(!pLmEqual(los[i].p,los[i2].p))
4001       {
4002   bw=TRUE;
4003   incr=FALSE;
4004       }
4005       else
4006       {
4007   if ((!incr) &&(step==1)) break;
4008       }
4009
4010
4011     }
4012     else
4013     {
4014
4015       step=si_min(i-i2,step);
4016       if (step==0) break;
4017       i2+=step;
4018       if(pLmEqual(los[i].p,los[i2].p)){
4019   if(step==1) break;
4020   else
4021   {
4022     bw=FALSE;
4023   }
4024       }
4025
4026     }
4027     if (incr)
4028       step*=2;
4029     else
4030     {
4031       if (step%2==1)
4032   step=(step+1)/2;
4033       else
4034   step/=2;
4035
4036     }
4037   }
4038   return i2;
4039}
4040static void canonicalize_region(red_object* los, int l, int u,slimgb_alg* c){
4041    assume(l<=u+1);
4042    int i;
4043    for(i=l;i<=u;i++){
4044        kBucketCanonicalize(los[i].bucket);
4045    }
4046
4047}
4048static void multi_reduction_find(red_object* los, int losl,slimgb_alg* c,int startf,find_erg & erg){
4049  kStrategy strat=c->strat;
4050
4051  assume(startf<=losl);
4052  assume((startf==losl-1)||(pLmCmp(los[startf].p,los[startf+1].p)==-1));
4053  int i=startf;
4054
4055  int j;
4056  while(i>=0){
4057    assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)<=0));
4058    assume(is_valid_ro(los[i]));
4059    assume((!(c->eliminationProblem))||(los[i].sugar>=c->pTotaldegree(los[i].p)));
4060    j=kFindDivisibleByInS_easy(strat,los[i]);
4061    if(j>=0){
4062
4063      erg.to_reduce_u=i;
4064      erg.reduce_by=j;
4065      erg.fromS=TRUE;
4066      int i2=fwbw(los,i);
4067      assume(pLmEqual(los[i].p,los[i2].p));
4068      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
4069      assume(i>=i2);
4070
4071
4072      erg.to_reduce_l=i2;
4073      assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
4074      canonicalize_region(los,erg.to_reduce_u+1,startf,c);
4075      return;
4076    }
4077    if (j<0){
4078
4079      //not reduceable, try to use this for reducing higher terms
4080      int i2=fwbw(los,i);
4081      assume(pLmEqual(los[i].p,los[i2].p));
4082      assume((i2==0)||(!pLmEqual(los[i2].p,los[i2-1].p)));
4083      assume(i>=i2);
4084      if(i2!=i){
4085
4086
4087  erg.to_reduce_u=i-1;
4088  erg.to_reduce_l=i2;
4089  erg.reduce_by=i;
4090  erg.fromS=FALSE;
4091  assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
4092  canonicalize_region(los,erg.to_reduce_u+1,startf,c);
4093  return;
4094      }
4095
4096      i--;
4097    }
4098  }
4099  erg.reduce_by=-1;//error code
4100  return;
4101}
4102
4103 //  nicht reduzierbare eintraege in ergebnisliste schreiben
4104//   nullen loeschen
4105//   while(finde_groessten leitterm reduzierbar(c,erg)){
4106
4107static int multi_reduction_clear_zeroes(red_object* los, int  losl, int l, int u)
4108{
4109
4110
4111  int deleted=0;
4112  int  i=l;
4113  int last=-1;
4114  while(i<=u)
4115  {
4116
4117    if(los[i].p==NULL){
4118      kBucketDestroy(&los[i].bucket);
4119//      delete los[i];//here we assume los are constructed with new
4120      //destroy resources, must be added here
4121     if (last>=0)
4122     {
4123       memmove(los+(int)(last+1-deleted),los+(last+1),sizeof(red_object)*(i-1-last));
4124     }
4125     last=i;
4126     deleted++;
4127    }
4128    i++;
4129  }
4130  if((last>=0)&&(last!=losl-1))
4131      memmove(los+(int)(last+1-deleted),los+last+1,sizeof(red_object)*(losl-1-last));
4132  return deleted;
4133
4134}
4135int search_red_object_pos(red_object* a, int top, red_object* key ){
4136
4137    int an = 0;
4138    int en= top;
4139    if (top==-1) return 0;
4140    if (pLmCmp(key->p,a[top].p)==1)
4141      return top+1;
4142    int i;
4143    loop
4144    {
4145      if (an >= en-1)
4146      {
4147        if (pLmCmp(key->p,a[an].p)==-1)
4148           return an;
4149        return en;
4150      }
4151      i=(an+en) / 2;
4152      if (pLmCmp(key->p,a[i].p)==-1)
4153        en=i;
4154      else
4155        an=i;
4156    }
4157
4158}
4159static void sort_region_down(red_object* los, int l, int u, slimgb_alg* c)
4160{
4161  int r_size=u-l+1;
4162  qsort(los+l,r_size,sizeof(red_object),red_object_better_gen);
4163  int i;
4164  int * new_indices=(int*) omalloc((r_size)*sizeof(int));
4165  int bound=0;
4166  BOOLEAN at_end=FALSE;
4167  for(i=l;i<=u;i++){
4168    if (!(at_end)){
4169      bound=new_indices[i-l]=bound+search_red_object_pos(los+bound,l-bound-1,los+i);
4170      if (bound==l) at_end=TRUE;
4171    }
4172    else{
4173      new_indices[i-l]=l;
4174    }
4175  }
4176  red_object* los_region=(red_object*) omalloc(sizeof(red_object)*(u-l+1));
4177  for (int i=0;i<r_size;i++){
4178    new_indices[i]+=i;
4179    los_region[i]=los[l+i];
4180    assume((i==0)||(new_indices[i]>new_indices[i-1]));
4181
4182  }
4183
4184  i=r_size-1;
4185  int j=u;
4186  int j2=l-1;
4187  while(i>=0){
4188    if (new_indices[i]==j){
4189      los[j]=los_region[i];
4190      i--;
4191      j--;
4192    } else{
4193      assume(new_indices[i]<j);
4194      los[j]=los[j2];
4195      assume(j2>=0);
4196      j2--;
4197      j--;
4198    }
4199  }
4200  omfree(los_region);
4201
4202  omfree(new_indices);
4203
4204}
4205
4206//assume that los is ordered ascending by leading term, all non zero
4207static void multi_reduction(red_object* los, int & losl, slimgb_alg* c)
4208{
4209  poly* delay=(poly*) omalloc(losl*sizeof(poly));
4210  int delay_s=0;
4211  //initialize;
4212  assume(c->strat->sl>=0);
4213  assume(losl>0);
4214  int i;
4215  wlen_type max_initial_quality=0;
4216
4217  for(i=0;i<losl;i++){
4218    los[i].sev=pGetShortExpVector(los[i].p);
4219//SetShortExpVector();
4220    los[i].p=kBucketGetLm(los[i].bucket);
4221    if (los[i].initial_quality>max_initial_quality)
4222        max_initial_quality=los[i].initial_quality;
4223    // else
4224//         Print("init2_qal=%lld;", los[i].initial_quality);
4225//     Print("initial_quality=%lld;",max_initial_quality);
4226  }
4227
4228  kStrategy strat=c->strat;
4229  int curr_pos=losl-1;
4230
4231
4232//  nicht reduzierbare eintrï¿œe in ergebnisliste schreiben
4233  // nullen loeschen
4234  while(curr_pos>=0){
4235    if ((c->use_noro_last_block)&&(lies_in_last_dp_block(los[curr_pos].p,c))){
4236        PrintS("L");
4237        int pn_noro=curr_pos+1;
4238        poly* p_noro=(poly*) omalloc(pn_noro*sizeof(poly));
4239        for(i=0;i<pn_noro;i++){
4240            int dummy_len;
4241            poly p;
4242            los[i].p=NULL;
4243            kBucketClear(los[i].bucket,&p,&dummy_len);
4244            p_noro[i]=p;
4245        }
4246
4247
4248        if (npPrimeM<255){
4249          noro_step<tgb_uint8>(p_noro,pn_noro,c);
4250        } else {
4251          if (npPrimeM<65000){
4252            noro_step<tgb_uint16>(p_noro,pn_noro,c);
4253          } else{
4254            noro_step<tgb_uint32>(p_noro,pn_noro,c);
4255          }
4256        }
4257        for(i=0;i<pn_noro;i++){
4258            los[i].p=p_noro[i];
4259            los[i].sev=pGetShortExpVector(los[i].p);
4260            //ignore quality
4261            kBucketInit(los[i].bucket,los[i].p,pLength(los[i].p));
4262
4263        }
4264        qsort(los,pn_noro,sizeof(red_object),red_object_better_gen);
4265        int deleted=multi_reduction_clear_zeroes(los, losl, pn_noro, curr_pos);
4266        losl -= deleted;
4267        curr_pos -= deleted;
4268        break;
4269    }
4270    find_erg erg;
4271
4272    multi_reduction_find(los, losl,c,curr_pos,erg);//last argument should be curr_pos
4273    if(erg.reduce_by<0) break;
4274
4275
4276
4277    erg.expand=NULL;
4278    int d=erg.to_reduce_u-erg.to_reduce_l+1;
4279
4280
4281    multi_reduction_lls_trick(los,losl,c,erg);
4282
4283
4284    int i;
4285    int len;
4286    //    wrp(los[erg.to_reduce_u].p);
4287    //PrintLn();
4288    multi_reduce_step(erg,los,c);
4289
4290
4291    if(!K_TEST_OPT_REDTHROUGH){
4292  for(i=erg.to_reduce_l;i<=erg.to_reduce_u;i++){
4293     if  (los[i].p!=NULL)  //the check (los[i].p!=NULL) might be invalid
4294     {
4295         //
4296         assume(los[i].initial_quality>0);
4297
4298               if(los[i].guess_quality(c)
4299                  >1.5*delay_factor*max_initial_quality){
4300                       if (TEST_OPT_PROT)
4301                           PrintS("v");
4302                       los[i].canonicalize();
4303                       if(los[i].guess_quality(c)
4304                           >delay_factor*max_initial_quality){
4305                               if (TEST_OPT_PROT)
4306                                   PrintS(".");
4307                               los[i].clear_to_poly();
4308                               //delay.push_back(los[i].p);
4309                               delay[delay_s]=los[i].p;
4310                               delay_s++;
4311
4312                               los[i].p=NULL;
4313
4314                      }
4315                  }
4316
4317            }
4318     }
4319  }
4320    int deleted=multi_reduction_clear_zeroes(los, losl, erg.to_reduce_l, erg.to_reduce_u);
4321    if(erg.fromS==FALSE)
4322      curr_pos=si_max(erg.to_reduce_u,erg.reduce_by);
4323    else
4324      curr_pos=erg.to_reduce_u;
4325    losl -= deleted;
4326    curr_pos -= deleted;
4327
4328    //Print("deleted %i \n",deleted);
4329    if ((TEST_V_UPTORADICAL) &&(!(erg.fromS)))
4330        sort_region_down(los,si_min(erg.to_reduce_l,erg.reduce_by),(si_max(erg.to_reduce_u,erg.reduce_by))-deleted,c);
4331    else
4332    sort_region_down(los, erg.to_reduce_l, erg.to_reduce_u-deleted, c);
4333
4334
4335    if(erg.expand)
4336    {
4337#ifdef FIND_DETERMINISTIC
4338      int i;
4339      for(i=0;c->expandS[i];i++);
4340      c->expandS=(poly*) omrealloc(c->expandS,(i+2)*sizeof(poly));
4341      c->expandS[i]=erg.expand;
4342      c->expandS[i+1]=NULL;
4343#else
4344      int ecart=0;
4345      if (c->eliminationProblem){
4346        ecart=c->pTotaldegree_full(erg.expand)-c->pTotaldegree(erg.expand);
4347      }
4348      add_to_reductors(c,erg.expand,erg.expand_length,ecart);
4349#endif
4350    }
4351
4352  }
4353
4354
4355  //sorted_pair_node** pairs=(sorted_pair_node**)
4356  //  omalloc(delay_s*sizeof(sorted_pair_node*));
4357  c->introduceDelayedPairs(delay,delay_s);
4358  /*
4359  for(i=0;i<delay_s;i++){
4360
4361      poly p=delay[i];
4362      //if (rPar(c->r)==0)
4363      simplify_poly(p,c->r);
4364      sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
4365      si->i=-1;
4366      si->j=-1;
4367       if (!rField_is_Zp(c->r)){
4368        if (!c->nc)
4369            p=redTailShort(p, c->strat);
4370        pCleardenom(p);
4371        pContent(p);
4372      }
4373      si->expected_length=pQuality(p,c,pLength(p));
4374      si->deg=pTotaldegree(p);
4375
4376      si->lcm_of_lm=p;
4377      pairs[i]=si;
4378  }
4379  qsort(pairs,delay_s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
4380  c->apairs=spn_merge(c->apairs,c->pair_top+1,pairs,delay_s,c);
4381  c->pair_top+=delay_s;*/
4382  omfree(delay);
4383  //omfree(pairs);
4384  return;
4385}
4386void red_object::flatten(){
4387  assume(p==kBucketGetLm(bucket));
4388}
4389void red_object::validate(){
4390  p=kBucketGetLm(bucket);
4391  if(p)
4392    sev=pGetShortExpVector(p);
4393}
4394int red_object::clear_to_poly(){
4395  flatten();
4396  int l;
4397  kBucketClear(bucket,&p,&l);
4398  return l;
4399}
4400
4401
4402
4403
4404
4405void reduction_step::reduce(red_object* r, int l, int u){}
4406void simple_reducer::do_reduce(red_object & ro)
4407{
4408  number coef;
4409#ifdef HAVE_PLURAL
4410  if (c->nc)
4411    nc_BucketPolyRed_Z(ro.bucket, p, &coef);
4412  else
4413#endif
4414    coef=kBucketPolyRed(ro.bucket,p,
4415       p_len,
4416       c->strat->kNoether);
4417  nDelete(&coef);
4418}
4419
4420
4421void simple_reducer::reduce(red_object* r, int l, int u){
4422  this->pre_reduce(r,l,u);
4423  int i;
4424//debug start
4425  int im;
4426
4427
4428  if(c->eliminationProblem){
4429    assume(p_LmEqual(r[l].p,r[u].p,c->r));
4430    /*int lm_deg=pTotaldegree(r[l].p);
4431    reducer_deg=lm_deg+pTotaldegree_full(p)-pTotaldegree(p);*/
4432  }
4433
4434  for(i=l;i<=u;i++){
4435
4436
4437
4438    this->do_reduce(r[i]);
4439    if (c->eliminationProblem){
4440        r[i].sugar=si_max(r[i].sugar,reducer_deg);
4441    }
4442  }
4443  for(i=l;i<=u;i++){
4444
4445    kBucketSimpleContent(r[i].bucket);
4446    r[i].validate();
4447    #ifdef TGB_DEBUG
4448    #endif
4449  }
4450}
4451reduction_step::~reduction_step(){}
4452simple_reducer::~simple_reducer(){
4453  if(fill_back!=NULL)
4454  {
4455    kBucketInit(fill_back,p,p_len);
4456  }
4457  fill_back=NULL;
4458
4459}
4460
4461void multi_reduce_step(find_erg & erg, red_object* r, slimgb_alg* c){
4462  static int id=0;
4463  id++;
4464  unsigned long sev;
4465    BOOLEAN lt_changed=FALSE;
4466  int rn=erg.reduce_by;
4467  poly red;
4468  int red_len;
4469  simple_reducer* pointer;
4470  BOOLEAN work_on_copy=FALSE;
4471  if(erg.fromS){
4472    red=c->strat->S[rn];
4473    red_len=c->strat->lenS[rn];
4474    assume(red_len==pLength(red));
4475  }
4476  else
4477  {
4478    r[rn].flatten();
4479    kBucketClear(r[rn].bucket,&red,&red_len);
4480
4481    if (!rField_is_Zp(c->r))
4482    {
4483      pCleardenom(red);//should be unnecessary
4484      pContent(red);
4485
4486
4487    }
4488    pNormalize(red);
4489    if (c->eliminationProblem){
4490        r[rn].sugar=c->pTotaldegree_full(red);
4491    }
4492
4493    if ((!(erg.fromS))&&(TEST_V_UPTORADICAL)){
4494
4495         if (polynomial_root(red,c->r))
4496            lt_changed=TRUE;
4497            sev=p_GetShortExpVector(red,c->r);}
4498    red_len=pLength(red);
4499  }
4500  if (((TEST_V_MODPSOLVSB)&&(red_len>1))||((c->nc)||(erg.to_reduce_u-erg.to_reduce_l>5))){
4501    work_on_copy=TRUE;
4502    // poly m=pOne();
4503    poly m=c->tmp_lm;
4504    pSetCoeff(m,nInit(1));
4505    pSetComp(m,0);
4506    for(int i=1;i<=pVariables;i++)
4507      pSetExp(m,i,(pGetExp(r[erg.to_reduce_l].p, i)-pGetExp(red,i)));
4508    pSetm(m);
4509    poly red_cp;
4510    #ifdef HAVE_PLURAL
4511    if (c->nc)
4512      red_cp = nc_mm_Mult_pp(m, red, c->r);
4513    else
4514    #endif
4515      red_cp=ppMult_mm(red,m);
4516    if(!erg.fromS){
4517      kBucketInit(r[rn].bucket,red,red_len);
4518    }
4519    //now reduce the copy
4520    //static poly redNF2 (poly h,slimgb_alg* c , int &len, number&  m,int n)
4521
4522    if (!c->nc)
4523      redTailShort(red_cp,c->strat);
4524    //number mul;
4525    // red_len--;
4526//     red_cp->next=redNF2(red_cp->next,c,red_len,mul,c->average_length);
4527//     pSetCoeff(red_cp,nMult(red_cp->coef,mul));
4528//     nDelete(&mul);
4529//     red_len++;
4530    red=red_cp;
4531    red_len=pLength(red);
4532    // pDelete(&m);
4533
4534  }
4535  int i;
4536
4537
4538
4539  assume(red_len==pLength(red));
4540
4541  int reducer_deg=0;
4542  if (c->eliminationProblem){
4543     int lm_deg=c->pTotaldegree(r[erg.to_reduce_l].p);
4544     int ecart;
4545     if (erg.fromS){
4546       ecart=c->strat->ecartS[erg.reduce_by];
4547     } else {
4548       ecart=c->pTotaldegree_full(red)-lm_deg;
4549     }
4550     reducer_deg=lm_deg+ecart;
4551  }
4552  pointer=new simple_reducer(red,red_len,reducer_deg,c);
4553
4554  if ((!work_on_copy) && (!erg.fromS))
4555    pointer->fill_back=r[rn].bucket;
4556  else
4557    pointer->fill_back=NULL;
4558  pointer->reduction_id=id;
4559  pointer->c=c;
4560
4561  pointer->reduce(r,erg.to_reduce_l, erg.to_reduce_u);
4562  if(work_on_copy) pDelete(&pointer->p);
4563  delete pointer;
4564  if (lt_changed){
4565    assume(!erg.fromS);
4566    r[erg.reduce_by].sev=sev;
4567  }
4568
4569};
4570
4571
4572
4573
4574void simple_reducer:: pre_reduce(red_object* r, int l, int u){}
4575
Note: See TracBrowser for help on using the repository browser.