source: git/Singular/kbuckets.cc @ 512a2b

spielwiese
Last change on this file since 512a2b was 512a2b, checked in by Olaf Bachmann <obachman@…>, 24 years ago
p_polys.h git-svn-id: file:///usr/local/Singular/svn/trunk@4606 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 22.1 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id: kbuckets.cc,v 1.17 2000-09-18 09:19:06 obachman Exp $ */
5
6#include "mod2.h"
7#include "tok.h"
8#include "structs.h"
9#include "omalloc.h"
10#include "polys.h"
11#include "febase.h"
12#include "kbuckets.h"
13#include "numbers.h"
14#include "p_Procs.h"
15
16static omBin kBucket_bin = omGetSpecBin(sizeof(kBucket));
17
18//////////////////////////////////////////////////////////////////////////
19///
20/// Some internal stuff
21///
22
23// returns ceil(log_4(l))
24inline unsigned int pLogLength(unsigned int l)
25{
26  unsigned int i = 0;
27
28  if (l == 0) return 0;
29  l--;
30  while ((l = (l >> 2))) i++;
31  return i+1;
32}
33
34// returns ceil(log_4(pLength(p)))
35inline unsigned int pLogLength(poly p)
36{
37  return pLogLength((unsigned int) pLength(p));
38}
39
40#if defined(PDEBUG) && ! defined(HAVE_PSEUDO_BUCKETS)
41
42void kbTest(kBucket_pt bucket, int i)
43{
44  if (_p_Test(bucket->buckets[i], bucket->bucket_ring, PDEBUG))
45  {
46    if (bucket->buckets_length[i] != pLength(bucket->buckets[i]))
47    {
48      dReportError("Bucket %d lengths difference should:%d has:%d",
49                   i, bucket->buckets_length[i], pLength(bucket->buckets[i]));
50    }
51    else if (i > 0 && (int) pLogLength(bucket->buckets_length[i]) > i)
52    {
53      dReportError("Bucket %d too long %d",
54                   i, bucket->buckets_length[i]);
55    }
56  }
57  if (i==0 && bucket->buckets_length[0] > 1)
58  {
59    dReportError("Bucket 0 too long");
60  }
61}
62
63
64void kbTests(kBucket_pt bucket)
65{
66  int i;
67  poly lm = bucket->buckets[0];
68
69  omCheckAddrBin(bucket, kBucket_bin);
70  kbTest(bucket, 0);
71  for (i=1; i<= (int) bucket->buckets_used; i++)
72  {
73    kbTest(bucket, i);
74    if (lm != NULL &&  bucket->buckets[i] != NULL
75        && pLmCmp(lm, bucket->buckets[i]) != 1)
76    {
77      dReportError("Bucket %d larger than lm", i);
78    }
79  }
80
81  for (; i<=MAX_BUCKET; i++)
82  {
83    if (bucket->buckets[i] != NULL || bucket->buckets_length[i] != 0)
84    {
85      dReportError("Bucket %d not zero", i);
86    }
87  }
88}
89
90#else
91
92#define kbTests(bucket) (NULL)
93#define kbTest(bucket, i) (NULL)
94
95#endif // PDEBUG
96
97//////////////////////////////////////////////////////////////////////////
98///
99/// Creation/Destruction of buckets
100///
101
102kBucket_pt kBucketCreate(ring bucket_ring)
103{
104  assume(bucket_ring != NULL);
105  kBucket_pt bucket = (kBucket_pt) omAlloc0Bin(kBucket_bin);
106  bucket->bucket_ring = bucket_ring;
107  return bucket;
108}
109
110void kBucketDestroy(kBucket_pt *bucket)
111{
112  omFreeBin(*bucket, kBucket_bin);
113  *bucket = NULL;
114}
115
116
117/////////////////////////////////////////////////////////////////////////////
118// Convertion from/to Bpolys
119//
120#ifndef HAVE_PSEUDO_BUCKETS
121
122inline void kBucketAdjustBucketsUsed(kBucket_pt bucket)
123{
124  while ( bucket->buckets_used > 0 &&
125          bucket->buckets[bucket->buckets_used] == NULL)
126    (bucket->buckets_used)--;
127}
128
129inline void kBucketMergeLm(kBucket_pt bucket)
130{
131  if (bucket->buckets[0] != NULL)
132  {
133    poly lm = bucket->buckets[0];
134    int i = 1;
135    int l = 4;
136    while ( bucket->buckets_length[i] >= l)
137    {
138      i++;
139      l = l << 2;
140    }
141    pNext(lm) = bucket->buckets[i];
142    bucket->buckets[i] = lm;
143    bucket->buckets_length[i]++;
144    assume(i <= bucket->buckets_used+1);
145    if (i > bucket->buckets_used)  bucket->buckets_used = i;
146    bucket->buckets[0] = NULL;
147    bucket->buckets_length[0] = 0;
148  }
149}
150
151
152static BOOLEAN kBucketIsCleared(kBucket_pt bucket)
153{
154  int i;
155
156  for (i = 0;i<=MAX_BUCKET;i++)
157  {
158    if (bucket->buckets[i] != NULL) return FALSE;
159    if (bucket->buckets_length[i] != 0) return FALSE;
160  }
161  return TRUE;
162}
163
164void kBucketInit(kBucket_pt bucket, poly lm, int length)
165{
166  int i;
167  assume(bucket != NULL);
168  assume(length <= 0 || length == pLength(lm));
169  assume(kBucketIsCleared(bucket));
170
171  if (length <= 0)
172    length = pLength(lm);
173
174  bucket->buckets[0] = lm;
175  if (length > 1)
176  {
177    unsigned int i = pLogLength(length-1);
178    bucket->buckets[i] = pNext(lm);
179    pNext(lm) = NULL;
180    bucket->buckets_length[i] = length-1;
181    bucket->buckets_used = i;
182    bucket->buckets_length[0] = 1;
183  }
184  else
185  {
186    bucket->buckets_used = 0;
187    bucket->buckets_length[0] = 0;
188  }
189}
190
191static int kBucketCanonicalize(kBucket_pt bucket);
192
193static void kBucketSetLm(kBucket_pt bucket)
194{
195  int j = 0;
196  poly lt;
197  BOOLEAN zero = FALSE;
198  assume(bucket->buckets[0] == NULL && bucket->buckets_length[0] == 0);
199
200  do
201  {
202    j = 0;
203    for (int i = 1; i<=bucket->buckets_used; i++)
204    {
205      if (bucket->buckets[i] != NULL)
206      {
207        int comp = (j == 0 ? 1 :
208                    p_LmCmp(bucket->buckets[i], bucket->buckets[j], bucket->bucket_ring));
209
210        if (comp > 0)
211        {
212          if (j > 0 &&
213              bucket->buckets[j] != NULL &&
214              nIsZero(pGetCoeff(bucket->buckets[j])))
215          {
216            p_DeleteLm(&(bucket->buckets[j]), bucket->bucket_ring);
217            (bucket->buckets_length[j])--;
218          }
219          j = i;
220        }
221        else if (comp == 0)
222        {
223          number tn = pGetCoeff(bucket->buckets[j]);
224          pSetCoeff0(bucket->buckets[j],
225                     nAdd(pGetCoeff(bucket->buckets[i]), tn));
226          nDelete(&tn);
227          p_DeleteLm(&(bucket->buckets[i]), bucket->bucket_ring);
228          (bucket->buckets_length[i])--;
229        }
230      }
231    }
232    if (j > 0 && nIsZero(pGetCoeff(bucket->buckets[j])))
233    {
234      p_DeleteLm(&(bucket->buckets[j]), bucket->bucket_ring);
235      (bucket->buckets_length[j])--;
236      j = -1;
237    }
238  }
239  while (j < 0);
240
241  if (j == 0)
242  {
243    return;
244  }
245
246  assume(bucket->buckets[j] != NULL);
247  lt = bucket->buckets[j];
248  bucket->buckets[j] = pNext(lt);
249  bucket->buckets_length[j]--;
250  pNext(lt) = NULL;
251  bucket->buckets[0] = lt;
252  bucket->buckets_length[0] = 1;
253
254  kBucketAdjustBucketsUsed(bucket);
255}
256
257const poly kBucketGetLm(kBucket_pt bucket)
258{
259  if (bucket->buckets[0] == NULL)
260    kBucketSetLm(bucket);
261  return bucket->buckets[0];
262}
263
264poly kBucketExtractLm(kBucket_pt bucket)
265{
266  poly lm;
267  if (bucket->buckets[0] == NULL) kBucketSetLm(bucket);
268  lm = bucket->buckets[0];
269  bucket->buckets[0] = NULL;
270  bucket->buckets_length[0] = 0;
271  return lm;
272}
273
274static int kBucketCanonicalize(kBucket_pt bucket)
275{
276  poly p = bucket->buckets[1];
277  poly lm;
278  int pl = bucket->buckets_length[1], i;
279  bucket->buckets[1] = NULL;
280  bucket->buckets_length[1] = 0;
281
282
283  for (i=2; i<=bucket->buckets_used; i++)
284  {
285    p = p_Add_q(p, bucket->buckets[i],
286                 pl, bucket->buckets_length[i], bucket->bucket_ring);
287    bucket->buckets[i] = NULL;
288    bucket->buckets_length[i] = 0;
289  }
290
291  lm = bucket->buckets[0];
292  if (lm != NULL)
293  {
294    pNext(lm) = p;
295    p = lm;
296    pl++;
297    bucket->buckets[0] = NULL;
298    bucket->buckets_length[0] = 0;
299  }
300  if (pl > 0)
301  {
302    i = pLogLength(pl);
303    bucket->buckets[i] = p;
304    bucket->buckets_length[i] = pl;
305  }
306  else
307  {
308    i = 0;
309  }
310  bucket->buckets_used = i;
311  assume(pLength(p) == (int) pl);
312  return i;
313}
314
315void kBucketClear(kBucket_pt bucket, poly *p, int *length)
316{
317  int i = kBucketCanonicalize(bucket);
318  if (i > 0)
319  {
320    *p = bucket->buckets[i];
321    *length = bucket->buckets_length[i];
322    bucket->buckets[i] = NULL;
323    bucket->buckets_length[i] = 0;
324    bucket->buckets_used = 0;
325  }
326  else
327  {
328    *p = NULL;
329    *length = 0;
330  }
331}
332
333#else // HAVE_PSEUDO_BUCKETS
334
335void kBucketInit(kBucket_pt bucket, poly lm, int length)
336{
337  int i;
338
339  assume(bucket != NULL);
340  assume(length <= 0 || length == pLength(lm));
341
342  bucket->p = lm;
343  if (length <= 0) bucket->l = pLength(lm);
344  else bucket->l = length;
345
346}
347
348const poly kBucketGetLm(kBucket_pt bucket)
349{
350  return bucket->p;
351}
352
353poly kBucketExtractLm(kBucket_pt bucket)
354{
355  poly lm = bucket->p;
356  assume(pLength(bucket->p) == bucket->l);
357  pIter(bucket->p);
358  (bucket->l)--;
359  pNext(lm) = NULL;
360  return lm;
361}
362
363void kBucketClear(kBucket_pt bucket, poly *p, int *length)
364{
365  assume(pLength(bucket->p) == bucket->l);
366  *p = bucket->p;
367  *length = bucket->l;
368  bucket->p = NULL;
369  bucket->l = 0;
370}
371
372#endif // ! HAVE_PSEUDO_BUCKETS
373
374//////////////////////////////////////////////////////////////////////////
375///
376/// Multiply Bucket by number ,i.e. Bpoly == n*Bpoly
377///
378void kBucket_Mult_n(kBucket_pt bucket, number n)
379{
380#ifndef HAVE_PSEUDO_BUCKETS
381  int i;
382
383  for (i=0; i<= bucket->buckets_used; i++)
384    if (bucket->buckets[i] != NULL)
385      bucket->buckets[i] = p_Mult_nn(bucket->buckets[i], n, bucket->bucket_ring);
386#else
387  bucket->p = p_Mult_nn(bucket->p, n, bucket->bucket_ring);
388#endif
389}
390
391//////////////////////////////////////////////////////////////////////////
392///
393/// Bpoly == Bpoly - m*p; where m is a monom
394/// Does not destroy p and m
395/// assume (*l <= 0 || pLength(p) == *l)
396void kBucket_Minus_m_Mult_p(kBucket_pt bucket, poly m, poly p, int *l,
397                            poly spNoether)
398{
399  assume(*l <= 0 || pLength(p) == *l);
400  int i, l1;
401  poly p1 = p;
402
403  if (*l <= 0)
404  {
405    l1 = pLength(p1);
406    *l = l1;
407  }
408  else
409    l1 = *l;
410
411  if (m == NULL || p == NULL) return;
412
413#ifndef HAVE_PSEUDO_BUCKETS
414  kBucketMergeLm(bucket);
415  kbTests(bucket);
416  i = pLogLength(l1);
417
418  if (i <= bucket->buckets_used && bucket->buckets[i] != NULL)
419  {
420    p1 = p_Minus_mm_Mult_qq(bucket->buckets[i], m, p1,
421                            bucket->buckets_length[i], l1, 
422                            spNoether,
423                            bucket->bucket_ring);
424    l1 = bucket->buckets_length[i];
425    bucket->buckets[i] = NULL;
426    bucket->buckets_length[i] = 0;
427    i = pLogLength(l1);
428    while (bucket->buckets[i] != NULL)
429    {
430      p1 = p_Add_q(p1, bucket->buckets[i],
431                     l1, bucket->buckets_length[i], bucket->bucket_ring);
432      bucket->buckets[i] = NULL;
433      bucket->buckets_length[i] = 0;
434      i = pLogLength(l1);
435    }
436  }
437  else
438  {
439    pSetCoeff0(m, nNeg(pGetCoeff(m)));
440    p1 = bucket->bucket_ring->p_Procs->pp_Mult_mm(p1, m, spNoether, bucket->bucket_ring);
441    pSetCoeff0(m, nNeg(pGetCoeff(m)));
442  }
443
444  bucket->buckets[i] = p1;
445  bucket->buckets_length[i]=l1;
446  if (i >= bucket->buckets_used)
447    bucket->buckets_used = i;
448  else
449    kBucketAdjustBucketsUsed(bucket);
450#else // HAVE_PSEUDO_BUCKETS
451  bucket->p = p_Minus_mm_Mult_qq(bucket->p, m,  p,
452                               bucket->l, l1,
453                               spNoether, bucket->bucket_ring);
454#endif
455  kbTests(bucket);
456}
457
458/////////////////////////////////////////////////////////////////////////////
459//
460// Extract all monomials from bucket with component comp
461// Return as a polynomial *p with length *l
462// In other words, afterwards
463// Bpoly == Bpoly - (poly consisting of all monomials with component comp)
464// and components of monomials of *p are all 0
465//
466void kBucketTakeOutComp(kBucket_pt bucket,
467                        Exponent_t comp,
468                        poly *r_p, int *l)
469{
470  poly p = NULL, q;
471  int i, lp = 0, lq;
472
473#ifndef HAVE_PSEUDO_BUCKETS
474  kBucketMergeLm(bucket);
475  for (i=1; i<=bucket->buckets_used; i++)
476  {
477    if (bucket->buckets[i] != NULL)
478    {
479      pTakeOutComp(&(bucket->buckets[i]), comp, &q, &lq);
480      if (q != NULL)
481      {
482        assume(pLength(q) == lq);
483        bucket->buckets_length[i] -= lq;
484        assume(pLength(bucket->buckets[i]) == bucket->buckets_length[i]);
485        p = p_Add_q(p, q, lp, lq, bucket->bucket_ring);
486      }
487    }
488  }
489  kBucketAdjustBucketsUsed(bucket);
490#else
491  pTakeOutComp(&(bucket->p), comp, &p, &lp);
492  (bucket->l) -= lp;
493#endif
494  *r_p = p;
495  *l = lp;
496
497  kbTests(bucket);
498}
499
500void kBucketDecrOrdTakeOutComp(kBucket_pt bucket,
501                               Exponent_t comp, Order_t order,
502                               poly *r_p, int *l)
503{
504  poly p = NULL, q;
505  int i, lp = 0, lq;
506
507#ifndef HAVE_PSEUDO_BUCKETS
508  kBucketMergeLm(bucket);
509  for (i=1; i<=bucket->buckets_used; i++)
510  {
511    if (bucket->buckets[i] != NULL)
512    {
513      pDecrOrdTakeOutComp(&(bucket->buckets[i]), comp, order, &q, &lq);
514      if (q != NULL)
515      {
516        bucket->buckets_length[i] -= lq;
517        p = p_Add_q(p, q, lp, lq, bucket->bucket_ring);
518      }
519    }
520  }
521  kBucketAdjustBucketsUsed(bucket);
522#else
523  pDecrOrdTakeOutComp(&(bucket->p), comp, order, &p, &lp);
524  (bucket->l) -= lp;
525#endif
526
527  *r_p = p;
528  *l = lp;
529}
530
531/////////////////////////////////////////////////////////////////////////////
532// Reduction of Bpoly with a given poly
533//
534
535extern int ksCheckCoeff(number *a, number *b);
536
537number kBucketPolyRed(kBucket_pt bucket,
538                      poly p1, int l1,
539                      poly spNoether)
540{
541  assume(p1 != NULL &&
542         pDivisibleBy(p1,  kBucketGetLm(bucket)));
543  assume(pLength(p1) == (int) l1);
544
545  poly a1 = pNext(p1), lm = kBucketExtractLm(bucket);
546  BOOLEAN reset_vec=FALSE;
547  number rn;
548
549  if(a1==NULL)
550  {
551    p_DeleteLm(&lm, bucket->bucket_ring);
552    return nInit(1);
553  }
554
555  if (! nIsOne(pGetCoeff(p1)))
556  {
557    number an = pGetCoeff(p1), bn = pGetCoeff(lm);
558    int ct = ksCheckCoeff(&an, &bn);
559    pSetCoeff(lm, bn);
560    if ((ct == 0) || (ct == 2)) kBucket_Mult_n(bucket, an);
561    rn = an;
562  }
563  else
564  {
565    rn = nInit(1);
566  }
567
568  if (pGetComp(p1) != pGetComp(lm))
569  {
570    pSetCompP(a1, pGetComp(lm));
571    reset_vec = TRUE;
572    pSetComp(lm, pGetComp(p1));
573    pSetm(lm);
574  }
575
576  pExpVectorSub(lm,p1);
577  l1--;
578
579  kBucket_Minus_m_Mult_p(bucket, lm, a1, &l1, spNoether);
580
581  p_DeleteLm(&lm, bucket->bucket_ring);
582  if (reset_vec) pSetCompP(a1, 0);
583  kbTests(bucket);
584  return rn;
585}
586
587/////////////////////////////////////////////////////////////////////////////
588// Reduction of a poly using buckets
589//
590
591#ifdef KB_USE_BUCKETS
592/*2
593*  reduction procedure for the homogeneous case
594*  and the case of a degree-ordering
595*/
596void kBucketRedHomog (LObject* h,kStrategy strat)
597{
598  if (strat->tl<0)
599  {
600    return;
601  }
602  int j = 0;
603  int k = 0;
604  poly lm;
605#ifdef KDEBUG
606  if (TEST_OPT_DEBUG)
607  {
608    PrintS("red:");
609    wrp(h->p);
610    PrintS(" ");
611  }
612#endif
613
614  kBucketInit(&(strat->bucket), h->p, h->length);
615
616  lm = kBucketGetLm(&(strat->bucket));
617
618  for(;;)
619  {
620#ifdef KDEBUG
621    if (TEST_OPT_DEBUG) Print("%d",j);
622#endif
623    j = 0;
624
625#ifdef KB_HAVE_SHORT_EVECTORS
626    unsigned long ev = ~ pGetShortExpVector(lm);
627#endif
628
629    if (strat->ak)
630    {
631      while(j <= strat->tl)
632      {
633        if (pLmDivisibleBy(strat->T[j].p,lm)) goto Found;
634        j++;
635      }
636    }
637    else
638    {
639      while(j <= strat->tl)
640      {
641#ifdef KB_HAVE_SHORT_EVECTORS
642        if ((strat->T[j].sev & ev) || ! pLmDivisibleByNoComp(strat->T[j].p,lm))
643          j++;
644        else
645          goto Found;
646#else
647        if (pLmDivisibleByNoComp(strat->T[j].p,lm)) goto Found;
648        j++;
649#endif
650      }
651    }
652
653    assume(j > strat->tl);
654    kBucketClear(&(strat->bucket), &(h->p), &(h->length));
655    if (TEST_OPT_INTSTRATEGY) pCleardenom((*h).p);
656    assume(pLength(h->p) == h->length);
657    return;
658
659    Found:
660    // now we found one to reduce with
661    assume(pDivisibleBy(strat->T[j].p,lm));
662
663#ifdef KDEBUG
664    if (TEST_OPT_DEBUG)
665    {
666        wrp(strat->T[j].p);
667    }
668#endif
669    assume(strat->T[j].length <= 0 ||
670           strat->T[j].length == pLength(strat->T[j].p));
671      /*- compute the s-polynomial -*/
672    if (strat->T[j].length <= 0)
673      strat->T[j].length = pLength(strat->T[j].p);
674
675    number up = kBucketPolyRed(&(strat->bucket),
676                               strat->T[j].p, strat->T[j].length,
677                               strat->kNoether);
678    nDelete(&up);
679    lm = kBucketGetLm(&(strat->bucket));
680
681    if (lm == NULL)
682    {
683      h->p = NULL;
684      h->length = 0;
685#ifdef KDEBUG
686      if (TEST_OPT_DEBUG) PrintS(" to 0\n");
687#endif
688      if (h->lcm!=NULL) pLmFree((*h).lcm);
689#ifdef KDEBUG
690      (*h).lcm=NULL;
691#endif
692      return;
693    }
694  }
695}
696
697
698/*2
699*  reduction procedure for the sugar-strategy (honey)
700* reduces h with elements from T choosing first possible
701* element in T with respect to the given ecart
702*/
703void kBucketRedHoney (LObject*  h,kStrategy strat)
704{
705  if (strat->tl<0)
706  {
707    return;
708  }
709
710  poly pi;
711  int i,j,at,reddeg,d,pass,ei, li;
712
713  pass = j = 0;
714  d = reddeg = pFDeg((*h).p)+(*h).ecart;
715#ifdef KDEBUG
716  if (TEST_OPT_DEBUG)
717  {
718    PrintS("red:");
719    wrp((*h).p);
720  }
721#endif
722  kTest(strat);
723
724  if (strat->fromT)
725  {
726    HeapPoly b_hp(h->heap), r_hp(h->p, h->length);
727    pCopyHeapPoly(&b_hp, &r_hp, FALSE);
728    kBucketInit(&(strat->bucket), b_hp.p, b_hp.length, h->heap);
729  }
730  else
731    kBucketInit(&(strat->bucket),  h->p, h->length, h->heap);
732
733  poly lm = kBucketGetLm(&(strat->bucket));
734
735  strat->fromT=FALSE;
736
737  while (1)
738  {
739    j = 0;
740#ifdef KB_HAVE_SHORT_EVECTORS
741    unsigned long ev = ~ pGetShortExpVector(lm);
742#endif
743
744    if (strat->ak)
745    {
746      while(j <= strat->tl)
747      {
748        if (pLmDivisibleBy(strat->T[j].p,lm)) goto Found;
749        j++;
750      }
751    }
752    else
753    {
754      while(j <= strat->tl)
755      {
756#ifdef KB_HAVE_SHORT_EVECTORS
757        if ((strat->T[j].sev & ev) || ! pLmDivisibleByNoComp(strat->T[j].p,lm))
758          j++;
759        else
760          goto Found;
761#else
762        if (pLmDivisibleByNoComp(strat->T[j].p,lm)) goto Found;
763        j++;
764#endif
765      }
766    }
767
768    assume(j > strat->tl);
769#ifdef KDEBUG
770    if (TEST_OPT_DEBUG) PrintLn();
771#endif
772    kBucketClear(&(strat->bucket), &(h->p), &(h->length));
773    if (TEST_OPT_INTSTRATEGY) pCleardenom(h->p);// also does a pContent
774    return;
775
776    Found:
777    assume(pDivisibleBy(strat->T[j].p,lm));
778#ifdef KDEBUG
779    if (TEST_OPT_DEBUG) Print(" T[%d]",j);
780#endif
781    pi = strat->T[j].p;
782    ei = strat->T[j].ecart;
783    li = strat->T[j].length;
784    if (li == 0)
785    {
786      strat->T[j].length = pLength(pi);
787      li =  strat->T[j].length;
788    }
789
790
791    /*
792     * the polynomial to reduce with (up to the moment) is;
793     * pi with ecart ei
794     */
795    i = j;
796    loop
797      {
798        /*- takes the first possible with respect to ecart -*/
799        i++;
800        if (i > strat->tl)
801          break;
802        if ((!TEST_OPT_REDBEST) && (ei <= (*h).ecart))
803          break;
804        if ((strat->T[i].ecart < ei) && pLmDivisibleBy(strat->T[i].p,lm))
805        {
806#ifdef KDEBUG
807          if (TEST_OPT_DEBUG) Print(" T[%d]",i);
808#endif
809          /*
810           * the polynomial to reduce with is now;
811           */
812          pi = strat->T[i].p;
813          ei = strat->T[i].ecart;
814          li = strat->T[i].length;
815          if (li == 0)
816          {
817            strat->T[i].length = pLength(pi);
818            li =  strat->T[i].length;
819          }
820        }
821      }
822
823    /*
824     * end of search: have to reduce with pi
825     */
826    if ((pass!=0) && (ei > (*h).ecart))
827    {
828      /*
829       * It is not possible to reduce h with smaller ecart;
830       * if possible h goes to the lazy-set L,i.e
831       * if its position in L would be not the last one
832       */
833      if (strat->Ll >= 0) /* L is not empty */
834      {
835        h->p = lm;
836        at = strat->posInL(strat->L,strat->Ll,*h,strat);
837        if(at <= strat->Ll)
838          /*- h will not become the next element to reduce -*/
839        {
840          kBucketClear(&(strat->bucket), &(h->p), &(h->length));
841          enterL(&strat->L,&strat->Ll,&strat->Lmax,*h,at);
842#ifdef KDEBUG
843          if (TEST_OPT_DEBUG) Print(" ecart too big: -> L%d\n",at);
844#endif
845          h->p = NULL;
846          h->length = 0;
847          h->heap = NULL;
848          return;
849        }
850      }
851    }
852#ifdef KDEBUG
853    if (TEST_OPT_DEBUG)
854    {
855      PrintS("\nwith ");
856      wrp(pi);
857    }
858#endif
859
860    // Poly Reduction
861    number up = kBucketPolyRed(&(strat->bucket),
862                               pi, li,
863                               strat->kNoether);
864    nDelete(&up);
865    lm = kBucketGetLm(&(strat->bucket));
866
867#ifdef KDEBUG
868    if (TEST_OPT_DEBUG)
869    {
870      PrintS(" to ");
871      wrp((*h).p);
872      PrintLn();
873    }
874#endif
875    if (lm == NULL)
876    {
877      if (h->lcm!=NULL) pLmFree((*h).lcm);
878#ifdef KDEBUG
879      (*h).lcm=NULL;
880#endif
881      h->p = NULL;
882      h->length = 0;
883      return;
884    }
885    /* compute the ecart */
886    if (ei <= (*h).ecart)
887      (*h).ecart = d-pFDeg(lm);
888    else
889      (*h).ecart = d-pFDeg(lm)+ei-(*h).ecart;
890    /*
891     * try to reduce the s-polynomial h
892     *test first whether h should go to the lazyset L
893     *-if the degree jumps
894     *-if the number of pre-defined reductions jumps
895     */
896    pass++;
897    d = pFDeg(lm)+(*h).ecart;
898    if ((strat->Ll >= 0) && ((d > reddeg) || (pass > strat->LazyPass)))
899    {
900      h->p = lm;
901      at = strat->posInL(strat->L,strat->Ll,*h,strat);
902      if (at <= strat->Ll)
903      {
904        kBucketClear(&(strat->bucket), &(h->p), &(h->length));
905        assume(pLength(h->p) == h->length);
906        /*test if h is already standardbasis element*/
907        i=strat->sl+1;
908        do
909        {
910          i--;
911          if (i<0)
912          {
913            // we actually should do a pCleardenom(h->p) here
914            return;
915          }
916        } while (!pLmDivisibleBy(strat->S[i], h->p));
917
918        // enter in Lazyset and return
919        enterL(&strat->L,&strat->Ll,&strat->Lmax,*h,at);
920#ifdef KDEBUG
921        if (TEST_OPT_DEBUG)
922          Print(" degree jumped: -> L%d\n",at);
923#endif
924        h->p = NULL;
925        h->length = 0;
926        h->heap = NULL;
927        return;
928      }
929    }
930    else if (TEST_OPT_PROT && (strat->Ll < 0) && (d > reddeg))
931    {
932      reddeg = d;
933      Print(".%d",d); mflush();
934    }
935  }
936}
937
938/*2
939*compute the normalform of the tail p->next of p
940*with respect to S
941*/
942void kBucketRedTail(LObject *hp, int pos, kStrategy strat)
943{
944  poly h, hn, redstart = NULL;
945  int length = 1, j;
946
947  if (strat->noTailReduction) return;
948  h = hp->p;
949  hn = pNext(h);
950  strat->redTailChange=FALSE;
951
952  while(hn != NULL)
953  {
954    j = 0;
955    while (j <= pos)
956    {
957      if (pDivisibleBy(strat->S[j], hn))
958      {
959        if (redstart == NULL)
960        {
961#ifdef KB_USE_HEAPS
962          int i;
963          if (hp->length <= 0) hp->length = pLength(hp->p);
964          h = hp->p;
965          for (i = 1, h = hp->p; i < length; i++, h = pNext(h));
966          hn = pNext(h);
967#endif
968          strat->redTailChange = TRUE;
969          redstart = h;
970
971          kBucketInit(&(strat->bucket), hn, -1, hp->heap);
972        }
973        number up = kBucketPolyRed(&(strat->bucket), strat->S[j],
974                                   pLength(strat->S[j]),
975                                   strat->kNoether);
976        nDelete(&up);
977        hn = kBucketGetLm(&(strat->bucket));
978        if (hn == NULL)
979        {
980          pNext(h) = NULL;
981          goto Finish;
982        }
983        j = 0;
984      }
985      else
986      {
987        j++;
988      }
989    }
990
991    if (redstart != NULL)
992    {
993      hn = kBucketExtractLm(&(strat->bucket));
994      pNext(h) = hn;
995      h = hn;
996      hn = kBucketGetLm(&(strat->bucket));
997    }
998    else
999    {
1000      pNext(h) = hn;
1001      h = hn;
1002      hn = pNext(h);
1003    }
1004    length++;
1005  }
1006
1007
1008  Finish:
1009  assume(pLength(hp->p) == length);
1010  hp->length = length;
1011}
1012
1013#endif // KB_USE_BUCKETS
1014
Note: See TracBrowser for help on using the repository browser.