1 //! \file tgb.cc
2 // multiple rings
3 // shorten_tails und dessen Aufrufe pruefen wlength!!!
4 /****************************************
5 * Computer Algebra System SINGULAR *
6 ****************************************/
7 /*
8 * ABSTRACT: slimgb and F4 implementation
9 */
10 //#include <vector>
11 //using namespace std;
12
13 ///@TODO: delay nur auf Sugarvergroesserung
14 ///@TODO: grade aus ecartS, setze dazu strat->honey; und nutze p.ecart
15 ///@TODO: no tail reductions in syz comp
16 #include "kernel/mod2.h"
17
18 #include "kernel/GBEngine/tgb.h"
19 #include "kernel/GBEngine/tgb_internal.h"
20 #include "kernel/GBEngine/tgbgauss.h"
21
22 #include "misc/options.h"
23 #include "kernel/digitech.h"
24 #include "polys/nc/nc.h"
25 #include "polys/nc/sca.h"
26 #include "polys/prCopy.h"
27
28 #include "coeffs/longrat.h" // nlQlogSize
29
30 #include <stdlib.h>
31 #include <stdio.h>
32 #include <queue>
33
34 #define BUCKETS_FOR_NORO_RED 1
35 #define SR_HDL(A) ((long)(A))
36 static const int bundle_size = 100;
37 static const int bundle_size_noro = 10000;
38 static const int delay_factor = 3;
39 #define ADD_LATER_SIZE 500
40 #if 1
41 STATIC_VAR omBin lm_bin = NULL;
42 static int add_to_reductors(slimgb_alg* c, poly h, int len, int ecart, BOOLEAN simplified=FALSE);
43 static void multi_reduction(red_object* los, int & losl, slimgb_alg* c);
44 static void multi_reduce_step(find_erg & erg, red_object* r, slimgb_alg* c);
45 static BOOLEAN extended_product_criterion(poly p1, poly gcd1, poly p2, poly gcd2, slimgb_alg* c);
46 static poly gcd_of_terms(poly p, ring r);
47 static int tgb_pair_better_gen(const void* ap,const void* bp);
48 static BOOLEAN pair_better(sorted_pair_node* a,sorted_pair_node* b, slimgb_alg* c=NULL);
49 static BOOLEAN state_is(calc_state state, const int & i, const int & j, slimgb_alg* c);
50 static void super_clean_top_of_pair_list(slimgb_alg* c);
51 static int simple_posInS (kStrategy strat, poly p,int len, wlen_type wlen);
52 static int* make_connections(int from, int to, poly bound, slimgb_alg* c);
53 static BOOLEAN has_t_rep(const int & arg_i, const int & arg_j, slimgb_alg* state);
54 static void shorten_tails(slimgb_alg* c, poly monom);
55 static poly redNF2 (poly h,slimgb_alg* c , int &len, number& m,int n=0);
56 static poly redNFTail (poly h,const int sl,kStrategy strat, int len);
57 static int bucket_guess(kBucket* bucket);
58
simplify_poly(poly p,ring r)59 static void simplify_poly (poly p, ring r)
60 {
61 assume (r == currRing);
62 if(TEST_OPT_INTSTRATEGY)
63 {
64 p_Cleardenom (p, r);
65 //includes p_Content(p,r);
66 }
67 else
68 pNorm (p);
69 }
70
71 //static const BOOLEAN up_to_radical=TRUE;
72
slim_nsize(number n,ring r)73 int slim_nsize (number n, ring r)
74 {
75 if(rField_is_Zp (r))
76 {
77 return 1;
78 }
79 if(rField_is_Q (r))
80 {
81 return nlQlogSize (n, r->cf);
82 }
83 else
84 {
85 return n_Size (n, r->cf);
86 }
87 }
88
monomial_root(poly m,ring r)89 static BOOLEAN monomial_root (poly m, ring r)
90 {
91 BOOLEAN changed = FALSE;
92 int i;
93 for(i = 1; i <= rVar (r); i++)
94 {
95 int e = p_GetExp (m, i, r);
96 if(e > 1)
97 {
98 p_SetExp (m, i, 1, r);
99 changed = TRUE;
100 }
101 }
102 if(changed)
103 {
104 p_Setm (m, r);
105 }
106 return changed;
107 }
108
polynomial_root(poly h,ring r)109 static BOOLEAN polynomial_root (poly h, ring r)
110 {
111 poly got = gcd_of_terms (h, r);
112 BOOLEAN changed = FALSE;
113 if((got != NULL) && (TEST_V_UPTORADICAL))
114 {
115 poly copy = p_Copy (got, r);
116 //p_wrp(got,c->r);
117 changed = monomial_root (got, r);
118 if(changed)
119 {
120 poly div_by = pMDivide (copy, got);
121 poly iter = h;
122 while(iter)
123 {
124 pExpVectorSub (iter, div_by);
125 pIter (iter);
126 }
127 p_Delete (&div_by, r);
128 if(TEST_OPT_PROT)
129 PrintS ("U");
130 }
131 p_Delete (©, r);
132 }
133 p_Delete (&got, r);
134 return changed;
135 }
136
p_Init_Special(const ring r)137 static inline poly p_Init_Special (const ring r)
138 {
139 return p_Init (r, lm_bin);
140 }
141
pOne_Special(const ring r=currRing)142 static inline poly pOne_Special (const ring r = currRing)
143 {
144 poly rc = p_Init_Special (r);
145 pSetCoeff0 (rc, n_Init (1, r->cf));
146 return rc;
147 }
148
149 // zum Initialiseren: in t_rep_gb plazieren:
150
151 #endif
152 #define LEN_VAR3
153 #define degbound(p) assume(pTotaldegree(p)<10)
154 //#define inDebug(p) assume((debug_Ideal==NULL)||(kNF(debug_Ideal,NULL,p,0,0)==0))
155
156 //die meisten Varianten stossen sich an coef_buckets
157
158 #ifdef LEN_VAR1
159 // erste Variante: Laenge: Anzahl der Monome
pSLength(poly p,int l)160 static inline int pSLength (poly p, int l)
161 {
162 return l;
163 }
164
kSBucketLength(kBucket * bucket,poly lm)165 static inline int kSBucketLength (kBucket * bucket, poly lm)
166 {
167 return bucket_guess (bucket);
168 }
169 #endif
170
171 #ifdef LEN_VAR2
172 // 2. Variante: Laenge: Platz fuer die Koeff.
pSLength(poly p,int l)173 int pSLength (poly p, int l)
174 {
175 int s = 0;
176 while(p != NULL)
177 {
178 s += nSize (pGetCoeff (p));
179 pIter (p);
180 }
181 return s;
182 }
183
kSBucketLength(kBucket * b,poly lm)184 int kSBucketLength (kBucket * b, poly lm)
185 {
186 int s = 0;
187 int i;
188 for(i = MAX_BUCKET; i >= 0; i--)
189 {
190 s += pSLength (b->buckets[i], 0);
191 }
192 return s;
193 }
194 #endif
195
196 #ifdef LEN_VAR3
pSLength(poly p,int l)197 static inline wlen_type pSLength (poly p, int l)
198 {
199 wlen_type c;
200 number coef = pGetCoeff (p);
201 if(rField_is_Q (currRing))
202 {
203 c = nlQlogSize (coef, currRing->cf);
204 }
205 else
206 c = nSize (coef);
207 if(!(TEST_V_COEFSTRAT))
208 {
209 return (wlen_type) c *(wlen_type) l /*pLength(p) */ ;
210 }
211 else
212 {
213 wlen_type res = l;
214 res *= c;
215 res *= c;
216 return res;
217 }
218 }
219
220 //! TODO CoefBuckets bercksichtigen
kSBucketLength(kBucket * b,poly lm=NULL)221 wlen_type kSBucketLength (kBucket * b, poly lm = NULL)
222 {
223 int s = 0;
224 wlen_type c;
225 number coef;
226 if(lm == NULL)
227 coef = pGetCoeff (kBucketGetLm (b));
228 //c=nSize(pGetCoeff(kBucketGetLm(b)));
229 else
230 coef = pGetCoeff (lm);
231 //c=nSize(pGetCoeff(lm));
232 if(rField_is_Q (currRing))
233 {
234 c = nlQlogSize (coef, currRing->cf);
235 }
236 else
237 c = nSize (coef);
238
239 int i;
240 for(i = b->buckets_used; i >= 0; i--)
241 {
242 assume ((b->buckets_length[i] == 0) || (b->buckets[i] != NULL));
243 s += b->buckets_length[i] /*pLength(b->buckets[i]) */ ;
244 }
245 #ifdef HAVE_COEF_BUCKETS
246 assume (b->buckets[0] == kBucketGetLm (b));
247 if(b->coef[0] != NULL)
248 {
249 if(rField_is_Q (currRing))
250 {
251 int modifier = nlQlogSize (pGetCoeff (b->coef[0]), currRing->cf);
252 c += modifier;
253 }
254 else
255 {
256 int modifier = nSize (pGetCoeff (b->coef[0]));
257 c *= modifier;
258 }
259 }
260 #endif
261 if(!(TEST_V_COEFSTRAT))
262 {
263 return s * c;
264 }
265 else
266 {
267 wlen_type res = s;
268 res *= c;
269 res *= c;
270 return res;
271 }
272 }
273 #endif
274 #ifdef LEN_VAR5
pSLength(poly p,int l)275 static inline wlen_type pSLength (poly p, int l)
276 {
277 int c;
278 number coef = pGetCoeff (p);
279 if(rField_is_Q (currRing))
280 {
281 c = nlQlogSize (coef, currRing->cf);
282 }
283 else
284 c = nSize (coef);
285 wlen_type erg = l;
286 erg *= c;
287 erg *= c;
288 //PrintS("lenvar 5");
289 assume (erg >= 0);
290 return erg; /*pLength(p) */ ;
291 }
292
293 //! TODO CoefBuckets beruecksichtigen
kSBucketLength(kBucket * b,poly lm=NULL)294 wlen_type kSBucketLength (kBucket * b, poly lm = NULL)
295 {
296 wlen_type s = 0;
297 int c;
298 number coef;
299 if(lm == NULL)
300 coef = pGetCoeff (kBucketGetLm (b));
301 //c=nSize(pGetCoeff(kBucketGetLm(b)));
302 else
303 coef = pGetCoeff (lm);
304 //c=nSize(pGetCoeff(lm));
305 if(rField_is_Q (currRing))
306 {
307 c = nlQlogSize (coef, currRing->cf);
308 }
309 else
310 c = nSize (coef);
311
312 int i;
313 for(i = b->buckets_used; i >= 0; i--)
314 {
315 assume ((b->buckets_length[i] == 0) || (b->buckets[i] != NULL));
316 s += b->buckets_length[i] /*pLength(b->buckets[i]) */ ;
317 }
318 #ifdef HAVE_COEF_BUCKETS
319 assume (b->buckets[0] == kBucketGetLm (b));
320 if(b->coef[0] != NULL)
321 {
322 if(rField_is_Q (currRing))
323 {
324 int modifier = nlQlogSize (pGetCoeff (b->coef[0]), currRing->cf);
325 c += modifier;
326 }
327 else
328 {
329 int modifier = nSize (pGetCoeff (b->coef[0]));
330 c *= modifier;
331 }
332 }
333 #endif
334 wlen_type erg = s;
335 erg *= c;
336 erg *= c;
337 return erg;
338 }
339 #endif
340
341 #ifdef LEN_VAR4
342 // 4.Variante: Laenge: Platz fuer Leitk * (1+Platz fuer andere Koeff.)
pSLength(poly p,int l)343 int pSLength (poly p, int l)
344 {
345 int s = 1;
346 int c = nSize (pGetCoeff (p));
347 pIter (p);
348 while(p != NULL)
349 {
350 s += nSize (pGetCoeff (p));
351 pIter (p);
352 }
353 return s * c;
354 }
355
kSBucketLength(kBucket * b)356 int kSBucketLength (kBucket * b)
357 {
358 int s = 1;
359 int c = nSize (pGetCoeff (kBucketGetLm (b)));
360 int i;
361 for(i = MAX_BUCKET; i > 0; i--)
362 {
363 if(b->buckets[i] == NULL)
364 continue;
365 s += pSLength (b->buckets[i], 0);
366 }
367 return s * c;
368 }
369 #endif
370 //BUG/TODO this stuff will fail on internal Schreyer orderings
elength_is_normal_length(poly p,slimgb_alg * c)371 static BOOLEAN elength_is_normal_length (poly p, slimgb_alg * c)
372 {
373 ring r = c->r;
374 if(p_GetComp (p, r) != 0)
375 return FALSE;
376 if(c->lastDpBlockStart <= (currRing->N))
377 {
378 int i;
379 for(i = 1; i < c->lastDpBlockStart; i++)
380 {
381 if(p_GetExp (p, i, r) != 0)
382 {
383 break;
384 }
385 }
386 if(i >= c->lastDpBlockStart)
387 {
388 //wrp(p);
389 //PrintS("\n");
390 return TRUE;
391 }
392 else
393 return FALSE;
394 }
395 else
396 return FALSE;
397 }
398
lies_in_last_dp_block(poly p,slimgb_alg * c)399 static BOOLEAN lies_in_last_dp_block (poly p, slimgb_alg * c)
400 {
401 ring r = c->r;
402 if(p_GetComp (p, r) != 0)
403 return FALSE;
404 if(c->lastDpBlockStart <= (currRing->N))
405 {
406 int i;
407 for(i = 1; i < c->lastDpBlockStart; i++)
408 {
409 if(p_GetExp (p, i, r) != 0)
410 {
411 break;
412 }
413 }
414 if(i >= c->lastDpBlockStart)
415 {
416 //wrp(p);
417 //PrintS("\n");
418 return TRUE;
419 }
420 else
421 return FALSE;
422 }
423 else
424 return FALSE;
425 }
426
get_last_dp_block_start(ring r)427 static int get_last_dp_block_start (ring r)
428 {
429 //ring r=c->r;
430 int last_block;
431
432 if(rRing_has_CompLastBlock (r))
433 {
434 last_block = rBlocks (r) - 3;
435 }
436 else
437 {
438 last_block = rBlocks (r) - 2;
439 }
440 assume (last_block >= 0);
441 if(r->order[last_block] == ringorder_dp)
442 return r->block0[last_block];
443 return (currRing->N) + 1;
444 }
445
do_pELength(poly p,slimgb_alg * c,int dlm=-1)446 static wlen_type do_pELength (poly p, slimgb_alg * c, int dlm = -1)
447 {
448 if(p == NULL)
449 return 0;
450 wlen_type s = 0;
451 poly pi = p;
452 if(dlm < 0)
453 {
454 dlm = c->pTotaldegree (p);
455 s = 1;
456 pi = p->next;
457 }
458
459 while(pi)
460 {
461 int d = c->pTotaldegree (pi);
462 if(d > dlm)
463 s += 1 + d - dlm;
464 else
465 ++s;
466 pi = pi->next;
467 }
468 return s;
469 }
470
pELength(poly p,slimgb_alg * c,ring)471 wlen_type pELength (poly p, slimgb_alg * c, ring /*r*/)
472 {
473 if(p == NULL)
474 return 0;
475 wlen_type s = 0;
476 poly pi = p;
477 int dlm;
478 dlm = c->pTotaldegree (p);
479 s = 1;
480 pi = p->next;
481
482 while(pi)
483 {
484 int d = c->pTotaldegree (pi);
485 if(d > dlm)
486 s += 1 + d - dlm;
487 else
488 ++s;
489 pi = pi->next;
490 }
491 return s;
492 }
493
kEBucketLength(kBucket * b,poly lm,slimgb_alg * ca)494 wlen_type kEBucketLength (kBucket * b, poly lm, slimgb_alg * ca)
495 {
496 wlen_type s = 0;
497 if(lm == NULL)
498 {
499 lm = kBucketGetLm (b);
500 }
501 if(lm == NULL)
502 return 0;
503 if(elength_is_normal_length (lm, ca))
504 {
505 return bucket_guess (b);
506 }
507 int d = ca->pTotaldegree (lm);
508 #if 0
509 assume (sugar >= d);
510 s = 1 + (bucket_guess (b) - 1) * (sugar - d + 1);
511 return s;
512 #else
513
514 //int d=pTotaldegree(lm,ca->r);
515 int i;
516 for(i = b->buckets_used; i >= 0; i--)
517 {
518 if(b->buckets[i] == NULL)
519 continue;
520
521 if((ca->pTotaldegree (b->buckets[i]) <= d)
522 && (elength_is_normal_length (b->buckets[i], ca)))
523 {
524 s += b->buckets_length[i];
525 }
526 else
527 {
528 s += do_pELength (b->buckets[i], ca, d);
529 }
530 }
531 return s;
532 #endif
533 }
534
pELength(poly p,slimgb_alg * c,int l)535 static inline int pELength (poly p, slimgb_alg * c, int l)
536 {
537 if(p == NULL)
538 return 0;
539 if((l > 0) && (elength_is_normal_length (p, c)))
540 return l;
541 return do_pELength (p, c);
542 }
543
pQuality(poly p,slimgb_alg * c,int l=-1)544 static inline wlen_type pQuality (poly p, slimgb_alg * c, int l = -1)
545 {
546 if(l < 0)
547 l = pLength (p);
548 if(c->isDifficultField)
549 {
550 if(c->eliminationProblem)
551 {
552 wlen_type cs;
553 number coef = pGetCoeff (p);
554 if(rField_is_Q (currRing))
555 {
556 cs = nlQlogSize (coef, currRing->cf);
557 }
558 else
559 cs = nSize (coef);
560 wlen_type erg = cs;
561 if(TEST_V_COEFSTRAT)
562 erg *= cs;
563 //erg*=cs;//for quadratic
564 erg *= pELength (p, c, l);
565 //FIXME: not quadratic coeff size
566 //return cs*pELength(p,c,l);
567 return erg;
568 }
569 //PrintS("I am here");
570 wlen_type r = pSLength (p, l);
571 assume (r >= 0);
572 return r;
573 }
574 if(c->eliminationProblem)
575 return pELength (p, c, l);
576 return l;
577 }
578
pTotaldegree_full(poly p)579 static inline int pTotaldegree_full (poly p)
580 {
581 int r = 0;
582 while(p)
583 {
584 int d = pTotaldegree (p);
585 r = si_max (r, d);
586 pIter (p);
587 }
588 return r;
589 }
590
guess_quality(slimgb_alg * c)591 wlen_type red_object::guess_quality (slimgb_alg * c)
592 {
593 //works at the moment only for lenvar 1, because in different
594 //case, you have to look on coefs
595 wlen_type s = 0;
596 if(c->isDifficultField)
597 {
598 //s=kSBucketLength(bucket,this->p);
599 if(c->eliminationProblem)
600 {
601 wlen_type cs;
602 number coef;
603
604 coef = pGetCoeff (kBucketGetLm (bucket));
605 //c=nSize(pGetCoeff(kBucketGetLm(b)));
606
607 //c=nSize(pGetCoeff(lm));
608 if(rField_is_Q (currRing))
609 {
610 cs = nlQlogSize (coef, currRing->cf);
611 }
612 else
613 cs = nSize (coef);
614 #ifdef HAVE_COEF_BUCKETS
615 if(bucket->coef[0] != NULL)
616 {
617 if(rField_is_Q (currRing))
618 {
619 int modifier = nlQlogSize (pGetCoeff (bucket->coef[0]), currRing->cf);
620 cs += modifier;
621 }
622 else
623 {
624 int modifier = nSize (pGetCoeff (bucket->coef[0]));
625 cs *= modifier;
626 }
627 }
628 #endif
629 //FIXME:not quadratic
630 wlen_type erg = kEBucketLength (this->bucket, this->p, c);
631 //erg*=cs;//for quadratic
632 erg *= cs;
633 if(TEST_V_COEFSTRAT)
634 erg *= cs;
635 //return cs*kEBucketLength(this->bucket,this->p,c);
636 return erg;
637 }
638 s = kSBucketLength (bucket, NULL);
639 }
640 else
641 {
642 if(c->eliminationProblem)
643 //if (false)
644 s = kEBucketLength (this->bucket, this->p, c);
645 else
646 s = bucket_guess (bucket);
647 }
648 return s;
649 }
650
651 #if 0 //currently unused
652 static void finalize_reduction_step (reduction_step * r)
653 {
654 delete r;
655 }
656 #endif
657 #if 0 //currently unused
658 static int LObject_better_gen (const void *ap, const void *bp)
659 {
660 LObject *a = *(LObject **) ap;
661 LObject *b = *(LObject **) bp;
662 return (pLmCmp (a->p, b->p));
663 }
664 #endif
red_object_better_gen(const void * ap,const void * bp)665 static int red_object_better_gen (const void *ap, const void *bp)
666 {
667 return (pLmCmp (((red_object *) ap)->p, ((red_object *) bp)->p));
668 }
669
670 #if 0 //currently unused
671 static int pLmCmp_func_inverted (const void *ap1, const void *ap2)
672 {
673 poly p1, p2;
674 p1 = *((poly *) ap1);
675 p2 = *((poly *) ap2);
676 return -pLmCmp (p1, p2);
677 }
678 #endif
679
tgb_pair_better_gen2(const void * ap,const void * bp)680 int tgb_pair_better_gen2 (const void *ap, const void *bp)
681 {
682 return (-tgb_pair_better_gen (ap, bp));
683 }
684
kFindDivisibleByInS_easy(kStrategy strat,const red_object & obj)685 int kFindDivisibleByInS_easy (kStrategy strat, const red_object & obj)
686 {
687 poly p = obj.p;
688 if ((strat->syzComp>0) && (pGetComp(p)>strat->syzComp)) return -1;
689 long not_sev = ~obj.sev;
690 for(int i = 0; i <= strat->sl; i++)
691 {
692 if(pLmShortDivisibleBy (strat->S[i], strat->sevS[i], p, not_sev))
693 return i;
694 }
695 return -1;
696 }
697
kFindDivisibleByInS_easy(kStrategy strat,poly p,long sev)698 int kFindDivisibleByInS_easy (kStrategy strat, poly p, long sev)
699 {
700 if ((strat->syzComp>0) && (pGetComp(p)>strat->syzComp)) return -1;
701 long not_sev = ~sev;
702 for(int i = 0; i <= strat->sl; i++)
703 {
704 if(pLmShortDivisibleBy (strat->S[i], strat->sevS[i], p, not_sev))
705 return i;
706 }
707 return -1;
708 }
709
710 static int
posInPairs(sorted_pair_node ** p,int pn,sorted_pair_node * qe,slimgb_alg * c,int an=0)711 posInPairs (sorted_pair_node ** p, int pn, sorted_pair_node * qe,
712 slimgb_alg * c, int an = 0)
713 {
714 if(pn == 0)
715 return 0;
716
717 int length = pn - 1;
718 int i;
719 //int an = 0;
720 int en = length;
721
722 if(pair_better (qe, p[en], c))
723 return length + 1;
724
725 while(1)
726 {
727 //if (an >= en-1)
728 if(en - 1 <= an)
729 {
730 if(pair_better (p[an], qe, c))
731 return an;
732 return en;
733 }
734 i = (an + en) / 2;
735 if(pair_better (p[i], qe, c))
736 en = i;
737 else
738 an = i;
739 }
740 }
741
ascending(int * i,int top)742 static BOOLEAN ascending (int *i, int top)
743 {
744 if(top < 1)
745 return TRUE;
746 if(i[top] < i[top - 1])
747 return FALSE;
748 return ascending (i, top - 1);
749 }
750
spn_merge(sorted_pair_node ** p,int pn,sorted_pair_node ** q,int qn,slimgb_alg * c)751 sorted_pair_node **spn_merge (sorted_pair_node ** p, int pn,
752 sorted_pair_node ** q, int qn, slimgb_alg * c)
753 {
754 int i;
755 int *a = (int *) omalloc (qn * sizeof (int));
756 // int mc;
757 // PrintS("Debug\n");
758 // for(mc=0;mc<qn;mc++)
759 // {
760 // wrp(q[mc]->lcm_of_lm);
761 // PrintS("\n");
762 // }
763 // PrintS("Debug they are in\n");
764 // for(mc=0;mc<pn;mc++)
765 // {
766 // wrp(p[mc]->lcm_of_lm);
767 // PrintS("\n");
768 // }
769 int lastpos = 0;
770 for(i = 0; i < qn; i++)
771 {
772 lastpos = posInPairs (p, pn, q[i], c, si_max (lastpos - 1, 0));
773 // cout<<lastpos<<"\n";
774 a[i] = lastpos;
775 }
776 if((pn + qn) > c->max_pairs)
777 {
778 p =
779 (sorted_pair_node **) omrealloc (p,
780 2 * (pn +
781 qn) *
782 sizeof (sorted_pair_node *));
783 c->max_pairs = 2 * (pn + qn);
784 }
785 for(i = qn - 1; i >= 0; i--)
786 {
787 size_t size;
788 if(qn - 1 > i)
789 size = (a[i + 1] - a[i]) * sizeof (sorted_pair_node *);
790 else
791 size = (pn - a[i]) * sizeof (sorted_pair_node *); //as indices begin with 0
792 memmove (p + a[i] + (1 + i), p + a[i], size);
793 p[a[i] + i] = q[i];
794 }
795 omFree (a);
796 return p;
797 }
798
799 static BOOLEAN
trivial_syzygie(int pos1,int pos2,poly bound,slimgb_alg * c)800 trivial_syzygie (int pos1, int pos2, poly bound, slimgb_alg * c)
801 {
802 poly p1 = c->S->m[pos1];
803 poly p2 = c->S->m[pos2];
804
805 if(pGetComp (p1) > 0 || pGetComp (p2) > 0)
806 return FALSE;
807 int i = 1;
808 poly m = NULL;
809 poly gcd1 = c->gcd_of_terms[pos1];
810 poly gcd2 = c->gcd_of_terms[pos2];
811
812 if((gcd1 != NULL) && (gcd2 != NULL))
813 {
814 gcd1->next = gcd2; //may ordered incorrect
815 m = gcd_of_terms (gcd1, c->r);
816 gcd1->next = NULL;
817 }
818 if(m == NULL)
819 {
820 loop
821 {
822 if(pGetExp (p1, i) + pGetExp (p2, i) > pGetExp (bound, i))
823 return FALSE;
824 if(i == (currRing->N))
825 {
826 //PrintS("trivial");
827 return TRUE;
828 }
829 i++;
830 }
831 }
832 else
833 {
834 loop
835 {
836 if(pGetExp (p1, i) - pGetExp (m, i) + pGetExp (p2, i) >
837 pGetExp (bound, i))
838 {
839 pDelete (&m);
840 return FALSE;
841 }
842 if(i == (currRing->N))
843 {
844 pDelete (&m);
845 //PrintS("trivial");
846 return TRUE;
847 }
848 i++;
849 }
850 }
851 }
852
853 //! returns position sets w as weight
find_best(red_object * r,int l,int u,wlen_type & w,slimgb_alg * c)854 int find_best (red_object * r, int l, int u, wlen_type & w, slimgb_alg * c)
855 {
856 int best = l;
857 int i;
858 w = r[l].guess_quality (c);
859 for(i = l + 1; i <= u; i++)
860 {
861 wlen_type w2 = r[i].guess_quality (c);
862 if(w2 < w)
863 {
864 w = w2;
865 best = i;
866 }
867 }
868 return best;
869 }
870
canonicalize()871 void red_object::canonicalize ()
872 {
873 kBucketCanonicalize (bucket);
874 }
875
good_has_t_rep(int i,int j,slimgb_alg * c)876 BOOLEAN good_has_t_rep (int i, int j, slimgb_alg * c)
877 {
878 assume (i >= 0);
879 assume (j >= 0);
880 if(has_t_rep (i, j, c))
881 return TRUE;
882 //poly lm=pOne();
883 assume (c->tmp_lm != NULL);
884 poly lm = c->tmp_lm;
885
886 pLcm (c->S->m[i], c->S->m[j], lm);
887 pSetm (lm);
888 assume (lm != NULL);
889 //int deciding_deg= pTotaldegree(lm);
890 int *i_con = make_connections (i, j, lm, c);
891 //p_Delete(&lm,c->r);
892
893 for(int n = 0; ((n < c->n) && (i_con[n] >= 0)); n++)
894 {
895 if(i_con[n] == j)
896 {
897 now_t_rep (i, j, c);
898 omFree (i_con);
899 return TRUE;
900 }
901 }
902 omFree (i_con);
903
904 return FALSE;
905 }
906
lenS_correct(kStrategy strat)907 BOOLEAN lenS_correct (kStrategy strat)
908 {
909 int i;
910 for(i = 0; i <= strat->sl; i++)
911 {
912 if(strat->lenS[i] != pLength (strat->S[i]))
913 return FALSE;
914 }
915 return TRUE;
916 }
917
918
cleanS(kStrategy strat,slimgb_alg * c)919 static void cleanS (kStrategy strat, slimgb_alg * c)
920 {
921 int i = 0;
922 LObject P;
923 while(i <= strat->sl)
924 {
925 P.p = strat->S[i];
926 P.sev = strat->sevS[i];
927 //int dummy=strat->sl;
928 //if(kFindDivisibleByInS(strat,&dummy,&P)!=i)
929 if(kFindDivisibleByInS_easy (strat, P.p, P.sev) != i)
930 {
931 deleteInS (i, strat);
932 //remember destroying poly
933 BOOLEAN found = FALSE;
934 int j;
935 for(j = 0; j < c->n; j++)
936 {
937 if(c->S->m[j] == P.p)
938 {
939 found = TRUE;
940 break;
941 }
942 }
943 if(!found)
944 pDelete (&P.p);
945 //remember additional reductors
946 }
947 else
948 i++;
949 }
950 }
951
bucket_guess(kBucket * bucket)952 static int bucket_guess (kBucket * bucket)
953 {
954 int sum = 0;
955 int i;
956 for(i = bucket->buckets_used; i >= 0; i--)
957 {
958 if(bucket->buckets[i])
959 sum += bucket->buckets_length[i];
960 }
961 return sum;
962 }
963
964 static int
add_to_reductors(slimgb_alg * c,poly h,int len,int ecart,BOOLEAN simplified)965 add_to_reductors (slimgb_alg * c, poly h, int len, int ecart,
966 BOOLEAN simplified)
967 {
968 //inDebug(h);
969 assume (lenS_correct (c->strat));
970 assume (len == pLength (h));
971 int i;
972 // if (c->isDifficultField)
973 // i=simple_posInS(c->strat,h,pSLength(h,len),c->isDifficultField);
974 // else
975 // i=simple_posInS(c->strat,h,len,c->isDifficultField);
976
977 LObject P;
978 memset (&P, 0, sizeof (P));
979 P.tailRing = c->r;
980 P.p = h; /*p_Copy(h,c->r); */
981 P.ecart = ecart;
982 P.FDeg = c->r->pFDeg (P.p, c->r);
983 if(!(simplified))
984 {
985 if(TEST_OPT_INTSTRATEGY)
986 {
987 p_Cleardenom (P.p, c->r); //includes p_Content(P.p,c->r );
988 }
989 else
990 pNorm (P.p);
991 //pNormalize (P.p);
992 }
993 wlen_type pq = pQuality (h, c, len);
994 i = simple_posInS (c->strat, h, len, pq);
995 c->strat->enterS (P, i, c->strat, -1);
996
997 c->strat->lenS[i] = len;
998 assume (pLength (c->strat->S[i]) == c->strat->lenS[i]);
999 if(c->strat->lenSw != NULL)
1000 c->strat->lenSw[i] = pq;
1001
1002 return i;
1003 }
1004
length_one_crit(slimgb_alg * c,int pos,int len)1005 static void length_one_crit (slimgb_alg * c, int pos, int len)
1006 {
1007 if(c->nc)
1008 return;
1009 if(len == 1)
1010 {
1011 int i;
1012 for(i = 0; i < pos; i++)
1013 {
1014 if(c->lengths[i] == 1)
1015 c->states[pos][i] = HASTREP;
1016 }
1017 for(i = pos + 1; i < c->n; i++)
1018 {
1019 if(c->lengths[i] == 1)
1020 c->states[i][pos] = HASTREP;
1021 }
1022 if(!c->nc)
1023 shorten_tails (c, c->S->m[pos]);
1024 }
1025 }
1026
move_forward_in_S(int old_pos,int new_pos,kStrategy strat)1027 static void move_forward_in_S (int old_pos, int new_pos, kStrategy strat)
1028 {
1029 assume (old_pos >= new_pos);
1030 poly p = strat->S[old_pos];
1031 int ecart = strat->ecartS[old_pos];
1032 long sev = strat->sevS[old_pos];
1033 int s_2_r = strat->S_2_R[old_pos];
1034 int length = strat->lenS[old_pos];
1035 assume (length == pLength (strat->S[old_pos]));
1036 wlen_type length_w;
1037 if(strat->lenSw != NULL)
1038 length_w = strat->lenSw[old_pos];
1039 int i;
1040 for(i = old_pos; i > new_pos; i--)
1041 {
1042 strat->S[i] = strat->S[i - 1];
1043 strat->ecartS[i] = strat->ecartS[i - 1];
1044 strat->sevS[i] = strat->sevS[i - 1];
1045 strat->S_2_R[i] = strat->S_2_R[i - 1];
1046 }
1047 if(strat->lenS != NULL)
1048 for(i = old_pos; i > new_pos; i--)
1049 strat->lenS[i] = strat->lenS[i - 1];
1050 if(strat->lenSw != NULL)
1051 for(i = old_pos; i > new_pos; i--)
1052 strat->lenSw[i] = strat->lenSw[i - 1];
1053
1054 strat->S[new_pos] = p;
1055 strat->ecartS[new_pos] = ecart;
1056 strat->sevS[new_pos] = sev;
1057 strat->S_2_R[new_pos] = s_2_r;
1058 strat->lenS[new_pos] = length;
1059 if(strat->lenSw != NULL)
1060 strat->lenSw[new_pos] = length_w;
1061 //assume(lenS_correct(strat));
1062 }
1063
move_backward_in_S(int old_pos,int new_pos,kStrategy strat)1064 static void move_backward_in_S (int old_pos, int new_pos, kStrategy strat)
1065 {
1066 assume (old_pos <= new_pos);
1067 poly p = strat->S[old_pos];
1068 int ecart = strat->ecartS[old_pos];
1069 long sev = strat->sevS[old_pos];
1070 int s_2_r = strat->S_2_R[old_pos];
1071 int length = strat->lenS[old_pos];
1072 assume (length == pLength (strat->S[old_pos]));
1073 wlen_type length_w;
1074 if(strat->lenSw != NULL)
1075 length_w = strat->lenSw[old_pos];
1076 int i;
1077 for(i = old_pos; i < new_pos; i++)
1078 {
1079 strat->S[i] = strat->S[i + 1];
1080 strat->ecartS[i] = strat->ecartS[i + 1];
1081 strat->sevS[i] = strat->sevS[i + 1];
1082 strat->S_2_R[i] = strat->S_2_R[i + 1];
1083 }
1084 if(strat->lenS != NULL)
1085 for(i = old_pos; i < new_pos; i++)
1086 strat->lenS[i] = strat->lenS[i + 1];
1087 if(strat->lenSw != NULL)
1088 for(i = old_pos; i < new_pos; i++)
1089 strat->lenSw[i] = strat->lenSw[i + 1];
1090
1091 strat->S[new_pos] = p;
1092 strat->ecartS[new_pos] = ecart;
1093 strat->sevS[new_pos] = sev;
1094 strat->S_2_R[new_pos] = s_2_r;
1095 strat->lenS[new_pos] = length;
1096 if(strat->lenSw != NULL)
1097 strat->lenSw[new_pos] = length_w;
1098 //assume(lenS_correct(strat));
1099 }
1100
make_connections(int from,int to,poly bound,slimgb_alg * c)1101 static int *make_connections (int from, int to, poly bound, slimgb_alg * c)
1102 {
1103 ideal I = c->S;
1104 int *cans = (int *) omAlloc (c->n * sizeof (int));
1105 int *connected = (int *) omAlloc (c->n * sizeof (int));
1106 cans[0] = to;
1107 int cans_length = 1;
1108 connected[0] = from;
1109 int last_cans_pos = -1;
1110 int connected_length = 1;
1111 long neg_bounds_short = ~p_GetShortExpVector (bound, c->r);
1112
1113 int not_yet_found = cans_length;
1114 int con_checked = 0;
1115 int pos;
1116
1117 while(TRUE)
1118 {
1119 if((con_checked < connected_length) && (not_yet_found > 0))
1120 {
1121 pos = connected[con_checked];
1122 for(int i = 0; i < cans_length; i++)
1123 {
1124 if(cans[i] < 0)
1125 continue;
1126 //FIXME: triv. syz. does not hold on noncommutative, check it for modules
1127 if((has_t_rep (pos, cans[i], c))
1128 || ((!rIsPluralRing (c->r))
1129 && (trivial_syzygie (pos, cans[i], bound, c))))
1130 {
1131 connected[connected_length] = cans[i];
1132 connected_length++;
1133 cans[i] = -1;
1134 --not_yet_found;
1135
1136 if(connected[connected_length - 1] == to)
1137 {
1138 if(connected_length < c->n)
1139 {
1140 connected[connected_length] = -1;
1141 }
1142 omFree (cans);
1143 return connected;
1144 }
1145 }
1146 }
1147 con_checked++;
1148 }
1149 else
1150 {
1151 for(last_cans_pos++; last_cans_pos <= c->n; last_cans_pos++)
1152 {
1153 if(last_cans_pos == c->n)
1154 {
1155 if(connected_length < c->n)
1156 {
1157 connected[connected_length] = -1;
1158 }
1159 omFree (cans);
1160 return connected;
1161 }
1162 if((last_cans_pos == from) || (last_cans_pos == to))
1163 continue;
1164 if(p_LmShortDivisibleBy
1165 (I->m[last_cans_pos], c->short_Exps[last_cans_pos], bound,
1166 neg_bounds_short, c->r))
1167 {
1168 cans[cans_length] = last_cans_pos;
1169 cans_length++;
1170 break;
1171 }
1172 }
1173 not_yet_found++;
1174 for(int i = 0; i < con_checked; i++)
1175 {
1176 if(has_t_rep (connected[i], last_cans_pos, c))
1177 {
1178 connected[connected_length] = last_cans_pos;
1179 connected_length++;
1180 cans[cans_length - 1] = -1;
1181 --not_yet_found;
1182 if(connected[connected_length - 1] == to)
1183 {
1184 if(connected_length < c->n)
1185 {
1186 connected[connected_length] = -1;
1187 }
1188 omFree (cans);
1189 return connected;
1190 }
1191 break;
1192 }
1193 }
1194 }
1195 }
1196 if(connected_length < c->n)
1197 {
1198 connected[connected_length] = -1;
1199 }
1200 omFree (cans);
1201 return connected;
1202 }
1203
1204 #ifdef HEAD_BIN
p_MoveHead(poly p,omBin b)1205 static inline poly p_MoveHead (poly p, omBin b)
1206 {
1207 poly np;
1208 omTypeAllocBin (poly, np, b);
1209 memmove (np, p, omSizeWOfBin(b) * sizeof (long));
1210 omFreeBinAddr (p);
1211 return np;
1212 }
1213 #endif
1214
replace_pair(int & i,int & j,slimgb_alg * c)1215 static void replace_pair (int &i, int &j, slimgb_alg * c)
1216 {
1217 if(i < 0)
1218 return;
1219 c->soon_free = NULL;
1220 int syz_deg;
1221 poly lm = pOne ();
1222
1223 pLcm (c->S->m[i], c->S->m[j], lm);
1224 pSetm (lm);
1225
1226 int *i_con = make_connections (i, j, lm, c);
1227
1228 for(int n = 0; ((n < c->n) && (i_con[n] >= 0)); n++)
1229 {
1230 if(i_con[n] == j)
1231 {
1232 now_t_rep (i, j, c);
1233 omFree (i_con);
1234 p_Delete (&lm, c->r);
1235 return;
1236 }
1237 }
1238
1239 int *j_con = make_connections (j, i, lm, c);
1240
1241 // if(c->n>1)
1242 // {
1243 // if (i_con[1]>=0)
1244 // i=i_con[1];
1245 // else
1246 // {
1247 // if (j_con[1]>=0)
1248 // j=j_con[1];
1249 // }
1250 // }
1251
1252 int sugar = syz_deg = c->pTotaldegree (lm);
1253
1254 p_Delete (&lm, c->r);
1255 if(c->T_deg_full) //Sugar
1256 {
1257 int t_i = c->T_deg_full[i] - c->T_deg[i];
1258 int t_j = c->T_deg_full[j] - c->T_deg[j];
1259 sugar += si_max (t_i, t_j);
1260 //Print("\n max: %d\n",max(t_i,t_j));
1261 }
1262
1263 for(int m = 0; ((m < c->n) && (i_con[m] >= 0)); m++)
1264 {
1265 if(c->T_deg_full != NULL)
1266 {
1267 int s1 = c->T_deg_full[i_con[m]] + syz_deg - c->T_deg[i_con[m]];
1268 if(s1 > sugar)
1269 continue;
1270 }
1271 if(c->weighted_lengths[i_con[m]] < c->weighted_lengths[i])
1272 i = i_con[m];
1273 }
1274 for(int m = 0; ((m < c->n) && (j_con[m] >= 0)); m++)
1275 {
1276 if(c->T_deg_full != NULL)
1277 {
1278 int s1 = c->T_deg_full[j_con[m]] + syz_deg - c->T_deg[j_con[m]];
1279 if(s1 > sugar)
1280 continue;
1281 }
1282 if(c->weighted_lengths[j_con[m]] < c->weighted_lengths[j])
1283 j = j_con[m];
1284 }
1285
1286 //can also try dependend search
1287 omFree (i_con);
1288 omFree (j_con);
1289 return;
1290 }
1291
add_later(poly p,const char * prot,slimgb_alg * c)1292 static void add_later (poly p, const char *prot, slimgb_alg * c)
1293 {
1294 int i = 0;
1295 //check, if it is already in the queue
1296
1297 while(c->add_later->m[i] != NULL)
1298 {
1299 if(p_LmEqual (c->add_later->m[i], p, c->r))
1300 return;
1301 i++;
1302 }
1303 if(TEST_OPT_PROT)
1304 PrintS (prot);
1305 c->add_later->m[i] = p;
1306 }
1307
simple_posInS(kStrategy strat,poly p,int len,wlen_type wlen)1308 static int simple_posInS (kStrategy strat, poly p, int len, wlen_type wlen)
1309 {
1310 if(strat->sl == -1)
1311 return 0;
1312 if(strat->lenSw)
1313 return pos_helper (strat, p, (wlen_type) wlen, (wlen_set) strat->lenSw,
1314 strat->S);
1315 return pos_helper (strat, p, len, strat->lenS, strat->S);
1316 }
1317
1318 /*2
1319 *if the leading term of p
1320 *divides the leading term of some S[i] it will be canceled
1321 */
1322 static inline void
clearS(poly p,unsigned long p_sev,int l,int * at,int * k,kStrategy strat)1323 clearS (poly p, unsigned long p_sev, int l, int *at, int *k, kStrategy strat)
1324 {
1325 assume (p_sev == pGetShortExpVector (p));
1326 if(!pLmShortDivisibleBy (p, p_sev, strat->S[*at], ~strat->sevS[*at]))
1327 return;
1328 if(l >= strat->lenS[*at])
1329 return;
1330 if(TEST_OPT_PROT)
1331 PrintS ("!");
1332 mflush ();
1333 //pDelete(&strat->S[*at]);
1334 deleteInS ((*at), strat);
1335 (*at)--;
1336 (*k)--;
1337 // assume(lenS_correct(strat));
1338 }
1339
iq_crit(const void * ap,const void * bp)1340 static int iq_crit (const void *ap, const void *bp)
1341 {
1342 sorted_pair_node *a = *((sorted_pair_node **) ap);
1343 sorted_pair_node *b = *((sorted_pair_node **) bp);
1344 assume (a->i > a->j);
1345 assume (b->i > b->j);
1346
1347 if(a->deg < b->deg)
1348 return -1;
1349 if(a->deg > b->deg)
1350 return 1;
1351 int comp = pLmCmp (a->lcm_of_lm, b->lcm_of_lm);
1352 if(comp != 0)
1353 return comp;
1354 if(a->expected_length < b->expected_length)
1355 return -1;
1356 if(a->expected_length > b->expected_length)
1357 return 1;
1358 if(a->j > b->j)
1359 return 1;
1360 if(a->j < b->j)
1361 return -1;
1362 return 0;
1363 }
1364
coeff_mult_size_estimate(int s1,int s2,ring r)1365 static wlen_type coeff_mult_size_estimate (int s1, int s2, ring r)
1366 {
1367 if(rField_is_Q (r))
1368 return s1 + s2;
1369 else
1370 return s1 * s2;
1371 }
1372
pair_weighted_length(int i,int j,slimgb_alg * c)1373 static wlen_type pair_weighted_length (int i, int j, slimgb_alg * c)
1374 {
1375 if((c->isDifficultField) && (c->eliminationProblem))
1376 {
1377 int c1 = slim_nsize (p_GetCoeff (c->S->m[i], c->r), c->r);
1378 int c2 = slim_nsize (p_GetCoeff (c->S->m[j], c->r), c->r);
1379 wlen_type el1 = c->weighted_lengths[i] / c1;
1380 assume (el1 != 0);
1381 assume (c->weighted_lengths[i] % c1 == 0);
1382 wlen_type el2 = c->weighted_lengths[j] / c2;
1383 assume (el2 != 0);
1384 //assume (c->weighted_lengths[j] % c2 == 0); // fails in Tst/Plural/dmod_lib.tst
1385 //should be * for function fields
1386 //return (c1+c2) * (el1+el2-2);
1387 wlen_type res = coeff_mult_size_estimate (c1, c2, c->r);
1388 res *= el1 + el2 - 2;
1389 return res;
1390
1391 }
1392 if(c->isDifficultField)
1393 {
1394 //int cs=slim_nsize(p_GetCoeff(c->S->m[i],c->r),c->r)+
1395 // slim_nsize(p_GetCoeff(c->S->m[j],c->r),c->r);
1396 if(!(TEST_V_COEFSTRAT))
1397 {
1398 wlen_type cs =
1399 coeff_mult_size_estimate (slim_nsize
1400 (p_GetCoeff (c->S->m[i], c->r), c->r),
1401 slim_nsize (p_GetCoeff (c->S->m[j], c->r),
1402 c->r), c->r);
1403 return (wlen_type) (c->lengths[i] + c->lengths[j] - 2) * (wlen_type) cs;
1404 }
1405 else
1406 {
1407
1408 wlen_type cs =
1409 coeff_mult_size_estimate (slim_nsize
1410 (p_GetCoeff (c->S->m[i], c->r), c->r),
1411 slim_nsize (p_GetCoeff (c->S->m[j], c->r),
1412 c->r), c->r);
1413 cs *= cs;
1414 return (wlen_type) (c->lengths[i] + c->lengths[j] - 2) * (wlen_type) cs;
1415 }
1416 }
1417 if(c->eliminationProblem)
1418 {
1419
1420 return (c->weighted_lengths[i] + c->weighted_lengths[j] - 2);
1421 }
1422 return c->lengths[i] + c->lengths[j] - 2;
1423
1424 }
1425
add_to_basis_ideal_quotient(poly h,slimgb_alg * c,int * ip)1426 sorted_pair_node **add_to_basis_ideal_quotient (poly h, slimgb_alg * c,
1427 int *ip)
1428 {
1429 p_Test (h, c->r);
1430 assume (h != NULL);
1431 poly got = gcd_of_terms (h, c->r);
1432 if((got != NULL) && (TEST_V_UPTORADICAL))
1433 {
1434 poly copy = p_Copy (got, c->r);
1435 //p_wrp(got,c->r);
1436 BOOLEAN changed = monomial_root (got, c->r);
1437 if(changed)
1438 {
1439 poly div_by = pMDivide (copy, got);
1440 poly iter = h;
1441 while(iter)
1442 {
1443 pExpVectorSub (iter, div_by);
1444 pIter (iter);
1445 }
1446 p_Delete (&div_by, c->r);
1447 PrintS ("U");
1448 }
1449 p_Delete (©, c->r);
1450 }
1451
1452 #define ENLARGE(pointer, type) pointer=(type*) omrealloc(pointer, c->array_lengths*sizeof(type))
1453
1454 #define ENLARGE_ALIGN(pointer, type) {if(pointer)\
1455 pointer=(type*)omReallocAligned(pointer, c->array_lengths*sizeof(type));\
1456 else pointer=(type*)omAllocAligned(c->array_lengths*sizeof(type));}
1457 // BOOLEAN corr=lenS_correct(c->strat);
1458 int sugar;
1459 int ecart = 0;
1460 ++(c->n);
1461 ++(c->S->ncols);
1462 int i, j;
1463 i = c->n - 1;
1464 sorted_pair_node **nodes =
1465 (sorted_pair_node **) omalloc (sizeof (sorted_pair_node *) * i);
1466 int spc = 0;
1467 if(c->n > c->array_lengths)
1468 {
1469 c->array_lengths = c->array_lengths * 2;
1470 assume (c->array_lengths >= c->n);
1471 ENLARGE (c->T_deg, int);
1472 ENLARGE_ALIGN (c->tmp_pair_lm, poly);
1473 ENLARGE_ALIGN (c->tmp_spn, sorted_pair_node *);
1474
1475 ENLARGE_ALIGN (c->short_Exps, long);
1476 ENLARGE (c->lengths, int);
1477 #ifndef HAVE_BOOST
1478 #ifndef USE_STDVECBOOL
1479
1480 ENLARGE_ALIGN (c->states, char *);
1481 #endif
1482 #endif
1483 ENLARGE_ALIGN (c->gcd_of_terms, poly);
1484 //if (c->weighted_lengths!=NULL) {
1485 ENLARGE_ALIGN (c->weighted_lengths, wlen_type);
1486 //}
1487 //ENLARGE_ALIGN(c->S->m,poly);
1488 }
1489 pEnlargeSet (&c->S->m, c->n - 1, 1);
1490 if(c->T_deg_full)
1491 ENLARGE (c->T_deg_full, int);
1492 sugar = c->T_deg[i] = c->pTotaldegree (h);
1493 if(c->T_deg_full)
1494 {
1495 sugar = c->T_deg_full[i] = c->pTotaldegree_full (h);
1496 ecart = sugar - c->T_deg[i];
1497 assume (ecart >= 0);
1498 }
1499 c->tmp_pair_lm[i] = pOne_Special (c->r);
1500
1501 c->tmp_spn[i] = (sorted_pair_node *) omAlloc (sizeof (sorted_pair_node));
1502
1503 c->lengths[i] = pLength (h);
1504
1505 //necessary for correct weighted length
1506
1507 if(TEST_OPT_INTSTRATEGY)
1508 {
1509 p_Cleardenom (h, c->r); //includes p_Content(h,c->r);
1510 }
1511 else
1512 pNorm (h);
1513 //pNormalize (h);
1514
1515 c->weighted_lengths[i] = pQuality (h, c, c->lengths[i]);
1516 c->gcd_of_terms[i] = got;
1517 #ifdef HAVE_BOOST
1518 c->states.push_back (dynamic_bitset <> (i));
1519
1520 #else
1521 #ifdef USE_STDVECBOOL
1522
1523 c->states.push_back (vector < bool > (i));
1524
1525 #else
1526 if(i > 0)
1527 c->states[i] = (char *) omAlloc (i * sizeof (char));
1528 else
1529 c->states[i] = NULL;
1530 #endif
1531 #endif
1532
1533 c->S->m[i] = h;
1534 c->short_Exps[i] = p_GetShortExpVector (h, c->r);
1535
1536 #undef ENLARGE
1537 #undef ENLARGE_ALIGN
1538 if(p_GetComp (h, currRing) <= c->syz_comp)
1539 {
1540 for(j = 0; j < i; j++)
1541 {
1542
1543
1544 #ifndef HAVE_BOOST
1545 c->states[i][j] = UNCALCULATED;
1546 #endif
1547 assume (p_LmDivisibleBy (c->S->m[i], c->S->m[j], c->r) ==
1548 p_LmShortDivisibleBy (c->S->m[i], c->short_Exps[i], c->S->m[j],
1549 ~(c->short_Exps[j]), c->r));
1550
1551 if(__p_GetComp (c->S->m[i], c->r) != __p_GetComp (c->S->m[j], c->r))
1552 {
1553 //c->states[i][j]=UNCALCULATED;
1554 //WARNUNG: be careful
1555 continue;
1556 }
1557 else if((!c->nc) && (c->lengths[i] == 1) && (c->lengths[j] == 1))
1558 {
1559 c->states[i][j] = HASTREP;
1560 }
1561 else if(((!c->nc) || (c->is_homog && rIsSCA (c->r)))
1562 && (pHasNotCF (c->S->m[i], c->S->m[j])))
1563 // else if ((!(c->nc)) && (pHasNotCF(c->S->m[i],c->S->m[j])))
1564 {
1565 c->easy_product_crit++;
1566 c->states[i][j] = HASTREP;
1567 continue;
1568 }
1569 else
1570 if(extended_product_criterion
1571 (c->S->m[i], c->gcd_of_terms[i], c->S->m[j], c->gcd_of_terms[j],
1572 c))
1573 {
1574 c->states[i][j] = HASTREP;
1575 c->extended_product_crit++;
1576 //PrintS("E");
1577 }
1578 // if (c->states[i][j]==UNCALCULATED)
1579 // {
1580
1581 if((TEST_V_FINDMONOM) && (!c->nc))
1582 {
1583 //PrintS("COMMU");
1584 // if (c->lengths[i]==c->lengths[j])
1585 // {
1586 // poly short_s=ksCreateShortSpoly(c->S->m[i],c->S->m[j],c->r);
1587 // if (short_s==NULL)
1588 // {
1589 // c->states[i][j]=HASTREP;
1590 // }
1591 // else
1592 // {
1593 // p_Delete(&short_s, currRing);
1594 // }
1595 // }
1596 if(c->lengths[i] + c->lengths[j] == 3)
1597 {
1598
1599
1600 poly short_s = ksCreateShortSpoly (c->S->m[i], c->S->m[j], c->r);
1601 if(short_s == NULL)
1602 {
1603 c->states[i][j] = HASTREP;
1604 }
1605 else
1606 {
1607 assume (pLength (short_s) == 1);
1608 if(TEST_V_UPTORADICAL)
1609 monomial_root (short_s, c->r);
1610 int iS = kFindDivisibleByInS_easy (c->strat, short_s,
1611 p_GetShortExpVector (short_s,
1612 c->r));
1613 if(iS < 0)
1614 {
1615 //PrintS("N");
1616 if(TRUE)
1617 {
1618 c->states[i][j] = HASTREP;
1619 add_later (short_s, "N", c);
1620 }
1621 else
1622 p_Delete (&short_s, currRing);
1623 }
1624 else
1625 {
1626 if(c->strat->lenS[iS] > 1)
1627 {
1628 //PrintS("O");
1629 if(TRUE)
1630 {
1631 c->states[i][j] = HASTREP;
1632 add_later (short_s, "O", c);
1633 }
1634 else
1635 p_Delete (&short_s, currRing);
1636 }
1637 else
1638 p_Delete (&short_s, currRing);
1639 c->states[i][j] = HASTREP;
1640 }
1641
1642
1643 }
1644 }
1645 }
1646 // if (short_s)
1647 // {
1648 assume (spc <= j);
1649 sorted_pair_node *s = c->tmp_spn[spc]; //(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1650 s->i = si_max (i, j);
1651 s->j = si_min (i, j);
1652 assume (s->j == j);
1653 s->expected_length = pair_weighted_length (i, j, c); //c->lengths[i]+c->lengths[j]-2;
1654
1655 poly lm = c->tmp_pair_lm[spc]; //=pOne_Special();
1656
1657 pLcm (c->S->m[i], c->S->m[j], lm);
1658 pSetm (lm);
1659 p_Test (lm, c->r);
1660 s->deg = c->pTotaldegree (lm);
1661
1662 if(c->T_deg_full) //Sugar
1663 {
1664 int t_i = c->T_deg_full[s->i] - c->T_deg[s->i];
1665 int t_j = c->T_deg_full[s->j] - c->T_deg[s->j];
1666 s->deg += si_max (t_i, t_j);
1667 //Print("\n max: %d\n",max(t_i,t_j));
1668 }
1669 p_Test (lm, c->r);
1670 s->lcm_of_lm = lm;
1671 // pDelete(&short_s);
1672 //assume(lm!=NULL);
1673 nodes[spc] = s;
1674 spc++;
1675
1676 // }
1677 //else
1678 //{
1679 //c->states[i][j]=HASTREP;
1680 //}
1681 }
1682 } //if syz_comp end
1683
1684 assume (spc <= i);
1685 //now ideal quotient crit
1686 qsort (nodes, spc, sizeof (sorted_pair_node *), iq_crit);
1687
1688 sorted_pair_node **nodes_final =
1689 (sorted_pair_node **) omalloc (sizeof (sorted_pair_node *) * (i+1));
1690 int spc_final = 0;
1691 j = 0;
1692 while(j < spc)
1693 {
1694 int lower = j;
1695 int upper;
1696 BOOLEAN has = FALSE;
1697 for(upper = lower + 1; upper < spc; upper++)
1698 {
1699 if(!pLmEqual (nodes[lower]->lcm_of_lm, nodes[upper]->lcm_of_lm))
1700 {
1701 break;
1702 }
1703 if(has_t_rep (nodes[upper]->i, nodes[upper]->j, c))
1704 has = TRUE;
1705 }
1706 upper = upper - 1;
1707 int z;
1708 assume (spc_final <= j);
1709 for(z = 0; z < spc_final; z++)
1710 {
1711 if(p_LmDivisibleBy
1712 (nodes_final[z]->lcm_of_lm, nodes[lower]->lcm_of_lm, c->r))
1713 {
1714 has = TRUE;
1715 break;
1716 }
1717 }
1718
1719 if(has)
1720 {
1721 for(; lower <= upper; lower++)
1722 {
1723 //free_sorted_pair_node(nodes[lower],c->r);
1724 //omfree(nodes[lower]);
1725 nodes[lower] = NULL;
1726 }
1727 j = upper + 1;
1728 continue;
1729 }
1730 else
1731 {
1732 p_Test (nodes[lower]->lcm_of_lm, c->r);
1733 nodes[lower]->lcm_of_lm = pCopy (nodes[lower]->lcm_of_lm);
1734 assume (__p_GetComp (c->S->m[nodes[lower]->i], c->r) ==
1735 __p_GetComp (c->S->m[nodes[lower]->j], c->r));
1736 nodes_final[spc_final] =
1737 (sorted_pair_node *) omAlloc (sizeof (sorted_pair_node));
1738
1739 *(nodes_final[spc_final++]) = *(nodes[lower]);
1740 //c->tmp_spn[nodes[lower]->j]=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
1741 nodes[lower] = NULL;
1742 for(lower = lower + 1; lower <= upper; lower++)
1743 {
1744 // free_sorted_pair_node(nodes[lower],c->r);
1745 //omfree(nodes[lower]);
1746 nodes[lower] = NULL;
1747 }
1748 j = upper + 1;
1749 continue;
1750 }
1751 }
1752
1753 // Print("i:%d,spc_final:%d",i,spc_final);
1754
1755 assume (spc_final <= spc);
1756 omFree (nodes);
1757 nodes = NULL;
1758
1759 add_to_reductors (c, h, c->lengths[c->n - 1], ecart, TRUE);
1760 //i=posInS(c->strat,c->strat->sl,h,0 ecart);
1761 if(!(c->nc))
1762 {
1763 if(c->lengths[c->n - 1] == 1)
1764 shorten_tails (c, c->S->m[c->n - 1]);
1765 }
1766 //you should really update c->lengths, c->strat->lenS, and the oder of polys in strat if you sort after lengths
1767
1768 //for(i=c->strat->sl; i>0;i--)
1769 // if(c->strat->lenS[i]<c->strat->lenS[i-1]) printf("fehler bei %d\n",i);
1770 if(c->Rcounter > 50)
1771 {
1772 c->Rcounter = 0;
1773 cleanS (c->strat, c);
1774 }
1775
1776 #ifdef HAVE_PLURAL
1777 // for SCA:
1778 // here write at the end of nodes_final[spc_final,...,spc_final+lmdeg-1]
1779 if(rIsSCA (c->r))
1780 {
1781 const poly pNext = pNext (h);
1782
1783 if(pNext != NULL)
1784 {
1785 // for additional polynomials
1786 const unsigned int m_iFirstAltVar = scaFirstAltVar (c->r);
1787 const unsigned int m_iLastAltVar = scaLastAltVar (c->r);
1788
1789 int N = // c->r->N;
1790 m_iLastAltVar - m_iFirstAltVar + 1; // should be enough
1791 // TODO: but we may also use got = gcd({m}_{m\in f}))!
1792
1793 poly *array_arg = (poly *) omalloc (N * sizeof (poly)); // !
1794 int j = 0;
1795
1796
1797 for(unsigned short v = m_iFirstAltVar; v <= m_iLastAltVar; v++)
1798 // for all x_v | Ann(lm(h))
1799 if(p_GetExp (h, v, c->r)) // TODO: use 'got' here!
1800 {
1801 assume (p_GetExp (h, v, c->r) == 1);
1802
1803 poly p = sca_pp_Mult_xi_pp (v, pNext, c->r); // x_v * h;
1804
1805 if(p != NULL) // if (x_v * h != 0)
1806 array_arg[j++] = p;
1807 } // for all x_v | Ann(lm(h))
1808
1809 c->introduceDelayedPairs (array_arg, j);
1810
1811 omFree (array_arg); // !!!
1812 }
1813 // PrintS("Saturation - done!!!\n");
1814 }
1815 #endif // if SCAlgebra
1816
1817
1818 if(!ip)
1819 {
1820 qsort (nodes_final, spc_final, sizeof (sorted_pair_node *),
1821 tgb_pair_better_gen2);
1822
1823
1824 c->apairs =
1825 spn_merge (c->apairs, c->pair_top + 1, nodes_final, spc_final, c);
1826 c->pair_top += spc_final;
1827 clean_top_of_pair_list (c);
1828 omFree (nodes_final);
1829 return NULL;
1830 }
1831 {
1832 *ip = spc_final;
1833 return nodes_final;
1834 }
1835 }
1836
redNF2(poly h,slimgb_alg * c,int & len,number & m,int n)1837 static poly redNF2 (poly h, slimgb_alg * c, int &len, number & m, int n)
1838 {
1839 m = nInit (1);
1840 if(h == NULL)
1841 return NULL;
1842
1843 assume (len == pLength (h));
1844 kStrategy strat = c->strat;
1845 if(0 > strat->sl)
1846 {
1847 return h;
1848 }
1849 int j;
1850
1851 LObject P (h);
1852 P.SetShortExpVector ();
1853 P.bucket = kBucketCreate (currRing);
1854 // BOOLEAN corr=lenS_correct(strat);
1855 kBucketInit (P.bucket, P.p, len /*pLength(P.p) */ );
1856 //wlen_set lenSw=(wlen_set) c->strat->lenS;
1857 //FIXME: plainly wrong
1858 //strat->lenS;
1859 //if (strat->lenSw!=NULL)
1860 // lenSw=strat->lenSw;
1861 //int max_pos=simple_posInS(strat,P.p);
1862 loop
1863 {
1864 //int dummy=strat->sl;
1865 j = kFindDivisibleByInS_easy (strat, P.p, P.sev);
1866 //j=kFindDivisibleByInS(strat,&dummy,&P);
1867 if((j >= 0) && ((!n) ||
1868 ((strat->lenS[j] <= n) &&
1869 ((strat->lenSw == NULL) || (strat->lenSw[j] <= n)))))
1870 {
1871 nNormalize (pGetCoeff (P.p));
1872 #ifdef KDEBUG
1873 if(TEST_OPT_DEBUG)
1874 {
1875 PrintS ("red:");
1876 wrp (h);
1877 PrintS (" with ");
1878 wrp (strat->S[j]);
1879 }
1880 #endif
1881
1882 number coef = kBucketPolyRed (P.bucket, strat->S[j],
1883 strat->lenS[j] /*pLength(strat->S[j]) */ ,
1884 strat->kNoether);
1885 number m2 = nMult (m, coef);
1886 nDelete (&m);
1887 m = m2;
1888 nDelete (&coef);
1889 h = kBucketGetLm (P.bucket);
1890
1891 if(h == NULL)
1892 {
1893 len = 0;
1894 kBucketDestroy (&P.bucket);
1895 return NULL;
1896 }
1897 P.p = h;
1898 P.t_p = NULL;
1899 P.SetShortExpVector ();
1900 #ifdef KDEBUG
1901 if(TEST_OPT_DEBUG)
1902 {
1903 PrintS ("\nto:");
1904 wrp (h);
1905 PrintLn ();
1906 }
1907 #endif
1908 }
1909 else
1910 {
1911 kBucketClear (P.bucket, &(P.p), &len);
1912 kBucketDestroy (&P.bucket);
1913 pNormalize (P.p);
1914 assume (len == (pLength (P.p)));
1915 return P.p;
1916 }
1917 }
1918 }
1919
redTailShort(poly h,kStrategy strat)1920 static poly redTailShort (poly h, kStrategy strat)
1921 {
1922 if(h == NULL)
1923 return NULL; //n_Init(1,currRing);
1924 if(TEST_V_MODPSOLVSB)
1925 {
1926 bit_reduce (pNext (h), strat->tailRing);
1927 }
1928 int i;
1929 int len = pLength (h);
1930 for(i = 0; i <= strat->sl; i++)
1931 {
1932 if((strat->lenS[i] > 2)
1933 || ((strat->lenSw != NULL) && (strat->lenSw[i] > 2)))
1934 break;
1935 }
1936 return (redNFTail (h, i - 1, strat, len));
1937 }
1938
line_of_extended_prod(int fixpos,slimgb_alg * c)1939 static void line_of_extended_prod (int fixpos, slimgb_alg * c)
1940 {
1941 if(c->gcd_of_terms[fixpos] == NULL)
1942 {
1943 c->gcd_of_terms[fixpos] = gcd_of_terms (c->S->m[fixpos], c->r);
1944 if(c->gcd_of_terms[fixpos])
1945 {
1946 int i;
1947 for(i = 0; i < fixpos; i++)
1948 if((c->states[fixpos][i] != HASTREP)
1949 &&
1950 (extended_product_criterion
1951 (c->S->m[fixpos], c->gcd_of_terms[fixpos], c->S->m[i],
1952 c->gcd_of_terms[i], c)))
1953 {
1954 c->states[fixpos][i] = HASTREP;
1955 c->extended_product_crit++;
1956 }
1957 for(i = fixpos + 1; i < c->n; i++)
1958 if((c->states[i][fixpos] != HASTREP)
1959 &&
1960 (extended_product_criterion
1961 (c->S->m[fixpos], c->gcd_of_terms[fixpos], c->S->m[i],
1962 c->gcd_of_terms[i], c)))
1963 {
1964 c->states[i][fixpos] = HASTREP;
1965 c->extended_product_crit++;
1966 }
1967 }
1968 }
1969 }
1970
c_S_element_changed_hook(int pos,slimgb_alg * c)1971 static void c_S_element_changed_hook (int pos, slimgb_alg * c)
1972 {
1973 length_one_crit (c, pos, c->lengths[pos]);
1974 if(!c->nc)
1975 line_of_extended_prod (pos, c);
1976 }
1977
1978 class poly_tree_node
1979 {
1980 public:
1981 poly p;
1982 poly_tree_node *l;
1983 poly_tree_node *r;
1984 int n;
poly_tree_node(int sn)1985 poly_tree_node (int sn):l (NULL), r (NULL), n (sn)
1986 {
1987 }
1988 };
1989 class exp_number_builder
1990 {
1991 public:
1992 poly_tree_node * top_level;
1993 int n;
1994 int get_n (poly p);
exp_number_builder()1995 exp_number_builder ():top_level (0), n (0)
1996 {
1997 }
1998 };
get_n(poly p)1999 int exp_number_builder::get_n (poly p)
2000 {
2001 poly_tree_node **node = &top_level;
2002 while(*node != NULL)
2003 {
2004 int c = pLmCmp (p, (*node)->p);
2005 if(c == 0)
2006 return (*node)->n;
2007 if(c == -1)
2008 node = &((*node)->r);
2009 else
2010 node = &((*node)->l);
2011 }
2012 (*node) = new poly_tree_node (n);
2013 n++;
2014 (*node)->p = pLmInit (p);
2015 return (*node)->n;
2016 }
2017
2018 //mac_polys exp are smaller iff they are greater by monomial ordering
2019 //corresponding to solving linear equations notation
2020
is_valid_ro(red_object & ro)2021 BOOLEAN is_valid_ro (red_object & ro)
2022 {
2023 red_object r2 = ro;
2024 ro.validate ();
2025 if((r2.p != ro.p) || (r2.sev != ro.sev))
2026 return FALSE;
2027 return TRUE;
2028 }
2029
terms_sort_crit(const void * a,const void * b)2030 int terms_sort_crit (const void *a, const void *b)
2031 {
2032 return -pLmCmp (*((poly *) a), *((poly *) b));
2033 }
2034
2035 #if 0 // currently unused
2036 static void unify_terms (poly * terms, int &sum)
2037 {
2038 if(sum == 0)
2039 return;
2040 int last = 0;
2041 int curr = 1;
2042 while(curr < sum)
2043 {
2044 if(!(pLmEqual (terms[curr], terms[last])))
2045 {
2046 terms[++last] = terms[curr];
2047 }
2048 ++curr;
2049 }
2050 sum = last + 1;
2051 }
2052 #endif
2053 #if 0 // currently unused
2054 static void
2055 export_mat (number * number_array, int pn, int tn, const char *format_str,
2056 int mat_nr)
2057 {
2058 char matname[20];
2059 sprintf (matname, format_str, mat_nr);
2060 FILE *out = fopen (matname, "w");
2061 int i, j;
2062 fprintf (out, "mat=[\n");
2063 for(i = 0; i < pn; i++)
2064 {
2065 fprintf (out, "[\n");
2066 for(j = 0; j < tn; j++)
2067 {
2068 if(j > 0)
2069 {
2070 fprintf (out, ", ");
2071 }
2072 fprintf (out, "%i", npInt (number_array[i * tn + j], currRing));
2073 }
2074 if(i < pn - 1)
2075 fprintf (out, "],\n");
2076 else
2077 fprintf (out, "],\n");
2078 }
2079 fprintf (out, "]\n");
2080 fclose (out);
2081 }
2082 #endif
2083 //typedef unsigned short number_type;
2084
2085
2086 #ifdef USE_NORO
2087 #ifndef NORO_CACHE
2088 static void
linalg_step_modp(poly * p,poly * p_out,int & pn,poly * terms,int tn,slimgb_alg * c)2089 linalg_step_modp (poly * p, poly * p_out, int &pn, poly * terms, int tn,
2090 slimgb_alg * c)
2091 {
2092 STATIC_VAR int export_n = 0;
2093 assume (terms[tn - 1] != NULL);
2094 assume (rField_is_Zp (c->r));
2095 //I don't do deletes, copies of number_types ...
2096 const number_type zero = 0; //npInit(0);
2097 int array_size = pn * tn;
2098 number_type *number_array =
2099 (number_type *) omalloc (pn * tn * sizeof (number_type));
2100 int i;
2101 for(i = 0; i < array_size; i++)
2102 {
2103 number_array[i] = zero;
2104 }
2105 for(i = 0; i < pn; i++)
2106 {
2107 poly h = p[i];
2108 //int base=tn*i;
2109 write_poly_to_row (number_array + tn * i, h, terms, tn, c->r);
2110
2111 }
2112 #if 0
2113 //export matrix
2114 export_mat (number_array, pn, tn, "mat%i.py", ++export_n);
2115 #endif
2116 int rank = pn;
2117 simplest_gauss_modp (number_array, rank, tn);
2118 int act_row = 0;
2119 int p_pos = 0;
2120 for(i = 0; i < pn; i++)
2121 {
2122 poly h = NULL;
2123 int j;
2124 int base = tn * i;
2125 number *row = number_array + base;
2126 h = row_to_poly (row, terms, tn, c->r);
2127
2128 if(h != NULL)
2129 {
2130 p_out[p_pos++] = h;
2131 }
2132 }
2133 pn = p_pos;
2134 //assert(p_pos==rank)
2135 while(p_pos < pn)
2136 {
2137 p_out[p_pos++] = NULL;
2138 }
2139 #if 0
2140 export_mat (number_array, pn, tn, "mat%i.py", ++export_n);
2141 #endif
2142 }
2143 #endif
2144 #endif
mass_add(poly * p,int pn,slimgb_alg * c)2145 static void mass_add (poly * p, int pn, slimgb_alg * c)
2146 {
2147 int j;
2148 int *ibuf = (int *) omalloc (pn * sizeof (int));
2149 sorted_pair_node ***sbuf =
2150 (sorted_pair_node ***) omalloc (pn * sizeof (sorted_pair_node **));
2151 for(j = 0; j < pn; j++)
2152 {
2153 p_Test (p[j], c->r);
2154 sbuf[j] = add_to_basis_ideal_quotient (p[j], c, ibuf + j);
2155 }
2156 int sum = 0;
2157 for(j = 0; j < pn; j++)
2158 {
2159 sum += ibuf[j];
2160 }
2161 sorted_pair_node **big_sbuf =
2162 (sorted_pair_node **) omalloc (sum * sizeof (sorted_pair_node *));
2163 int partsum = 0;
2164 for(j = 0; j < pn; j++)
2165 {
2166 memmove (big_sbuf + partsum, sbuf[j],
2167 ibuf[j] * sizeof (sorted_pair_node *));
2168 omFree (sbuf[j]);
2169 partsum += ibuf[j];
2170 }
2171
2172 qsort (big_sbuf, sum, sizeof (sorted_pair_node *), tgb_pair_better_gen2);
2173 c->apairs = spn_merge (c->apairs, c->pair_top + 1, big_sbuf, sum, c);
2174 c->pair_top += sum;
2175 clean_top_of_pair_list (c);
2176 omFree (big_sbuf);
2177 omFree (sbuf);
2178 omFree (ibuf);
2179 //omfree(buf);
2180 #ifdef TGB_DEBUG
2181 int z;
2182 for(z = 1; z <= c->pair_top; z++)
2183 {
2184 assume (pair_better (c->apairs[z], c->apairs[z - 1], c));
2185 }
2186 #endif
2187
2188 }
2189
2190 #ifdef NORO_CACHE
2191 #ifndef NORO_NON_POLY
evaluateRows()2192 void NoroCache::evaluateRows ()
2193 {
2194 //after that can evaluate placeholders
2195 int i;
2196 buffer = (number *) omAlloc (nIrreducibleMonomials * sizeof (number));
2197 for(i = 0; i < root.branches_len; i++)
2198 {
2199 evaluateRows (1, root.branches[i]);
2200 }
2201 omFree (buffer);
2202 buffer = NULL;
2203 }
2204
evaluateRows(int level,NoroCacheNode * node)2205 void NoroCache::evaluateRows (int level, NoroCacheNode * node)
2206 {
2207 assume (level >= 0);
2208 if(node == NULL)
2209 return;
2210 if(level < (currRing->N))
2211 {
2212 int i, sum;
2213 for(i = 0; i < node->branches_len; i++)
2214 {
2215 evaluateRows (level + 1, node->branches[i]);
2216 }
2217 }
2218 else
2219 {
2220 DataNoroCacheNode *dn = (DataNoroCacheNode *) node;
2221 if(dn->value_len != backLinkCode)
2222 {
2223 poly p = dn->value_poly;
2224 #ifndef NORO_SPARSE_ROWS_PRE
2225 dn->row = new DenseRow ();
2226 DenseRow *row = dn->row;
2227 memset (buffer, 0, sizeof (number) * nIrreducibleMonomials);
2228
2229 if(p == NULL)
2230 {
2231 row->array = NULL;
2232 row->begin = 0;
2233 row->end = 0;
2234 return;
2235 }
2236 int i = 0;
2237 int idx;
2238 number *a = buffer;
2239 while(p)
2240 {
2241 DataNoroCacheNode *ref = getCacheReference (p);
2242
2243 idx = ref->term_index;
2244 assume (idx >= 0);
2245 a[idx] = p_GetCoeff (p, currRing);
2246 if(i == 0)
2247 row->begin = idx;
2248 i++;
2249 pIter (p);
2250 }
2251 row->end = idx + 1;
2252 assume (row->end > row->begin);
2253 int len = row->end - row->begin;
2254 row->array = (number *) omalloc ((len) * sizeof (number));
2255 memcpy (row->array, a + row->begin, len * sizeof (number));
2256 #else
2257 assume (dn->value_len == pLength (dn->value_poly));
2258 dn->row = new SparseRow (dn->value_len);
2259 SparseRow *row = dn->row;
2260 int i = 0;
2261 while(p)
2262 {
2263 DataNoroCacheNode *ref = getCacheReference (p);
2264
2265 int idx = ref->term_index;
2266 assume (idx >= 0);
2267 row->idx_array[i] = idx;
2268 row->coef_array[i] = p_GetCoeff (p, currRing);
2269 i++;
2270 pIter (p);
2271 }
2272 if(i != dn->value_len)
2273 {
2274 PrintS ("F4 calc wrong, as poly len was wrong\n");
2275 }
2276 assume (i == dn->value_len);
2277 #endif
2278 }
2279 }
2280 }
2281
2282 void
evaluatePlaceHolder(number * row,std::vector<NoroPlaceHolder> & place_holders)2283 NoroCache::evaluatePlaceHolder (number * row,
2284 std::vector < NoroPlaceHolder >
2285 &place_holders)
2286 {
2287 int i;
2288 int s = place_holders.size ();
2289
2290 if (currRing->cf-ch<=NV_MAX_PRIME)
2291 {
2292 for(i = 0; i < s; i++)
2293 {
2294 DataNoroCacheNode *ref = place_holders[i].ref;
2295 number coef = place_holders[i].coef;
2296 if(ref->value_len == backLinkCode)
2297 {
2298 row[ref->term_index] = npAddM (row[ref->term_index], coef);
2299 }
2300 else
2301 {
2302 #ifndef NORO_SPARSE_ROWS_PRE
2303 DenseRow *ref_row = ref->row;
2304 if(ref_row == NULL)
2305 continue;
2306 number *ref_begin = ref_row->array;
2307 number *ref_end = ref_row->array + (ref_row->end - ref_row->begin);
2308 number *my_pos = row + ref_row->begin;
2309 //TODO npisOne distinction
2310 if(!(npIsOne (coef)))
2311 {
2312 while(ref_begin != ref_end)
2313 {
2314 *my_pos = npAddM (*my_pos, npMult (coef, *ref_begin));
2315 ++ref_begin;
2316 ++my_pos;
2317 }
2318 }
2319 else
2320 {
2321 while(ref_begin != ref_end)
2322 {
2323
2324 *my_pos = npAddM (*my_pos, *ref_begin);
2325 ++ref_begin;
2326 ++my_pos;
2327 }
2328 }
2329 #else
2330 SparseRow *ref_row = ref->row;
2331 if(ref_row == NULL)
2332 continue;
2333 int n = ref_row->len;
2334 int j;
2335 int *idx_array = ref_row->idx_array;
2336 number *coef_array = ref_row->coef_array;
2337 if(!(npIsOne (coef)))
2338 {
2339 for(j = 0; j < n; j++)
2340 {
2341 int idx = idx_array[j];
2342 number ref_coef = coef_array[j];
2343 row[idx] = npAddM (row[idx], npMult (coef, ref_coef));
2344 }
2345 }
2346 else
2347 {
2348 for(j = 0; j < n; j++)
2349 {
2350 int idx = idx_array[j];
2351 number ref_coef = coef_array[j];
2352 row[idx] = npAddM (row[idx], ref_coef);
2353 }
2354 }
2355 #endif
2356 }
2357 }
2358 }
2359 else /*ch >NV_MAX_PRIME */
2360 {
2361 for(i = 0; i < s; i++)
2362 {
2363 DataNoroCacheNode *ref = place_holders[i].ref;
2364 number coef = place_holders[i].coef;
2365 if(ref->value_len == backLinkCode)
2366 {
2367 row[ref->term_index] = npAddM (row[ref->term_index], coef);
2368 }
2369 else
2370 {
2371 #ifndef NORO_SPARSE_ROWS_PRE
2372 DenseRow *ref_row = ref->row;
2373 if(ref_row == NULL)
2374 continue;
2375 number *ref_begin = ref_row->array;
2376 number *ref_end = ref_row->array + (ref_row->end - ref_row->begin);
2377 number *my_pos = row + ref_row->begin;
2378 //TODO npisOne distinction
2379 if(!(npIsOne (coef)))
2380 {
2381 while(ref_begin != ref_end)
2382 {
2383 *my_pos = npAddM (*my_pos, nvMult (coef, *ref_begin));
2384 ++ref_begin;
2385 ++my_pos;
2386 }
2387 }
2388 else
2389 {
2390 while(ref_begin != ref_end)
2391 {
2392 *my_pos = npAddM (*my_pos, *ref_begin);
2393 ++ref_begin;
2394 ++my_pos;
2395 }
2396 }
2397 #else
2398 SparseRow *ref_row = ref->row;
2399 if(ref_row == NULL)
2400 continue;
2401 int n = ref_row->len;
2402 int j;
2403 int *idx_array = ref_row->idx_array;
2404 number *coef_array = ref_row->coef_array;
2405 if(!(npIsOne (coef)))
2406 {
2407 for(j = 0; j < n; j++)
2408 {
2409 int idx = idx_array[j];
2410 number ref_coef = coef_array[j];
2411 row[idx] = npAddM (row[idx], nvMult (coef, ref_coef));
2412 }
2413 }
2414 else
2415 {
2416 for(j = 0; j < n; j++)
2417 {
2418 int idx = idx_array[j];
2419 number ref_coef = coef_array[j];
2420 row[idx] = npAddM (row[idx], ref_coef);
2421 }
2422 }
2423 #endif
2424 }
2425 }
2426 }
2427 }
2428 #endif
2429
2430 //poly noro_red_non_unique(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2431
2432 #ifndef NORO_NON_POLY
2433 MonRedRes
noro_red_mon(poly t,BOOLEAN force_unique,NoroCache * cache,slimgb_alg * c)2434 noro_red_mon (poly t, BOOLEAN force_unique, NoroCache * cache, slimgb_alg * c)
2435 {
2436 MonRedRes res_holder;
2437
2438 //wrp(t);
2439 res_holder.changed = TRUE;
2440 if(force_unique)
2441 {
2442 DataNoroCacheNode *ref = cache->getCacheReference (t);
2443 if(ref != NULL)
2444 {
2445 res_holder.len = ref->value_len;
2446 if(res_holder.len == NoroCache::backLinkCode)
2447 {
2448 res_holder.len = 1;
2449 }
2450 res_holder.coef = p_GetCoeff (t, c->r);
2451 res_holder.p = ref->value_poly;
2452 res_holder.ref = ref;
2453 res_holder.onlyBorrowed = TRUE;
2454 res_holder.changed = TRUE;
2455 p_Delete (&t, c->r);
2456 return res_holder;
2457 }
2458 }
2459 else
2460 {
2461 BOOLEAN succ;
2462 poly cache_lookup = cache->lookup (t, succ, res_holder.len); //don't own this yet
2463 if(succ)
2464 {
2465 if(cache_lookup == t)
2466 {
2467 //know they are equal
2468 //res_holder.len=1;
2469
2470 res_holder.changed = FALSE;
2471 res_holder.p = t;
2472 res_holder.coef = npInit (1);
2473
2474 res_holder.onlyBorrowed = FALSE;
2475 return res_holder;
2476 }
2477
2478 res_holder.coef = p_GetCoeff (t, c->r);
2479 p_Delete (&t, c->r);
2480
2481 res_holder.p = cache_lookup;
2482
2483 res_holder.onlyBorrowed = TRUE;
2484 return res_holder;
2485
2486 }
2487 }
2488
2489 unsigned long sev = p_GetShortExpVector (t, currRing);
2490 int i = kFindDivisibleByInS_easy (c->strat, t, sev);
2491 if(i >= 0)
2492 {
2493 number coef_bak = p_GetCoeff (t, c->r);
2494
2495 p_SetCoeff (t, npInit (1), c->r);
2496 assume (npIsOne (p_GetCoeff (c->strat->S[i], c->r)));
2497 number coefstrat = p_GetCoeff (c->strat->S[i], c->r);
2498
2499 //poly t_copy_mon=p_Copy(t,c->r);
2500 poly exp_diff = cache->temp_term;
2501 p_ExpVectorDiff (exp_diff, t, c->strat->S[i], c->r);
2502 p_SetCoeff (exp_diff, npNeg (nInvers (coefstrat)), c->r);
2503 // nInvers may be npInvers or nvInvers
2504 p_Setm (exp_diff, c->r);
2505 assume (c->strat->S[i] != NULL);
2506 //poly t_to_del=t;
2507 poly res;
2508 res = pp_Mult_mm (pNext (c->strat->S[i]), exp_diff, c->r);
2509
2510 res_holder.len = c->strat->lenS[i] - 1;
2511 res = noro_red_non_unique (res, res_holder.len, cache, c);
2512
2513 DataNoroCacheNode *ref = cache->insert (t, res, res_holder.len);
2514 p_Delete (&t, c->r);
2515 //p_Delete(&t_copy_mon,c->r);
2516 //res=pMult_nn(res,coef_bak);
2517 res_holder.changed = TRUE;
2518 res_holder.p = res;
2519 res_holder.coef = coef_bak;
2520 res_holder.onlyBorrowed = TRUE;
2521 res_holder.ref = ref;
2522 return res_holder;
2523 }
2524 else
2525 {
2526 number coef_bak = p_GetCoeff (t, c->r);
2527 number one = npInit (1);
2528 p_SetCoeff (t, one, c->r);
2529 res_holder.len = 1;
2530 if(!(force_unique))
2531 {
2532 res_holder.ref = cache->insert (t, t, res_holder.len);
2533 p_SetCoeff (t, coef_bak, c->r);
2534 //return t;
2535
2536 //we need distinction
2537 res_holder.changed = FALSE;
2538 res_holder.p = t;
2539
2540 res_holder.coef = npInit (1);
2541 res_holder.onlyBorrowed = FALSE;
2542 return res_holder;
2543 }
2544 else
2545 {
2546 res_holder.ref = cache->insertAndTransferOwnerShip (t, c->r);
2547 res_holder.coef = coef_bak;
2548 res_holder.onlyBorrowed = TRUE;
2549 res_holder.changed = FALSE;
2550 res_holder.p = t;
2551 return res_holder;
2552 }
2553 }
2554
2555 }
2556 #endif
2557 //SparseRow* noro_red_to_non_poly(poly p, int &len, NoroCache* cache,slimgb_alg* c);
2558 #ifndef NORO_NON_POLY
2559 //len input and out: Idea: reverse addition
noro_red_non_unique(poly p,int & len,NoroCache * cache,slimgb_alg * c)2560 poly noro_red_non_unique (poly p, int &len, NoroCache * cache, slimgb_alg * c)
2561 {
2562 assume (len == pLength (p));
2563 poly orig_p = p;
2564 if(p == NULL)
2565 {
2566 len = 0;
2567 return NULL;
2568 }
2569 kBucket_pt bucket = kBucketCreate (currRing);
2570 kBucketInit (bucket, NULL, 0);
2571 poly unchanged_head = NULL;
2572 poly unchanged_tail = NULL;
2573 int unchanged_size = 0;
2574
2575 while(p)
2576 {
2577 poly t = p;
2578 pIter (p);
2579 pNext (t) = NULL;
2580 #ifndef SING_NDEBUG
2581 number coef_debug = p_GetCoeff (t, currRing);
2582 #endif
2583 MonRedRes red = noro_red_mon (t, FALSE, cache, c);
2584 if((!(red.changed)) && (!(red.onlyBorrowed)))
2585 {
2586 unchanged_size++;
2587 assume (npIsOne (red.coef));
2588 assume (p_GetCoeff (red.p, currRing) == coef_debug);
2589 if(unchanged_head)
2590 {
2591 pNext (unchanged_tail) = red.p;
2592 pIter (unchanged_tail);
2593 }
2594 else
2595 {
2596 unchanged_tail = red.p;
2597 unchanged_head = red.p;
2598 }
2599 }
2600 else
2601 {
2602 assume (red.len == pLength (red.p));
2603 if(red.onlyBorrowed)
2604 {
2605 if(npIsOne (red.coef))
2606 {
2607 t = p_Copy (red.p, currRing);
2608 }
2609 else
2610 t = __pp_Mult_nn (red.p, red.coef, currRing);
2611 }
2612 else
2613 {
2614 if(npIsOne (red.coef))
2615 t = red.p;
2616 else
2617 t = __p_Mult_nn (red.p, red.coef, currRing);
2618 }
2619 kBucket_Add_q (bucket, t, &red.len);
2620 }
2621 }
2622 poly res = NULL;
2623 len = 0;
2624 kBucket_Add_q (bucket, unchanged_head, &unchanged_size);
2625 kBucketClear (bucket, &res, &len);
2626 kBucketDestroy (&bucket);
2627 return res;
2628 }
2629 #endif
2630 #ifdef NORO_SPARSE_ROWS_PRE
2631 //len input and out: Idea: reverse addition
2632
2633 /*template <class number_type> SparseRow<number_type>* noro_red_to_non_poly(poly p, int &len, NoroCache<number_type>* cache,slimgb_alg* c)
2634 * {
2635 if (n_GetChar(currRing->cf)<255)
2636 {
2637 return noro_red_to_non_poly_t<tgb_uint8>(p,len,cache,c);
2638 }
2639 else
2640 {
2641 if (n_GetChar(currRing->cf)<65000)
2642 {
2643 return noro_red_to_non_poly_t<tgb_uint16>(p,len,cache,c);
2644 }
2645 else
2646 {
2647 return noro_red_to_non_poly_t<tgb_uint32>(p,len,cache,c);
2648 }
2649 }
2650 }*/
2651 #endif
2652 //len input and out: Idea: reverse addition
2653 #ifndef NORO_NON_POLY
noro_red(poly p,int & len,NoroCache * cache,slimgb_alg * c)2654 std::vector < NoroPlaceHolder > noro_red (poly p, int &len, NoroCache * cache,
2655 slimgb_alg * c)
2656 {
2657 std::vector < NoroPlaceHolder > res;
2658 while(p)
2659 {
2660 poly t = p;
2661 pIter (p);
2662 pNext (t) = NULL;
2663
2664 MonRedRes red = noro_red_mon (t, TRUE, cache, c);
2665 assume (red.onlyBorrowed);
2666 assume (red.coef);
2667 assume (red.ref);
2668 NoroPlaceHolder h;
2669 h.ref = red.ref;
2670 h.coef = red.coef;
2671 assume (!((h.ref->value_poly == NULL) && (h.ref->value_len != 0)));
2672 if(h.ref->value_poly)
2673 res.push_back (h);
2674 }
2675 return res;
2676 }
2677 #endif
2678
2679 #endif
2680 #ifdef USE_NORO
2681 #ifndef NORO_CACHE
noro_step(poly * p,int & pn,slimgb_alg * c)2682 void noro_step (poly * p, int &pn, slimgb_alg * c)
2683 {
2684 poly *reduced = (poly *) omalloc (pn * sizeof (poly));
2685 int j;
2686 int *reduced_len = (int *) omalloc (pn * sizeof (int));
2687 int reduced_c = 0;
2688 //if (TEST_OPT_PROT)
2689 // PrintS("reduced system:\n");
2690 #ifdef NORO_CACHE
2691 NoroCache cache;
2692 #endif
2693 for(j = 0; j < pn; j++)
2694 {
2695
2696 poly h = p[j];
2697 int h_len = pLength (h);
2698
2699 number coef;
2700 #ifndef NORO_CACHE
2701 h = redNF2 (p_Copy (h, c->r), c, h_len, coef, 0);
2702 #else
2703 h = noro_red (p_Copy (h, c->r), h_len, &cache, c);
2704 assume (pLength (h) == h_len);
2705 #endif
2706 if(h != NULL)
2707 {
2708 #ifndef NORO_CACHE
2709
2710 h = redNFTail (h, c->strat->sl, c->strat, h_len);
2711 h_len = pLength (h);
2712 #endif
2713 reduced[reduced_c] = h;
2714 reduced_len[reduced_c] = h_len;
2715 reduced_c++;
2716 if(TEST_OPT_PROT)
2717 Print ("%d ", h_len);
2718 }
2719 }
2720 int reduced_sum = 0;
2721 for(j = 0; j < reduced_c; j++)
2722 {
2723 reduced_sum += reduced_len[j];
2724 }
2725 poly *terms = (poly *) omalloc (reduced_sum * sizeof (poly));
2726 int tc = 0;
2727 for(j = 0; j < reduced_c; j++)
2728 {
2729 poly h = reduced[j];
2730
2731 while(h != NULL)
2732 {
2733 terms[tc++] = h;
2734 pIter (h);
2735 assume (tc <= reduced_sum);
2736 }
2737 }
2738 assume (tc == reduced_sum);
2739 qsort (terms, reduced_sum, sizeof (poly), terms_sort_crit);
2740 int nterms = reduced_sum;
2741 //if (TEST_OPT_PROT)
2742 //Print("orig estimation:%i\n",reduced_sum);
2743 unify_terms (terms, nterms);
2744 //if (TEST_OPT_PROT)
2745 // Print("actual number of columns:%i\n",nterms);
2746 int rank = reduced_c;
2747 linalg_step_modp (reduced, p, rank, terms, nterms, c);
2748 omFree (terms);
2749
2750 pn = rank;
2751 omFree (reduced);
2752
2753 if(TEST_OPT_PROT)
2754 PrintS ("\n");
2755 }
2756 #else
2757
2758 #endif
2759 #endif
go_on(slimgb_alg * c)2760 static void go_on (slimgb_alg * c)
2761 {
2762 //set limit of 1000 for multireductions, at the moment for
2763 //programming reasons
2764 #ifdef USE_NORO
2765 //Print("module rank%d\n",c->S->rank);
2766 const BOOLEAN use_noro = c->use_noro;
2767 #else
2768 const BOOLEAN use_noro = FALSE;
2769 #endif
2770 int i = 0;
2771 c->average_length = 0;
2772 for(i = 0; i < c->n; i++)
2773 {
2774 c->average_length += c->lengths[i];
2775 }
2776 c->average_length = c->average_length / c->n;
2777 i = 0;
2778 int max_pairs = bundle_size;
2779
2780 #ifdef USE_NORO
2781 if((use_noro) || (c->use_noro_last_block))
2782 max_pairs = bundle_size_noro;
2783 #endif
2784 poly *p = (poly *) omAlloc ((max_pairs + 1) * sizeof (poly)); //nullterminated
2785
2786 int curr_deg = -1;
2787 while(i < max_pairs)
2788 {
2789 sorted_pair_node *s = top_pair (c); //here is actually chain criterium done
2790
2791 if(!s)
2792 break;
2793
2794 if(curr_deg >= 0)
2795 {
2796 if(s->deg > curr_deg)
2797 break;
2798 }
2799
2800 else
2801 curr_deg = s->deg;
2802 quick_pop_pair (c);
2803 if(s->i >= 0)
2804 {
2805 //be careful replace_pair use createShortSpoly which is not noncommutative
2806 now_t_rep (s->i, s->j, c);
2807 replace_pair (s->i, s->j, c);
2808
2809 if(s->i == s->j)
2810 {
2811 free_sorted_pair_node (s, c->r);
2812 continue;
2813 }
2814 now_t_rep (s->i, s->j, c);
2815 }
2816 poly h;
2817 if(s->i >= 0)
2818 {
2819 #ifdef HAVE_PLURAL
2820 if(c->nc)
2821 {
2822 h = nc_CreateSpoly (c->S->m[s->i], c->S->m[s->j] /*, NULL */ , c->r);
2823
2824 if(h != NULL)
2825 p_Cleardenom (h, c->r);
2826 }
2827 else
2828 #endif
2829 h = ksOldCreateSpoly (c->S->m[s->i], c->S->m[s->j], NULL, c->r);
2830 p_Test (h, c->r);
2831 }
2832 else
2833 {
2834 h = s->lcm_of_lm;
2835 p_Test (h, c->r);
2836 }
2837 // if(s->i>=0)
2838 // now_t_rep(s->j,s->i,c);
2839 number coef;
2840 int mlen = pLength (h);
2841 p_Test (h, c->r);
2842 if((!c->nc) & (!(use_noro)))
2843 {
2844 h = redNF2 (h, c, mlen, coef, 2);
2845 redTailShort (h, c->strat);
2846 nDelete (&coef);
2847 }
2848 p_Test (h, c->r);
2849 free_sorted_pair_node (s, c->r);
2850 if(!h)
2851 continue;
2852
2853 if(TEST_OPT_IDLIFT
2854 && p_GetComp(h, currRing) > c->syz_comp)
2855 {
2856 pDelete(&h);
2857 continue;
2858 }
2859
2860 p[i] = h;
2861 i++;
2862 }
2863 p[i] = NULL;
2864 // pre_comp(p,i,c);
2865 if(i == 0)
2866 {
2867 omFree (p);
2868 return;
2869 }
2870 #ifdef TGB_RESORT_PAIRS
2871 c->replaced = new bool[c->n];
2872 c->used_b = FALSE;
2873 #endif
2874
2875 c->normal_forms += i;
2876 int j;
2877 #ifdef USE_NORO
2878 //if ((!(c->nc))&&(rField_is_Zp(c->r)))
2879 //{
2880 if(use_noro)
2881 {
2882 int pn = i;
2883 if(pn == 0)
2884 {
2885 omFree (p);
2886 return;
2887 }
2888 {
2889 if(n_GetChar(currRing->cf) < 255)
2890 {
2891 noro_step < tgb_uint8 > (p, pn, c);
2892 }
2893 else
2894 {
2895 if(n_GetChar(currRing->cf) < 65000)
2896 {
2897 noro_step < tgb_uint16 > (p, pn, c);
2898 }
2899 else
2900 {
2901 noro_step < tgb_uint32 > (p, pn, c);
2902 }
2903 }
2904 }
2905
2906 //if (TEST_OPT_PROT)
2907 //{
2908 // Print("reported rank:%i\n",pn);
2909 //}
2910 mass_add (p, pn, c);
2911 omFree (p);
2912 return;
2913 /*if (TEST_OPT_PROT)
2914 for(j=0;j<pn;j++)
2915 {
2916 p_wrp(p[j],c->r);
2917 } */
2918 }
2919 #endif
2920 red_object *buf = (red_object *) omAlloc (i * sizeof (red_object)); /*i>0*/
2921 for(j = 0; j < i; j++)
2922 {
2923 p_Test (p[j], c->r);
2924 buf[j].p = p[j];
2925 buf[j].sev = pGetShortExpVector (p[j]);
2926 buf[j].bucket = kBucketCreate (currRing);
2927 p_Test (p[j], c->r);
2928 int len = pLength (p[j]);
2929 kBucketInit (buf[j].bucket, buf[j].p, len);
2930 buf[j].initial_quality = buf[j].guess_quality (c);
2931 assume (buf[j].initial_quality >= 0);
2932 }
2933 omFree (p);
2934 qsort (buf, i, sizeof (red_object), red_object_better_gen);
2935 // Print("\ncurr_deg:%i\n",curr_deg);
2936 if(TEST_OPT_PROT)
2937 {
2938 Print ("%dM[%d,", curr_deg, i);
2939 }
2940
2941 multi_reduction (buf, i, c);
2942 #ifdef TGB_RESORT_PAIRS
2943 if(c->used_b)
2944 {
2945 if(TEST_OPT_PROT)
2946 PrintS ("B");
2947 int e;
2948 for(e = 0; e <= c->pair_top; e++)
2949 {
2950 if(c->apairs[e]->i < 0)
2951 continue;
2952 assume (c->apairs[e]->j >= 0);
2953 if((c->replaced[c->apairs[e]->i]) || (c->replaced[c->apairs[e]->j]))
2954 {
2955 sorted_pair_node *s = c->apairs[e];
2956 s->expected_length = pair_weighted_length (s->i, s->j, c);
2957 }
2958 }
2959 qsort (c->apairs, c->pair_top + 1, sizeof (sorted_pair_node *),
2960 tgb_pair_better_gen2);
2961 }
2962 #endif
2963 #ifdef TGB_DEBUG
2964 {
2965 int k;
2966 for(k = 0; k < i; k++)
2967 {
2968 assume (kFindDivisibleByInS_easy (c->strat, buf[k]) < 0);
2969 int k2;
2970 for(k2 = 0; k2 < i; k2++)
2971 {
2972 if(k == k2)
2973 continue;
2974 assume ((!(p_LmDivisibleBy (buf[k].p, buf[k2].p, c->r)))
2975 || (wrp (buf[k].p), Print (" k %d k2 %d ", k, k2),
2976 wrp (buf[k2].p), FALSE));
2977 }
2978 }
2979 }
2980 #endif
2981 //resort S
2982
2983 if(TEST_OPT_PROT)
2984 Print ("%i]", i);
2985
2986 poly *add_those = (poly *) omalloc (i * sizeof (poly));
2987 for(j = 0; j < i; j++)
2988 {
2989 int len;
2990 poly p;
2991 buf[j].flatten ();
2992 kBucketClear (buf[j].bucket, &p, &len);
2993 kBucketDestroy (&buf[j].bucket);
2994 p_Test (p, c->r);
2995 //if (!c->nc) {
2996 if((c->tailReductions) || (lies_in_last_dp_block (p, c)))
2997 {
2998 p = redNFTail (p, c->strat->sl, c->strat, 0);
2999 }
3000 else
3001 {
3002 p = redTailShort (p, c->strat);
3003 }
3004 //}
3005 p_Test (p, c->r);
3006 add_those[j] = p;
3007
3008 //sbuf[j]=add_to_basis(p,-1,-1,c,ibuf+j);
3009 }
3010 mass_add (add_those, i, c);
3011 omFree (add_those);
3012 omFree (buf);
3013
3014 if(TEST_OPT_PROT)
3015 Print ("(%d)", c->pair_top + 1);
3016 //TODO: implement that while(!(idIs0(c->add_later)))
3017 #ifdef TGB_RESORT_PAIRS
3018 delete c->replaced;
3019 c->replaced = NULL;
3020 c->used_b = FALSE;
3021 #endif
3022 return;
3023 }
3024
3025 #ifdef REDTAIL_S
3026
redNFTail(poly h,const int sl,kStrategy strat,int len)3027 static poly redNFTail (poly h, const int sl, kStrategy strat, int len)
3028 {
3029 if(h == NULL)
3030 return NULL;
3031 pTest (h);
3032 if(0 > sl)
3033 return h;
3034 if(pNext (h) == NULL)
3035 return h;
3036 BOOLEAN nc = rIsPluralRing (currRing);
3037
3038 int j;
3039 poly res = h;
3040 poly act = res;
3041 LObject P (pNext (h));
3042 pNext (res) = NULL;
3043 P.bucket = kBucketCreate (currRing);
3044 len--;
3045 h = P.p;
3046 if(len <= 0)
3047 len = pLength (h);
3048 kBucketInit (P.bucket, h /*P.p */ , len /*pLength(P.p) */ );
3049 pTest (h);
3050 loop
3051 {
3052 P.p = h;
3053 P.t_p = NULL;
3054 P.SetShortExpVector ();
3055 loop
3056 {
3057 //int dummy=strat->sl;
3058 j = kFindDivisibleByInS_easy (strat, P.p, P.sev); //kFindDivisibleByInS(strat,&dummy,&P);
3059 if(j >= 0)
3060 {
3061 #ifdef REDTAIL_PROT
3062 PrintS ("r");
3063 #endif
3064 nNormalize (pGetCoeff (P.p));
3065 #ifdef KDEBUG
3066 if(TEST_OPT_DEBUG)
3067 {
3068 PrintS ("red tail:");
3069 wrp (h);
3070 PrintS (" with ");
3071 wrp (strat->S[j]);
3072 }
3073 #endif
3074 number coef;
3075 pTest (strat->S[j]);
3076 #ifdef HAVE_PLURAL
3077 if(nc)
3078 {
3079 nc_kBucketPolyRed_Z (P.bucket, strat->S[j], &coef);
3080 }
3081 else
3082 #endif
3083 coef = kBucketPolyRed (P.bucket, strat->S[j],
3084 strat->lenS[j] /*pLength(strat->S[j]) */ ,
3085 strat->kNoether);
3086 res=__p_Mult_nn (res, coef, currRing);
3087 nDelete (&coef);
3088 h = kBucketGetLm (P.bucket);
3089 if(h == NULL)
3090 {
3091 #ifdef REDTAIL_PROT
3092 PrintS (" ");
3093 #endif
3094 kBucketDestroy (&P.bucket);
3095 return res;
3096 }
3097 P.p = h;
3098 P.t_p = NULL;
3099 P.SetShortExpVector ();
3100 #ifdef KDEBUG
3101 if(TEST_OPT_DEBUG)
3102 {
3103 PrintS ("\nto tail:");
3104 wrp (h);
3105 PrintLn ();
3106 }
3107 #endif
3108 }
3109 else
3110 {
3111 #ifdef REDTAIL_PROT
3112 PrintS ("n");
3113 #endif
3114 break;
3115 }
3116 } /* end loop current mon */
3117 // poly tmp=pHead(h /*kBucketGetLm(P.bucket)*/);
3118 //act->next=tmp;pIter(act);
3119 act->next = kBucketExtractLm (P.bucket);
3120 pIter (act);
3121 h = kBucketGetLm (P.bucket);
3122 if(h == NULL)
3123 {
3124 #ifdef REDTAIL_PROT
3125 PrintS (" ");
3126 #endif
3127 kBucketDestroy (&P.bucket);
3128 return res;
3129 }
3130 pTest (h);
3131 }
3132 }
3133 #endif
3134
3135
3136 //try to fill, return FALSE iff queue is empty
3137
3138 //transfers ownership of m to mat
init_with_mac_poly(tgb_sparse_matrix * mat,int row,mac_poly m)3139 void init_with_mac_poly (tgb_sparse_matrix * mat, int row, mac_poly m)
3140 {
3141 assume (mat->mp[row] == NULL);
3142 mat->mp[row] = m;
3143 #ifdef TGB_DEBUG
3144 mac_poly r = m;
3145 while(r)
3146 {
3147 assume (r->exp < mat->columns);
3148 r = r->next;
3149 }
3150 #endif
3151 }
3152
3153 poly
free_row_to_poly(tgb_sparse_matrix * mat,int row,poly * monoms,int monom_index)3154 free_row_to_poly (tgb_sparse_matrix * mat, int row, poly * monoms,
3155 int monom_index)
3156 {
3157 poly p = NULL;
3158 poly *set_this = &p;
3159 mac_poly r = mat->mp[row];
3160 mat->mp[row] = NULL;
3161 while(r)
3162 {
3163 (*set_this) = pLmInit (monoms[monom_index - 1 - r->exp]);
3164 pSetCoeff ((*set_this), r->coef);
3165 set_this = &((*set_this)->next);
3166 mac_poly old = r;
3167 r = r->next;
3168 delete old;
3169
3170 }
3171 return p;
3172 }
3173
poly_crit(const void * ap1,const void * ap2)3174 static int poly_crit (const void *ap1, const void *ap2)
3175 {
3176 poly p1, p2;
3177 p1 = *((poly *) ap1);
3178 p2 = *((poly *) ap2);
3179
3180 int c = pLmCmp (p1, p2);
3181 if(c != 0)
3182 return c;
3183 int l1 = pLength (p1);
3184 int l2 = pLength (p2);
3185 if(l1 < l2)
3186 return -1;
3187 if(l1 > l2)
3188 return 1;
3189 return 0;
3190 }
3191
introduceDelayedPairs(poly * pa,int s)3192 void slimgb_alg::introduceDelayedPairs (poly * pa, int s)
3193 {
3194 if(s == 0)
3195 return;
3196 sorted_pair_node **si_array =
3197 (sorted_pair_node **) omAlloc (s * sizeof (sorted_pair_node *));
3198
3199 for(int i = 0; i < s; i++)
3200 {
3201 sorted_pair_node *si =
3202 (sorted_pair_node *) omAlloc (sizeof (sorted_pair_node));
3203 si->i = -1;
3204 si->j = -2;
3205 poly p = pa[i];
3206 simplify_poly (p, r);
3207 si->expected_length = pQuality (p, this, pLength (p));
3208 p_Test (p, r);
3209 si->deg = this->pTotaldegree_full (p);
3210 /*if (!rField_is_Zp(r))
3211 {
3212 p_Content(p,r);
3213 p_Cleardenom(p,r);
3214 } */
3215
3216 si->lcm_of_lm = p;
3217
3218 // c->apairs[n-1-i]=si;
3219 si_array[i] = si;
3220 }
3221
3222 qsort (si_array, s, sizeof (sorted_pair_node *), tgb_pair_better_gen2);
3223 apairs = spn_merge (apairs, pair_top + 1, si_array, s, this);
3224 pair_top += s;
3225 omFree (si_array);
3226 }
3227
slimgb_alg(ideal I,int syz_comp,BOOLEAN F4,int deg_pos)3228 slimgb_alg::slimgb_alg (ideal I, int syz_comp, BOOLEAN F4, int deg_pos)
3229 {
3230 this->deg_pos = deg_pos;
3231 lastCleanedDeg = -1;
3232 completed = FALSE;
3233 this->syz_comp = syz_comp;
3234 r = currRing;
3235 nc = rIsPluralRing (r);
3236 this->lastDpBlockStart = get_last_dp_block_start (r);
3237 //Print("last dp Block start: %i\n", this->lastDpBlockStart);
3238 is_homog = TRUE;
3239 {
3240 int hzz;
3241 for(hzz = 0; hzz < IDELEMS (I); hzz++)
3242 {
3243 assume (I->m[hzz] != NULL);
3244 int d = this->pTotaldegree (I->m[hzz]);
3245 poly t = I->m[hzz]->next;
3246 while(t)
3247 {
3248 if(d != this->pTotaldegree (t))
3249 {
3250 is_homog = FALSE;
3251 break;
3252 }
3253 t = t->next;
3254 }
3255 if(!(is_homog))
3256 break;
3257 }
3258 }
3259 eliminationProblem = ((!(is_homog)) && ((currRing->pLexOrder) || (I->rank > 1)));
3260 tailReductions = ((is_homog) || ((TEST_OPT_REDTAIL) && (!(I->rank > 1))));
3261 // Print("is homog:%d",c->is_homog);
3262 void *h;
3263 int i;
3264 to_destroy = NULL;
3265 easy_product_crit = 0;
3266 extended_product_crit = 0;
3267 if(rField_is_Zp (r))
3268 isDifficultField = FALSE;
3269 else
3270 isDifficultField = TRUE;
3271 //not fully correct
3272 //(rChar()==0);
3273 F4_mode = F4;
3274
3275 reduction_steps = 0;
3276 last_index = -1;
3277
3278 F = NULL;
3279 F_minus = NULL;
3280
3281 Rcounter = 0;
3282
3283 soon_free = NULL;
3284
3285 tmp_lm = pOne ();
3286
3287 normal_forms = 0;
3288 current_degree = 1;
3289
3290 max_pairs = 5 * IDELEMS (I);
3291
3292 apairs =
3293 (sorted_pair_node **) omAlloc (sizeof (sorted_pair_node *) * max_pairs);
3294 pair_top = -1;
3295
3296 int n = IDELEMS (I);
3297 array_lengths = n;
3298
3299
3300 i = 0;
3301 this->n = 0;
3302 T_deg = (int *) omAlloc (n * sizeof (int));
3303 if(eliminationProblem)
3304 T_deg_full = (int *) omAlloc (n * sizeof (int));
3305 else
3306 T_deg_full = NULL;
3307 tmp_pair_lm = (poly *) omAlloc (n * sizeof (poly));
3308 tmp_spn = (sorted_pair_node **) omAlloc (n * sizeof (sorted_pair_node *));
3309 lm_bin = omGetSpecBin (POLYSIZE + (r->ExpL_Size) * sizeof (long));
3310 #ifdef HEAD_BIN
3311 HeadBin = omGetSpecBin (POLYSIZE + (currRing->ExpL_Size) * sizeof (long));
3312 #endif
3313 /* omUnGetSpecBin(&(c->HeadBin)); */
3314 #ifndef HAVE_BOOST
3315 #ifdef USE_STDVECBOOL
3316 #else
3317 h = omAlloc (n * sizeof (char *));
3318
3319 states = (char **) h;
3320 #endif
3321 #endif
3322 h = omAlloc (n * sizeof (int));
3323 lengths = (int *) h;
3324 weighted_lengths = (wlen_type *) omAllocAligned (n * sizeof (wlen_type));
3325 gcd_of_terms = (poly *) omAlloc (n * sizeof (poly));
3326
3327 short_Exps = (long *) omAlloc (n * sizeof (long));
3328 if(F4_mode)
3329 S = idInit (n, I->rank);
3330 else
3331 S = idInit (1, I->rank);
3332 strat = new skStrategy;
3333 if(eliminationProblem)
3334 strat->honey = TRUE;
3335 strat->syzComp = syz_comp;
3336 initBuchMoraCrit (strat);
3337 initBuchMoraPos (strat);
3338 strat->initEcart = initEcartBBA;
3339 strat->tailRing = r;
3340 strat->enterS = enterSBba;
3341 strat->sl = -1;
3342 i = n;
3343 i = 1; //some strange bug else
3344 /* initS(c->S,NULL,c->strat); */
3345 /* intS start: */
3346 // i=((i+IDELEMS(c->S)+15)/16)*16;
3347 strat->ecartS = (intset) omAlloc (i * sizeof (int)); /*initec(i); */
3348 strat->sevS = (unsigned long *) omAlloc0 (i * sizeof (unsigned long));
3349 /*initsevS(i); */
3350 strat->S_2_R = (int *) omAlloc0 (i * sizeof (int)); /*initS_2_R(i); */
3351 strat->fromQ = NULL;
3352 strat->Shdl = idInit (1, 1);
3353 strat->S = strat->Shdl->m;
3354 strat->lenS = (int *) omAlloc0 (i * sizeof (int));
3355 if((isDifficultField) || (eliminationProblem))
3356 strat->lenSw = (wlen_type *) omAlloc0 (i * sizeof (wlen_type));
3357 else
3358 strat->lenSw = NULL;
3359 assume (n > 0);
3360 add_to_basis_ideal_quotient (I->m[0], this, NULL);
3361
3362 assume (strat->sl == IDELEMS (strat->Shdl) - 1);
3363 if(!(F4_mode))
3364 {
3365 poly *array_arg = I->m;
3366 array_arg++;
3367 introduceDelayedPairs (array_arg, n - 1);
3368 /*
3369 for (i=1;i<n;i++)//the 1 is wanted, because first element is added to basis
3370 {
3371 // add_to_basis(I->m[i],-1,-1,c);
3372 si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
3373 si->i=-1;
3374 si->j=-2;
3375 si->expected_length=pQuality(I->m[i],this,pLength(I->m[i]));
3376 si->deg=pTotaldegree(I->m[i]);
3377 if (!rField_is_Zp(r))
3378 {
3379 p_Cleardenom(I->m[i], r);
3380 }
3381 si->lcm_of_lm=I->m[i];
3382
3383 // c->apairs[n-1-i]=si;
3384 apairs[n-i-1]=si;
3385 ++(pair_top);
3386 } */
3387 }
3388 else
3389 {
3390 for(i = 1; i < n; i++) //the 1 is wanted, because first element is added to basis
3391 add_to_basis_ideal_quotient (I->m[i], this, NULL);
3392 }
3393 for(i = 0; i < IDELEMS (I); i++)
3394 {
3395 I->m[i] = NULL;
3396 }
3397 idDelete (&I);
3398 add_later = idInit (ADD_LATER_SIZE, S->rank);
3399 #ifdef USE_NORO
3400 use_noro = ((!(nc)) && (S->rank <= 1) && (rField_is_Zp (r))
3401 && (!(eliminationProblem)) && (n_GetChar(currRing->cf) <= NV_MAX_PRIME));
3402 use_noro_last_block = false;
3403 if((!(use_noro)) && (lastDpBlockStart <= (currRing->N)))
3404 {
3405 use_noro_last_block = ((!(nc)) && (S->rank <= 1) && (rField_is_Zp (r))
3406 && (n_GetChar(currRing->cf) <= NV_MAX_PRIME));
3407 }
3408 #else
3409 use_noro = false;
3410 use_noro_last_block = false;
3411 #endif
3412 //Print("NORO last block %i",use_noro_last_block);
3413 memset (add_later->m, 0, ADD_LATER_SIZE * sizeof (poly));
3414 }
3415
~slimgb_alg()3416 slimgb_alg::~slimgb_alg ()
3417 {
3418
3419 if(!(completed))
3420 {
3421 poly *add = (poly *) omAlloc ((pair_top + 2) * sizeof (poly));
3422 int piter;
3423 int pos = 0;
3424 for(piter = 0; piter <= pair_top; piter++)
3425 {
3426 sorted_pair_node *s = apairs[piter];
3427 if(s->i < 0)
3428 {
3429 //delayed element
3430 if(s->lcm_of_lm != NULL)
3431 {
3432 add[pos] = s->lcm_of_lm;
3433 pos++;
3434 }
3435 }
3436 free_sorted_pair_node (s, r);
3437 apairs[piter] = NULL;
3438 }
3439 pair_top = -1;
3440 add[pos] = NULL;
3441 pos = 0;
3442 while(add[pos] != NULL)
3443 {
3444 add_to_basis_ideal_quotient (add[pos], this, NULL);
3445 pos++;
3446 }
3447 for(piter = 0; piter <= pair_top; piter++)
3448 {
3449 sorted_pair_node *s = apairs[piter];
3450 assume (s->i >= 0);
3451 free_sorted_pair_node (s, r);
3452 apairs[piter] = NULL;
3453 }
3454 pair_top = -1;
3455 }
3456 id_Delete (&add_later, r);
3457 int i, j;
3458 slimgb_alg *c = this;
3459 while(c->to_destroy)
3460 {
3461 pDelete (&(c->to_destroy->p));
3462 poly_list_node *old = c->to_destroy;
3463 c->to_destroy = c->to_destroy->next;
3464 omFree (old);
3465 }
3466 while(c->F)
3467 {
3468 for(i = 0; i < c->F->size; i++)
3469 {
3470 pDelete (&(c->F->mp[i].m));
3471 }
3472 omFree (c->F->mp);
3473 c->F->mp = NULL;
3474 mp_array_list *old = c->F;
3475 c->F = c->F->next;
3476 omFree (old);
3477 }
3478 while(c->F_minus)
3479 {
3480 for(i = 0; i < c->F_minus->size; i++)
3481 {
3482 pDelete (&(c->F_minus->p[i]));
3483 }
3484 omFree (c->F_minus->p);
3485 c->F_minus->p = NULL;
3486 poly_array_list *old = c->F_minus;
3487 c->F_minus = c->F_minus->next;
3488 omFree (old);
3489 }
3490 #ifndef HAVE_BOOST
3491 #ifndef USE_STDVECBOOL
3492 for(int z = 1 /* zero length at 0 */ ; z < c->n; z++)
3493 {
3494 omFree (c->states[z]);
3495 }
3496 omFree (c->states);
3497 #endif
3498 #endif
3499
3500 omFree (c->lengths);
3501 omFree (c->weighted_lengths);
3502 for(int z = 0; z < c->n; z++)
3503 {
3504 pDelete (&c->tmp_pair_lm[z]);
3505 omFree (c->tmp_spn[z]);
3506 }
3507 omFree (c->tmp_pair_lm);
3508 omFree (c->tmp_spn);
3509
3510 omFree (c->T_deg);
3511 omfree (c->T_deg_full); /*c->T_deg_full my be NULL*/
3512
3513 omFree (c->strat->ecartS);
3514 omFree (c->strat->sevS);
3515 // initsevS(i);
3516 omFree (c->strat->S_2_R);
3517
3518
3519 omFree (c->strat->lenS);
3520
3521 if(c->strat->lenSw)
3522 omFree (c->strat->lenSw);
3523
3524 for(i = 0; i < c->n; i++)
3525 {
3526 if(c->gcd_of_terms[i])
3527 pDelete (&(c->gcd_of_terms[i]));
3528 }
3529 omFree (c->gcd_of_terms);
3530
3531 omFree (c->apairs);
3532 if(TEST_OPT_PROT)
3533 {
3534 //Print("calculated %d NFs\n",c->normal_forms);
3535 Print ("\nNF:%i product criterion:%i, ext_product criterion:%i \n",
3536 c->normal_forms, c->easy_product_crit, c->extended_product_crit);
3537 }
3538
3539 for(i = 0; i <= c->strat->sl; i++)
3540 {
3541 if(!c->strat->S[i])
3542 continue;
3543 BOOLEAN found = FALSE;
3544 for(j = 0; j < c->n; j++)
3545 {
3546 if(c->S->m[j] == c->strat->S[i])
3547 {
3548 found = TRUE;
3549 break;
3550 }
3551 }
3552 if(!found)
3553 pDelete (&c->strat->S[i]);
3554 }
3555 // for(i=0;i<c->n;i++)
3556 // {
3557 // if (c->rep[i]!=i)
3558 // {
3559 // // for(j=0;j<=c->strat->sl;j++)
3560 // {
3561 // // if(c->strat->S[j]==c->S->m[i])
3562 // {
3563 // // c->strat->S[j]=NULL;
3564 // // break;
3565 // // }
3566 // // }
3567 // // PrintS("R_delete");
3568 // pDelete(&c->S->m[i]);
3569 // }
3570 // }
3571
3572 if(completed)
3573 {
3574 for(i = 0; i < c->n; i++)
3575 {
3576 assume (c->S->m[i] != NULL);
3577 if(p_GetComp (c->S->m[i], currRing) > this->syz_comp)
3578 continue;
3579 for(j = 0; j < c->n; j++)
3580 {
3581 if((c->S->m[j] == NULL) || (i == j))
3582 continue;
3583 assume (p_LmShortDivisibleBy (c->S->m[j], c->short_Exps[j],
3584 c->S->m[i], ~c->short_Exps[i],
3585 c->r) == p_LmDivisibleBy (c->S->m[j],
3586 c->S->m[i],
3587 c->r));
3588 if(p_LmShortDivisibleBy (c->S->m[j], c->short_Exps[j],
3589 c->S->m[i], ~c->short_Exps[i], c->r))
3590 {
3591 pDelete (&c->S->m[i]);
3592 break;
3593 }
3594 }
3595 }
3596 }
3597 omFree (c->short_Exps);
3598
3599 ideal I = c->S;
3600 IDELEMS (I) = c->n;
3601 idSkipZeroes (I);
3602 for(i = 0; i <= c->strat->sl; i++)
3603 c->strat->S[i] = NULL;
3604 id_Delete (&c->strat->Shdl, c->r);
3605 pDelete (&c->tmp_lm);
3606 omUnGetSpecBin (&lm_bin);
3607 delete c->strat;
3608 }
3609
t_rep_gb(const ring r,ideal arg_I,int syz_comp,BOOLEAN F4_mode)3610 ideal t_rep_gb (const ring r, ideal arg_I, int syz_comp, BOOLEAN F4_mode)
3611 {
3612 assume (r == currRing);
3613 ring orig_ring = r;
3614 int pos;
3615 ring new_ring = rAssure_TDeg (orig_ring, pos);
3616 ideal s_h;
3617 if(orig_ring != new_ring)
3618 {
3619 rChangeCurrRing (new_ring);
3620 s_h = idrCopyR_NoSort (arg_I, orig_ring, new_ring);
3621 /*int i;
3622 for(i=0;i<IDELEMS(s_h);i++)
3623 {
3624 poly p=s_h->m[i];
3625 while(p)
3626 {
3627 p_Setm(p,new_ring);
3628 pIter(p);
3629 }
3630 } */
3631 }
3632 else
3633 {
3634 s_h = id_Copy (arg_I, orig_ring);
3635 }
3636 idTest (s_h);
3637
3638 ideal s_result = do_t_rep_gb (new_ring, s_h, syz_comp, F4_mode, pos);
3639 ideal result;
3640 if(orig_ring != new_ring)
3641 {
3642 idTest (s_result);
3643 rChangeCurrRing (orig_ring);
3644 result = idrMoveR_NoSort (s_result, new_ring, orig_ring);
3645
3646 idTest (result);
3647 //rChangeCurrRing(new_ring);
3648 rDelete(new_ring);
3649 //rChangeCurrRing(orig_ring);
3650 }
3651 else
3652 result = s_result;
3653 idTest (result);
3654 return result;
3655 }
3656
3657 ideal
do_t_rep_gb(ring,ideal arg_I,int syz_comp,BOOLEAN F4_mode,int deg_pos)3658 do_t_rep_gb (ring /*r*/, ideal arg_I, int syz_comp, BOOLEAN F4_mode, int deg_pos)
3659 {
3660 // Print("QlogSize(0) %d, QlogSize(1) %d,QlogSize(-2) %d, QlogSize(5) %d\n", QlogSize(nlInit(0)),QlogSize(nlInit(1)),QlogSize(nlInit(-2)),QlogSize(nlInit(5)));
3661
3662 if(TEST_OPT_PROT)
3663 if(F4_mode)
3664 PrintS ("F4 Modus\n");
3665
3666 //debug_Ideal=arg_debug_Ideal;
3667 //if (debug_Ideal) PrintS("DebugIdeal received\n");
3668 // Print("Idelems %i \n----------\n",IDELEMS(arg_I));
3669 ideal I = arg_I;
3670 id_Compactify (I,currRing);
3671 if(idIs0 (I))
3672 return I;
3673 int i;
3674 for(i = 0; i < IDELEMS (I); i++)
3675 {
3676 assume (I->m[i] != NULL);
3677 simplify_poly (I->m[i], currRing);
3678 }
3679
3680 qsort (I->m, IDELEMS (I), sizeof (poly), poly_crit);
3681 //Print("Idelems %i \n----------\n",IDELEMS(I));
3682 //slimgb_alg* c=(slimgb_alg*) omalloc(sizeof(slimgb_alg));
3683 //int syz_comp=arg_I->rank;
3684 slimgb_alg *c = new slimgb_alg (I, syz_comp, F4_mode, deg_pos);
3685
3686 while((c->pair_top >= 0)
3687 && ((!(TEST_OPT_DEGBOUND))
3688 || (c->apairs[c->pair_top]->deg <= Kstd1_deg)))
3689 {
3690 #ifdef HAVE_F4
3691 if(F4_mode)
3692 go_on_F4 (c);
3693 else
3694 #endif
3695 go_on (c);
3696 }
3697 if(c->pair_top < 0)
3698 c->completed = TRUE;
3699 I = c->S;
3700 delete c;
3701 if(TEST_OPT_REDSB)
3702 {
3703 ideal erg = kInterRed (I, NULL);
3704 assume (I != erg);
3705 id_Delete (&I, currRing);
3706 return erg;
3707 }
3708 //qsort(I->m, IDELEMS(I),sizeof(poly),pLmCmp_func);
3709 assume (I->rank >= id_RankFreeModule (I,currRing));
3710 return (I);
3711 }
3712
now_t_rep(const int & arg_i,const int & arg_j,slimgb_alg * c)3713 void now_t_rep (const int &arg_i, const int &arg_j, slimgb_alg * c)
3714 {
3715 int i, j;
3716 if(arg_i == arg_j)
3717 {
3718 return;
3719 }
3720 if(arg_i > arg_j)
3721 {
3722 i = arg_j;
3723 j = arg_i;
3724 }
3725 else
3726 {
3727 i = arg_i;
3728 j = arg_j;
3729 }
3730 c->states[j][i] = HASTREP;
3731 }
3732
3733 static BOOLEAN
has_t_rep(const int & arg_i,const int & arg_j,slimgb_alg * state)3734 has_t_rep (const int &arg_i, const int &arg_j, slimgb_alg * state)
3735 {
3736 assume (0 <= arg_i);
3737 assume (0 <= arg_j);
3738 assume (arg_i < state->n);
3739 assume (arg_j < state->n);
3740 if(arg_i == arg_j)
3741 {
3742 return (TRUE);
3743 }
3744 if(arg_i > arg_j)
3745 {
3746 return (state->states[arg_i][arg_j] == HASTREP);
3747 }
3748 else
3749 {
3750 return (state->states[arg_j][arg_i] == HASTREP);
3751 }
3752 }
3753
3754 #if 0 // unused
3755 static int pLcmDeg (poly a, poly b)
3756 {
3757 int i;
3758 int n = 0;
3759 for(i = (currRing->N); i; i--)
3760 {
3761 n += si_max (pGetExp (a, i), pGetExp (b, i));
3762 }
3763 return n;
3764 }
3765 #endif
3766
shorten_tails(slimgb_alg * c,poly monom)3767 static void shorten_tails (slimgb_alg * c, poly monom)
3768 {
3769 return;
3770 // BOOLEAN corr=lenS_correct(c->strat);
3771 for(int i = 0; i < c->n; i++)
3772 {
3773 //enter tail
3774
3775 if(c->S->m[i] == NULL)
3776 continue;
3777 poly tail = c->S->m[i]->next;
3778 poly prev = c->S->m[i];
3779 BOOLEAN did_something = FALSE;
3780 while((tail != NULL) && (pLmCmp (tail, monom) >= 0))
3781 {
3782 if(p_LmDivisibleBy (monom, tail, c->r))
3783 {
3784 did_something = TRUE;
3785 prev->next = tail->next;
3786 tail->next = NULL;
3787 p_Delete (&tail, c->r);
3788 tail = prev;
3789 //PrintS("Shortened");
3790 c->lengths[i]--;
3791 }
3792 prev = tail;
3793 tail = tail->next;
3794 }
3795 if(did_something)
3796 {
3797 int new_pos;
3798 wlen_type q;
3799 q = pQuality (c->S->m[i], c, c->lengths[i]);
3800 new_pos = simple_posInS (c->strat, c->S->m[i], c->lengths[i], q);
3801
3802 int old_pos = -1;
3803 //assume new_pos<old_pos
3804 for(int z = 0; z <= c->strat->sl; z++)
3805 {
3806 if(c->strat->S[z] == c->S->m[i])
3807 {
3808 old_pos = z;
3809 break;
3810 }
3811 }
3812 if(old_pos == -1)
3813 for(int z = new_pos - 1; z >= 0; z--)
3814 {
3815 if(c->strat->S[z] == c->S->m[i])
3816 {
3817 old_pos = z;
3818 break;
3819 }
3820 }
3821 assume (old_pos >= 0);
3822 assume (new_pos <= old_pos);
3823 assume (pLength (c->strat->S[old_pos]) == c->lengths[i]);
3824 c->strat->lenS[old_pos] = c->lengths[i];
3825 if(c->strat->lenSw)
3826 c->strat->lenSw[old_pos] = q;
3827 if(new_pos < old_pos)
3828 move_forward_in_S (old_pos, new_pos, c->strat);
3829 length_one_crit (c, i, c->lengths[i]);
3830 }
3831 }
3832 }
3833
3834 #if 0 // currently unused
3835 static sorted_pair_node *pop_pair (slimgb_alg * c)
3836 {
3837 clean_top_of_pair_list (c);
3838
3839 if(c->pair_top < 0)
3840 return NULL;
3841 else
3842 return (c->apairs[c->pair_top--]);
3843 }
3844 #endif
3845
cleanDegs(int lower,int upper)3846 void slimgb_alg::cleanDegs (int lower, int upper)
3847 {
3848 assume (is_homog);
3849 int deg;
3850 if(TEST_OPT_PROT)
3851 {
3852 PrintS ("C");
3853 }
3854 for(deg = lower; deg <= upper; deg++)
3855 {
3856 int i;
3857 for(i = 0; i < n; i++)
3858 {
3859 if(T_deg[i] == deg)
3860 {
3861 poly h;
3862 h = S->m[i];
3863 h = redNFTail (h, strat->sl, strat, lengths[i]);
3864 if(TEST_OPT_INTSTRATEGY)
3865 {
3866 p_Cleardenom (h, r); //includes p_Content(h,r);
3867 }
3868 else
3869 pNorm (h);
3870 //TODO:GCD of TERMS
3871 poly got =::gcd_of_terms (h, r);
3872 p_Delete (&gcd_of_terms[i], r);
3873 gcd_of_terms[i] = got;
3874 int len = pLength (h);
3875 wlen_type wlen = pQuality (h, this, len);
3876 if(weighted_lengths)
3877 weighted_lengths[i] = wlen;
3878 lengths[i] = len;
3879 assume (h == S->m[i]);
3880 int j;
3881 for(j = 0; j <= strat->sl; j++)
3882 {
3883 if(h == strat->S[j])
3884 {
3885 int new_pos = simple_posInS (strat, h, len, wlen);
3886 if(strat->lenS)
3887 {
3888 strat->lenS[j] = len;
3889 }
3890 if(strat->lenSw)
3891 {
3892 strat->lenSw[j] = wlen;
3893 }
3894 if(new_pos < j)
3895 {
3896 move_forward_in_S (j, new_pos, strat);
3897 }
3898 else
3899 {
3900 if(new_pos > j)
3901 new_pos = new_pos - 1; //is identical with one element
3902 if(new_pos > j)
3903 move_backward_in_S (j, new_pos, strat);
3904 }
3905 break;
3906 }
3907 }
3908 }
3909 }
3910 }
3911 {
3912 int i, j;
3913 for(i = 0; i < this->n; i++)
3914 {
3915 for(j = 0; j < i; j++)
3916 {
3917 if(T_deg[i] + T_deg[j] <= upper)
3918 {
3919 now_t_rep (i, j, this);
3920 }
3921 }
3922 }
3923 }
3924 //TODO resort and update strat->S,strat->lenSw
3925 //TODO mark pairs
3926 }
3927
top_pair(slimgb_alg * c)3928 sorted_pair_node *top_pair (slimgb_alg * c)
3929 {
3930 while(c->pair_top >= 0)
3931 {
3932 super_clean_top_of_pair_list (c); //yeah, I know, it's odd that I use a different proc here
3933 if((c->is_homog) && (c->pair_top >= 0)
3934 && (c->apairs[c->pair_top]->deg >= c->lastCleanedDeg + 2))
3935 {
3936 int upper = c->apairs[c->pair_top]->deg - 1;
3937 c->cleanDegs (c->lastCleanedDeg + 1, upper);
3938 c->lastCleanedDeg = upper;
3939 }
3940 else
3941 {
3942 break;
3943 }
3944 }
3945
3946 if(c->pair_top < 0)
3947 return NULL;
3948 else
3949 return (c->apairs[c->pair_top]);
3950 }
3951
quick_pop_pair(slimgb_alg * c)3952 sorted_pair_node *quick_pop_pair (slimgb_alg * c)
3953 {
3954 if(c->pair_top < 0)
3955 return NULL;
3956 else
3957 return (c->apairs[c->pair_top--]);
3958 }
3959
super_clean_top_of_pair_list(slimgb_alg * c)3960 static void super_clean_top_of_pair_list (slimgb_alg * c)
3961 {
3962 while((c->pair_top >= 0)
3963 && (c->apairs[c->pair_top]->i >= 0)
3964 &&
3965 (good_has_t_rep
3966 (c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i, c)))
3967 {
3968 free_sorted_pair_node (c->apairs[c->pair_top], c->r);
3969 c->pair_top--;
3970 }
3971 }
3972
clean_top_of_pair_list(slimgb_alg * c)3973 void clean_top_of_pair_list (slimgb_alg * c)
3974 {
3975 while((c->pair_top >= 0) && (c->apairs[c->pair_top]->i >= 0)
3976 &&
3977 (!state_is
3978 (UNCALCULATED, c->apairs[c->pair_top]->j, c->apairs[c->pair_top]->i,
3979 c)))
3980 {
3981 free_sorted_pair_node (c->apairs[c->pair_top], c->r);
3982 c->pair_top--;
3983 }
3984 }
3985
3986 static BOOLEAN
state_is(calc_state state,const int & arg_i,const int & arg_j,slimgb_alg * c)3987 state_is (calc_state state, const int &arg_i, const int &arg_j,
3988 slimgb_alg * c)
3989 {
3990 assume (0 <= arg_i);
3991 assume (0 <= arg_j);
3992 assume (arg_i < c->n);
3993 assume (arg_j < c->n);
3994 if(arg_i == arg_j)
3995 {
3996 return (TRUE);
3997 }
3998 if(arg_i > arg_j)
3999 {
4000 return (c->states[arg_i][arg_j] == state);
4001 }
4002 else
4003 return (c->states[arg_j][arg_i] == state);
4004 }
4005
free_sorted_pair_node(sorted_pair_node * s,const ring r)4006 void free_sorted_pair_node (sorted_pair_node * s, const ring r)
4007 {
4008 if(s->i >= 0)
4009 p_Delete (&s->lcm_of_lm, r);
4010 omFree (s);
4011 }
4012
4013 static BOOLEAN
pair_better(sorted_pair_node * a,sorted_pair_node * b,slimgb_alg *)4014 pair_better (sorted_pair_node * a, sorted_pair_node * b, slimgb_alg * /*c*/)
4015 {
4016 if(a->deg < b->deg)
4017 return TRUE;
4018 if(a->deg > b->deg)
4019 return FALSE;
4020
4021 int comp = pLmCmp (a->lcm_of_lm, b->lcm_of_lm);
4022 if(comp == 1)
4023 return FALSE;
4024 if(-1 == comp)
4025 return TRUE;
4026 if(a->expected_length < b->expected_length)
4027 return TRUE;
4028 if(a->expected_length > b->expected_length)
4029 return FALSE;
4030 if(a->i + a->j < b->i + b->j)
4031 return TRUE;
4032 if(a->i + a->j > b->i + b->j)
4033 return FALSE;
4034 if(a->i < b->i)
4035 return TRUE;
4036 if(a->i > b->i)
4037 return FALSE;
4038 return TRUE;
4039 }
4040
tgb_pair_better_gen(const void * ap,const void * bp)4041 static int tgb_pair_better_gen (const void *ap, const void *bp)
4042 {
4043 sorted_pair_node *a = *((sorted_pair_node **) ap);
4044 sorted_pair_node *b = *((sorted_pair_node **) bp);
4045 assume ((a->i > a->j) || (a->i < 0));
4046 assume ((b->i > b->j) || (b->i < 0));
4047 if(a->deg < b->deg)
4048 return -1;
4049 if(a->deg > b->deg)
4050 return 1;
4051
4052 int comp = pLmCmp (a->lcm_of_lm, b->lcm_of_lm);
4053
4054 if(comp == 1)
4055 return 1;
4056 if(-1 == comp)
4057 return -1;
4058 if(a->expected_length < b->expected_length)
4059 return -1;
4060 if(a->expected_length > b->expected_length)
4061 return 1;
4062 if(a->i + a->j < b->i + b->j)
4063 return -1;
4064 if(a->i + a->j > b->i + b->j)
4065 return 1;
4066 if(a->i < b->i)
4067 return -1;
4068 if(a->i > b->i)
4069 return 1;
4070 return 0;
4071 }
4072
gcd_of_terms(poly p,ring r)4073 static poly gcd_of_terms (poly p, ring r)
4074 {
4075 int max_g_0 = 0;
4076 assume (p != NULL);
4077 int i;
4078 poly m = pOne ();
4079 poly t;
4080 for(i = (currRing->N); i; i--)
4081 {
4082 pSetExp (m, i, pGetExp (p, i));
4083 if(max_g_0 == 0)
4084 if(pGetExp (m, i) > 0)
4085 max_g_0 = i;
4086 }
4087
4088 t = p->next;
4089 while(t != NULL)
4090 {
4091 if(max_g_0 == 0)
4092 break;
4093 for(i = max_g_0; i; i--)
4094 {
4095 pSetExp (m, i, si_min (pGetExp (t, i), pGetExp (m, i)));
4096 if(max_g_0 == i)
4097 if(pGetExp (m, i) == 0)
4098 max_g_0 = 0;
4099 if((max_g_0 == 0) && (pGetExp (m, i) > 0))
4100 {
4101 max_g_0 = i;
4102 }
4103 }
4104 t = t->next;
4105 }
4106 p_Setm (m, r);
4107 if(max_g_0 > 0)
4108 return m;
4109 pDelete (&m);
4110 return NULL;
4111 }
4112
pHasNotCFExtended(poly p1,poly p2,poly m)4113 static inline BOOLEAN pHasNotCFExtended (poly p1, poly p2, poly m)
4114 {
4115
4116 if(pGetComp (p1) > 0 || pGetComp (p2) > 0)
4117 return FALSE;
4118 int i = 1;
4119 loop
4120 {
4121 if((pGetExp (p1, i) - pGetExp (m, i) > 0)
4122 && (pGetExp (p2, i) - pGetExp (m, i) > 0))
4123 return FALSE;
4124 if(i == (currRing->N))
4125 return TRUE;
4126 i++;
4127 }
4128 }
4129
4130 //for impl reasons may return false if the the normal product criterion matches
4131 static inline BOOLEAN
extended_product_criterion(poly p1,poly gcd1,poly p2,poly gcd2,slimgb_alg * c)4132 extended_product_criterion (poly p1, poly gcd1, poly p2, poly gcd2,
4133 slimgb_alg * c)
4134 {
4135 if(c->nc)
4136 return FALSE;
4137 if(gcd1 == NULL)
4138 return FALSE;
4139 if(gcd2 == NULL)
4140 return FALSE;
4141 gcd1->next = gcd2; //may ordered incorrect
4142 poly m = gcd_of_terms (gcd1, c->r);
4143 gcd1->next = NULL;
4144 if(m == NULL)
4145 return FALSE;
4146
4147 BOOLEAN erg = pHasNotCFExtended (p1, p2, m);
4148 pDelete (&m);
4149 return erg;
4150 }
4151
4152 #if 0 //currently unused
4153 static poly kBucketGcd (kBucket * b, ring r)
4154 {
4155 int s = 0;
4156 int i;
4157 poly m, n;
4158 BOOLEAN initialized = FALSE;
4159 for(i = MAX_BUCKET - 1; i >= 0; i--)
4160 {
4161 if(b->buckets[i] != NULL)
4162 {
4163 if(!initialized)
4164 {
4165 m = gcd_of_terms (b->buckets[i], r);
4166 initialized = TRUE;
4167 if(m == NULL)
4168 return NULL;
4169 }
4170 else
4171 {
4172 n = gcd_of_terms (b->buckets[i], r);
4173 if(n == NULL)
4174 {
4175 pDelete (&m);
4176 return NULL;
4177 }
4178 n->next = m;
4179 poly t = gcd_of_terms (n, r);
4180 n->next = NULL;
4181 pDelete (&m);
4182 pDelete (&n);
4183 m = t;
4184 if(m == NULL)
4185 return NULL;
4186
4187 }
4188 }
4189 }
4190 return m;
4191 }
4192 #endif
4193
quality_of_pos_in_strat_S(int pos,slimgb_alg * c)4194 static inline wlen_type quality_of_pos_in_strat_S (int pos, slimgb_alg * c)
4195 {
4196 if(c->strat->lenSw != NULL)
4197 return c->strat->lenSw[pos];
4198 return c->strat->lenS[pos];
4199 }
4200
4201 #ifdef HAVE_PLURAL
4202 static inline wlen_type
quality_of_pos_in_strat_S_mult_high(int pos,poly high,slimgb_alg * c)4203 quality_of_pos_in_strat_S_mult_high (int pos, poly high, slimgb_alg * c)
4204 //meant only for nc
4205 {
4206 poly m = pOne ();
4207 pExpVectorDiff (m, high, c->strat->S[pos]);
4208 poly product = nc_mm_Mult_pp (m, c->strat->S[pos], c->r);
4209 wlen_type erg = pQuality (product, c);
4210 pDelete (&m);
4211 pDelete (&product);
4212 return erg;
4213 }
4214 #endif
4215
4216 static void
multi_reduction_lls_trick(red_object * los,int,slimgb_alg * c,find_erg & erg)4217 multi_reduction_lls_trick (red_object * los, int /*losl*/, slimgb_alg * c,
4218 find_erg & erg)
4219 {
4220 erg.expand = NULL;
4221 BOOLEAN swap_roles; //from reduce_by, to_reduce_u if fromS
4222 if(erg.fromS)
4223 {
4224 if(pLmEqual (c->strat->S[erg.reduce_by], los[erg.to_reduce_u].p))
4225 {
4226 wlen_type quality_a = quality_of_pos_in_strat_S (erg.reduce_by, c);
4227 int best = erg.to_reduce_u + 1;
4228 /*
4229 for (i=erg.to_reduce_u;i>=erg.to_reduce_l;i--)
4230 {
4231 int qc=los[i].guess_quality(c);
4232 if (qc<quality_a)
4233 {
4234 best=i;
4235 quality_a=qc;
4236 }
4237 }
4238 if(best!=erg.to_reduce_u+1)
4239 {*/
4240 wlen_type qc;
4241 best = find_best (los, erg.to_reduce_l, erg.to_reduce_u, qc, c);
4242 if(qc < quality_a)
4243 {
4244 los[best].flatten ();
4245 int b_pos = kBucketCanonicalize (los[best].bucket);
4246 los[best].p = los[best].bucket->buckets[b_pos];
4247 qc = pQuality (los[best].bucket->buckets[b_pos], c);
4248 if(qc < quality_a)
4249 {
4250 red_object h = los[erg.to_reduce_u];
4251 los[erg.to_reduce_u] = los[best];
4252 los[best] = h;
4253 swap_roles = TRUE;
4254 }
4255 else
4256 swap_roles = FALSE;
4257 }
4258 else
4259 {
4260 swap_roles = FALSE;
4261 }
4262 }
4263 else
4264 {
4265 if(erg.to_reduce_u > erg.to_reduce_l)
4266 {
4267 wlen_type quality_a = quality_of_pos_in_strat_S (erg.reduce_by, c);
4268 #ifdef HAVE_PLURAL
4269 if((c->nc) && (!(rIsSCA (c->r))))
4270 quality_a =
4271 quality_of_pos_in_strat_S_mult_high (erg.reduce_by,
4272 los[erg.to_reduce_u].p, c);
4273 #endif
4274 int best = erg.to_reduce_u + 1;
4275 wlen_type qc;
4276 best = find_best (los, erg.to_reduce_l, erg.to_reduce_u, qc, c);
4277 assume (qc == los[best].guess_quality (c));
4278 if(qc < quality_a)
4279 {
4280 los[best].flatten ();
4281 int b_pos = kBucketCanonicalize (los[best].bucket);
4282 los[best].p = los[best].bucket->buckets[b_pos];
4283 qc = pQuality (los[best].bucket->buckets[b_pos], c);
4284 //(best!=erg.to_reduce_u+1)
4285 if(qc < quality_a)
4286 {
4287 red_object h = los[erg.to_reduce_u];
4288 los[erg.to_reduce_u] = los[best];
4289 los[best] = h;
4290 erg.reduce_by = erg.to_reduce_u;
4291 erg.fromS = FALSE;
4292 erg.to_reduce_u--;
4293 }
4294 }
4295 }
4296 else
4297 {
4298 assume (erg.to_reduce_u == erg.to_reduce_l);
4299 wlen_type quality_a = quality_of_pos_in_strat_S (erg.reduce_by, c);
4300 wlen_type qc = los[erg.to_reduce_u].guess_quality (c);
4301 if(qc < 0)
4302 PrintS ("Wrong wlen_type");
4303 if(qc < quality_a)
4304 {
4305 int best = erg.to_reduce_u;
4306 los[best].flatten ();
4307 int b_pos = kBucketCanonicalize (los[best].bucket);
4308 los[best].p = los[best].bucket->buckets[b_pos];
4309 qc = pQuality (los[best].bucket->buckets[b_pos], c);
4310 assume (qc >= 0);
4311 if(qc < quality_a)
4312 {
4313 BOOLEAN exp = FALSE;
4314 if(qc <= 2)
4315 {
4316 //Print("\n qc is %lld \n",qc);
4317 exp = TRUE;
4318 }
4319 else
4320 {
4321 if(qc < quality_a / 2)
4322 exp = TRUE;
4323 else if(erg.reduce_by < c->n / 4)
4324 exp = TRUE;
4325 }
4326 if(exp)
4327 {
4328 poly clear_into;
4329 los[erg.to_reduce_u].flatten ();
4330 kBucketClear (los[erg.to_reduce_u].bucket, &clear_into,
4331 &erg.expand_length);
4332 erg.expand = pCopy (clear_into);
4333 kBucketInit (los[erg.to_reduce_u].bucket, clear_into,
4334 erg.expand_length);
4335 if(TEST_OPT_PROT)
4336 PrintS ("e");
4337 }
4338 }
4339 }
4340 }
4341
4342 swap_roles = FALSE;
4343 return;
4344 }
4345 }
4346 else
4347 {
4348 if(erg.reduce_by > erg.to_reduce_u)
4349 {
4350 //then lm(rb)>= lm(tru) so =
4351 assume (erg.reduce_by == erg.to_reduce_u + 1);
4352 int best = erg.reduce_by;
4353 wlen_type quality_a = los[erg.reduce_by].guess_quality (c);
4354 wlen_type qc;
4355 best = find_best (los, erg.to_reduce_l, erg.to_reduce_u, qc, c);
4356
4357 if(qc < quality_a)
4358 {
4359 red_object h = los[erg.reduce_by];
4360 los[erg.reduce_by] = los[best];
4361 los[best] = h;
4362 }
4363 swap_roles = FALSE;
4364 return;
4365 }
4366 else
4367 {
4368 assume (!pLmEqual (los[erg.reduce_by].p, los[erg.to_reduce_l].p));
4369 assume (erg.to_reduce_u == erg.to_reduce_l);
4370 //further assume, that reduce_by is the above all other polys
4371 //with same leading term
4372 int il = erg.reduce_by;
4373 wlen_type quality_a = los[erg.reduce_by].guess_quality (c);
4374 wlen_type qc;
4375 while((il > 0) && pLmEqual (los[il - 1].p, los[il].p))
4376 {
4377 il--;
4378 qc = los[il].guess_quality (c);
4379 if(qc < quality_a)
4380 {
4381 quality_a = qc;
4382 erg.reduce_by = il;
4383 }
4384 }
4385 swap_roles = FALSE;
4386 }
4387 }
4388 if(swap_roles)
4389 {
4390 if(TEST_OPT_PROT)
4391 PrintS ("b");
4392 poly clear_into;
4393 int new_length;
4394 int bp = erg.to_reduce_u; //bucket_positon
4395 //kBucketClear(los[bp].bucket,&clear_into,&new_length);
4396 new_length = los[bp].clear_to_poly ();
4397 clear_into = los[bp].p;
4398 poly p = c->strat->S[erg.reduce_by];
4399 int j = erg.reduce_by;
4400 int old_length = c->strat->lenS[j]; // in view of S
4401 los[bp].p = p;
4402 kBucketInit (los[bp].bucket, p, old_length);
4403 wlen_type qal = pQuality (clear_into, c, new_length);
4404 int pos_in_c = -1;
4405 int z;
4406 int new_pos;
4407 new_pos = simple_posInS (c->strat, clear_into, new_length, qal);
4408 assume (new_pos <= j);
4409 for(z = c->n; z; z--)
4410 {
4411 if(p == c->S->m[z - 1])
4412 {
4413 pos_in_c = z - 1;
4414 break;
4415 }
4416 }
4417
4418 int tdeg_full = -1;
4419 int tdeg = -1;
4420 if(pos_in_c >= 0)
4421 {
4422 #ifdef TGB_RESORT_PAIRS
4423 c->used_b = TRUE;
4424 c->replaced[pos_in_c] = TRUE;
4425 #endif
4426 tdeg = c->T_deg[pos_in_c];
4427 c->S->m[pos_in_c] = clear_into;
4428 c->lengths[pos_in_c] = new_length;
4429 c->weighted_lengths[pos_in_c] = qal;
4430 if(c->gcd_of_terms[pos_in_c] == NULL)
4431 c->gcd_of_terms[pos_in_c] = gcd_of_terms (clear_into, c->r);
4432 if(c->T_deg_full)
4433 tdeg_full = c->T_deg_full[pos_in_c] =
4434 c->pTotaldegree_full (clear_into);
4435 else
4436 tdeg_full = tdeg;
4437 c_S_element_changed_hook (pos_in_c, c);
4438 }
4439 else
4440 {
4441 if(c->eliminationProblem)
4442 {
4443 tdeg_full = c->pTotaldegree_full (clear_into);
4444 tdeg = c->pTotaldegree (clear_into);
4445 }
4446 }
4447 c->strat->S[j] = clear_into;
4448 c->strat->lenS[j] = new_length;
4449
4450 assume (pLength (clear_into) == new_length);
4451 if(c->strat->lenSw != NULL)
4452 c->strat->lenSw[j] = qal;
4453 if(TEST_OPT_INTSTRATEGY)
4454 {
4455 p_Cleardenom (clear_into, c->r); //should be unnecessary
4456 //includes p_Content(clear_into, c->r);
4457 }
4458 else
4459 pNorm (clear_into);
4460 #ifdef FIND_DETERMINISTIC
4461 erg.reduce_by = j;
4462 //resort later see diploma thesis, find_in_S must be deterministic
4463 //during multireduction if spolys are only in the span of the
4464 //input polys
4465 #else
4466 if(new_pos < j)
4467 {
4468 if(c->strat->honey)
4469 c->strat->ecartS[j] = tdeg_full - tdeg;
4470 move_forward_in_S (j, new_pos, c->strat);
4471 erg.reduce_by = new_pos;
4472 }
4473 #endif
4474 }
4475 }
4476
fwbw(red_object * los,int i)4477 static int fwbw (red_object * los, int i)
4478 {
4479 int i2 = i;
4480 int step = 1;
4481
4482 BOOLEAN bw = FALSE;
4483 BOOLEAN incr = TRUE;
4484
4485 while(1)
4486 {
4487 if(!bw)
4488 {
4489 step = si_min (i2, step);
4490 if(step == 0)
4491 break;
4492 i2 -= step;
4493
4494 if(!pLmEqual (los[i].p, los[i2].p))
4495 {
4496 bw = TRUE;
4497 incr = FALSE;
4498 }
4499 else
4500 {
4501 if((!incr) && (step == 1))
4502 break;
4503 }
4504 }
4505 else
4506 {
4507 step = si_min (i - i2, step);
4508 if(step == 0)
4509 break;
4510 i2 += step;
4511 if(pLmEqual (los[i].p, los[i2].p))
4512 {
4513 if(step == 1)
4514 break;
4515 else
4516 {
4517 bw = FALSE;
4518 }
4519 }
4520 }
4521 if(incr)
4522 step *= 2;
4523 else
4524 {
4525 if(step % 2 == 1)
4526 step = (step + 1) / 2;
4527 else
4528 step /= 2;
4529 }
4530 }
4531 return i2;
4532 }
4533
4534 static void
canonicalize_region(red_object * los,int l,int u,slimgb_alg *)4535 canonicalize_region (red_object * los, int l, int u, slimgb_alg * /*c*/)
4536 {
4537 assume (l <= u + 1);
4538 int i;
4539 for(i = l; i <= u; i++)
4540 {
4541 kBucketCanonicalize (los[i].bucket);
4542 }
4543 }
4544
4545 #ifdef SING_NDEBUG
4546 static void
multi_reduction_find(red_object * los,int,slimgb_alg * c,int startf,find_erg & erg)4547 multi_reduction_find (red_object * los, int /*losl*/, slimgb_alg * c, int startf,
4548 find_erg & erg)
4549 #else
4550 static void
4551 multi_reduction_find (red_object * los, int losl, slimgb_alg * c, int startf,
4552 find_erg & erg)
4553 #endif
4554 {
4555 kStrategy strat = c->strat;
4556
4557 #ifndef SING_NDEBUG
4558 assume (startf <= losl);
4559 assume ((startf == losl - 1)
4560 || (pLmCmp (los[startf].p, los[startf + 1].p) == -1));
4561 #endif
4562 int i = startf;
4563
4564 int j;
4565 while(i >= 0)
4566 {
4567 #ifndef SING_NDEBUG
4568 assume ((i == losl - 1) || (pLmCmp (los[i].p, los[i + 1].p) <= 0));
4569 #endif
4570 assume (is_valid_ro (los[i]));
4571 j = kFindDivisibleByInS_easy (strat, los[i]);
4572 if(j >= 0)
4573 {
4574 erg.to_reduce_u = i;
4575 erg.reduce_by = j;
4576 erg.fromS = TRUE;
4577 int i2 = fwbw (los, i);
4578 assume (pLmEqual (los[i].p, los[i2].p));
4579 assume ((i2 == 0) || (!pLmEqual (los[i2].p, los[i2 - 1].p)));
4580 assume (i >= i2);
4581
4582 erg.to_reduce_l = i2;
4583 #ifndef SING_NDEBUG
4584 assume ((i == losl - 1) || (pLmCmp (los[i].p, los[i + 1].p) == -1));
4585 #endif
4586 canonicalize_region (los, erg.to_reduce_u + 1, startf, c);
4587 return;
4588 }
4589 else /*if(j < 0)*/
4590 {
4591 //not reduceable, try to use this for reducing higher terms
4592 int i2 = fwbw (los, i);
4593 assume (pLmEqual (los[i].p, los[i2].p));
4594 assume ((i2 == 0) || (!pLmEqual (los[i2].p, los[i2 - 1].p)));
4595 assume (i >= i2);
4596 if(i2 != i)
4597 {
4598 erg.to_reduce_u = i - 1;
4599 erg.to_reduce_l = i2;
4600 erg.reduce_by = i;
4601 erg.fromS = FALSE;
4602 #ifndef SING_NDEBUG
4603 assume ((i == losl - 1) || (pLmCmp (los[i].p, los[i + 1].p) == -1));
4604 #endif
4605 canonicalize_region (los, erg.to_reduce_u + 1, startf, c);
4606 return;
4607 }
4608 i--;
4609 }
4610 }
4611 erg.reduce_by = -1; //error code
4612 return;
4613 }
4614
4615 // nicht reduzierbare eintrage in ergnisliste schreiben
4616 // nullen loeschen
4617 // while(finde_groessten leitterm reduzierbar(c,erg))
4618 // {
4619
4620 static int
multi_reduction_clear_zeroes(red_object * los,int losl,int l,int u)4621 multi_reduction_clear_zeroes (red_object * los, int losl, int l, int u)
4622 {
4623 int deleted = 0;
4624 int i = l;
4625 int last = -1;
4626 while(i <= u)
4627 {
4628 if(los[i].p == NULL)
4629 {
4630 kBucketDestroy (&los[i].bucket);
4631 // delete los[i];//here we assume los are constructed with new
4632 //destroy resources, must be added here
4633 if(last >= 0)
4634 {
4635 memmove (los + (int) (last + 1 - deleted), los + (last + 1),
4636 sizeof (red_object) * (i - 1 - last));
4637 }
4638 last = i;
4639 deleted++;
4640 }
4641 i++;
4642 }
4643 if((last >= 0) && (last != losl - 1))
4644 memmove (los + (int) (last + 1 - deleted), los + last + 1,
4645 sizeof (red_object) * (losl - 1 - last));
4646 return deleted;
4647 }
4648
search_red_object_pos(red_object * a,int top,red_object * key)4649 int search_red_object_pos (red_object * a, int top, red_object * key)
4650 {
4651 int an = 0;
4652 int en = top;
4653 if(top == -1)
4654 return 0;
4655 if(pLmCmp (key->p, a[top].p) == 1)
4656 return top + 1;
4657 int i;
4658 loop
4659 {
4660 if(an >= en - 1)
4661 {
4662 if(pLmCmp (key->p, a[an].p) == -1)
4663 return an;
4664 return en;
4665 }
4666 i = (an + en) / 2;
4667 if(pLmCmp (key->p, a[i].p) == -1)
4668 en = i;
4669 else
4670 an = i;
4671 }
4672 }
4673
sort_region_down(red_object * los,int l,int u,slimgb_alg *)4674 static void sort_region_down (red_object * los, int l, int u, slimgb_alg * /*c*/)
4675 {
4676 int r_size = u - l + 1;
4677 qsort (los + l, r_size, sizeof (red_object), red_object_better_gen);
4678 int i;
4679 int *new_indices = (int *) omalloc ((r_size) * sizeof (int));
4680 int bound = 0;
4681 BOOLEAN at_end = FALSE;
4682 for(i = l; i <= u; i++)
4683 {
4684 if(!(at_end))
4685 {
4686 bound = new_indices[i - l] =
4687 bound + search_red_object_pos (los + bound, l - bound - 1, los + i);
4688 if(bound == l)
4689 at_end = TRUE;
4690 }
4691 else
4692 {
4693 new_indices[i - l] = l;
4694 }
4695 }
4696 red_object *los_region =
4697 (red_object *) omalloc (sizeof (red_object) * (u - l + 1));
4698 for(int i = 0; i < r_size; i++)
4699 {
4700 new_indices[i] += i;
4701 los_region[i] = los[l + i];
4702 assume ((i == 0) || (new_indices[i] > new_indices[i - 1]));
4703 }
4704
4705 i = r_size - 1;
4706 int j = u;
4707 int j2 = l - 1;
4708 while(i >= 0)
4709 {
4710 if(new_indices[i] == j)
4711 {
4712 los[j] = los_region[i];
4713 i--;
4714 j--;
4715 }
4716 else
4717 {
4718 assume (new_indices[i] < j);
4719 los[j] = los[j2];
4720 assume (j2 >= 0);
4721 j2--;
4722 j--;
4723 }
4724 }
4725 omFree (los_region);
4726 omFree (new_indices);
4727 }
4728
4729 //assume that los is ordered ascending by leading term, all non zero
multi_reduction(red_object * los,int & losl,slimgb_alg * c)4730 static void multi_reduction (red_object * los, int &losl, slimgb_alg * c)
4731 {
4732 poly *delay = (poly *) omAlloc (losl * sizeof (poly));
4733 int delay_s = 0;
4734 //initialize;
4735 assume (c->strat->sl >= 0);
4736 assume (losl > 0);
4737 int i;
4738 wlen_type max_initial_quality = 0;
4739
4740 for(i = 0; i < losl; i++)
4741 {
4742 los[i].sev = pGetShortExpVector (los[i].p);
4743 //SetShortExpVector();
4744 los[i].p = kBucketGetLm (los[i].bucket);
4745 if(los[i].initial_quality > max_initial_quality)
4746 max_initial_quality = los[i].initial_quality;
4747 // else
4748 // Print("init2_qal=%lld;", los[i].initial_quality);
4749 // Print("initial_quality=%lld;",max_initial_quality);
4750 }
4751
4752 int curr_pos = losl - 1;
4753
4754 // nicht reduzierbare eintrage in ergnisliste schreiben
4755 // nullen loeschen
4756 while(curr_pos >= 0)
4757 {
4758 if((c->use_noro_last_block)
4759 && (lies_in_last_dp_block (los[curr_pos].p, c)))
4760 {
4761 int pn_noro = curr_pos + 1;
4762 poly *p_noro = (poly *) omAlloc (pn_noro * sizeof (poly));
4763 for(i = 0; i < pn_noro; i++)
4764 {
4765 int dummy_len;
4766 poly p;
4767 los[i].p = NULL;
4768 kBucketClear (los[i].bucket, &p, &dummy_len);
4769 p_noro[i] = p;
4770 }
4771 if(n_GetChar(currRing->cf) < 255)
4772 {
4773 noro_step < tgb_uint8 > (p_noro, pn_noro, c);
4774 }
4775 else
4776 {
4777 if(n_GetChar(currRing->cf) < 65000)
4778 {
4779 noro_step < tgb_uint16 > (p_noro, pn_noro, c);
4780 }
4781 else
4782 {
4783 noro_step < tgb_uint32 > (p_noro, pn_noro, c);
4784 }
4785 }
4786 for(i = 0; i < pn_noro; i++)
4787 {
4788 los[i].p = p_noro[i];
4789 los[i].sev = pGetShortExpVector (los[i].p);
4790 //ignore quality
4791 kBucketInit (los[i].bucket, los[i].p, pLength (los[i].p));
4792 }
4793 qsort (los, pn_noro, sizeof (red_object), red_object_better_gen);
4794 int deleted =
4795 multi_reduction_clear_zeroes (los, losl, pn_noro, curr_pos);
4796 losl -= deleted;
4797 curr_pos -= deleted;
4798 break;
4799 }
4800 find_erg erg;
4801
4802 multi_reduction_find (los, losl, c, curr_pos, erg); //last argument should be curr_pos
4803 if(erg.reduce_by < 0)
4804 break;
4805
4806 erg.expand = NULL;
4807
4808 multi_reduction_lls_trick (los, losl, c, erg);
4809
4810 int i;
4811 // wrp(los[erg.to_reduce_u].p);
4812 //PrintLn();
4813 multi_reduce_step (erg, los, c);
4814
4815
4816 if(!TEST_OPT_REDTHROUGH)
4817 {
4818 for(i = erg.to_reduce_l; i <= erg.to_reduce_u; i++)
4819 {
4820 if(los[i].p != NULL) //the check (los[i].p!=NULL) might be invalid
4821 {
4822 //
4823 assume (los[i].initial_quality > 0);
4824 if(los[i].guess_quality (c)
4825 > 1.5 * delay_factor * max_initial_quality)
4826 {
4827 if(TEST_OPT_PROT)
4828 PrintS ("v");
4829 los[i].canonicalize ();
4830 if(los[i].guess_quality (c) > delay_factor * max_initial_quality)
4831 {
4832 if(TEST_OPT_PROT)
4833 PrintS (".");
4834 los[i].clear_to_poly ();
4835 //delay.push_back(los[i].p);
4836 delay[delay_s] = los[i].p;
4837 delay_s++;
4838 los[i].p = NULL;
4839 }
4840 }
4841 }
4842 }
4843 }
4844 int deleted = multi_reduction_clear_zeroes (los, losl, erg.to_reduce_l,
4845 erg.to_reduce_u);
4846 if(erg.fromS == FALSE)
4847 curr_pos = si_max (erg.to_reduce_u, erg.reduce_by);
4848 else
4849 curr_pos = erg.to_reduce_u;
4850 losl -= deleted;
4851 curr_pos -= deleted;
4852
4853 //Print("deleted %i \n",deleted);
4854 if((TEST_V_UPTORADICAL) && (!(erg.fromS)))
4855 sort_region_down (los, si_min (erg.to_reduce_l, erg.reduce_by),
4856 (si_max (erg.to_reduce_u, erg.reduce_by)) - deleted,
4857 c);
4858 else
4859 sort_region_down (los, erg.to_reduce_l, erg.to_reduce_u - deleted, c);
4860
4861 if(erg.expand)
4862 {
4863 #ifdef FIND_DETERMINISTIC
4864 int i;
4865 for(i = 0; c->expandS[i]; i++) ;
4866 c->expandS = (poly *) omrealloc (c->expandS, (i + 2) * sizeof (poly));
4867 c->expandS[i] = erg.expand;
4868 c->expandS[i + 1] = NULL;
4869 #else
4870 int ecart = 0;
4871 if(c->eliminationProblem)
4872 {
4873 ecart =
4874 c->pTotaldegree_full (erg.expand) - c->pTotaldegree (erg.expand);
4875 }
4876 add_to_reductors (c, erg.expand, erg.expand_length, ecart);
4877 #endif
4878 }
4879 }
4880
4881 //sorted_pair_node** pairs=(sorted_pair_node**)
4882 // omalloc(delay_s*sizeof(sorted_pair_node*));
4883 c->introduceDelayedPairs (delay, delay_s);
4884 /*
4885 for(i=0;i<delay_s;i++)
4886 {
4887 poly p=delay[i];
4888 //if (rPar(c->r)==0)
4889 simplify_poly(p,c->r);
4890 sorted_pair_node* si=(sorted_pair_node*) omalloc(sizeof(sorted_pair_node));
4891 si->i=-1;
4892 si->j=-1;
4893 if (!rField_is_Zp(c->r))
4894 {
4895 if (!c->nc)
4896 p=redTailShort(p, c->strat);
4897 p_Cleardenom(p, c->r);
4898 p_Content(p, c->r);
4899 }
4900 si->expected_length=pQuality(p,c,pLength(p));
4901 si->deg=pTotaldegree(p);
4902 si->lcm_of_lm=p;
4903 pairs[i]=si;
4904 }
4905 qsort(pairs,delay_s,sizeof(sorted_pair_node*),tgb_pair_better_gen2);
4906 c->apairs=spn_merge(c->apairs,c->pair_top+1,pairs,delay_s,c);
4907 c->pair_top+=delay_s; */
4908 omFree (delay);
4909 //omfree(pairs);
4910 return;
4911 }
4912
flatten()4913 void red_object::flatten ()
4914 {
4915 assume (p == kBucketGetLm (bucket));
4916 }
4917
validate()4918 void red_object::validate ()
4919 {
4920 p = kBucketGetLm (bucket);
4921 if(p)
4922 sev = pGetShortExpVector (p);
4923 }
4924
clear_to_poly()4925 int red_object::clear_to_poly ()
4926 {
4927 flatten ();
4928 int l;
4929 kBucketClear (bucket, &p, &l);
4930 return l;
4931 }
4932
reduce(red_object *,int,int)4933 void reduction_step::reduce (red_object * /*r*/, int /*l*/, int /*u*/)
4934 {
4935 }
4936
do_reduce(red_object & ro)4937 void simple_reducer::do_reduce (red_object & ro)
4938 {
4939 number coef;
4940 #ifdef HAVE_PLURAL
4941 if(c->nc)
4942 nc_kBucketPolyRed_Z (ro.bucket, p, &coef);
4943 else
4944 #endif
4945 coef = kBucketPolyRed (ro.bucket, p, p_len, c->strat->kNoether);
4946 nDelete (&coef);
4947 }
4948
reduce(red_object * r,int l,int u)4949 void simple_reducer::reduce (red_object * r, int l, int u)
4950 {
4951 this->pre_reduce (r, l, u);
4952 int i;
4953 //debug start
4954
4955 if(c->eliminationProblem)
4956 {
4957 assume (p_LmEqual (r[l].p, r[u].p, c->r));
4958 /*int lm_deg=pTotaldegree(r[l].p);
4959 reducer_deg=lm_deg+pTotaldegree_full(p)-pTotaldegree(p); */
4960 }
4961
4962 for(i = l; i <= u; i++)
4963 {
4964 this->do_reduce (r[i]);
4965 }
4966 for(i = l; i <= u; i++)
4967 {
4968 kBucketSimpleContent (r[i].bucket);
4969 r[i].validate ();
4970 }
4971 }
4972
~reduction_step()4973 reduction_step::~reduction_step ()
4974 {
4975 }
4976
~simple_reducer()4977 simple_reducer::~simple_reducer ()
4978 {
4979 if(fill_back != NULL)
4980 {
4981 kBucketInit (fill_back, p, p_len);
4982 }
4983 fill_back = NULL;
4984 }
4985
multi_reduce_step(find_erg & erg,red_object * r,slimgb_alg * c)4986 void multi_reduce_step (find_erg & erg, red_object * r, slimgb_alg * c)
4987 {
4988 STATIC_VAR int id = 0;
4989 id++;
4990 unsigned long sev;
4991 BOOLEAN lt_changed = FALSE;
4992 int rn = erg.reduce_by;
4993 poly red;
4994 int red_len;
4995 simple_reducer *pointer;
4996 BOOLEAN work_on_copy = FALSE;
4997 if(erg.fromS)
4998 {
4999 red = c->strat->S[rn];
5000 red_len = c->strat->lenS[rn];
5001 assume (red_len == pLength (red));
5002 }
5003 else
5004 {
5005 r[rn].flatten ();
5006 kBucketClear (r[rn].bucket, &red, &red_len);
5007
5008 if(TEST_OPT_INTSTRATEGY)
5009 {
5010 p_Cleardenom (red, c->r); //should be unnecessary
5011 //includes p_Content(red, c->r);
5012 }
5013 //pNormalize (red);
5014
5015 if((!(erg.fromS)) && (TEST_V_UPTORADICAL))
5016 {
5017 if(polynomial_root (red, c->r))
5018 lt_changed = TRUE;
5019 sev = p_GetShortExpVector (red, c->r);
5020 }
5021 red_len = pLength (red);
5022 }
5023 if(((TEST_V_MODPSOLVSB) && (red_len > 1))
5024 || ((c->nc) || (erg.to_reduce_u - erg.to_reduce_l > 5)))
5025 {
5026 work_on_copy = TRUE;
5027 // poly m=pOne();
5028 poly m = c->tmp_lm;
5029 pSetCoeff (m, nInit (1));
5030 pSetComp (m, 0);
5031 for(int i = 1; i <= (currRing->N); i++)
5032 pSetExp (m, i, (pGetExp (r[erg.to_reduce_l].p, i) - pGetExp (red, i)));
5033 pSetm (m);
5034 poly red_cp;
5035 #ifdef HAVE_PLURAL
5036 if(c->nc)
5037 red_cp = nc_mm_Mult_pp (m, red, c->r);
5038 else
5039 #endif
5040 red_cp = ppMult_mm (red, m);
5041 if(!erg.fromS)
5042 {
5043 kBucketInit (r[rn].bucket, red, red_len);
5044 }
5045 //now reduce the copy
5046 //static poly redNF2 (poly h,slimgb_alg* c , int &len, number& m,int n)
5047
5048 if(!c->nc)
5049 redTailShort (red_cp, c->strat);
5050 //number mul;
5051 // red_len--;
5052 // red_cp->next=redNF2(red_cp->next,c,red_len,mul,c->average_length);
5053 // pSetCoeff(red_cp,nMult(red_cp->coef,mul));
5054 // nDelete(&mul);
5055 // red_len++;
5056 red = red_cp;
5057 red_len = pLength (red);
5058 // pDelete(&m);
5059 }
5060
5061 assume (red_len == pLength (red));
5062
5063 int reducer_deg = 0;
5064 if(c->eliminationProblem)
5065 {
5066 int lm_deg = c->pTotaldegree (r[erg.to_reduce_l].p);
5067 int ecart;
5068 if(erg.fromS)
5069 {
5070 ecart = c->strat->ecartS[erg.reduce_by];
5071 }
5072 else
5073 {
5074 ecart = c->pTotaldegree_full (red) - lm_deg;
5075 }
5076 reducer_deg = lm_deg + ecart;
5077 }
5078 pointer = new simple_reducer (red, red_len, reducer_deg, c);
5079
5080 if((!work_on_copy) && (!erg.fromS))
5081 pointer->fill_back = r[rn].bucket;
5082 else
5083 pointer->fill_back = NULL;
5084 pointer->reduction_id = id;
5085 pointer->c = c;
5086
5087 pointer->reduce (r, erg.to_reduce_l, erg.to_reduce_u);
5088 if(work_on_copy)
5089 pDelete (&pointer->p);
5090 delete pointer;
5091 if(lt_changed)
5092 {
5093 assume (!erg.fromS);
5094 r[erg.reduce_by].sev = sev;
5095 }
5096 }
5097
pre_reduce(red_object *,int,int)5098 void simple_reducer::pre_reduce (red_object * /*r*/, int /*l*/, int /*u*/)
5099 {
5100 }
5101
5102 #if 0
5103 template int pos_helper<int, int*>(skStrategy*, spolyrec*, int, int*, spolyrec**);
5104 template int pos_helper<long, long*>(skStrategy*, spolyrec*, long, long*, spolyrec**);
5105
5106 template void noro_step<unsigned char>(spolyrec**, int&, slimgb_alg*);
5107 template void noro_step<unsigned int>(spolyrec**, int&, slimgb_alg*);
5108 template void noro_step<unsigned short>(spolyrec**, int&, slimgb_alg*);
5109
5110
5111 template int term_nodes_sort_crit<unsigned char>(void const*, void const*);
5112 template int term_nodes_sort_crit<unsigned int>(void const*, void const*);
5113 template int term_nodes_sort_crit<unsigned short>(void const*, void const*);
5114
5115 template spolyrec* row_to_poly<unsigned char>(unsigned char*, spolyrec**, int, ip_sring*);
5116 template spolyrec* row_to_poly<unsigned int>(unsigned int*, spolyrec**, int, ip_sring*);
5117 template spolyrec* row_to_poly<unsigned short>(unsigned short*, spolyrec**, int, ip_sring*);
5118
5119 template void simplest_gauss_modp<unsigned char>(unsigned char*, int, int);
5120 template void simplest_gauss_modp<unsigned int>(unsigned int*, int, int);
5121 template void simplest_gauss_modp<unsigned short>(unsigned short*, int, int);
5122
5123
5124 template int modP_lastIndexRow<unsigned char>(unsigned char*, int);
5125 template int modP_lastIndexRow<unsigned int>(unsigned int*, int);
5126 template int modP_lastIndexRow<unsigned short>(unsigned short*, int);
5127
5128 template SparseRow<unsigned char>* noro_red_to_non_poly_t<unsigned char>(spolyrec*, int&, NoroCache<unsigned char>*, slimgb_alg*);
5129 template SparseRow<unsigned int>* noro_red_to_non_poly_t<unsigned int>(spolyrec*, int&, NoroCache<unsigned int>*, slimgb_alg*);
5130 template SparseRow<unsigned short>* noro_red_to_non_poly_t<unsigned short>(spolyrec*, int&, NoroCache<unsigned short>*, slimgb_alg*);
5131
5132
5133 template MonRedResNP<unsigned char> noro_red_mon_to_non_poly<unsigned char>(spolyrec*, NoroCache<unsigned char>*, slimgb_alg*);
5134 template MonRedResNP<unsigned int> noro_red_mon_to_non_poly<unsigned int>(spolyrec*, NoroCache<unsigned int>*, slimgb_alg*);
5135 template MonRedResNP<unsigned short> noro_red_mon_to_non_poly<unsigned short>(spolyrec*, NoroCache<unsigned short>*, slimgb_alg*);
5136
5137 template SparseRow<unsigned char>* noro_red_to_non_poly_dense<unsigned char>(MonRedResNP<unsigned char>*, int, NoroCache<unsigned char>*);
5138 template SparseRow<unsigned char>* noro_red_to_non_poly_sparse<unsigned char>(MonRedResNP<unsigned char>*, int, NoroCache<unsigned char>*);
5139 template SparseRow<unsigned int>* noro_red_to_non_poly_dense<unsigned int>(MonRedResNP<unsigned int>*, int, NoroCache<unsigned int>*);
5140 template SparseRow<unsigned int>* noro_red_to_non_poly_sparse<unsigned int>(MonRedResNP<unsigned int>*, int, NoroCache<unsigned int>*);
5141 template SparseRow<unsigned short>* noro_red_to_non_poly_dense<unsigned short>(MonRedResNP<unsigned short>*, int, NoroCache<unsigned short>*);
5142 template SparseRow<unsigned short>* noro_red_to_non_poly_sparse<unsigned short>(MonRedResNP<unsigned short>*, int, NoroCache<unsigned short>*);
5143
5144
5145
5146 template class DataNoroCacheNode<unsigned char>;
5147 template class DataNoroCacheNode<unsigned int>;
5148 template class DataNoroCacheNode<unsigned short>;
5149
5150 template class NoroCache<unsigned char>;
5151 template class NoroCache<unsigned int>;
5152 template class NoroCache<unsigned short>;
5153
5154
5155
5156 template void add_coef_times_dense<unsigned char>(unsigned char*, int, unsigned char const*, int, snumber*);
5157 template void add_coef_times_dense<unsigned int>(unsigned int*, int, unsigned int const*, int, snumber*);
5158 template void add_coef_times_dense<unsigned short>(unsigned short*, int, unsigned short const*, int, snumber*);
5159 template void add_coef_times_sparse<unsigned char>(unsigned char*, int, SparseRow<unsigned char>*, snumber*);
5160 template void add_coef_times_sparse<unsigned int>(unsigned int*, int, SparseRow<unsigned int>*, snumber*);
5161 template void add_coef_times_sparse<unsigned short>(unsigned short*, int, SparseRow<unsigned short>*, snumber*);
5162 template void add_dense<unsigned char>(unsigned char*, int, unsigned char const*, int);
5163 template void add_dense<unsigned int>(unsigned int*, int, unsigned int const*, int);
5164 template void add_dense<unsigned short>(unsigned short*, int, unsigned short const*, int);
5165 template void add_sparse<unsigned char>(unsigned char*, int, SparseRow<unsigned char>*);
5166 template void add_sparse<unsigned int>(unsigned int*, int, SparseRow<unsigned int>*);
5167 template void add_sparse<unsigned short>(unsigned short*, int, SparseRow<unsigned short>*);
5168
5169
5170 template void sub_dense<unsigned char>(unsigned char*, int, unsigned char const*, int);
5171 template void sub_dense<unsigned int>(unsigned int*, int, unsigned int const*, int);
5172 template void sub_dense<unsigned short>(unsigned short*, int, unsigned short const*, int);
5173 template void sub_sparse<unsigned char>(unsigned char*, int, SparseRow<unsigned char>*);
5174 template void sub_sparse<unsigned int>(unsigned int*, int, SparseRow<unsigned int>*);
5175 template void sub_sparse<unsigned short>(unsigned short*, int, SparseRow<unsigned short>*);
5176 template void write_coef_idx_to_buffer_dense<unsigned char>(CoefIdx<unsigned char>*, int&, unsigned char*, int);
5177 template void write_coef_idx_to_buffer_dense<unsigned int>(CoefIdx<unsigned int>*, int&, unsigned int*, int);
5178 template void write_coef_idx_to_buffer_dense<unsigned short>(CoefIdx<unsigned short>*, int&, unsigned short*, int);
5179 template void write_coef_idx_to_buffer<unsigned char>(CoefIdx<unsigned char>*, int&, int*, unsigned char*, int);
5180 template void write_coef_idx_to_buffer<unsigned int>(CoefIdx<unsigned int>*, int&, int*, unsigned int*, int);
5181 template void write_coef_idx_to_buffer<unsigned short>(CoefIdx<unsigned short>*, int&, int*, unsigned short*, int);
5182 template void write_coef_times_xx_idx_to_buffer_dense<unsigned char>(CoefIdx<unsigned char>*, int&, unsigned char*, int, snumber*);
5183 template void write_coef_times_xx_idx_to_buffer_dense<unsigned int>(CoefIdx<unsigned int>*, int&, unsigned int*, int, snumber*);
5184 template void write_coef_times_xx_idx_to_buffer_dense<unsigned short>(CoefIdx<unsigned short>*, int&, unsigned short*, int, snumber*);
5185 template void write_coef_times_xx_idx_to_buffer<unsigned char>(CoefIdx<unsigned char>*, int&, int*, unsigned char*, int, snumber*);
5186 template void write_coef_times_xx_idx_to_buffer<unsigned int>(CoefIdx<unsigned int>*, int&, int*, unsigned int*, int, snumber*);
5187 template void write_coef_times_xx_idx_to_buffer<unsigned short>(CoefIdx<unsigned short>*, int&, int*, unsigned short*, int, snumber*);
5188 template void write_minus_coef_idx_to_buffer_dense<unsigned char>(CoefIdx<unsigned char>*, int&, unsigned char*, int);
5189 template void write_minus_coef_idx_to_buffer_dense<unsigned int>(CoefIdx<unsigned int>*, int&, unsigned int*, int);
5190 template void write_minus_coef_idx_to_buffer_dense<unsigned short>(CoefIdx<unsigned short>*, int&, unsigned short*, int);
5191 template void write_minus_coef_idx_to_buffer<unsigned char>(CoefIdx<unsigned char>*, int&, int*, unsigned char*, int);
5192 template void write_minus_coef_idx_to_buffer<unsigned int>(CoefIdx<unsigned int>*, int&, int*, unsigned int*, int);
5193 template void write_minus_coef_idx_to_buffer<unsigned short>(CoefIdx<unsigned short>*, int&, int*, unsigned short*, int);
5194
5195
5196 template class std::vector<DataNoroCacheNode<unsigned char>*>;
5197 template class std::vector<DataNoroCacheNode<unsigned int>*>;
5198 template class std::vector<DataNoroCacheNode<unsigned short>*>;
5199 template class std::vector<PolySimple>;
5200
5201 template void std::sort( CoefIdx<unsigned char>* , CoefIdx<unsigned char>* );
5202 template void std::sort( CoefIdx<unsigned int>* , CoefIdx<unsigned int>* );
5203 template void std::sort( CoefIdx<unsigned short>*, CoefIdx<unsigned short>* );
5204
5205 template void std::sort_heap<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5206 template void std::sort_heap<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5207 template void std::sort_heap<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5208
5209 template void std::make_heap<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5210 template void std::make_heap<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5211 template void std::make_heap<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5212 #endif
5213
5214 #if 0
5215 template void std::__final_insertion_sort<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5216 template void std::__final_insertion_sort<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5217 template void std::__final_insertion_sort<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5218
5219 template void std::__introsort_loop<CoefIdx<unsigned char>*, long>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*, long);
5220 template void std::__introsort_loop<CoefIdx<unsigned int>*, long>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*, long);
5221 template void std::__introsort_loop<CoefIdx<unsigned short>*, long>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*, long);
5222
5223 template void std::__heap_select<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5224 template void std::__heap_select<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5225 template void std::__heap_select<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5226
5227 template void std::__insertion_sort<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5228 template void std::__insertion_sort<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5229 template void std::__insertion_sort<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5230
5231 template void std::__move_median_first<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*, CoefIdx<unsigned char>*);
5232 template void std::__move_median_first<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*, CoefIdx<unsigned int>*);
5233 template void std::__move_median_first<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*, CoefIdx<unsigned short>*);
5234
5235 template void std::__unguarded_linear_insert<CoefIdx<unsigned char>*>(CoefIdx<unsigned char>*);
5236 template void std::__unguarded_linear_insert<CoefIdx<unsigned int>*>(CoefIdx<unsigned int>*);
5237 template void std::__unguarded_linear_insert<CoefIdx<unsigned short>*>(CoefIdx<unsigned short>*);
5238
5239 template void std::__adjust_heap<CoefIdx<unsigned char>*, long, CoefIdx<unsigned char> >(CoefIdx<unsigned char>*, long, long, CoefIdx<unsigned char>);
5240 template void std::__adjust_heap<CoefIdx<unsigned int>*, long, CoefIdx<unsigned int> >(CoefIdx<unsigned int>*, long, long, CoefIdx<unsigned int>);
5241 template void std::__adjust_heap<CoefIdx<unsigned short>*, long, CoefIdx<unsigned short> >(CoefIdx<unsigned short>*, long, long, CoefIdx<unsigned short>);
5242
5243 template void std::__push_heap<CoefIdx<unsigned char>*, long, CoefIdx<unsigned char> >(CoefIdx<unsigned char>*, long, long, CoefIdx<unsigned char>);
5244 template void std::__push_heap<CoefIdx<unsigned int>*, long, CoefIdx<unsigned int> >(CoefIdx<unsigned int>*, long, long, CoefIdx<unsigned int>);
5245 template void std::__push_heap<CoefIdx<unsigned short>*, long, CoefIdx<unsigned short> >(CoefIdx<unsigned short>*, long, long, CoefIdx<unsigned short>);
5246
5247 template CoefIdx<unsigned char>* std::__unguarded_partition<CoefIdx<unsigned char>*, CoefIdx<unsigned char> >(CoefIdx<unsigned char>*, CoefIdx<unsigned char>*, CoefIdx<unsigned char> const&);
5248 template CoefIdx<unsigned int>* std::__unguarded_partition<CoefIdx<unsigned int>*, CoefIdx<unsigned int> >(CoefIdx<unsigned int>*, CoefIdx<unsigned int>*, CoefIdx<unsigned int> const&);
5249 template CoefIdx<unsigned short>* std::__unguarded_partition<CoefIdx<unsigned short>*, CoefIdx<unsigned short> >(CoefIdx<unsigned short>*, CoefIdx<unsigned short>*, CoefIdx<unsigned short> const&);
5250
5251 #endif
5252
5253