source: git/Singular/kbuckets.cc @ c232af

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