source: git/Singular/tgb.cc @ 1e89ba

fieker-DuValspielwiese
Last change on this file since 1e89ba was 1e89ba, checked in by Michael Brickenstein <bricken@…>, 21 years ago
*bricken: PAR_N constant git-svn-id: file:///usr/local/Singular/svn/trunk@6582 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 56.3 KB
Line 
1
2// #define OM_CHECK 3
3// #define OM_TRACK 5
4// #define OM_KEEP  1
5// TODO:
6//       deg -> poly_crit
7//       multiple rings
8//       shorten_tails und dessen Aufrufe pruefen wlength!!!
9//       calculating with formal sums
10//       try to create spolys as formal sums
11
12#include "tgb.h"
13#define OM_KEEP 0
14#define LEN_VAR1
15
16#ifdef LEN_VAR1
17// erste Variante: Laenge: Anzahl der Monome
18int pSLength(poly p, int l) {
19  return l; }
20int kSBucketLength(kBucket* bucket) {return bucket_guess(bucket);}
21#endif
22
23#ifdef LEN_VAR2
24// 2. Variante: Laenge: Platz fuer die Koeff.
25int pSLength(poly p,int l)
26{
27  int s=0;
28  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
29  return s;
30}
31int kSBucketLength(kBucket* b)
32{
33  int s=0;
34  int i;
35  for (i=MAX_BUCKET;i>=0;i--)
36  {
37    s+=pSLength(b->buckets[i],0);
38  }
39  return s;
40}
41#endif
42
43#ifdef LEN_VAR3
44// 3.Variante: Laenge: Platz fuer Leitk * Monomanzahl
45int pSLength(poly p,int l)
46{
47  int c=nSize(pGetCoeff(p));
48  return c*l /*pLength(p)*/;
49}
50int kSBucketLength(kBucket* b)
51{
52  int s=0;
53  int c=nSize(pGetCoeff(kBucketGetLm(b)))+1;
54  int i;
55  for (i=MAX_BUCKET;i>0;i--)
56  {
57    s+=b->buckets_length[i] /*pLength(b->buckets[i])*/;
58  }
59  return s*c;
60}
61#endif
62
63#ifdef LEN_VAR4
64// 4.Variante: Laenge: Platz fuer Leitk * (1+Platz fuer andere Koeff.)
65int pSLength(poly p, int l)
66{
67  int s=1;
68  int c=nSize(pGetCoeff(p));
69  pIter(p);
70  while (p!=NULL) { s+=nSize(pGetCoeff(p));pIter(p); }
71  return s*c;
72}
73int kSBucketLength(kBucket* b)
74{
75  int s=1;
76  int c=nSize(pGetCoeff(kBucketGetLm(b)));
77  int i;
78  for (i=MAX_BUCKET;i>0;i--)
79  {
80    if(b->buckets[i]==NULL) continue;
81    s+=pSLength(b->buckets[i],0);
82  }
83  return s*c;
84}
85#endif
86
87static int LObject_better_gen(const void* ap, const void* bp)
88{
89  LObject* a=*(LObject**)ap;
90  LObject* b=*(LObject**)bp;
91  return(pLmCmp(a->p,b->p));
92}
93static int red_object_better_gen(const void* ap, const void* bp)
94{
95
96
97  return(pLmCmp(((red_object*) ap)->p,((red_object*) bp)->p));
98}
99static int pair_better_gen2(const void* ap,const void* bp){
100  return(-pair_better_gen(ap,bp));
101}
102static int kFindDivisibleByInS_easy(kStrategy strat,const red_object & obj){
103  int i;
104  long not_sev=~obj.sev;
105  poly p=obj.p;
106  for(i=0;i<=strat->sl;i++){
107    if (pLmShortDivisibleBy(strat->S[i],strat->sevS[i],p,not_sev))
108      return i;
109  }
110  return -1;
111}
112static int posInPairs (sorted_pair_node**  p, int pn, sorted_pair_node* qe,calc_dat* c,int an=0)
113{
114  if(pn==0) return 0;
115
116  int length=pn-1;
117  int i;
118  //int an = 0;
119  int en= length;
120
121  if (pair_better(qe,p[en],c))
122    return length+1;
123
124  while(1)
125    {
126      //if (an >= en-1)
127      if(en-1<=an)
128      {
129        if (pair_better(p[an],qe,c)) return an;
130        return en;
131      }
132      i=(an+en) / 2;
133        if (pair_better(p[i],qe,c))
134          en=i;
135      else an=i;
136    }
137}
138static BOOLEAN  ascending(int* i,int top){
139  if(top<1) return TRUE;
140  if(i[top]<i[top-1]) return FALSE;
141  return ascending(i,top-1);
142}
143
144sorted_pair_node**  merge(sorted_pair_node** p, int pn,sorted_pair_node **q, int qn,calc_dat* c){
145  int i;
146  int* a= (int*) omalloc(qn*sizeof(int));
147//   int mc;
148//   PrintS("Debug\n");
149//   for(mc=0;mc<qn;mc++)
150// {
151
152//     wrp(q[mc]->lcm_of_lm);
153//     PrintS("\n");
154// }
155//    PrintS("Debug they are in\n");
156//   for(mc=0;mc<pn;mc++)
157// {
158
159//     wrp(p[mc]->lcm_of_lm);
160//     PrintS("\n");
161// }
162  int lastpos=0;
163  for(i=0;i<qn;i++){
164    lastpos=posInPairs(p,pn,q[i],c, max(lastpos-1,0));
165    //   cout<<lastpos<<"\n";
166    a[i]=lastpos;
167
168  }
169  if((pn+qn)>c->max_pairs){
170    p=(sorted_pair_node**) omrealloc(p,2*(pn+qn)*sizeof(sorted_pair_node*));
171    c->max_pairs=2*(pn+qn);
172  }
173  for(i=qn-1;i>=0;i--){
174    size_t size;
175    if(qn-1>i)
176      size=(a[i+1]-a[i])*sizeof(sorted_pair_node*);
177    else
178      size=(pn-a[i])*sizeof(sorted_pair_node*); //as indices begin with 0
179    memmove (p+a[i]+(1+i), p+a[i], size);
180    p[a[i]+i]=q[i];
181  }
182  omfree(a);
183  return p;
184}
185
186
187static BOOLEAN trivial_syzygie(int pos1,int pos2,poly bound,calc_dat* c){
188
189
190  poly p1=c->S->m[pos1];
191  poly p2=c->S->m[pos2];
192  ring r=c->r;
193 
194
195  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
196    return FALSE;
197  int i = 1;
198  poly m=NULL;
199  poly gcd1=c->gcd_of_terms[pos1];
200  poly gcd2=c->gcd_of_terms[pos2];
201 
202  if((gcd1!=NULL) && (gcd2!=NULL)) 
203    {
204      gcd1->next=gcd2; //may ordered incorrect
205      poly m=gcd_of_terms(gcd1,c->r);
206      gcd1->next=NULL;
207     
208    } 
209
210  if (m==NULL) 
211  {
212     loop
213      {
214        if (pGetExp(p1, i)+ pGetExp(p2, i) > pGetExp(bound,i))   return FALSE;
215        if (i == pVariables){
216          PrintS("trivial");
217          return TRUE;
218        }
219        i++;
220      }
221  }
222  else 
223  {
224    loop
225      {
226        if (pGetExp(p1, i)-pGetExp(m,i) + pGetExp(p2, i) > pGetExp(bound,i))   return FALSE;
227        if (i == pVariables){
228          pDelete(&m);
229          PrintS("trivial");
230          return TRUE;
231        }
232        i++;
233      }
234  }
235
236 
237
238 
239}
240
241BOOLEAN good_has_t_rep(int i, int j,calc_dat* c){
242  assume(i>=0);
243    assume(j>=0);
244  if (has_t_rep(i,j,c)) return TRUE;
245  poly lm=pOne();
246
247  pLcm(c->S->m[i], c->S->m[j], lm);
248  pSetm(lm);
249  assume(lm!=NULL);
250  int deciding_deg= pTotaldegree(lm);
251  int* i_con =make_connections(i,j,lm,c);
252  p_Delete(&lm,c->r);
253
254
255  for (int n=0;((n<c->n) && (i_con[n]>=0));n++){
256    if (i_con[n]==j){
257      now_t_rep(i,j,c);
258      omfree(i_con);
259
260      return TRUE;
261    }
262  }
263  omfree(i_con);
264
265  return FALSE;
266}
267BOOLEAN lenS_correct(kStrategy strat){
268  int i;
269  for(i=0;i<=strat->sl;i++){
270    if (strat->lenS[i]!=pLength(strat->S[i]))
271      return FALSE;
272  }
273  return TRUE;
274}
275
276static void notice_miss(int i, int j, calc_dat* c){
277  PrintS("-");
278 
279}
280
281static void cleanS(kStrategy strat){
282  int i=0;
283  LObject P;
284  while(i<=strat->sl){
285    P.p=strat->S[i];
286    P.sev=strat->sevS[i];
287    if(kFindDivisibleByInS(strat->S,strat->sevS,strat->sl,&P)!=i){
288      deleteInS(i,strat);
289      //remember destroying poly
290    }
291    else i++;
292  }
293}
294static int bucket_guess(kBucket* bucket){
295  int sum=0;
296  int i;
297  for (i=MAX_BUCKET;i>=0;i--){
298    sum+=bucket->buckets_length[i];
299  }
300  return sum;
301}
302
303
304
305
306
307
308static int add_to_reductors(calc_dat* c, poly h, int len){
309  assume(lenS_correct(c->strat));
310 
311  int i;
312  if (c->is_char0)
313       i=simple_posInS(c->strat,h,pSLength(h,len),c->is_char0);
314  else
315    i=simple_posInS(c->strat,h,len,c->is_char0);
316 
317  LObject P; memset(&P,0,sizeof(P));
318  P.tailRing=c->r;
319  P.p=h; /*p_Copy(h,c->r);*/
320  P.FDeg=pFDeg(P.p,c->r);
321  if (!rField_is_Zp(c->r)){ 
322    pCleardenom(P.p);
323    pContent(P.p); //is a duplicate call, but belongs here
324  }
325 
326  else                     
327    pNorm(P.p);
328 
329 
330  c->strat->enterS(P,i,c->strat);
331 
332 
333
334  c->strat->lenS[i]=len;
335 
336  if(c->strat->lenSw)
337    c->strat->lenSw[i]=pSLength(P.p,len);
338  return i;
339 
340}
341static void length_one_crit(calc_dat* c, int pos, int len)
342{
343  if (len==1)
344  {
345    int i;
346    for ( i=0;i<pos;i++)
347    {
348      if (c->lengths[i]==1)
349        c->states[pos][i]=HASTREP;
350    }
351    for ( i=pos+1;i<c->n;i++){
352      if (c->lengths[i]==1)
353        c->states[i][pos]=HASTREP;
354    }
355    shorten_tails(c,c->S->m[pos]);
356  }
357}
358static sorted_pair_node* find_next_pair2(calc_dat* c, BOOLEAN go_higher){
359  clean_top_of_pair_list(c);
360  sorted_pair_node* s=pop_pair(c);
361
362
363  return s;
364}
365
366static void move_forward_in_S(int old_pos, int new_pos,kStrategy strat, BOOLEAN is_char0)
367{
368  assume(old_pos>=new_pos);
369  poly p=strat->S[old_pos];
370  int ecart=strat->ecartS[old_pos];
371  long sev=strat->sevS[old_pos];
372  int s_2_r=strat->S_2_R[old_pos];
373  int length=strat->lenS[old_pos];
374  int length_w;
375  if(is_char0)
376    length_w=strat->lenSw[old_pos];
377  int i;
378  for (i=old_pos; i>new_pos; i--)
379  {
380    strat->S[i] = strat->S[i-1];
381    strat->ecartS[i] = strat->ecartS[i-1];
382    strat->sevS[i] = strat->sevS[i-1];
383    strat->S_2_R[i] = strat->S_2_R[i-1];
384  }
385  if (strat->lenS!=NULL)
386    for (i=old_pos; i>new_pos; i--)
387      strat->lenS[i] = strat->lenS[i-1];
388  if (strat->lenSw!=NULL)
389    for (i=old_pos; i>new_pos; i--)
390      strat->lenSw[i] = strat->lenSw[i-1];
391
392  strat->S[new_pos]=p;
393  strat->ecartS[new_pos]=ecart;
394  strat->sevS[new_pos]=sev;
395  strat->S_2_R[new_pos]=s_2_r;
396  strat->lenS[new_pos]=length;
397  if(is_char0)
398    strat->lenSw[new_pos]=length_w;
399  //assume(lenS_correct(strat));
400}
401static void replace_pair(int & i, int & j, calc_dat* c)
402{
403  c->soon_free=NULL;
404  int curr_deg;
405  poly lm=pOne();
406
407  pLcm(c->S->m[i], c->S->m[j], lm);
408  pSetm(lm);
409  int deciding_deg= pTotaldegree(lm);
410  int* i_con =make_connections(i,j,lm,c);
411  int z=0;
412
413  for (int n=0;((n<c->n) && (i_con[n]>=0));n++){
414    if (i_con[n]==j){
415      now_t_rep(i,j,c);
416      omfree(i_con);
417      p_Delete(&lm,c->r);
418      return;
419    }
420  }
421
422  int* j_con =make_connections(j,lm,c);
423  i= i_con[0];
424  j=j_con[0];
425  if(c->n>1){
426    if (i_con[1]>=0)
427      i=i_con[1];
428    else {
429      if (j_con[1]>=0)
430        j=j_con[1];
431    }
432  }
433  pLcm(c->S->m[i], c->S->m[j], lm);
434  pSetm(lm);
435  poly short_s;
436  curr_deg=pTotaldegree(lm);
437  int_pair_node* last=NULL;
438
439  for (int n=0;((n<c->n) && (j_con[n]>=0));n++){
440    for (int m=0;((m<c->n) && (i_con[m]>=0));m++){
441      pLcm(c->S->m[i_con[m]], c->S->m[j_con[n]], lm);
442      pSetm(lm);
443      if (pTotaldegree(lm)>=deciding_deg)
444      {
445        soon_t_rep(i_con[m],j_con[n],c);
446        int_pair_node* h= (int_pair_node*)omalloc(sizeof(int_pair_node));
447        if (last!=NULL)
448          last->next=h;
449        else
450          c->soon_free=h;
451        h->next=NULL;
452        h->a=i_con[m];
453        h->b=j_con[n];
454        last=h;
455      }
456      //      if ((comp_deg<curr_deg)
457      //  ||
458      //  ((comp_deg==curr_deg) &&
459      short_s=ksCreateShortSpoly(c->S->m[i_con[m]],c->S->m[j_con[n]],c->r);
460      if (short_s==NULL) {
461        i=i_con[m];
462        j=j_con[n];
463        now_t_rep(i_con[m],j_con[n],c);
464        p_Delete(&lm,c->r);
465        omfree(i_con);
466        omfree(j_con);
467
468        return;
469      }
470#ifdef QUICK_SPOLY_TEST
471      for (int dz=0;dz<=c->n;dz++){
472        if (dz==c->n) {
473          //have found not head reducing pair
474          i=i_con[m];
475          j=j_con[n];
476          p_Delete(&short_s,c->r);
477          p_Delete(&lm,c->r);
478          omfree(i_con);
479          omfree(j_con);
480
481          return;
482        }
483        if (p_LmDivisibleBy(c->S->m[dz],short_s,c->r)) break;
484      }
485#endif
486      int comp_deg(pTotaldegree(short_s));
487      p_Delete(&short_s,c->r);
488      if ((comp_deg<curr_deg))
489         
490      {
491        curr_deg=comp_deg;
492        i=i_con[m];
493        j=j_con[n];
494      }
495    }
496  }
497  p_Delete(&lm,c->r);
498  omfree(i_con);
499  omfree(j_con);
500  return;
501}
502
503
504static int* make_connections(int from, poly bound, calc_dat* c)
505{
506  ideal I=c->S;
507  int s=pTotaldegree(bound);
508  int* cans=(int*) omalloc(c->n*sizeof(int));
509  int* connected=(int*) omalloc(c->n*sizeof(int));
510  int cans_length=0;
511  connected[0]=from;
512  int connected_length=1;
513  long neg_bounds_short= ~p_GetShortExpVector(bound,c->r);
514  for (int i=0;i<c->n;i++){
515    if (c->T_deg[i]>s) continue;
516    if (i!=from){
517      if(p_LmShortDivisibleBy(I->m[i],c->short_Exps[i],bound,neg_bounds_short,c->r)){
518        cans[cans_length]=i;
519        cans_length++;
520      }
521    }
522  }
523  int not_yet_found=cans_length;
524  int con_checked=0;
525  int pos;
526  while((not_yet_found>0) && (con_checked<connected_length)){
527    pos=connected[con_checked];
528    for(int i=0;i<cans_length;i++){
529      if (cans[i]<0) continue;
530      if (has_t_rep(pos,cans[i],c))
531      {
532        connected[connected_length]=cans[i];
533        connected_length++;
534        cans[i]=-1;
535        --not_yet_found;
536      }
537    }
538    con_checked++;
539  }
540  if (connected_length<c->n){
541    connected[connected_length]=-1;
542  }
543  omfree(cans);
544  return connected;
545}
546static int* make_connections(int from, int to, poly bound, calc_dat* c)
547{
548  ideal I=c->S;
549  int s=pTotaldegree(bound);
550  int* cans=(int*) omalloc(c->n*sizeof(int));
551  int* connected=(int*) omalloc(c->n*sizeof(int));
552  cans[0]=to;
553  int cans_length=1;
554  connected[0]=from;
555  int last_cans_pos=-1;
556  int connected_length=1;
557  long neg_bounds_short= ~p_GetShortExpVector(bound,c->r);
558
559  int not_yet_found=cans_length;
560  int con_checked=0;
561  int pos;
562  BOOLEAN can_find_more=TRUE;
563  while(((not_yet_found>0) && (con_checked<connected_length))||can_find_more){
564    if ((con_checked<connected_length)&& (not_yet_found>0)){
565      pos=connected[con_checked];
566      for(int i=0;i<cans_length;i++){
567        if (cans[i]<0) continue;
568        if (has_t_rep(pos,cans[i],c))//||(trivial_syzygie(pos,cans[i],bound,c))
569{
570
571          connected[connected_length]=cans[i];
572          connected_length++;
573          cans[i]=-1;
574          --not_yet_found;
575
576          if (connected[connected_length-1]==to){
577            if (connected_length<c->n){
578              connected[connected_length]=-1;
579            }
580            omfree(cans);
581            return connected;
582          }
583        }
584      }
585      con_checked++;
586    }
587    else
588    {
589      for(last_cans_pos++;last_cans_pos<=c->n;last_cans_pos++){
590        if (last_cans_pos==c->n){
591          if (connected_length<c->n){
592            connected[connected_length]=-1;
593          }
594          omfree(cans);
595          return connected;
596        }
597        if ((last_cans_pos==from)||(last_cans_pos==to))
598          continue;
599        if(p_LmShortDivisibleBy(I->m[last_cans_pos],c->short_Exps[last_cans_pos],bound,neg_bounds_short,c->r)){
600          cans[cans_length]=last_cans_pos;
601          cans_length++;
602          break;
603        }
604      }
605      not_yet_found++;
606      for (int i=0;i<con_checked;i++){
607        if (has_t_rep(connected[i],last_cans_pos,c)){
608
609          connected[connected_length]=last_cans_pos;
610          connected_length++;
611          cans[cans_length-1]=-1;
612
613          --not_yet_found;
614          if (connected[connected_length-1]==to){
615            if (connected_length<c->n){
616              connected[connected_length]=-1;
617            }
618
619            omfree(cans);
620            return connected;
621          }
622          break;
623        }
624      }
625    }
626  }
627  if (connected_length<c->n){
628    connected[connected_length]=-1;
629  }
630
631  omfree(cans);
632  return connected;
633}
634#ifdef HEAD_BIN
635static inline poly p_MoveHead(poly p, omBin b)
636{
637  poly np;
638  omTypeAllocBin(poly, np, b);
639  memmove(np, p, b->sizeW*sizeof(long));
640  omFreeBinAddr(p);
641  return np;
642}
643#endif
644
645
646
647//very important: ILM
648
649//len should be weighted length in char 0
650static int simple_posInS (kStrategy strat, poly p,int len, BOOLEAN is_char0)
651{
652
653
654  if(strat->sl==-1) return 0;
655  polyset set=strat->S;
656  intset setL=strat->lenS;
657  if (is_char0) setL=strat->lenSw;
658  int length=strat->sl;
659  int i;
660  int an = 0;
661  int en= length;
662
663  if ((len>setL[length])
664      || ((len==setL[length]) && (pLmCmp(set[length],p)== -1)))
665    return length+1;
666
667  loop
668  {
669    if (an >= en-1)
670    {
671      if ((len<setL[an])
672          || ((len==setL[an]) && (pLmCmp(set[an],p) == 1))) return an;
673      return en;
674    }
675    i=(an+en) / 2;
676    if ((len<setL[i])
677        || ((len==setL[i]) && (pLmCmp(set[i],p) == 1))) en=i;
678    //else if ((len>setL[i])
679    //|| ((len==setL[i]) && (pLmCmp(set[i],p) == -1))) an=i;
680    else an=i;
681  }
682}
683/*2
684 *if the leading term of p
685 *divides the leading term of some S[i] it will be canceled
686 */
687static inline void clearS (poly p, unsigned long p_sev,int l, int* at, int* k,
688                           kStrategy strat)
689{
690  assume(p_sev == pGetShortExpVector(p));
691  if (!pLmShortDivisibleBy(p,p_sev, strat->S[*at], ~ strat->sevS[*at])) return;
692  if (l>=strat->lenS[*at]) return;
693  PrintS("!");mflush();
694  //pDelete(&strat->S[*at]);
695  deleteInS((*at),strat);
696  (*at)--;
697  (*k)--;
698//  assume(lenS_correct(strat));
699}
700static sorted_pair_node** add_to_basis(poly h, int i_pos, int j_pos,calc_dat* c, int* ip)
701{
702
703  assume(h!=NULL);
704//  BOOLEAN corr=lenS_correct(c->strat);
705  BOOLEAN R_found=FALSE;
706  void* hp;
707
708  ++(c->n);
709  ++(c->S->ncols);
710  int i,j;
711  i=c->n-1;
712  sorted_pair_node** nodes=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*i);
713  int spc=0;
714  c->T_deg=(int*) omrealloc(c->T_deg,c->n*sizeof(int));
715  c->T_deg[i]=pTotaldegree(h);
716  hp=omrealloc(c->rep, c->n *sizeof(int));
717  if (hp!=NULL){
718    c->rep=(int*) hp;
719  } else {
720    exit(1);
721  }
722  c->short_Exps=(long *) omrealloc(c->short_Exps ,c->n*sizeof(long));
723
724  hp=omrealloc(c->lengths, c->n *sizeof(int));
725  if (hp!=NULL){
726    c->lengths=(int*) hp;
727  } else {
728    exit(1);
729  }
730  c->lengths[i]=pLength(h);
731  hp=omrealloc(c->states, c->n * sizeof(char*));
732 
733    c->states=(char**) hp;
734  c->gcd_of_terms=(poly*) omrealloc(c->gcd_of_terms, c->n *sizeof(poly));
735  c->gcd_of_terms[i]=gcd_of_terms(h,c->r);
736  c->rep[i]=i;
737  hp=omalloc(i*sizeof(char));
738  if (hp!=NULL){
739    c->states[i]=(char*) hp;
740  } else {
741    exit(1);
742  }
743  hp=omrealloc(c->S->m,c->n*sizeof(poly));
744  if (hp!=NULL){
745    c->S->m=(poly*) hp;
746  } else {
747    exit(1);
748  }
749  c->S->m[i]=h;
750  c->short_Exps[i]=p_GetShortExpVector(h,c->r);
751  for (j=0;j<i;j++){
752    if (c->rep[j]==j){
753      //check product criterion
754
755      c->states[i][j]=UNCALCULATED;
756
757      //lies I[i] under I[j] ?
758      if(p_LmShortDivisibleBy(c->S->m[i],c->short_Exps[i],c->S->m[j],~(c->short_Exps[j]),c->r)){
759        c->rep[j]=i;
760       
761        PrintS("R"); R_found=TRUE;
762
763        c->Rcounter++;
764        if((i_pos>=0) && (j_pos>=0)){
765       
766        }
767        for(int z=0;z<j;z++){
768          if(c->rep[z]!=z) continue;
769          if (c->states[j][z]==UNCALCULATED){
770            c->states[j][z]=UNIMPORTANT;
771          }
772        }
773        for(int z=j+1;z<i;z++){
774          if(c->rep[z]!=z) continue;
775          if (c->states[z][j]==UNCALCULATED){
776            c->states[z][j]=UNIMPORTANT;
777          }
778        }
779      }
780    }
781    else {
782      c->states[i][j]=UNIMPORTANT;
783    }
784    if ((c->lengths[i]==1) && (c->lengths[j]==1))
785      c->states[i][j]=HASTREP;
786    else if (pHasNotCF(c->S->m[i],c->S->m[j])){
787      c->easy_product_crit++;
788      c->states[i][j]=HASTREP;
789    }
790                        else if(extended_product_criterion(c->S->m[i],c->gcd_of_terms[i],c->S->m[j],c->gcd_of_terms[j],c)){
791                                        c->states[i][j]=HASTREP;
792                                        c->extended_product_crit++;
793                                        //PrintS("E");
794                        }
795    if (c->states[i][j]==UNCALCULATED){
796
797     
798      poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
799      if (short_s)
800      {
801        sorted_pair_node* s=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
802        s->i=max(i,j);
803        s->j=min(i,j);
804        s->expected_length=c->lengths[i]+c->lengths[j]-2;
805        s->deg=pTotaldegree(short_s);
806        poly lm=pOne();
807
808        pLcm(c->S->m[i], c->S->m[j], lm);
809        pSetm(lm);
810        s->lcm_of_lm=lm;
811          pDelete(&short_s);
812        //assume(lm!=NULL);
813        nodes[spc]=s;
814        spc++;
815      }
816      else
817      {
818        c->states[i][j]=HASTREP;
819      }
820    }
821  }
822
823  add_to_reductors(c, h, c->lengths[c->n-1]);
824  //i=posInS(c->strat,c->strat->sl,h,0 /*ecart*/);
825
826  if (c->lengths[c->n-1]==1)
827    shorten_tails(c,c->S->m[c->n-1]);
828  //you should really update c->lengths, c->strat->lenS, and the oder of polys in strat if you sort after lengths
829
830  //for(i=c->strat->sl; i>0;i--)
831  //  if(c->strat->lenS[i]<c->strat->lenS[i-1]) printf("fehler bei %d\n",i);
832  if (c->Rcounter>50) {
833    c->Rcounter=0;
834    cleanS(c->strat);
835  }
836  if(!ip){
837    qsort(nodes,spc,sizeof(sorted_pair_node*),pair_better_gen2);
838 
839   
840    c->apairs=merge(c->apairs,c->pair_top+1,nodes,spc,c);
841    c->pair_top+=spc;
842    clean_top_of_pair_list(c);
843    omfree(nodes);
844    return NULL;
845  }
846  {
847    *ip=spc;
848    return nodes;
849  }
850
851 
852
853}
854#if 0
855static poly redNF (poly h,kStrategy strat)
856{
857  int j = 0;
858  int z = 3;
859  unsigned long not_sev;
860
861  if (0 > strat->sl)
862  {
863    return h;
864  }
865  not_sev = ~ pGetShortExpVector(h);
866  loop
867    {
868      if (pLmShortDivisibleBy(strat->S[j], strat->sevS[j], h, not_sev))
869      {
870        //if (strat->interpt) test_int_std(strat->kIdeal);
871        /*- compute the s-polynomial -*/
872#ifdef KDEBUG
873        if (TEST_OPT_DEBUG)
874        {
875          PrintS("red:");
876          wrp(h);
877          PrintS(" with ");
878          wrp(strat->S[j]);
879        }
880#endif
881        h = ksOldSpolyRed(strat->S[j],h,strat->kNoether);
882#ifdef KDEBUG
883        if (TEST_OPT_DEBUG)
884        {
885          PrintS("\nto:");
886          wrp(h);
887          PrintLn();
888        }
889#endif
890        if (h == NULL) return NULL;
891        z++;
892        if (z>=10)
893        {
894          z=0;
895          pNormalize(h);
896        }
897        /*- try to reduce the s-polynomial -*/
898        j = 0;
899        not_sev = ~ pGetShortExpVector(h);
900      }
901      else
902      {
903        if (j >= strat->sl) return h;
904        j++;
905      }
906    }
907}
908#else
909
910static poly redNF2 (poly h,calc_dat* c , int &len)
911{
912  len=0;
913  if (h==NULL) return NULL;
914
915  len=pLength(h);
916  kStrategy strat=c->strat;
917  if (0 > strat->sl)
918  {
919    return h;
920  }
921  int j;
922  int len_upper_bound=len;
923  LObject P(h);
924  P.SetShortExpVector();
925  P.bucket = kBucketCreate(currRing);
926  // BOOLEAN corr=lenS_correct(strat);
927  kBucketInit(P.bucket,P.p,len /*pLength(P.p)*/);
928  //int max_pos=simple_posInS(strat,P.p);
929  loop
930    {
931//       if (corr){
932
933//      corr=lenS_correct(strat);
934//      if(!corr){
935//        PrintS("korupt");
936//      }
937//       }
938      int compare_bound;
939      compare_bound=bucket_guess(P.bucket);
940      len_upper_bound=min(compare_bound,len_upper_bound);
941      j=kFindDivisibleByInS(strat->S,strat->sevS,strat->sl,&P);
942      if (j>=0)
943      {
944        poly sec_copy=NULL;
945        //pseudo code
946        BOOLEAN must_expand=FALSE;
947        BOOLEAN must_replace_in_basis=(len_upper_bound<strat->lenS[j]);//first test
948        if (must_replace_in_basis)
949        {
950          //second test
951          if (pLmEqual(P.p,strat->S[j]))
952          {
953            PrintS("b");
954            sec_copy=kBucketClear(P.bucket);
955            sec_copy=redTailShort(sec_copy, strat);
956            kBucketInit(P.bucket,pCopy(sec_copy),pLength(sec_copy));
957          }
958          else
959          {
960            must_replace_in_basis=FALSE;
961            if ((len_upper_bound==1)
962                ||(len_upper_bound==2)
963                ||(len_upper_bound<strat->lenS[j]/2))
964            {
965              PrintS("e");
966              int dummy_len;
967              kBucketClear(P.bucket,&sec_copy,&dummy_len);
968              kBucketInit(P.bucket,pCopy(sec_copy),dummy_len
969                                                  /*pLength(sec_copy)*/);
970              must_expand=TRUE;
971            }
972          }
973        }
974//        must_expand=FALSE;
975//        must_replace_in_basis=FALSE;
976        nNormalize(pGetCoeff(P.p));
977#ifdef KDEBUG
978        if (TEST_OPT_DEBUG)
979        {
980          PrintS("red:");
981          wrp(h);
982          PrintS(" with ");
983          wrp(strat->S[j]);
984        }
985#endif
986        len_upper_bound=len_upper_bound+strat->lenS[j]-2;
987        number coef=kBucketPolyRed(P.bucket,strat->S[j],
988                                   strat->lenS[j]/*pLength(strat->S[j])*/,
989                                   strat->kNoether);
990        nDelete(&coef);
991        h = kBucketGetLm(P.bucket);
992
993        if (must_replace_in_basis){
994          int pos_in_c=-1;
995          poly p=strat->S[j];
996          int z;
997
998          int new_length=pLength(sec_copy);
999          Print("%i",strat->lenS[j]-new_length);
1000          len_upper_bound=new_length +strat->lenS[j]-2;//old entries length
1001          int new_pos;
1002          if(c->is_char0)
1003            new_pos=simple_posInS(c->strat,sec_copy,
1004                                  pSLength(sec_copy,new_length),
1005                                  c->is_char0);//hac
1006          else
1007            new_pos=simple_posInS(c->strat,sec_copy,new_length,c->is_char0);//hack
1008          assume(new_pos<=j);
1009//          p=NULL;
1010          for (z=c->n;z;z--)
1011          {
1012            if(p==c->S->m[z-1])
1013            {
1014
1015
1016              pos_in_c=z-1;
1017
1018              break;
1019            }
1020          }
1021          if (z<=0){
1022            //not in c->S
1023            //LEAVE
1024            deleteInS(j,c->strat);
1025
1026            add_to_reductors(c,sec_copy,pLength(sec_copy));
1027          }
1028          else {
1029//shorten_tails may alter position (not the length, even not by recursion in GLOBAL case)
1030
1031            strat->S[j]=sec_copy;
1032            c->strat->lenS[j]=new_length;
1033            pDelete(&p);
1034
1035            //        replace_quietly(c,j,sec_copy);
1036            // have to do many additional things for consistency
1037            {
1038
1039
1040
1041
1042              int old_pos=j;
1043              new_pos=min(old_pos, new_pos);
1044              assume(new_pos<=old_pos);
1045
1046
1047              c->strat->lenS[old_pos]=new_length;
1048              if(c->strat->lenSw)
1049                c->strat->lenSw[old_pos]=pSLength(sec_copy,new_length);
1050              int i=0;
1051              for(i=new_pos;i<old_pos;i++){
1052                if (strat->lenS[i]<=new_length)
1053                  new_pos++;
1054                else
1055                  break;
1056              }
1057              if (new_pos<old_pos)
1058                move_forward_in_S(old_pos,new_pos,c->strat, c->is_char0);
1059
1060              c->S->m[pos_in_c]=sec_copy;
1061
1062              c->lengths[pos_in_c]=new_length;
1063              if (new_length==1)
1064              {
1065                int i;
1066                for ( i=0;i<pos_in_c;i++)
1067                {
1068                  if (c->lengths[i]==1)
1069                    c->states[pos_in_c][i]=HASTREP;
1070                }
1071                for ( i=z;i<c->n;i++){
1072                  if (c->lengths[i]==1)
1073                    c->states[i][pos_in_c]=HASTREP;
1074                }
1075                shorten_tails(c,sec_copy);
1076              }
1077            }
1078          }
1079        }
1080        if(must_expand){
1081
1082          add_to_reductors(c,sec_copy,pLength(sec_copy));
1083        }
1084        if (h==NULL) return NULL;
1085        P.p=h;
1086        P.t_p=NULL;
1087        P.SetShortExpVector();
1088#ifdef KDEBUG
1089        if (TEST_OPT_DEBUG)
1090        {
1091          PrintS("\nto:");
1092          wrp(h);
1093          PrintLn();
1094        }
1095#endif
1096      }
1097      else
1098      {
1099        kBucketClear(P.bucket,&(P.p),&len);
1100        kBucketDestroy(&P.bucket);
1101        pNormalize(P.p);
1102        return P.p;
1103      }
1104    }
1105}
1106
1107
1108static poly redTailShort(poly h, kStrategy strat){
1109
1110  int sl=strat->sl;
1111  int i;
1112  int len=pLength(h);
1113  for(i=0;i<=strat->sl;i++){
1114    if(strat->lenS[i]>2)
1115      break;
1116  }
1117  return(redNFTail(h,i-1,strat, len));
1118}
1119
1120static void line_of_extended_prod(int fixpos,calc_dat* c){
1121    if (c->gcd_of_terms[fixpos]==NULL)
1122  {
1123    c->gcd_of_terms[fixpos]=gcd_of_terms(c->S->m[fixpos],c->r);
1124    if (c->gcd_of_terms[fixpos])
1125    {
1126      int i;
1127      for(i=0;i<fixpos;i++)
1128        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)))
1129{
1130          c->states[fixpos][i]=HASTREP;
1131          c->extended_product_crit++;
1132}     
1133      for(i=fixpos+1;i<c->n;i++)
1134        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)))
1135        {        c->states[i][fixpos]=HASTREP;
1136        c->extended_product_crit++;
1137        }
1138    }
1139  }
1140}
1141static void c_S_element_changed_hook(int pos, calc_dat* c){
1142  length_one_crit(c,pos, c->lengths[pos]);
1143  line_of_extended_prod(pos,c);
1144}
1145
1146static void go_on (calc_dat* c){
1147  //set limit of 1000 for multireductions, at the moment for
1148  //programming reasons
1149  int i=0;
1150  red_object* buf=(red_object*) omalloc(PAR_N*sizeof(red_object));
1151  int curr_deg=-1;
1152  while(i<PAR_N){
1153    sorted_pair_node* s=top_pair(c);
1154    if (!s) break;
1155    if(curr_deg>=0){
1156      if (s->deg >curr_deg) break;
1157    }
1158
1159    else curr_deg=s->deg;
1160    quick_pop_pair(c);
1161    if(s->i>=0){
1162      //replace_pair(s->i,s->j,c);
1163    if(s->i==s->j) {
1164      free_sorted_pair_node(s,c->r);
1165      continue;
1166        }
1167    }
1168    poly h;
1169    if(s->i>=0)
1170      h=ksOldCreateSpoly(c->S->m[s->i], c->S->m[s->j], NULL, c->r);
1171    else
1172      h=s->lcm_of_lm;
1173    if(s->i>=0)
1174      now_t_rep(s->j,s->i,c);
1175    free_sorted_pair_node(s,c->r);
1176    if(!h) continue;
1177    int len=pLength(h);
1178    buf[i].p=h;
1179    buf[i].sev=pGetShortExpVector(h);
1180    buf[i].sum=NULL;
1181    buf[i].bucket = kBucketCreate(currRing);
1182    kBucketInit(buf[i].bucket,buf[i].p,len);
1183    i++;
1184  }
1185  c->normal_forms+=i;
1186  qsort(buf,i,sizeof(red_object),red_object_better_gen);
1187//    Print("\ncurr_deg:%i\n",curr_deg);
1188  Print("M[%i, ",i); 
1189  multi_reduction(buf, i, c);
1190  Print("%i]",i);
1191  int j;
1192 //  for(j=0;j<i;j++){
1193//     if(buf[j].p==NULL) PrintS("\n ZERO ALERT \n");
1194//     int z;
1195//      for(z=0;z<j;z++){
1196//       if (pLmEqual(buf[z].p, buf[j].p))
1197//      PrintS("\n Critical Warning!!!! \n");
1198     
1199//     }
1200//   }
1201  int* ibuf=(int*) omalloc(i*sizeof(int));
1202  sorted_pair_node*** sbuf=(sorted_pair_node***) omalloc(i*sizeof(sorted_pair_node**));
1203  for(j=0;j<i;j++){
1204 
1205 
1206    int len;
1207    poly p;
1208    kBucketClear(buf[j].bucket,&p, &len);
1209    kBucketDestroy(&buf[j].bucket);
1210    // delete buf[j];
1211    //remember to free res here
1212    p=redTailShort(p, c->strat);
1213    sbuf[j]=add_to_basis(p,-1,-1,c,ibuf+j);
1214   
1215  }
1216  int sum=0;
1217  for(j=0;j<i;j++){
1218    sum+=ibuf[j];
1219  }
1220  sorted_pair_node** big_sbuf=(sorted_pair_node**) omalloc(sum*sizeof(sorted_pair_node*));
1221  int partsum=0;
1222  for(j=0;j<i;j++){
1223    memmove(big_sbuf+partsum, sbuf[j],ibuf[j]*sizeof(sorted_pair_node*));
1224    omfree(sbuf[j]);
1225    partsum+=ibuf[j];
1226  }
1227
1228  qsort(big_sbuf,sum,sizeof(sorted_pair_node*),pair_better_gen2);
1229  c->apairs=merge(c->apairs,c->pair_top+1,big_sbuf,sum,c);
1230  c->pair_top+=sum;
1231  clean_top_of_pair_list(c);
1232  omfree(big_sbuf);
1233  omfree(sbuf);
1234  omfree(ibuf);
1235  omfree(buf);
1236  return;
1237}
1238
1239static poly redNF (poly h,kStrategy strat, int &len)
1240{
1241  len=0;
1242  if (h==NULL) return NULL;
1243  int j;
1244
1245  len=pLength(h);
1246  if (0 > strat->sl)
1247  {
1248    return h;
1249  }
1250  LObject P(h);
1251  P.SetShortExpVector();
1252  P.bucket = kBucketCreate(currRing);
1253  kBucketInit(P.bucket,P.p,len /*pLength(P.p)*/);
1254  //int max_pos=simple_posInS(strat,P.p);
1255  loop
1256    {
1257      j=kFindDivisibleByInS(strat->S,strat->sevS,strat->sl,&P);
1258      if (j>=0)
1259      {
1260        nNormalize(pGetCoeff(P.p));
1261#ifdef KDEBUG
1262        if (TEST_OPT_DEBUG)
1263        {
1264          PrintS("red:");
1265          wrp(h);
1266          PrintS(" with ");
1267          wrp(strat->S[j]);
1268        }
1269#endif
1270        number coef=kBucketPolyRed(P.bucket,strat->S[j],
1271                                   strat->lenS[j]/*pLength(strat->S[j])*/,
1272                                   strat->kNoether);
1273        nDelete(&coef);
1274        h = kBucketGetLm(P.bucket);
1275        if (h==NULL) return NULL;
1276        P.p=h;
1277        P.t_p=NULL;
1278        P.SetShortExpVector();
1279#ifdef KDEBUG
1280        if (TEST_OPT_DEBUG)
1281        {
1282          PrintS("\nto:");
1283          wrp(h);
1284          PrintLn();
1285        }
1286#endif
1287      }
1288      else
1289      {
1290        kBucketClear(P.bucket,&(P.p),&len);
1291        kBucketDestroy(&P.bucket);
1292        pNormalize(P.p);
1293        return P.p;
1294      }
1295    }
1296}
1297#endif
1298#ifdef REDTAIL_S
1299
1300static poly redNFTail (poly h,const int sl,kStrategy strat, int len)
1301{
1302  if (h==NULL) return NULL;
1303  pTest(h);
1304  if (0 > sl)
1305    return h;
1306  if (pNext(h)==NULL) return h;
1307
1308  int j;
1309  poly res=h;
1310  poly act=res;
1311  LObject P(pNext(h));
1312  pNext(res)=NULL;
1313  P.bucket = kBucketCreate(currRing);
1314  len--;
1315  h=P.p;
1316  if (len <=0) len=pLength(h);
1317  kBucketInit(P.bucket,h /*P.p*/,len /*pLength(P.p)*/);
1318  pTest(h);
1319  loop
1320    {
1321      P.p=h;
1322      P.t_p=NULL;
1323      P.SetShortExpVector();
1324      loop
1325        {
1326          j=kFindDivisibleByInS(strat->S,strat->sevS,sl,&P);
1327          if (j>=0)
1328          {
1329#ifdef REDTAIL_PROT
1330            PrintS("r");
1331#endif
1332            nNormalize(pGetCoeff(P.p));
1333#ifdef KDEBUG
1334            if (TEST_OPT_DEBUG)
1335            {
1336              PrintS("red tail:");
1337              wrp(h);
1338              PrintS(" with ");
1339              wrp(strat->S[j]);
1340            }
1341#endif
1342            number coef;
1343            pTest(strat->S[j]);
1344            coef=kBucketPolyRed(P.bucket,strat->S[j],
1345                                strat->lenS[j]/*pLength(strat->S[j])*/,strat->kNoether);
1346            pMult_nn(res,coef);
1347            nDelete(&coef);
1348            h = kBucketGetLm(P.bucket);
1349            pTest(h);
1350            if (h==NULL)
1351            {
1352#ifdef REDTAIL_PROT
1353              PrintS(" ");
1354#endif
1355              return res;
1356            }
1357            pTest(h);
1358            P.p=h;
1359            P.t_p=NULL;
1360            P.SetShortExpVector();
1361#ifdef KDEBUG
1362            if (TEST_OPT_DEBUG)
1363            {
1364              PrintS("\nto tail:");
1365              wrp(h);
1366              PrintLn();
1367            }
1368#endif
1369          }
1370          else
1371          {
1372#ifdef REDTAIL_PROT
1373            PrintS("n");
1374#endif
1375            break;
1376          }
1377        } /* end loop current mon */
1378      poly tmp=pHead(h /*kBucketGetLm(P.bucket)*/);
1379      act->next=tmp;pIter(act);
1380      poly tmp2=pHead(h);
1381      pNeg(tmp2);
1382      int ltmp2=1;
1383      pTest(tmp2);
1384      kBucket_Add_q(P.bucket, tmp2, &ltmp2);
1385
1386      h = kBucketGetLm(P.bucket);
1387      if (h==NULL)
1388      {
1389#ifdef REDTAIL_PROT
1390        PrintS(" ");
1391#endif
1392        return res;
1393      }
1394      pTest(h);
1395    }
1396}
1397#endif
1398
1399static void do_this_spoly_stuff(int i,int j,calc_dat* c){
1400  poly f=c->S->m[i];
1401  poly g=c->S->m[j];
1402  poly h=ksOldCreateSpoly(f, g, NULL, c->r);
1403  poly hr=NULL;
1404#ifdef FULLREDUCTIONS
1405  if (h!=NULL)
1406  {
1407    int len;
1408
1409    hr=redNF2(h,c,len);
1410//      hr=redNF(h,c->strat,len);
1411
1412    if (hr!=NULL)
1413#ifdef REDTAIL_S
1414      hr = redNFTail(hr,c->strat->sl,c->strat,len);
1415#else
1416    hr = redtailBba(hr,c->strat->sl,c->strat);
1417#endif
1418
1419  }
1420#else
1421  if (h!=NULL)
1422  {
1423    int len;
1424    hr=redNF2(h,c,len);
1425  }
1426#endif
1427  c->normal_forms++;
1428  if (hr==NULL)
1429  {
1430    notice_miss(i, j, c);
1431
1432  }
1433  else
1434  {
1435
1436#ifdef HEAD_BIN
1437    hr=p_MoveHead(hr,c->HeadBin);
1438#endif
1439    add_to_basis(hr, i,j,c);
1440  }
1441}
1442//try to fill, return FALSE iff queue is empty
1443
1444static int poly_crit(const void* ap1, const void* ap2){
1445  poly p1,p2;
1446  p1=*((poly*) ap1);
1447  p2=*((poly*)ap2);
1448
1449  int c=pLmCmp(p1,p2);
1450  if (c !=0) return c;
1451  int l1=pLength(p1);
1452  int l2=pLength(p2);
1453  if (l1<l2) return -1;
1454  if (l1>l2) return 1;
1455  return 0;
1456}
1457ideal t_rep_gb(ring r,ideal arg_I){
1458   Print("Idelems %i \n----------\n",IDELEMS(arg_I));
1459  ideal I=idCompactify(arg_I);
1460  qsort(I->m,IDELEMS(I),sizeof(poly),poly_crit);
1461  Print("Idelems %i \n----------\n",IDELEMS(I));
1462  calc_dat* c=(calc_dat*) omalloc(sizeof(calc_dat));
1463  c->r=currRing;
1464  void* h;
1465  poly hp;
1466  int i,j;
1467  c->easy_product_crit=0;
1468  c->extended_product_crit=0;
1469  c->is_char0=(rChar()==0);
1470  c->reduction_steps=0;
1471  c->last_index=-1;
1472
1473
1474
1475  c->Rcounter=0;
1476
1477  c->soon_free=NULL;
1478
1479
1480  c->normal_forms=0;
1481  c->current_degree=1;
1482 
1483  c->max_pairs=5*I->idelems();
1484 
1485  c->apairs=(sorted_pair_node**) omalloc(sizeof(sorted_pair_node*)*c->max_pairs);
1486  c->pair_top=-1;
1487  int n=I->idelems();
1488  for (i=0;i<n;i++){
1489    wrp(I->m[i]);
1490    PrintS("\n");
1491  }
1492    i=0;
1493  c->n=0;
1494  c->T_deg=(int*) omalloc(n*sizeof(int));
1495 
1496#ifdef HEAD_BIN
1497  c->HeadBin=omGetSpecBin(POLYSIZE + (currRing->ExpL_Size)*sizeof(long));
1498#endif
1499  /* omUnGetSpecBin(&(c->HeadBin)); */
1500  h=omalloc(n*sizeof(char*));
1501  c->states=(char**) h;
1502  h=omalloc(n*sizeof(int));
1503  c->lengths=(int*) h;
1504  h=omalloc(n*sizeof(int));
1505        c->gcd_of_terms=(poly*) omalloc(n*sizeof(poly));
1506  c->rep=(int*) h;
1507  c->short_Exps=(long*) omalloc(n*sizeof(long));
1508  c->S=idInit(n,1);
1509  c->strat=new skStrategy;
1510  c->strat->syzComp = 0;
1511  initBuchMoraCrit(c->strat);
1512  initBuchMoraPos(c->strat);
1513  c->strat->initEcart = initEcartBBA;
1514  c->strat->enterS = enterSBba;
1515  c->strat->sl = -1;
1516  i=n;
1517  /* initS(c->S,NULL,c->strat); */
1518/* intS start: */
1519  i=((i+IDELEMS(c->S)+15)/16)*16;
1520  c->strat->ecartS=(intset)omAlloc(i*sizeof(int)); /*initec(i);*/
1521  c->strat->sevS=(unsigned long*)omAlloc0(i*sizeof(unsigned long));
1522  /*initsevS(i);*/
1523  c->strat->S_2_R=(int*)omAlloc0(i*sizeof(int));/*initS_2_R(i);*/
1524  c->strat->fromQ=NULL;
1525  c->strat->Shdl=idInit(1,1);
1526  c->strat->S=c->strat->Shdl->m;
1527  c->strat->lenS=(int*)omAlloc0(i*sizeof(int));
1528  if(c->is_char0)
1529    c->strat->lenSw=(int*)omAlloc0(i*sizeof(int));
1530  else
1531    c->strat->lenSw=NULL;
1532  sorted_pair_node* si;
1533  assume(n>0);
1534  add_to_basis(I->m[0],-1,-1,c);
1535
1536  assume(c->strat->sl==c->strat->Shdl->idelems()-1);
1537
1538  for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
1539   {
1540//     add_to_basis(I->m[i],-1,-1,c);
1541     si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1542      si->i=-1;
1543      si->j=-1;
1544      si->expected_length=pLength(I->m[i]);
1545      si->deg=pTotaldegree(I->m[i]);
1546      si->lcm_of_lm=I->m[i];
1547
1548//      c->apairs[n-1-i]=si;
1549      c->apairs[n-i-1]=si;
1550      ++(c->pair_top);
1551   }
1552 
1553
1554  while(c->pair_top>=0)
1555    go_on(c);
1556
1557  for(int z=0;z<c->n;z++){
1558    omfree(c->states[z]);
1559  }
1560  omfree(c->states);
1561  omfree(c->lengths);
1562
1563
1564  omfree(c->short_Exps);
1565  omfree(c->T_deg);
1566
1567     omFree(c->strat->ecartS);
1568     omFree(c->strat->sevS);
1569//   /*initsevS(i);*/
1570   omFree(c->strat->S_2_R);
1571   
1572
1573  omFree(c->strat->lenS);
1574
1575   if(c->strat->lenSw)  omFree(c->strat->lenSw);
1576
1577
1578
1579
1580  for(i=0;i<c->n;i++){
1581    if(c->gcd_of_terms[i])
1582      pDelete(&(c->gcd_of_terms[i]));
1583  }
1584  omfree(c->gcd_of_terms);
1585
1586  omfree(c->apairs);
1587  printf("calculated %d NFs\n",c->normal_forms);
1588  printf("applied %i product crit, %i extended_product crit \n", c->easy_product_crit, c->extended_product_crit);
1589  int deleted_form_c_s=0;
1590
1591  for(i=0;i<c->n;i++){
1592    if (c->rep[i]!=i){
1593      for(j=0;j<=c->strat->sl;j++){
1594        if(c->strat->S[j]==c->S->m[i]){
1595          c->strat->S[j]=NULL;
1596          break;
1597        }
1598      }
1599      PrintS("R_delete");
1600      pDelete(&c->S->m[i]);
1601    }
1602  }
1603  for(i=0;i<=c->strat->sl;i++){
1604    if (!c->strat->S[i]) continue;
1605    BOOLEAN found=FALSE;
1606    for(j=0;j<c->n;j++){
1607      if (c->S->m[j]==c->strat->S[i]){
1608        found=TRUE;
1609        break;
1610      }
1611    }
1612    if(!found) pDelete(&c->strat->S[i]);
1613  }
1614  omfree(c->rep);
1615  I=c->S;
1616 
1617  IDELEMS(I)=c->n;
1618
1619  idSkipZeroes(c->S);
1620  for(i=0;i<=c->strat->sl;i++)
1621    c->strat->S[i]=NULL;
1622  id_Delete(&c->strat->Shdl,c->r);
1623  delete c->strat;
1624  omfree(c);
1625
1626  return(I);
1627}
1628static void now_t_rep(const int & arg_i, const int & arg_j, calc_dat* c){
1629  int i,j;
1630  if (arg_i==arg_j){
1631    return;
1632  }
1633  if (arg_i>arg_j){
1634    i=arg_j;
1635    j=arg_i;
1636  } else {
1637    i=arg_i;
1638    j=arg_j;
1639  }
1640  c->states[j][i]=HASTREP;
1641}
1642static void soon_t_rep(const int& arg_i, const int& arg_j, calc_dat* c)
1643{
1644  assume(0<=arg_i);
1645  assume(0<=arg_j);
1646  assume(arg_i<c->n);
1647  assume(arg_j<c->n);
1648  int i,j;
1649  if (arg_i==arg_j){
1650    return;
1651  }
1652  if (arg_i>arg_j){
1653    i=arg_j;
1654    j=arg_i;
1655  } else {
1656    i=arg_i;
1657    j=arg_j;
1658  }
1659  if (!
1660      (c->states[j][i]==HASTREP))
1661    c->states[j][i]=SOONTREP;
1662}
1663static BOOLEAN has_t_rep(const int & arg_i, const  int & arg_j, calc_dat* state){
1664  assume(0<=arg_i);
1665  assume(0<=arg_j);
1666  assume(arg_i<state->n);
1667  assume(arg_j<state->n);
1668  if (arg_i==arg_j)
1669  {
1670    return (TRUE);
1671  }
1672  if (arg_i>arg_j)
1673  {
1674    return (state->states[arg_i][arg_j]==HASTREP);
1675  } else
1676  {
1677    return (state->states[arg_j][arg_i]==HASTREP);
1678  }
1679}
1680static int pLcmDeg(poly a, poly b)
1681{
1682  int i;
1683  int n=0;
1684  for (i=pVariables; i; i--)
1685  {
1686    n+=max( pGetExp(a,i), pGetExp(b,i));
1687  }
1688  return n;
1689
1690}
1691
1692
1693
1694static void shorten_tails(calc_dat* c, poly monom)
1695{
1696  return;
1697// BOOLEAN corr=lenS_correct(c->strat);
1698  for(int i=0;i<c->n;i++)
1699  {
1700    //enter tail
1701    if (c->rep[i]!=i) continue;
1702    if (c->S->m[i]==NULL) continue;
1703    poly tail=c->S->m[i]->next;
1704    poly prev=c->S->m[i];
1705    BOOLEAN did_something=FALSE;
1706    while((tail!=NULL)&& (pLmCmp(tail, monom)>=0))
1707    {
1708      if (p_LmDivisibleBy(monom,tail,c->r))
1709      {
1710        did_something=TRUE;
1711        prev->next=tail->next;
1712        tail->next=NULL;
1713        p_Delete(& tail,c->r);
1714        tail=prev;
1715        //PrintS("Shortened");
1716        c->lengths[i]--;
1717      }
1718      prev=tail;
1719      tail=tail->next;
1720    }
1721    if (did_something)
1722    {
1723      int new_pos;
1724      int q;
1725      q=quality(c->S->m[i],c->lengths[i],c);
1726      new_pos=simple_posInS(c->strat,c->S->m[i],q,c->is_char0);
1727
1728      int old_pos=-1;
1729      //assume new_pos<old_pos
1730      for (int z=0;z<=c->strat->sl;z++)
1731      {
1732        if (c->strat->S[z]==c->S->m[i])
1733        {
1734          old_pos=z;
1735          break;
1736        }
1737      }
1738      if (old_pos== -1)
1739        for (int z=new_pos-1;z>=0;z--)
1740        {
1741          if (c->strat->S[z]==c->S->m[i])
1742          {
1743            old_pos=z;
1744            break;
1745          }
1746        }
1747      assume(old_pos>=0);
1748      assume(new_pos<=old_pos);
1749      assume(pLength(c->strat->S[old_pos])==c->lengths[i]);
1750      c->strat->lenS[old_pos]=c->lengths[i];
1751      if (c->strat->lenSw)
1752        c->strat->lenSw[old_pos]=q;
1753
1754      if (new_pos<old_pos)
1755        move_forward_in_S(old_pos,new_pos,c->strat, c->is_char0);
1756
1757      length_one_crit(c,i,c->lengths[i]);
1758    }
1759  }
1760}
1761static sorted_pair_node* pop_pair(calc_dat* c){
1762  clean_top_of_pair_list(c);
1763
1764  if(c->pair_top<0) return NULL;
1765  else return (c->apairs[c->pair_top--]);
1766}
1767static sorted_pair_node* top_pair(calc_dat* c){
1768  super_clean_top_of_pair_list(c);//yeah, I know, it's odd that I use a different proc here
1769
1770  if(c->pair_top<0) return NULL;
1771  else return (c->apairs[c->pair_top]);
1772}
1773static sorted_pair_node* quick_pop_pair(calc_dat* c){
1774  if(c->pair_top<0) return NULL;
1775  else return (c->apairs[c->pair_top--]);
1776}
1777static BOOLEAN no_pairs(calc_dat* c){
1778  clean_top_of_pair_list(c);
1779  return (c->pair_top==-1);
1780}
1781
1782
1783static void super_clean_top_of_pair_list(calc_dat* c){
1784  while((c->pair_top>=0) && (c->apairs[c->pair_top]->i>=0) && (good_has_t_rep(c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,c))){
1785
1786    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
1787    c->pair_top--;
1788
1789  }
1790}
1791static void clean_top_of_pair_list(calc_dat* c){
1792  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))){
1793
1794    free_sorted_pair_node(c->apairs[c->pair_top],c->r);
1795    c->pair_top--;
1796
1797  }
1798}
1799static BOOLEAN state_is(calc_state state, const int & arg_i, const  int & arg_j, calc_dat* c){
1800  assume(0<=arg_i);
1801  assume(0<=arg_j);
1802  assume(arg_i<c->n);
1803  assume(arg_j<c->n);
1804  if (arg_i==arg_j)
1805  {
1806    return (TRUE);
1807  }
1808  if (arg_i>arg_j)
1809  {
1810    return (c->states[arg_i][arg_j]==state);
1811  }
1812  else return(c->states[arg_j][arg_i]==state);
1813}
1814static void free_sorted_pair_node(sorted_pair_node* s, ring r){
1815  if (s->i>=0)
1816    p_Delete(&s->lcm_of_lm,r);
1817  omfree(s);
1818}
1819static BOOLEAN pair_better(sorted_pair_node* a,sorted_pair_node* b, calc_dat* c){
1820  if (a->deg<b->deg) return TRUE;
1821  if (a->deg>b->deg) return FALSE;
1822
1823  if (a->expected_length<b->expected_length) return TRUE;
1824  if (a->expected_length>b->expected_length) return FALSE;
1825  int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
1826  if (comp==1) return FALSE;
1827  if (-1==comp) return TRUE;
1828  if (a->i<b->i) return TRUE;
1829  if (a->j<b->j) return TRUE;
1830  return FALSE;
1831}
1832
1833static int pair_better_gen(const void* ap,const void* bp){
1834
1835  sorted_pair_node* a=*((sorted_pair_node**)ap);
1836  sorted_pair_node* b=*((sorted_pair_node**)bp);
1837  assume(a->i>a->j);
1838  assume(b->i>b->j);
1839  if (a->deg<b->deg) return -1;
1840  if (a->deg>b->deg) return 1;
1841
1842
1843  if (a->expected_length<b->expected_length) return -1;
1844  if (a->expected_length>b->expected_length) return 1;
1845 int comp=pLmCmp(a->lcm_of_lm, b->lcm_of_lm);
1846 
1847  if (comp==1) return 1;
1848  if (-1==comp) return -1;
1849  if (a->i<b->i) return -1;
1850  if (a->j<b->j) return -1;
1851  return 0;
1852}
1853
1854
1855static poly gcd_of_terms(poly p, ring r){
1856  int max_g_0=0;
1857  assume(p!=NULL);
1858  int i;
1859  poly m=pOne();
1860  poly t;
1861  for (i=pVariables; i; i--)
1862  {
1863      pSetExp(m,i, pGetExp(p,i));
1864      if (max_g_0==0)
1865        if (pGetExp(m,i)>0)
1866          max_g_0=i;
1867  }
1868 
1869  t=p->next;
1870  while (t!=NULL){
1871   
1872    if (max_g_0==0) break;
1873    for (i=pVariables; i; i--)
1874    {
1875      pSetExp(m,i, min(pGetExp(t,i),pGetExp(m,i)));
1876      if (max_g_0==i)
1877        if (pGetExp(m,i)==0)
1878          max_g_0=0;
1879      if ((max_g_0==0) && (pGetExp(m,i)>0)){
1880        max_g_0=i;
1881      }
1882    }
1883                t=t->next;
1884  }
1885        for (i=pVariables;i;i--)
1886        {
1887                if(pGetExp(m,i)>0)
1888                        return m;
1889  }
1890        pDelete(&m);
1891  return NULL;
1892}
1893static BOOLEAN pHasNotCFExtended(poly p1, poly p2, poly m)
1894{
1895
1896  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
1897    return FALSE;
1898  int i = 1;
1899  loop
1900  {
1901    if ((pGetExp(p1, i)-pGetExp(m,i) >0) && (pGetExp(p2, i) -pGetExp(m,i)> 0))   return FALSE;
1902    if (i == pVariables)                                return TRUE;
1903    i++;
1904  }
1905}
1906
1907
1908//for impl reasons may return false if the the normal product criterion matches
1909static BOOLEAN extended_product_criterion(poly p1, poly gcd1, poly p2, poly gcd2, calc_dat* c){
1910  if(gcd1==NULL) return FALSE;
1911        if(gcd2==NULL) return FALSE;
1912        gcd1->next=gcd2; //may ordered incorrect
1913        poly m=gcd_of_terms(gcd1,c->r);
1914        gcd1->next=NULL;
1915        if (m==NULL) return FALSE;
1916
1917        BOOLEAN erg=pHasNotCFExtended(p1,p2,m);
1918        pDelete(&m);
1919        return erg;
1920}
1921static poly kBucketGcd(kBucket* b, ring r)
1922{
1923  int s=0;
1924  int i;
1925  poly m, n;
1926  BOOLEAN initialized=FALSE;
1927  for (i=MAX_BUCKET-1;i>=0;i--)
1928  { 
1929    if (b->buckets[i]!=NULL){
1930      if (!initialized){
1931        m=gcd_of_terms(b->buckets[i],r);
1932        initialized=TRUE;
1933        if (m==NULL) return NULL;
1934      }
1935      else
1936        {
1937          n=gcd_of_terms(b->buckets[i],r);
1938          if (n==NULL) {
1939            pDelete(&m);
1940            return NULL;   
1941          }
1942          n->next=m;
1943          poly t=gcd_of_terms(n,r);
1944          n->next=NULL;
1945          pDelete(&m);
1946          pDelete(&n);
1947          m=t;
1948          if (m==NULL) return NULL;
1949         
1950        }
1951    }
1952  }
1953  return m;
1954}
1955
1956
1957struct find_erg{
1958  poly expand;
1959  int expand_length;
1960  int to_reduce_u;
1961  int to_reduce_l;
1962  int reduce_by;//index of reductor
1963  BOOLEAN fromS;//else from los
1964
1965};
1966static int guess_quality(const red_object & p, calc_dat* c){
1967  //looks only on bucket
1968  if (c->is_char0) return kSBucketLength(p.bucket);
1969  return (bucket_guess(p.bucket));
1970}
1971static int quality_of_pos_in_strat_S(int pos, calc_dat* c){
1972  if (c->is_char0) return c->strat->lenSw[pos];
1973  return c->strat->lenS[pos];
1974}
1975static int quality(poly p, int len, calc_dat* c){
1976  if (c->is_char0) return pSLength(p,len);
1977  return pLength(p);
1978}
1979static void multi_reduction_lls_trick(red_object* los, int losl,calc_dat* c,find_erg & erg){
1980  erg.expand=NULL;
1981  BOOLEAN swap_roles; //from reduce_by, to_reduce_u if fromS
1982  if(erg.fromS){
1983    if(pLmEqual(c->strat->S[erg.reduce_by],los[erg.to_reduce_u].p))
1984    {
1985      int i;
1986      int quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
1987      int best=erg.to_reduce_u+1;
1988      for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
1989        int qc=los[i].guess_quality(c);
1990        if (qc<quality_a){
1991          best=i;
1992          quality_a=qc;
1993        }
1994      }
1995      if(best!=erg.to_reduce_u+1){
1996        red_object h=los[erg.to_reduce_u];
1997        los[erg.to_reduce_u]=los[best];
1998        los[best]=h;
1999        swap_roles=TRUE;
2000      }
2001      else{
2002       
2003        swap_roles=FALSE;
2004      }
2005 
2006    }
2007      else
2008    {
2009      if (erg.to_reduce_u>erg.to_reduce_l){
2010
2011        int i;
2012        int quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
2013        int best=erg.to_reduce_u+1;
2014        for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
2015          int qc=los[i].guess_quality(c);
2016          if (qc<quality_a){
2017            best=i;
2018            quality_a=qc;
2019          }
2020        }
2021        if(best!=erg.to_reduce_u+1){
2022          red_object h=los[erg.to_reduce_l];
2023          los[erg.to_reduce_l]=los[best];
2024          los[best]=h;
2025          erg.reduce_by=erg.to_reduce_l;
2026          erg.fromS=FALSE;
2027          erg.to_reduce_l++;
2028         
2029        }
2030      }
2031      else 
2032      {
2033        assume(erg.to_reduce_u==erg.to_reduce_l);
2034        int quality_a=quality_of_pos_in_strat_S(erg.reduce_by,c);
2035        int qc=los[erg.to_reduce_u].guess_quality(c);
2036        if(qc<quality_a){
2037          BOOLEAN exp=FALSE;
2038          if(qc<=2)
2039            exp=TRUE;
2040          else {
2041            if (qc<quality_a/2)
2042              exp=TRUE;
2043            else
2044              if(erg.reduce_by<c->n/4)
2045                exp=TRUE;
2046          }
2047          if (exp){
2048            poly clear_into;
2049           
2050            kBucketClear(los[erg.to_reduce_u].bucket,&clear_into,&erg.expand_length);
2051            erg.expand=pCopy(clear_into);
2052            kBucketInit(los[erg.to_reduce_u].bucket,clear_into,erg.expand_length);
2053            PrintS("e");
2054           
2055          }
2056        }
2057
2058       
2059      }
2060     
2061      swap_roles=FALSE;
2062      return;
2063      }
2064   
2065  }
2066  else{
2067    if(erg.reduce_by>erg.to_reduce_u){
2068      //then lm(rb)>= lm(tru) so =
2069      assume(erg.reduce_by==erg.to_reduce_u+1);
2070      int best=erg.reduce_by;
2071      int quality_a=los[erg.reduce_by].guess_quality(c);
2072      int i;
2073        for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--){
2074          int qc=los[i].guess_quality(c);
2075          if (qc<quality_a){
2076            best=i;
2077            quality_a=qc;
2078          }
2079        }
2080        if(best!=erg.reduce_by){
2081          red_object h=los[erg.reduce_by];
2082          los[erg.reduce_by]=los[best];
2083          los[best]=h;
2084        }
2085        swap_roles=FALSE;
2086        return;
2087       
2088         
2089    }
2090    else
2091    {
2092      assume(!pLmEqual(los[erg.reduce_by].p,los[erg.to_reduce_l].p));
2093      //further assume, that reduce_by is the above all other polys
2094      //with same leading term
2095      int il=erg.reduce_by;
2096      int quality_a =los[erg.reduce_by].guess_quality(c);
2097      int qc;
2098      while((il>0) && pLmEqual(los[il-1].p,los[il].p)){
2099        il--;
2100        qc=los[il].guess_quality(c);
2101        if (qc<quality_a){
2102          quality_a=qc;
2103          erg.reduce_by=il;
2104        }
2105      }
2106      swap_roles=FALSE;
2107    }
2108 
2109  }
2110  if(swap_roles){
2111    PrintS("b");
2112    poly clear_into;
2113    int dummy_len;
2114    int new_length;
2115    int bp=erg.to_reduce_u;//bucket_positon
2116    //kBucketClear(los[bp].bucket,&clear_into,&new_length);
2117    new_length=los[bp].clear_to_poly();
2118    clear_into=los[bp].p;
2119    poly p=c->strat->S[erg.reduce_by];
2120    int j=erg.reduce_by;
2121    int old_length=c->strat->lenS[j];// in view of S
2122    los[bp].p=p;
2123    kBucketInit(los[bp].bucket,p,old_length);
2124    int qal=quality(clear_into,new_length,c);
2125    int pos_in_c=-1;   
2126    int z;
2127    int new_pos;
2128    new_pos=simple_posInS(c->strat,clear_into,qal,c->is_char0);
2129    assume(new_pos<=j);
2130    for (z=c->n;z;z--)
2131    {
2132      if(p==c->S->m[z-1])
2133      {
2134        pos_in_c=z-1;
2135        break;
2136      }
2137    }
2138    if(pos_in_c>=0)
2139    {
2140      c->S->m[pos_in_c]=clear_into;
2141      c->lengths[pos_in_c]=new_length;
2142      c_S_element_changed_hook(pos_in_c,c);
2143    }
2144    c->strat->S[j]=clear_into;
2145    c->strat->lenS[j]=new_length;
2146    if(c->strat->lenSw)
2147      c->strat->lenS[j]=qal;
2148    if(c->is_char0)
2149    {
2150      pContent(clear_into);
2151      pCleardenom(clear_into);
2152    }
2153  else                     
2154    pNorm(clear_into);
2155    if (new_pos<j){
2156      move_forward_in_S(j,new_pos,c->strat,c->is_char0);
2157      erg.reduce_by=new_pos;
2158    }
2159  }
2160}
2161static void multi_reduction_find(red_object* los, int losl,calc_dat* c,int startf,find_erg & erg){
2162  kStrategy strat=c->strat;
2163
2164  assume(startf<=losl);
2165  assume((startf==losl-1)||(pLmCmp(los[startf].p,los[startf+1].p)==-1));
2166  int i=startf;
2167 
2168  int j;
2169  while(i>=0){
2170    assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)<=0));
2171    j=kFindDivisibleByInS_easy(strat,los[i]);
2172    if(j>=0){
2173     
2174      erg.to_reduce_u=i;
2175      erg.reduce_by=j;
2176      erg.fromS=TRUE;
2177      int i2;
2178      for(i2=i-1;i2>=0;i2--){
2179        if(!pLmEqual(los[i].p,los[i2].p))
2180          break;
2181      }
2182      erg.to_reduce_l=i2+1;
2183      assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2184      return;
2185    }
2186    if (j<0){
2187     
2188      //not reduceable, try to use this for reducing higher terms
2189      int i2;
2190      i2=i;
2191      while((i2>0)&&(pLmEqual(los[i].p,los[i2-1].p)))
2192        i2--;
2193      if(i2!=i){
2194       
2195        erg.to_reduce_u=i-1;
2196        erg.to_reduce_l=i2;
2197        erg.reduce_by=i;
2198        erg.fromS=FALSE;
2199        assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2200        return;
2201      }
2202 
2203      for (i2=i+1;i2<losl;i2++){
2204        if (p_LmShortDivisibleBy(los[i].p,los[i].sev,los[i2].p,~los[i2].sev,
2205                                c->r)){
2206          int i3=i2;
2207          while((i3+1<losl) && (pLmEqual(los[i2].p, los[i3+1].p)))
2208            i3++;
2209          erg.to_reduce_u=i3;
2210          erg.to_reduce_l=i2;
2211          erg.reduce_by=i;
2212          erg.fromS=FALSE;
2213          assume((i==losl-1)||(pLmCmp(los[i].p,los[i+1].p)==-1));
2214          return;
2215        }
2216//      else {assume(!p_LmDivisibleBy(los[i].p, los[i2].p,c->r));}
2217      }
2218
2219      i--;
2220    }
2221  }
2222  erg.reduce_by=-1;//error code
2223  return;
2224}
2225
2226 //  nicht reduzierbare eintraege in ergebnisliste schreiben
2227//   nullen loeschen
2228//   while(finde_groessten leitterm reduzierbar(c,erg)){
2229 
2230static int multi_reduction_clear_zeroes(red_object* los, int  losl, int l, int u)
2231{
2232
2233
2234  int deleted=0;
2235  int  i=l;
2236  int last=-1;
2237  while(i<=u)
2238  {
2239   
2240    if(los[i].p==NULL){
2241      kBucketDestroy(&los[i].bucket);
2242//      delete los[i];//here we assume los are constructed with new
2243      //destroy resources, must be added here   
2244     if (last>=0)
2245     {
2246       memmove(los+(int)(last+1-deleted),los+(last+1),sizeof(red_object)*(i-1-last));
2247     }
2248     last=i;
2249     deleted++;
2250    }
2251    i++;
2252  }
2253  if((last>=0)&&(last!=losl-1))
2254      memmove(los+(int)(last+1-deleted),los+last+1,sizeof(red_object)*(losl-1-last));
2255  return deleted;
2256 
2257}
2258
2259static void sort_region_down(red_object* los, int l, int u, calc_dat* c)
2260{
2261  qsort(los+l,u-l+1,sizeof(red_object),red_object_better_gen);
2262  int i;
2263
2264  for(i=l;i<=u;i++)
2265  {
2266    BOOLEAN moved=FALSE;
2267    int j;
2268    for(j=i;j;j--)
2269    {
2270      if(pLmCmp(los[j].p,los[j-1].p)==-1){
2271        red_object h=los[j];
2272        los[j]=los[j-1];
2273        los[j-1]=h;
2274        moved=TRUE;
2275      }
2276      else break;
2277    }
2278    if(!moved) return;
2279  }
2280}
2281
2282//assume that los is ordered ascending by leading term, all non zero
2283static void multi_reduction(red_object* los, int & losl, calc_dat* c)
2284{
2285 
2286  //initialize;
2287  assume(c->strat->sl>=0);
2288  assume(losl>0);
2289  int i;
2290  for(i=0;i<losl;i++){
2291    los[i].sev=pGetShortExpVector(los[i].p);
2292//SetShortExpVector();
2293    los[i].p=kBucketGetLm(los[i].bucket);
2294  }
2295
2296  kStrategy strat=c->strat;
2297  int curr_pos=losl-1;
2298
2299
2300//  nicht reduzierbare einträge in ergebnisliste schreiben
2301  // nullen loeschen
2302  while(curr_pos>=0){
2303    find_erg erg;
2304    multi_reduction_find(los, losl,c,curr_pos,erg);//last argument should be curr_pos
2305   //  PrintS("\n erg:\n");
2306//     Print("upper:%i\n",erg.to_reduce_u);
2307//     Print("lower:%i\n",erg.to_reduce_l);
2308//     Print("reduce_by:%i\n",erg.reduce_by);
2309//     Print("fromS:%i\n",erg.fromS);
2310    if(erg.reduce_by<0) break;
2311    multi_reduction_lls_trick(los,losl,c,erg);
2312    //erweitern? muß noch implementiert werden
2313    int i;
2314    int len;
2315    poly reductor;
2316    if(erg.fromS){
2317      reductor=strat->S[erg.reduce_by];
2318      len=strat->lenS[erg.reduce_by];
2319     
2320    }
2321    else 
2322    {
2323      //bucket aufloesen reduzieren, neu füllen
2324     
2325   
2326      int bn=kBucketCanonicalize(los[erg.reduce_by].bucket);
2327      reductor=los[erg.reduce_by].bucket->buckets[bn];
2328      len=los[erg.reduce_by].bucket->buckets_length[bn];
2329      if(c->is_char0)
2330        pContent(reductor);
2331 
2332    }
2333    for(i=erg.to_reduce_l;i<=erg.to_reduce_u;i++)
2334    {
2335   
2336      assume((erg.fromS)||(i!=erg.reduce_by));
2337      assume(reductor!=NULL);
2338       number coef=kBucketPolyRed(los[i].bucket,reductor,
2339                                  len,
2340                                  strat->kNoether);
2341       nDelete(&coef);
2342       los[i].p = kBucketGetLm(los[i].bucket);
2343       if(los[i].p!=NULL)
2344         if((i>0)&&(los[i-1].p!=NULL)&&(pLmEqual(los[i-1].p,los[i].p)))
2345             los[i].sev=los[i-1].sev;
2346         else
2347           los[i].sev=pGetShortExpVector(los[i].p);
2348       //better would be first sorting before sev
2349    }
2350 
2351                 
2352    int deleted=multi_reduction_clear_zeroes(los, losl, erg.to_reduce_l, erg.to_reduce_u);
2353    curr_pos=erg.to_reduce_u;
2354    losl -= deleted;
2355    curr_pos -= deleted;
2356
2357    //Print("deleted %i \n",deleted);
2358    sort_region_down(los, erg.to_reduce_l, erg.to_reduce_u-deleted, c);
2359//   sort_region_down(los, 0, losl-1, c);
2360    //  qsort(los,losl,sizeof(red_object),red_object_better_gen);
2361    if(erg.expand)
2362      add_to_reductors(c,erg.expand,erg.expand_length);
2363  }
2364  return;
2365}
2366void red_object::flatten(){
2367  if (sum!=NULL)
2368  {
2369
2370 
2371    if(kBucketGetLm(sum->ac->bucket)!=NULL){
2372      number mult_my=n_Mult(sum->c_my,sum->ac->multiplied,currRing);
2373      poly add_this;
2374      if(!nIsOne(mult_my))
2375        kBucket_Mult_n(bucket,mult_my);
2376      int len;
2377      poly clear_into;
2378      kBucketClear(sum->ac->bucket,&clear_into,&len);
2379      if(sum->ac->counter>1){
2380        add_this=pCopy(clear_into);
2381        kBucketInit(bucket,clear_into,len);
2382      }
2383      else
2384        add_this=clear_into;
2385      pMult_nn(add_this, sum->c_ac);
2386      nDelete(&sum->c_ac);
2387      nDelete(&sum->c_my);
2388      nDelete(&mult_my);
2389      delete sum;
2390      kBucket_Add_q(bucket,add_this, &len);
2391      sum->ac->decrease_counter();
2392     
2393    }
2394  }
2395}
2396void red_object::validate(){
2397  if(sum!=NULL)
2398  {
2399    poly lm=kBucketGetLm(bucket);
2400    poly lm_ac=kBucketGetLm(sum->ac->bucket);
2401    if ((lm_ac==NULL)||((lm!=NULL) && (pLmCmp(lm,lm_ac)!=-1))){
2402      flatten();
2403      p=kBucketGetLm(bucket);
2404    } 
2405    else
2406    {
2407      p=lm_ac;
2408    }
2409   
2410  }
2411  else
2412    p=kBucketGetLm(bucket);
2413
2414}
2415int red_object::clear_to_poly(){
2416  flatten();
2417  int l;
2418  kBucketClear(bucket,&p,&l);
2419  return l;
2420}
2421void red_object::reduction_step(int reduction_id, poly reductor_full, int full_len, poly reductor_part, reduction_accumulator* join_to, calc_dat* c)
2422{
2423  //we have to add support later for building new sums at this points, this involves a change in the interface
2424  if(this->sum==NULL)
2425    kBucketPolyRed(this->bucket,reductor_full,
2426                   full_len,
2427                   c->strat->kNoether);
2428  else 
2429  {
2430    assume(sum->ac!=NULL);
2431    if(sum->ac->last_reduction_id!=reduction_id){
2432     
2433
2434
2435     
2436     
2437      number n1=kBucketPolyRed(sum->ac->bucket,reductor_full, full_len, c->strat->kNoether);
2438      number n2=nMult(n1,sum->ac->multiplied);
2439      nDelete(&sum->ac->multiplied);
2440      nDelete(&n1);
2441      sum->ac->multiplied=n2;
2442    }
2443      //reduce and adjust multiplied
2444      sum->ac->last_reduction_id=reduction_id;
2445     
2446  }
2447   
2448     
2449}
2450 
2451
2452void red_object::adjust_coefs(number c_r, number c_ac_r){
2453  assume(this->sum!=NULL);
2454  number n1=nMult(sum->c_my, c_ac_r);
2455  number n2=nMult(sum->c_ac,c_r);
2456  nDelete(&sum->c_my);
2457  nDelete(&sum->c_ac);
2458 
2459  int ct = ksCheckCoeff(&n1, &n2);
2460  sum->c_my=n1;
2461  sum->c_ac=nNeg(n2);
2462  nDelete(&n2);
2463 
2464
2465}
2466int red_object::guess_quality(calc_dat* c){
2467    //works at the moment only for lenvar 1, because in different
2468    //case, you have to look on coefs
2469    int s=0;
2470    if (c->is_char0)
2471      s=kSBucketLength(bucket);
2472    else 
2473      s=bucket_guess(bucket);
2474    if (sum!=NULL){
2475      if (c->is_char0)
2476      s+=kSBucketLength(sum->ac->bucket);
2477    else 
2478      s+=bucket_guess(sum->ac->bucket);
2479    }
2480    return s;
2481}
Note: See TracBrowser for help on using the repository browser.