1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4
5 #include "misc/auxiliary.h"
6 #include "misc/options.h"
7
8 #include "polys/monomials/p_polys.h"
9 #include "coeffs/coeffs.h"
10 #include "coeffs/numbers.h"
11 #include "polys/monomials/ring.h"
12 #include "polys/kbuckets.h"
13
14 #ifdef HAVE_SHIFTBBA
15 #include "polys/shiftop.h"
16 #endif
17
18 #ifdef HAVE_COEF_BUCKETS
19 #define USE_COEF_BUCKETS
20 #endif
21
22 #ifdef USE_COEF_BUCKETS
23 #ifdef HAVE_RINGS_OLD
24 #define MULTIPLY_BUCKET(B,I) do \
25 { if (B->coef[I]!=NULL) \
26 { \
27 assume(p_IsConstant(B->Coef[i],B->bucket->ring)); \
28 B->buckets[I]=p_Mult_q(B->buckets[I],B->coef[I],B->bucket_ring); \
29 B->coef[I]=NULL; \
30 } \
31 } while(0) \
32 if (rField_is_Ring(B->bucket_ring)) B->buckets_length[i] = pLength(B->buckets[i]);
33 #else
34 #define MULTIPLY_BUCKET(B,I) do \
35 { if (B->coef[I]!=NULL) \
36 { \
37 B->buckets[I]=p_Mult_q(B->buckets[I],B->coef[I],B->bucket_ring); \
38 B->coef[I]=NULL; \
39 } \
40 } while(0)
41 #endif
42 #else
43 #define MULTIPLY_BUCKET(B,I)
44 #endif
45 STATIC_VAR omBin kBucket_bin = omGetSpecBin(sizeof(kBucket));
46 #ifdef USE_COEF_BUCKETS
47 STATIC_VAR int coef_start=1;
48 #endif
49 //////////////////////////////////////////////////////////////////////////
50 ///
51 /// Some internal stuff
52 ///
53
54 // https://graphics.stanford.edu/~seander/bithacks.html#IntegerLog
55 #ifndef BUCKET_TWO_BASE
LOG4(int v)56 static inline int LOG4(int v)
57 {
58 const unsigned int b[] = {0x2, 0xC, 0xF0, 0xFF00, 0xFFFF0000};
59 const unsigned int S[] = {1, 2, 4, 8, 16};
60
61 unsigned int r = 0; // result of log4(v) will go here
62 if (v & b[4]) { v >>= S[4]; r |= S[3]; }
63 if (v & b[3]) { v >>= S[3]; r |= S[2]; }
64 if (v & b[2]) { v >>= S[2]; r |= S[1]; }
65 if (v & b[1]) { v >>= S[1]; r |= S[0]; }
66 return (int)r;
67 }
68 #endif
69
70 // returns ceil(log_4(l))
pLogLength(unsigned int l)71 static inline unsigned int pLogLength(unsigned int l)
72 {
73 unsigned int i = 0;
74
75 if (l == 0) return 0;
76 l--;
77 #ifdef BUCKET_TWO_BASE
78 i=SI_LOG2(l);
79 #else
80 i=LOG4(l);
81 #endif
82 return i+1;
83 }
84
85 // returns ceil(log_4(pLength(p)))
pLogLength(poly p)86 static inline unsigned int pLogLength(poly p)
87 {
88 return pLogLength((unsigned int) pLength(p));
89 }
90
91 #ifdef KDEBUG
92
93 #ifndef HAVE_PSEUDO_BUCKETS
kbTest_i(kBucket_pt bucket,int i)94 BOOLEAN kbTest_i(kBucket_pt bucket, int i)
95 {//sBucketSortMerge
96 #ifdef USE_COEF_BUCKETS
97 assume(bucket->coef[0]==NULL);
98 if ((bucket->coef[i]!=NULL) && (bucket->buckets[i]==NULL))
99 {
100 dReportError("Bucket %d coef not NULL", i);
101 }
102 if (bucket->coef[i]!=NULL)
103 {
104 assume(bucket->buckets[i]!=NULL);
105 p_Test(bucket->coef[i],bucket->bucket_ring);
106 }
107 #endif
108 pFalseReturn(p_Test(bucket->buckets[i], bucket->bucket_ring));
109 if ((unsigned)bucket->buckets_length[i] != pLength(bucket->buckets[i]))
110 {
111 dReportError("Bucket %d lengths difference should:%d has:%d",
112 i, bucket->buckets_length[i], pLength(bucket->buckets[i]));
113 }
114 else if (i > 0 && (int) pLogLength(bucket->buckets_length[i]) > i)
115 {
116 dReportError("Bucket %d too long %d",
117 i, bucket->buckets_length[i]);
118 }
119 if (i==0 && bucket->buckets_length[0] > 1)
120 {
121 dReportError("Bucket 0 too long");
122 }
123 return TRUE;
124 }
125
126
kbTest(kBucket_pt bucket)127 BOOLEAN kbTest(kBucket_pt bucket)
128 {
129 #ifdef HAVE_COEF_BUCKETS
130 assume(bucket->coef[0]==NULL);
131 #endif
132 int i;
133 poly lm = bucket->buckets[0];
134
135 omCheckAddrBin(bucket, kBucket_bin);
136 assume(bucket->buckets_used <= MAX_BUCKET);
137 if (! kbTest_i(bucket, 0)) return FALSE;
138 for (i=1; i<= (int) bucket->buckets_used; i++)
139 {
140 if (!kbTest_i(bucket, i)) return FALSE;
141 if (lm != NULL && bucket->buckets[i] != NULL
142 && p_LmCmp(lm, bucket->buckets[i], bucket->bucket_ring) != 1)
143 {
144 dReportError("Bucket %d larger or equal than lm", i);
145 if (p_LmCmp(lm, bucket->buckets[i], bucket->bucket_ring) ==0)
146 dReportError("Bucket %d equal to lm", i);
147 return FALSE;
148 }
149 if (!p_Test(bucket->buckets[i],bucket->bucket_ring))
150 {
151 dReportError("Bucket %d is not =0(4)", i);
152 return FALSE;
153 }
154 }
155
156 for (; i<=MAX_BUCKET; i++)
157 {
158 if (bucket->buckets[i] != NULL || bucket->buckets_length[i] != 0)
159 {
160 dReportError("Bucket %d not zero", i);
161 return FALSE;
162 }
163 }
164 for(i=0;i<=MAX_BUCKET;i++)
165 {
166 if (bucket->buckets[i]!=NULL)
167 {
168 int j;
169 for(j=i+1;j<=MAX_BUCKET;j++)
170 {
171 if (bucket->buckets[j]==bucket->buckets[i])
172 {
173 dReportError("Bucket %d %d equal", i,j);
174 return FALSE;
175 }
176 }
177 }
178 #ifdef HAVE_COEF_BUCKETS
179 if (bucket->coef[i]!=NULL)
180 {
181 int j;
182 for(j=i+1;j<=MAX_BUCKET;j++)
183 {
184 if (bucket->coef[j]==bucket->coef[i])
185 {
186 dReportError("internal coef %d %d equal", i,j);
187 return FALSE;
188 }
189 }
190 }
191 #endif
192 }
193 return TRUE;
194 }
195
196 #else // HAVE_PSEUDO_BUCKETS
kbTest(kBucket_pt bucket)197 BOOLEAN kbTest(kBucket_pt bucket)
198 {
199 return TRUE;
200 }
201 #endif // ! HAVE_PSEUDO_BUCKETS
202 #endif // KDEBUG
203
204 //////////////////////////////////////////////////////////////////////////
205 ///
206 /// Creation/Destruction of buckets
207 ///
208
kBucketCreate(const ring bucket_ring)209 kBucket_pt kBucketCreate(const ring bucket_ring)
210 {
211 assume(bucket_ring != NULL);
212 kBucket_pt bucket = (kBucket_pt) omAlloc0Bin(kBucket_bin);
213 bucket->bucket_ring = bucket_ring;
214 return bucket;
215 }
kBucketDestroy(kBucket_pt * bucket_pt)216 void kBucketDestroy(kBucket_pt *bucket_pt)
217 {
218 omFreeBin(*bucket_pt, kBucket_bin);
219 *bucket_pt = NULL;
220 }
221
222
kBucketDeleteAndDestroy(kBucket_pt * bucket_pt)223 void kBucketDeleteAndDestroy(kBucket_pt *bucket_pt)
224 {
225 kBucket_pt bucket = *bucket_pt;
226 kbTest(bucket);
227 int i;
228 for (i=0; i<= bucket->buckets_used; i++)
229 {
230 p_Delete(&(bucket->buckets[i]), bucket->bucket_ring);
231 #ifdef USE_COEF_BUCKETS
232 p_Delete(&(bucket->coef[i]), bucket->bucket_ring);
233 #endif
234 }
235 omFreeBin(bucket, kBucket_bin);
236 *bucket_pt = NULL;
237 }
238
239 /////////////////////////////////////////////////////////////////////////////
240 // Convertion from/to Bpolys
241 //
242 #ifndef HAVE_PSEUDO_BUCKETS
243
kBucketMergeLm(kBucket_pt bucket)244 inline void kBucketMergeLm(kBucket_pt bucket)
245 {
246 kbTest(bucket);
247 if (bucket->buckets[0] != NULL)
248 {
249 poly lm = bucket->buckets[0];
250 int i = 1;
251 #ifdef BUCKET_TWO_BASE
252 int l = 2;
253 while ( bucket->buckets_length[i] >= l)
254 {
255 i++;
256 l = l << 1;
257 }
258 #else
259 int l = 4;
260 while ( bucket->buckets_length[i] >= l)
261 {
262 i++;
263 l = l << 2;
264 }
265 #endif
266 #ifndef USE_COEF_BUCKETS
267 MULTIPLY_BUCKET(bucket,i);
268 pNext(lm) = bucket->buckets[i];
269 bucket->buckets[i] = lm;
270 bucket->buckets_length[i]++;
271 assume(i <= bucket->buckets_used+1);
272 if (i > bucket->buckets_used) bucket->buckets_used = i;
273 bucket->buckets[0] = NULL;
274 bucket->buckets_length[0] = 0;
275 kbTest(bucket);
276 #else
277 if (i > bucket->buckets_used) bucket->buckets_used = i;
278 assume(i!=0);
279 if (bucket->buckets[i]!=NULL)
280 {
281 MULTIPLY_BUCKET(bucket,i);
282 pNext(lm) = bucket->buckets[i];
283 bucket->buckets[i] = lm;
284 bucket->buckets_length[i]++;
285 assume(i <= bucket->buckets_used+1);
286 }
287 else
288 {
289 #if 1
290 assume(bucket->buckets[i]==NULL);
291 assume(bucket->coef[0]==NULL);
292 assume(pLength(lm)==1);
293 assume(pNext(lm)==NULL);
294 number coef=p_GetCoeff(lm,bucket->bucket_ring);
295 //WARNING: not thread_safe
296 p_SetCoeff0(lm, n_Init(1,bucket->bucket_ring), bucket->bucket_ring);
297 bucket->buckets[i]=lm;
298 bucket->buckets_length[i]=1;
299 bucket->coef[i]=p_NSet(n_Copy(coef,bucket->bucket_ring),bucket->bucket_ring);
300
301 bucket->buckets[i]=lm;
302 bucket->buckets_length[i]=1;
303 #else
304 MULTIPLY_BUCKET(bucket,i);
305 pNext(lm) = bucket->buckets[i];
306 bucket->buckets[i] = lm;
307 bucket->buckets_length[i]++;
308 assume(i <= bucket->buckets_used+1);
309 #endif
310 }
311 bucket->buckets[0]=NULL;
312 bucket->buckets_length[0] = 0;
313 bucket->coef[0]=NULL;
314 kbTest(bucket);
315 #endif
316 }
317
318 }
319
kBucketIsCleared(kBucket_pt bucket)320 BOOLEAN kBucketIsCleared(kBucket_pt bucket)
321 {
322 int i;
323
324 for (i = 0;i<=MAX_BUCKET;i++)
325 {
326 if (bucket->buckets[i] != NULL) return FALSE;
327 #ifdef HAVE_COEF_BUCKETS
328 if (bucket->coef[i] != NULL) return FALSE;
329 #endif
330 if (bucket->buckets_length[i] != 0) return FALSE;
331 }
332 return TRUE;
333 }
334
kBucketInit(kBucket_pt bucket,poly lm,int length)335 void kBucketInit(kBucket_pt bucket, poly lm, int length)
336 {
337 //assume(false);
338 assume(bucket != NULL);
339 assume(length <= 0 || (unsigned)length == pLength(lm));
340 assume(kBucketIsCleared(bucket));
341
342 if (lm == NULL) return;
343
344 if (length <= 0)
345 length = pLength(lm);
346
347 bucket->buckets[0] = lm;
348 #ifdef HAVE_COEF_BUCKETS
349 assume(bucket->coef[0]==NULL);
350 #endif
351 #ifdef USE_COEF_BUCKETS
352 bucket->coef[0]=NULL;
353 #endif
354 if (lm!=NULL)
355 bucket->buckets_length[0] = 1;
356 else
357 bucket->buckets_length[0]= 0;
358 if (length > 1)
359 {
360 unsigned int i = pLogLength(length-1);
361 bucket->buckets[i] = pNext(lm);
362 pNext(lm) = NULL;
363 bucket->buckets_length[i] = length-1;
364 bucket->buckets_used = i;
365 }
366 else
367 {
368 bucket->buckets_used = 0;
369 }
370 }
371
kBucketCanonicalize(kBucket_pt bucket)372 int kBucketCanonicalize(kBucket_pt bucket)
373 {
374 #ifndef HAVE_PSEUDO_BUCKETS
375 assume(bucket->buckets_used<=MAX_BUCKET);
376 MULTIPLY_BUCKET(bucket,1);
377 kbTest(bucket);
378 poly p = bucket->buckets[1];
379 poly lm;
380 int pl = bucket->buckets_length[1];//, i;
381 int i;
382 bucket->buckets[1] = NULL;
383 bucket->buckets_length[1] = 0;
384 #ifdef USE_COEF_BUCKETS
385 assume(bucket->coef[1]==NULL);
386 #endif
387 ring r=bucket->bucket_ring;
388
389
390 for (i=1; i<=bucket->buckets_used; i++)
391 {
392 #ifdef USE_COEF_BUCKETS
393 if (bucket->coef[i]!=NULL)
394 {
395 assume(bucket->buckets[i]!=NULL);
396 p = p_Plus_mm_Mult_qq(p, bucket->coef[i], bucket->buckets[i],
397 pl, bucket->buckets_length[i], r);
398 p_Delete(&bucket->coef[i],r);
399 p_Delete(&bucket->buckets[i],r);
400 }
401 else
402 p = p_Add_q(p, bucket->buckets[i],
403 pl, bucket->buckets_length[i], r);
404 #else
405 p = p_Add_q(p, bucket->buckets[i],
406 pl, bucket->buckets_length[i], r);
407 #endif
408 if (i==1) continue;
409 bucket->buckets[i] = NULL;
410 bucket->buckets_length[i] = 0;
411 }
412 #ifdef HAVE_COEF_BUCKETS
413 assume(bucket->coef[0]==NULL);
414 #endif
415 lm = bucket->buckets[0];
416 if (lm != NULL)
417 {
418 pNext(lm) = p;
419 p = lm;
420 pl++;
421 bucket->buckets[0] = NULL;
422 bucket->buckets_length[0] = 0;
423 }
424 if (pl > 0)
425 {
426 i = pLogLength(pl);
427 bucket->buckets[i] = p;
428 bucket->buckets_length[i] = pl;
429 }
430 else
431 {
432 i = 0;
433 }
434 bucket->buckets_used = i;
435 assume(bucket->buckets_used <= MAX_BUCKET);
436 #ifdef USE_COEF_BUCKETS
437 assume(bucket->coef[0]==NULL);
438 assume(bucket->coef[i]==NULL);
439 #endif
440 assume(pLength(p) == (unsigned)pl);
441 //if (TEST_OPT_PROT) { Print("C(%d)",pl); }
442 kbTest(bucket);
443 return i;
444 #endif
445 }
446
kBucketNormalize(kBucket_pt bucket)447 void kBucketNormalize(kBucket_pt bucket)
448 {
449 #ifdef HAVE_PSEUDO_BUCKETS
450 p_Normalize(bucket->p,bucket->bucket_ring);
451 #else
452 MULTIPLY_BUCKET(bucket,1);
453 for (int i=0; i<=bucket->buckets_used; i++)
454 {
455 p_Normalize(bucket->buckets[i],bucket->bucket_ring);
456 }
457 #endif
458 }
459
kBucketClear(kBucket_pt bucket,poly * p,int * length)460 void kBucketClear(kBucket_pt bucket, poly *p, int *length)
461 {
462 int i = kBucketCanonicalize(bucket);
463 if (i > 0)
464 {
465 #ifdef USE_COEF_BUCKETS
466 MULTIPLY_BUCKET(bucket,i);
467 //bucket->coef[i]=NULL;
468 #endif
469 *p = bucket->buckets[i];
470 *length = bucket->buckets_length[i];
471 bucket->buckets[i] = NULL;
472 bucket->buckets_length[i] = 0;
473 bucket->buckets_used = 0;
474
475 }
476 else
477 {
478 *p = NULL;
479 *length = 0;
480 }
481 }
482
kBucketSetLm(kBucket_pt bucket,poly lm)483 void kBucketSetLm(kBucket_pt bucket, poly lm)
484 {
485 kBucketMergeLm(bucket);
486 pNext(lm) = NULL;
487 bucket->buckets[0] = lm;
488 bucket->buckets_length[0] = 1;
489 }
490
491 #else // HAVE_PSEUDO_BUCKETS
492
kBucketInit(kBucket_pt bucket,poly lm,int length)493 void kBucketInit(kBucket_pt bucket, poly lm, int length)
494 {
495 int i;
496
497 assume(bucket != NULL);
498 assume(length <= 0 || length == pLength(lm));
499
500 bucket->p = lm;
501 if (length <= 0) bucket->l = pLength(lm);
502 else bucket->l = length;
503
504 }
505
kBucketGetLm(kBucket_pt bucket)506 const poly kBucketGetLm(kBucket_pt bucket)
507 {
508 return bucket->p;
509 }
510
kBucketExtractLm(kBucket_pt bucket)511 poly kBucketExtractLm(kBucket_pt bucket)
512 {
513 poly lm = bucket->p;
514 assume(pLength(bucket->p) == bucket->l);
515 pIter(bucket->p);
516 (bucket->l)--;
517 pNext(lm) = NULL;
518 return lm;
519 }
520
kBucketClear(kBucket_pt bucket,poly * p,int * length)521 void kBucketClear(kBucket_pt bucket, poly *p, int *length)
522 {
523 assume(pLength(bucket->p) == bucket->l);
524 *p = bucket->p;
525 *length = bucket->l;
526 bucket->p = NULL;
527 bucket->l = 0;
528 }
529
530 #endif // ! HAVE_PSEUDO_BUCKETS
531 //////////////////////////////////////////////////////////////////////////
532 ///
533 /// For changing the ring of the Bpoly to new_tailBin
534 ///
kBucketShallowCopyDelete(kBucket_pt bucket,ring new_tailRing,omBin new_tailBin,pShallowCopyDeleteProc p_shallow_copy_delete)535 void kBucketShallowCopyDelete(kBucket_pt bucket,
536 ring new_tailRing, omBin new_tailBin,
537 pShallowCopyDeleteProc p_shallow_copy_delete)
538 {
539 #ifndef HAVE_PSEUDO_BUCKETS
540 int i;
541
542 kBucketCanonicalize(bucket);
543 for (i=0; i<= bucket->buckets_used; i++)
544 if (bucket->buckets[i] != NULL)
545 {
546 MULTIPLY_BUCKET(bucket,i);
547 bucket->buckets[i] = p_shallow_copy_delete(bucket->buckets[i],
548 bucket->bucket_ring,
549 new_tailRing,
550 new_tailBin);
551 }
552 #else
553 bucket->p = p_shallow_copy_delete(p,
554 bucket_ring,
555 new_tailRing,
556 new_tailBin);
557 #endif
558 bucket->bucket_ring = new_tailRing;
559 }
560
561 //////////////////////////////////////////////////////////////////////////
562 ///
563 /// Bucket number i from bucket is out of length sync, resync
564 ///
kBucketAdjust(kBucket_pt bucket,int i)565 void kBucketAdjust(kBucket_pt bucket, int i) {
566
567 MULTIPLY_BUCKET(bucket,i);
568
569 int l1 = bucket->buckets_length[i];
570 poly p1 = bucket->buckets[i];
571 bucket->buckets[i] = NULL;
572 bucket->buckets_length[i] = 0;
573 i = pLogLength(l1);
574
575 while (bucket->buckets[i] != NULL)
576 {
577 //kbTest(bucket);
578 MULTIPLY_BUCKET(bucket,i);
579 p1 = p_Add_q(p1, bucket->buckets[i],
580 l1, bucket->buckets_length[i], bucket->bucket_ring);
581 bucket->buckets[i] = NULL;
582 bucket->buckets_length[i] = 0;
583 i = pLogLength(l1);
584 }
585
586 bucket->buckets[i] = p1;
587 bucket->buckets_length[i]=l1;
588 if (i >= bucket->buckets_used)
589 bucket->buckets_used = i;
590 else
591 kBucketAdjustBucketsUsed(bucket);
592 }
593
594 //////////////////////////////////////////////////////////////////////////
595 ///
596 /// Multiply Bucket by number ,i.e. Bpoly == n*Bpoly
597 ///
kBucket_Mult_n(kBucket_pt bucket,number n)598 void kBucket_Mult_n(kBucket_pt bucket, number n)
599 {
600 #ifndef HAVE_PSEUDO_BUCKETS
601 kbTest(bucket);
602 ring r=bucket->bucket_ring;
603 int i;
604
605 for (i=0; i<= bucket->buckets_used; i++)
606 {
607 if (bucket->buckets[i] != NULL)
608 {
609 #ifdef USE_COEF_BUCKETS
610 if (i<coef_start)
611 bucket->buckets[i] = __p_Mult_nn(bucket->buckets[i], n, r);
612 /* Frank Seelisch on March 11, 2010:
613 This looks a bit strange: The following "if" is indented
614 like the previous line of code. But coded as it is,
615 it should actually be two spaces less indented.
616 Question: Should the following "if" also only be
617 performed when "(i<coef_start)" is true?
618 For the time being, I leave it as it is. */
619 if (rField_is_Ring(r) && !(rField_is_Domain(r)))
620 {
621 bucket->buckets_length[i] = pLength(bucket->buckets[i]);
622 kBucketAdjust(bucket, i);
623 }
624 else
625 if (bucket->coef[i]!=NULL)
626 {
627 bucket->coef[i] = __p_Mult_nn(bucket->coef[i],n,r);
628 }
629 else
630 {
631 bucket->coef[i] = p_NSet(n_Copy(n,r),r);
632 }
633 #else
634 bucket->buckets[i] = __p_Mult_nn(bucket->buckets[i], n, r);
635 #endif
636 }
637 }
638 if (rField_is_Ring(r) && !(rField_is_Domain(r)))
639 {
640 for (i=0; i<= bucket->buckets_used; i++)
641 {
642 if (bucket->buckets[i] != NULL)
643 {
644 bucket->buckets_length[i] = pLength(bucket->buckets[i]);
645 kBucketAdjust(bucket, i);
646 }
647 }
648 }
649 kbTest(bucket);
650 #else
651 bucket->p = __p_Mult_nn(bucket->p, n, bucket->bucket_ring);
652 #endif
653 }
654
655
656 //////////////////////////////////////////////////////////////////////////
657 ///
658 /// Add to Bucket a poly ,i.e. Bpoly == q+Bpoly
659 ///
kBucket_Add_q(kBucket_pt bucket,poly q,int * l)660 void kBucket_Add_q(kBucket_pt bucket, poly q, int *l)
661 {
662 if (q == NULL) return;
663 assume(*l <= 0 || pLength(q) == *l);
664
665 int i, l1;
666 ring r = bucket->bucket_ring;
667
668 if (*l <= 0)
669 {
670 l1 = pLength(q);
671 *l = l1;
672 }
673 else
674 l1 = *l;
675
676 kBucketMergeLm(bucket);
677 kbTest(bucket);
678 i = pLogLength(l1);
679
680 while (bucket->buckets[i] != NULL)
681 {
682 //MULTIPLY_BUCKET(bucket,i);
683 #ifdef USE_COEF_BUCKETS
684 if (bucket->coef[i]!=NULL)
685 {
686 q = p_Plus_mm_Mult_qq(q, bucket->coef[i], bucket->buckets[i],
687 l1, bucket->buckets_length[i], r);
688 p_Delete(&bucket->coef[i],r);
689 p_Delete(&bucket->buckets[i],r);
690 }
691 else
692 q = p_Add_q(q, bucket->buckets[i],
693 l1, bucket->buckets_length[i], r);
694 #else
695 q = p_Add_q(q, bucket->buckets[i],
696 l1, bucket->buckets_length[i], r);
697 #endif
698 bucket->buckets[i] = NULL;
699 bucket->buckets_length[i] = 0;
700 i = pLogLength(l1);
701 assume(i<= MAX_BUCKET);
702 assume(bucket->buckets_used<= MAX_BUCKET);
703 }
704
705 kbTest(bucket);
706 bucket->buckets[i] = q;
707 bucket->buckets_length[i]=l1;
708 if (i >= bucket->buckets_used)
709 bucket->buckets_used = i;
710 else
711 kBucketAdjustBucketsUsed(bucket);
712 kbTest(bucket);
713 }
714
715
716
717 //////////////////////////////////////////////////////////////////////////
718 ///
719 /// Bpoly == Bpoly - m*p; where m is a monom
720 /// Does not destroy p and m
721 /// assume (*l <= 0 || pLength(p) == *l)
kBucket_Minus_m_Mult_p(kBucket_pt bucket,poly m,poly p,int * l,poly spNoether)722 void kBucket_Minus_m_Mult_p(kBucket_pt bucket, poly m, poly p, int *l,
723 poly spNoether)
724 {
725 assume(*l <= 0 || pLength(p) == *l);
726 int i, l1;
727 poly p1 = p;
728 ring r = bucket->bucket_ring;
729
730 if (*l <= 0)
731 {
732 l1 = pLength(p1);
733 *l = l1;
734 }
735 else
736 l1 = *l;
737
738 if (m == NULL || p == NULL) return;
739
740 #ifndef HAVE_PSEUDO_BUCKETS
741 kBucketMergeLm(bucket);
742 kbTest(bucket);
743 i = pLogLength(l1);
744
745 #if defined(HAVE_PLURAL)
746 if ((rField_is_Ring(r) && !(rField_is_Domain(r)))
747 ||(rIsPluralRing(r)))
748 {
749 pSetCoeff0(m, n_InpNeg(pGetCoeff(m),r->cf));
750 p1=r->p_Procs->pp_mm_Mult(p,m,r);
751 pSetCoeff0(m, n_InpNeg(pGetCoeff(m),r->cf));
752 l1=pLength(p1);
753 i = pLogLength(l1);
754 }
755 else
756 #endif
757 {
758 if ((i <= bucket->buckets_used) && (bucket->buckets[i] != NULL))
759 {
760 assume(pLength(bucket->buckets[i])==(unsigned)bucket->buckets_length[i]);
761 //#ifdef USE_COEF_BUCKETS
762 // if(bucket->coef[i]!=NULL)
763 // {
764 // poly mult=p_Mult_mm(bucket->coef[i],m,r);
765 // bucket->coef[i]=NULL;
766 // p1 = p_Minus_mm_Mult_qq(bucket->buckets[i], mult, p1,
767 // bucket->buckets_length[i], l1,
768 // spNoether, r);
769 // }
770 // else
771 //#endif
772 MULTIPLY_BUCKET(bucket,i);
773 p1 = p_Minus_mm_Mult_qq(bucket->buckets[i], m, p1,
774 bucket->buckets_length[i], l1,
775 spNoether, r);
776 l1 = bucket->buckets_length[i];
777 bucket->buckets[i] = NULL;
778 bucket->buckets_length[i] = 0;
779 i = pLogLength(l1);
780 }
781 else
782 {
783 pSetCoeff0(m, n_InpNeg(pGetCoeff(m),r->cf));
784 if (spNoether != NULL)
785 {
786 l1 = -1;
787 p1 = r->p_Procs->pp_Mult_mm_Noether(p1, m, spNoether, l1, r);
788 i = pLogLength(l1);
789 }
790 else
791 {
792 p1 = r->p_Procs->pp_mm_Mult(p1, m, r);
793 }
794 pSetCoeff0(m, n_InpNeg(pGetCoeff(m),r->cf));
795 }
796 }
797
798 while (bucket->buckets[i] != NULL)
799 {
800 //kbTest(bucket);
801 MULTIPLY_BUCKET(bucket,i);
802 p1 = p_Add_q(p1, bucket->buckets[i],
803 l1, bucket->buckets_length[i], r);
804 bucket->buckets[i] = NULL;
805 bucket->buckets_length[i] = 0;
806 i = pLogLength(l1);
807 }
808
809 bucket->buckets[i] = p1;
810 bucket->buckets_length[i]=l1;
811 if (i >= bucket->buckets_used)
812 bucket->buckets_used = i;
813 else
814 kBucketAdjustBucketsUsed(bucket);
815 #else // HAVE_PSEUDO_BUCKETS
816 bucket->p = p_Minus_mm_Mult_qq(bucket->p, m, p,
817 bucket->l, l1,
818 spNoether, r);
819 #endif
820 }
821
822 //////////////////////////////////////////////////////////////////////////
823 ///
824 /// Bpoly == Bpoly + m*p; where m is a monom
825 /// Does not destroy p and m
826 /// assume (l <= 0 || pLength(p) == l)
kBucket_Plus_mm_Mult_pp(kBucket_pt bucket,poly m,poly p,int l)827 void kBucket_Plus_mm_Mult_pp(kBucket_pt bucket, poly m, poly p, int l)
828 {
829 assume((!rIsPluralRing(bucket->bucket_ring))||p_IsConstant(m, bucket->bucket_ring));
830 assume(l <= 0 || pLength(p) == (unsigned)l);
831 int i, l1;
832 poly p1 = p;
833 ring r = bucket->bucket_ring;
834
835 if (m == NULL || p == NULL) return;
836
837 if (l <= 0)
838 {
839 l1 = pLength(p1);
840 l = l1;
841 }
842 else
843 l1 = l;
844
845 kBucketMergeLm(bucket);
846 kbTest(bucket);
847 i = pLogLength(l1);
848 #ifdef USE_COEF_BUCKETS
849 number n=n_Init(1,r->cf);
850 #endif
851 if (i <= bucket->buckets_used && bucket->buckets[i] != NULL)
852 {
853 //if (FALSE){
854 #ifdef USE_COEF_BUCKETS
855 if ((bucket->coef[i]!=NULL) &&(i>=coef_start))
856 {
857 number orig_coef=p_GetCoeff(bucket->coef[i],r);
858 //we take ownership:
859 p_SetCoeff0(bucket->coef[i],n_Init(0,r),r);
860 number add_coef=n_Copy(p_GetCoeff(m,r),r);
861 number gcd=n_Gcd(add_coef, orig_coef,r);
862
863 if (!(n_IsOne(gcd,r)))
864 {
865 number orig_coef2=n_ExactDiv(orig_coef,gcd,r);
866 number add_coef2=n_ExactDiv(add_coef, gcd,r);
867 n_Delete(&orig_coef,r);
868 n_Delete(&add_coef,r);
869 orig_coef=orig_coef2;
870 add_coef=add_coef2;
871
872 //p_Mult_nn(bucket->buckets[i], orig_coef,r);
873 n_Delete(&n,r);
874 n=gcd;
875 }
876
877 //assume(n_IsOne(n,r));
878 number backup=p_GetCoeff(m,r);
879
880 p_SetCoeff0(m,add_coef,r);
881 bucket->buckets[i]=__p_Mult_nn(bucket->buckets[i],orig_coef,r);
882
883 n_Delete(&orig_coef,r);
884 p_Delete(&bucket->coef[i],r);
885
886 p1 = p_Plus_mm_Mult_qq(bucket->buckets[i], m, p1,
887 bucket->buckets_length[i], l1, r);
888 l1=bucket->buckets_length[i];
889 bucket->buckets[i]=NULL;
890 bucket->buckets_length[i] = 0;
891 i = pLogLength(l1);
892 assume(l1==pLength(p1));
893
894 p_SetCoeff(m,backup,r); //deletes add_coef
895 }
896 else
897 #endif
898 {
899 MULTIPLY_BUCKET(bucket,i);
900 p1 = p_Plus_mm_Mult_qq(bucket->buckets[i], m, p1,
901 bucket->buckets_length[i], l1, r);
902 l1 = bucket->buckets_length[i];
903 bucket->buckets[i] = NULL;
904 bucket->buckets_length[i] = 0;
905 i = pLogLength(l1);
906 }
907 }
908 else
909 {
910 #ifdef USE_COEF_BUCKETS
911 number swap_n=p_GetCoeff(m,r);
912
913 assume(n_IsOne(n,r));
914 p_SetCoeff0(m,n,r);
915 n=swap_n;
916 //p_SetCoeff0(n, swap_n, r);
917 //p_GetCoeff0(n, swap_n,r);
918 #endif
919 p1 = r->p_Procs->pp_Mult_mm(p1, m, r);
920 #ifdef USE_COEF_BUCKETS
921 //m may not be changed
922 p_SetCoeff(m,n_Copy(n,r),r);
923 #endif
924 }
925
926 while ((bucket->buckets[i] != NULL) && (p1!=NULL))
927 {
928 assume(i!=0);
929 #ifdef USE_COEF_BUCKETS
930 if ((bucket->coef[i]!=NULL) &&(i>=coef_start))
931 {
932 number orig_coef=p_GetCoeff(bucket->coef[i],r);
933 //we take ownership:
934 p_SetCoeff0(bucket->coef[i],n_Init(0,r),r);
935 number add_coef=n_Copy(n,r);
936 number gcd=n_Gcd(add_coef, orig_coef,r);
937
938 if (!(n_IsOne(gcd,r)))
939 {
940 number orig_coef2=n_ExactDiv(orig_coef,gcd,r);
941 number add_coef2=n_ExactDiv(add_coef, gcd,r);
942 n_Delete(&orig_coef,r);
943 n_Delete(&n,r);
944 n_Delete(&add_coef,r);
945 orig_coef=orig_coef2;
946 add_coef=add_coef2;
947 //p_Mult_nn(bucket->buckets[i], orig_coef,r);
948 n=gcd;
949 }
950 //assume(n_IsOne(n,r));
951 bucket->buckets[i]=__p_Mult_nn(bucket->buckets[i],orig_coef,r);
952 p1=__p_Mult_nn(p1,add_coef,r);
953
954 p1 = p_Add_q(p1, bucket->buckets[i],r);
955 l1=pLength(p1);
956
957 bucket->buckets[i]=NULL;
958 n_Delete(&orig_coef,r);
959 p_Delete(&bucket->coef[i],r);
960 //l1=bucket->buckets_length[i];
961 assume(l1==pLength(p1));
962 }
963 else
964 #endif
965 {
966 //don't do that, pull out gcd
967 #ifdef USE_COEF_BUCKETS
968 if(!(n_IsOne(n,r)))
969 {
970 p1=__p_Mult_nn(p1, n, r);
971 n_Delete(&n,r);
972 n=n_Init(1,r);
973 }
974 #endif
975 MULTIPLY_BUCKET(bucket,i);
976 p1 = p_Add_q(p1, bucket->buckets[i],
977 l1, bucket->buckets_length[i], r);
978 bucket->buckets[i] = NULL;
979 bucket->buckets_length[i] = 0;
980 }
981 i = pLogLength(l1);
982 }
983
984 bucket->buckets[i] = p1;
985 #ifdef USE_COEF_BUCKETS
986 assume(bucket->coef[i]==NULL);
987
988 if (!(n_IsOne(n,r)))
989 {
990 bucket->coef[i]=p_NSet(n,r);
991 }
992 else
993 {
994 bucket->coef[i]=NULL;
995 n_Delete(&n,r);
996 }
997
998 if (p1==NULL)
999 p_Delete(&bucket->coef[i],r);
1000 #endif
1001 bucket->buckets_length[i]=l1;
1002 if (i > bucket->buckets_used)
1003 bucket->buckets_used = i;
1004 else
1005 kBucketAdjustBucketsUsed(bucket);
1006
1007 kbTest(bucket);
1008 }
1009
kBucket_ExtractLarger(kBucket_pt bucket,poly q,poly append)1010 poly kBucket_ExtractLarger(kBucket_pt bucket, poly q, poly append)
1011 {
1012 if (q == NULL) return append;
1013 poly lm;
1014 loop
1015 {
1016 lm = kBucketGetLm(bucket);
1017 if (lm == NULL) return append;
1018 if (p_LmCmp(lm, q, bucket->bucket_ring) == 1)
1019 {
1020 lm = kBucketExtractLm(bucket);
1021 pNext(append) = lm;
1022 pIter(append);
1023 }
1024 else
1025 {
1026 return append;
1027 }
1028 }
1029 }
1030
1031 /////////////////////////////////////////////////////////////////////////////
1032 //
1033 // Extract all monomials from bucket with component comp
1034 // Return as a polynomial *p with length *l
1035 // In other words, afterwards
1036 // Bpoly = Bpoly - (poly consisting of all monomials with component comp)
1037 // and components of monomials of *p are all 0
1038 //
1039
1040 // Hmm... for now I'm too lazy to implement those independent of currRing
1041 // But better declare it extern than including polys.h
1042 extern void p_TakeOutComp(poly *p, long comp, poly *q, int *lq, const ring r);
1043
kBucketTakeOutComp(kBucket_pt bucket,long comp,poly * r_p,int * l)1044 void kBucketTakeOutComp(kBucket_pt bucket,
1045 long comp,
1046 poly *r_p, int *l)
1047 {
1048 poly p = NULL, q;
1049 int i, lp = 0, lq;
1050
1051 #ifndef HAVE_PSEUDO_BUCKETS
1052 kBucketMergeLm(bucket);
1053 for (i=1; i<=bucket->buckets_used; i++)
1054 {
1055 if (bucket->buckets[i] != NULL)
1056 {
1057 MULTIPLY_BUCKET(bucket,i);
1058 p_TakeOutComp(&(bucket->buckets[i]), comp, &q, &lq, bucket->bucket_ring);
1059 if (q != NULL)
1060 {
1061 assume(pLength(q) == (unsigned)lq);
1062 bucket->buckets_length[i] -= lq;
1063 assume(pLength(bucket->buckets[i]) == (unsigned)bucket->buckets_length[i]);
1064 p = p_Add_q(p, q, lp, lq, bucket->bucket_ring);
1065 }
1066 }
1067 }
1068 kBucketAdjustBucketsUsed(bucket);
1069 #else
1070 p_TakeOutComp(&(bucket->p), comp, &p, &lp,bucket->bucket_ring);
1071 (bucket->l) -= lp;
1072 #endif
1073 *r_p = p;
1074 *l = lp;
1075
1076 kbTest(bucket);
1077 }
1078
1079 /////////////////////////////////////////////////////////////////////////////
1080 // Reduction of Bpoly with a given poly
1081 //
1082
1083 extern int ksCheckCoeff(number *a, number *b);
1084
kBucketPolyRed(kBucket_pt bucket,poly p1,int l1,poly spNoether)1085 number kBucketPolyRed(kBucket_pt bucket,
1086 poly p1, int l1,
1087 poly spNoether)
1088 {
1089 ring r=bucket->bucket_ring;
1090 assume((!rIsPluralRing(r))||p_LmEqual(p1,kBucketGetLm(bucket), r));
1091 assume(p1 != NULL &&
1092 p_DivisibleBy(p1, kBucketGetLm(bucket), r));
1093 assume(pLength(p1) == (unsigned) l1);
1094
1095 poly a1 = pNext(p1), lm = kBucketExtractLm(bucket);
1096 BOOLEAN reset_vec=FALSE;
1097 number rn;
1098
1099 /* we shall reduce bucket=bn*lm+... by p1=an*t+a1 where t=lm(p1)
1100 and an,bn shall be defined further down only if lc(p1)!=1
1101 we already know: an|bn and t|lm */
1102 if(a1==NULL)
1103 {
1104 p_LmDelete(&lm, r);
1105 return n_Init(1,r->cf);
1106 }
1107
1108 if (! n_IsOne(pGetCoeff(p1),r->cf))
1109 {
1110 number an = pGetCoeff(p1), bn = pGetCoeff(lm);
1111 //StringSetS("##### an = "); nWrite(an); PrintS(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
1112 //StringSetS("##### bn = "); nWrite(bn); PrintS(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
1113 /* ksCheckCoeff: divide out gcd from an and bn: */
1114 int ct = ksCheckCoeff(&an, &bn,r->cf);
1115 /* the previous command returns ct=0 or ct=2 iff an!=1
1116 note: an is now 1 or -1 */
1117
1118 /* setup factor for p1 which cancels leading terms */
1119 p_SetCoeff(lm, bn, r);
1120 if ((ct == 0) || (ct == 2))
1121 {
1122 /* next line used to be here before but is WRONG:
1123 kBucket_Mult_n(bucket, an);
1124 its use would result in a wrong sign for the tail of bucket
1125 in the reduction */
1126
1127 /* correct factor for cancelation by changing sign if an=-1 */
1128 if (rField_is_Ring(r))
1129 lm = __p_Mult_nn(lm, an, r);
1130 else
1131 kBucket_Mult_n(bucket, an);
1132 }
1133 rn = an;
1134 }
1135 else
1136 {
1137 rn = n_Init(1,r->cf);
1138 }
1139
1140 if (p_GetComp(p1, r) != p_GetComp(lm, r))
1141 {
1142 p_SetCompP(a1, p_GetComp(lm, r), r);
1143 reset_vec = TRUE;
1144 p_SetComp(lm, p_GetComp(p1, r), r);
1145 p_Setm(lm, r);
1146 }
1147
1148 p_ExpVectorSub(lm, p1, r);
1149 l1--;
1150
1151 assume((unsigned)l1==pLength(a1));
1152
1153 #ifdef HAVE_SHIFTBBA
1154 poly lmRight;
1155 if (r->isLPring) {
1156 int firstBlock = p_mFirstVblock(p1, r);
1157 k_SplitFrame(lm, lmRight, si_max(firstBlock, 1), r);
1158 }
1159 #endif
1160 #if 0
1161 BOOLEAN backuped=FALSE;
1162 number coef;
1163 //@Viktor, don't ignore coefficients on monomials
1164 if(l1==1) {
1165
1166 //if (rField_is_Q(r)) {
1167 //avoid this for function fields, as gcds are expensive at the moment
1168
1169
1170 coef=p_GetCoeff(a1,r);
1171 lm=p_Mult_nn(lm, coef, r);
1172 p_SetCoeff0(a1, n_Init(1,r), r);
1173 backuped=TRUE;
1174 //WARNING: not thread_safe
1175 //deletes coef as side effect
1176 //}
1177 }
1178 #endif
1179
1180 #ifdef HAVE_SHIFTBBA
1181 if (r->isLPring)
1182 {
1183 kBucket_Minus_m_Mult_p(bucket, lm, r->p_Procs->pp_Mult_mm(a1, lmRight, r), &l1, spNoether);
1184 }
1185 else
1186 #endif
1187 {
1188 kBucket_Minus_m_Mult_p(bucket, lm, a1, &l1, spNoether);
1189 }
1190
1191 #if 0
1192 if (backuped)
1193 p_SetCoeff0(a1,coef,r);
1194 #endif
1195
1196 p_LmDelete(&lm, r);
1197 #ifdef HAVE_SHIFTBBA
1198 if (r->isLPring)
1199 {
1200 p_LmDelete(&lmRight, r);
1201 }
1202 #endif
1203 if (reset_vec) p_SetCompP(a1, 0, r);
1204 kbTest(bucket);
1205 return rn;
1206 }
1207
1208 #ifndef USE_COEF_BUCKETS
kBucketSimpleContent(kBucket_pt bucket)1209 void kBucketSimpleContent(kBucket_pt bucket)
1210 {
1211 if (bucket->buckets[0]==NULL) return;
1212
1213 ring r=bucket->bucket_ring;
1214 if (rField_is_Ring(r)) return;
1215
1216 coeffs cf=r->cf;
1217 if (cf->cfSubringGcd==ndGcd) /* trivial gcd*/ return;
1218
1219 number nn=pGetCoeff(bucket->buckets[0]);
1220 //if ((bucket->buckets_used==0)
1221 //&&(!n_IsOne(nn,cf)))
1222 //{
1223 // if (TEST_OPT_PROT) PrintS("@");
1224 // p_SetCoeff(bucket->buckets[0],n_Init(1,cf),r);
1225 // return;
1226 //}
1227
1228 if (n_Size(nn,cf)<2) return;
1229
1230 //kBucketAdjustBucketsUsed(bucket);
1231 number coef=n_Copy(nn,cf);
1232 // find an initial guess of a gcd
1233 for (int i=1; i<=bucket->buckets_used;i++)
1234 {
1235 if (bucket->buckets[i]!=NULL)
1236 {
1237 number t=p_InitContent(bucket->buckets[i],r);
1238 if (n_Size(t,cf)<2)
1239 {
1240 n_Delete(&t,cf);
1241 n_Delete(&coef,cf);
1242 return;
1243 }
1244 number t2=n_SubringGcd(coef,t,cf);
1245 n_Delete(&t,cf);
1246 n_Delete(&coef,cf);
1247 coef=t2;
1248 if (n_Size(coef,cf)<2) { n_Delete(&coef,cf);return;}
1249 }
1250 }
1251 // find the gcd
1252 for (int i=0; i<=bucket->buckets_used;i++)
1253 {
1254 if (bucket->buckets[i]!=NULL)
1255 {
1256 poly p=bucket->buckets[i];
1257 while(p!=NULL)
1258 {
1259 number t=n_SubringGcd(coef,pGetCoeff(p),cf);
1260 if (n_Size(t,cf)<2)
1261 {
1262 n_Delete(&t,cf);
1263 n_Delete(&coef,cf);
1264 return;
1265 }
1266 pIter(p);
1267 }
1268 }
1269 }
1270 // divided by the gcd
1271 if (TEST_OPT_PROT) PrintS("@");
1272 for (int i=bucket->buckets_used;i>=0;i--)
1273 {
1274 if (bucket->buckets[i]!=NULL)
1275 {
1276 poly p=bucket->buckets[i];
1277 while(p!=NULL)
1278 {
1279 number d = n_ExactDiv(pGetCoeff(p),coef,cf);
1280 p_SetCoeff(p,d,r);
1281 pIter(p);
1282 }
1283 }
1284 }
1285 n_Delete(&coef,cf);
1286 }
1287 #else
nIsPseudoUnit(number n,ring r)1288 static BOOLEAN nIsPseudoUnit(number n, ring r)
1289 {
1290 if (rField_is_Zp(r))
1291 return TRUE;
1292
1293 if (rParameter(r)==NULL)
1294 {
1295 return (n_Size(n,r->cf)==1);
1296 }
1297 //if (r->parameter!=NULL)
1298 return (n_IsOne(n,r->cf) || n_IsMOne(n,r->cf));
1299 }
1300
kBucketSimpleContent(kBucket_pt bucket)1301 void kBucketSimpleContent(kBucket_pt bucket)
1302 {
1303 ring r=bucket->bucket_ring;
1304 int i;
1305 //PrintS("HHHHHHHHHHHHH");
1306 for (i=0;i<=MAX_BUCKET;i++)
1307 {
1308 //if ((bucket->buckets[i]!=NULL) && (bucket->coef[i]!=NULL))
1309 // PrintS("H2H2H2");
1310 if (i==0)
1311 {
1312 assume(bucket->buckets[i]==NULL);
1313 }
1314 if ((bucket->buckets[i]!=NULL) && (bucket->coef[i]==NULL))
1315 return;
1316 }
1317 for (i=0;i<=MAX_BUCKET;i++)
1318 {
1319 //if ((bucket->buckets[i]!=NULL) && (bucket->coef[i]!=NULL))
1320 // PrintS("H2H2H2");
1321 if (i==0)
1322 {
1323 assume(bucket->buckets[i]==NULL);
1324 }
1325 if ((bucket->buckets[i]!=NULL)
1326 && (nIsPseudoUnit(p_GetCoeff(bucket->coef[i],r),r)))
1327 return;
1328 }
1329 //return;
1330
1331 number coef=n_Init(0,r);
1332 //ATTENTION: will not work correct for GB over ring
1333 //if (TEST_OPT_PROT)
1334 // PrintS("CCCCCCCCCCCCC");
1335 for (i=MAX_BUCKET;i>=0;i--)
1336 {
1337 if (i==0)
1338 {
1339 assume(bucket->buckets[i]==NULL);
1340 }
1341 if (bucket->buckets[i]!=NULL)
1342 {
1343 assume(bucket->coef[i]!=NULL);
1344 assume(!(n_IsZero(pGetCoeff(bucket->coef[i]),r)));
1345
1346 //in this way it should crash on programming errors, yeah
1347 number temp=n_Gcd(coef, pGetCoeff(bucket->coef[i]),r);
1348 n_Delete(&coef,r );
1349 coef=temp;
1350 if (nIsPseudoUnit(coef,r))
1351 {
1352 n_Delete(&coef,r);
1353 return;
1354 }
1355 assume(!(n_IsZero(coef,r)));
1356 }
1357 }
1358 if (n_IsZero(coef,r))
1359 {
1360 n_Delete(&coef,r);
1361 return;
1362 }
1363 if (TEST_OPT_PROT)
1364 PrintS("S");
1365 for(i=0;i<=MAX_BUCKET;i++)
1366 {
1367 if (bucket->buckets[i]!=NULL)
1368 {
1369 assume(!(n_IsZero(coef,r)));
1370 assume(bucket->coef[i]!=NULL);
1371 number lc=p_GetCoeff(bucket->coef[i],r);
1372 p_SetCoeff(bucket->coef[i], n_ExactDiv(lc,coef,r),r);
1373 assume(!(n_IsZero(p_GetCoeff(bucket->coef[i],r),r)));
1374 }
1375 }
1376 n_Delete(&coef,r);
1377 }
1378 #endif
1379
1380
kBucketExtractLmOfBucket(kBucket_pt bucket,int i)1381 poly kBucketExtractLmOfBucket(kBucket_pt bucket, int i)
1382 {
1383 assume(bucket->buckets[i]!=NULL);
1384
1385 poly p=bucket->buckets[i];
1386 bucket->buckets_length[i]--;
1387 #ifdef USE_COEF_BUCKETS
1388 ring r=bucket->bucket_ring;
1389 if (bucket->coef[i]!=NULL)
1390 {
1391 poly next=pNext(p);
1392 if (next==NULL)
1393 {
1394 MULTIPLY_BUCKET(bucket,i);
1395 p=bucket->buckets[i];
1396 bucket->buckets[i]=NULL;
1397 return p;
1398 }
1399 else
1400 {
1401 bucket->buckets[i]=next;
1402 number c=p_GetCoeff(bucket->coef[i],r);
1403 pNext(p)=NULL;
1404 p=__p_Mult_nn(p,c,r);
1405 assume(p!=NULL);
1406 return p;
1407 }
1408 }
1409 else
1410 #endif
1411 {
1412 bucket->buckets[i]=pNext(bucket->buckets[i]);
1413 pNext(p)=NULL;
1414 assume(p!=NULL);
1415 return p;
1416 }
1417 }
1418
1419 /*
1420 * input - output: a, b
1421 * returns:
1422 * a := a/gcd(a,b), b := b/gcd(a,b)
1423 * and return value
1424 * 0 -> a != 1, b != 1
1425 * 1 -> a == 1, b != 1
1426 * 2 -> a != 1, b == 1
1427 * 3 -> a == 1, b == 1
1428 * this value is used to control the spolys
1429 */
ksCheckCoeff(number * a,number * b,const coeffs r)1430 int ksCheckCoeff(number *a, number *b, const coeffs r)
1431 {
1432 int c = 0;
1433 number an = *a, bn = *b;
1434 n_Test(an,r);
1435 n_Test(bn,r);
1436
1437 number cn = n_SubringGcd(an, bn, r);
1438
1439 if(n_IsOne(cn, r))
1440 {
1441 an = n_Copy(an, r);
1442 bn = n_Copy(bn, r);
1443 }
1444 else
1445 {
1446 an = n_ExactDiv(an, cn, r); n_Normalize(an,r);
1447 bn = n_ExactDiv(bn, cn, r); n_Normalize(bn,r);
1448 }
1449 n_Delete(&cn, r);
1450 if (n_IsOne(an, r))
1451 {
1452 c = 1;
1453 }
1454 if (n_IsOne(bn, r))
1455 {
1456 c += 2;
1457 }
1458 *a = an;
1459 *b = bn;
1460 return c;
1461 }
1462
1463