source: git/kernel/GBEngine/kspoly.cc @ 1427c5c

spielwiese
Last change on this file since 1427c5c was 1427c5c, checked in by Karim Abou Zeid <karim23697@…>, 4 years ago
Fix bug in ksReducePolyGCD() for Letterplace. The commutative team should probably also look into this.
  • Property mode set to 100644
File size: 44.8 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5*  ABSTRACT -  Routines for Spoly creation and reductions
6*/
7
8// #define PDEBUG 2
9
10
11
12#include "kernel/mod2.h"
13#include "misc/options.h"
14#include "kernel/GBEngine/kutil.h"
15#include "coeffs/numbers.h"
16#include "polys/monomials/p_polys.h"
17#include "polys/templates/p_Procs.h"
18#include "polys/nc/nc.h"
19#ifdef HAVE_RINGS
20#include "kernel/polys.h"
21#endif
22#ifdef HAVE_SHIFTBBA
23#include "polys/shiftop.h"
24#endif
25
26#ifdef KDEBUG
27VAR int red_count = 0;
28VAR int create_count = 0;
29// define this if reductions are reported on TEST_OPT_DEBUG
30#define TEST_OPT_DEBUG_RED
31#endif
32
33/***************************************************************
34 *
35 * Reduces PR with PW
36 * Assumes PR != NULL, PW != NULL, Lm(PW) divides Lm(PR)
37 *
38 * returns 0: okay
39 *         1: tailRing changed
40 *         -1: cannot change tailRing
41 *         2: cannot change tailRing: strat==NULL
42 *
43 ***************************************************************/
44int ksReducePolyZ(LObject* PR,
45                 TObject* PW,
46                 poly spNoether,
47                 number *coef,
48                 kStrategy strat)
49{
50#ifdef KDEBUG
51  red_count++;
52#ifdef TEST_OPT_DEBUG_RED
53//  if (TEST_OPT_DEBUG)
54//  {
55//    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
56//    PW->wrp();
57//    //printf("\necart(PR)-ecart(PW): %i\n",PR->ecart-PW->ecart);
58//    //pWrite(PR->p);
59//  }
60#endif
61#endif
62  int ret = 0;
63  ring tailRing = PR->tailRing;
64  kTest_L(PR,tailRing);
65  kTest_T(PW);
66
67  poly p1 = PR->GetLmTailRing();   // p2 | p1
68  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
69  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
70  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
71  p_CheckPolyRing(p1, tailRing);
72  p_CheckPolyRing(p2, tailRing);
73
74  pAssume1(p2 != NULL && p1 != NULL &&
75           p_DivisibleBy(p2,  p1, tailRing));
76
77  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
78           (p_GetComp(p2, tailRing) == 0 &&
79            p_MaxComp(pNext(p2),tailRing) == 0));
80
81#ifdef HAVE_PLURAL
82  if (rIsPluralRing(currRing))
83  {
84    // for the time being: we know currRing==strat->tailRing
85    // no exp-bound checking needed
86    // (only needed if exp-bound(tailring)<exp-b(currRing))
87    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
88    else
89    {
90      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
91      assume(_p != NULL);
92      nc_PolyPolyRed(_p, p2,coef, currRing);
93      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
94      PR->pLength=0; // usually not used, GetpLength re-computes it if needed
95    }
96    return 0;
97  }
98#endif
99
100  if (t2==NULL)           // Divisor is just one term, therefore it will
101  {                       // just cancel the leading term
102    // adjust lead coefficient if needed
103    if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
104    {
105      number bn = pGetCoeff(lm);
106      number an = pGetCoeff(p2);
107      int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
108      p_SetCoeff(lm, bn, tailRing);
109      if ((ct == 0) || (ct == 2))
110      PR->Tail_Mult_nn(an);
111      if (coef != NULL) *coef = an;
112      else n_Delete(&an, tailRing->cf);
113    }
114    PR->LmDeleteAndIter();
115    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
116    return 0;
117  }
118
119  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
120
121  //if (tailRing != currRing)
122  {
123    // check that reduction does not violate exp bound
124    while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
125    {
126      // undo changes of lm
127      p_ExpVectorAdd(lm, p2, tailRing);
128      if (strat == NULL) return 2;
129      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
130      tailRing = strat->tailRing;
131      p1 = PR->GetLmTailRing();
132      p2 = PW->GetLmTailRing();
133      t2 = pNext(p2);
134      lm = p1;
135      p_ExpVectorSub(lm, p2, tailRing);
136      ret = 1;
137    }
138  }
139
140#ifdef HAVE_SHIFTBBA
141  poly lmRight;
142  if (tailRing->isLPring)
143  {
144    assume(PR->shift == 0);
145    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
146    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
147  }
148#endif
149
150  // take care of coef buisness
151  if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
152  {
153    number bn = pGetCoeff(lm);
154    number an = pGetCoeff(p2);
155    int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
156    p_SetCoeff(lm, bn, tailRing);
157    if ((ct == 0) || (ct == 2))
158      PR->Tail_Mult_nn(an);
159    if (coef != NULL) *coef = an;
160    else n_Delete(&an, tailRing->cf);
161  }
162  else
163  {
164    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
165  }
166
167
168  // and finally,
169#ifdef HAVE_SHIFTBBA
170  if (tailRing->isLPring)
171  {
172    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
173  }
174  else
175#endif
176  {
177    PR->Tail_Minus_mm_Mult_qq(lm, t2, pLength(t2) /*PW->GetpLength() - 1*/, spNoether);
178  }
179  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
180  PR->LmDeleteAndIter();
181
182  return ret;
183}
184
185int ksReducePoly(LObject* PR,
186                 TObject* PW,
187                 poly spNoether,
188                 number *coef,
189                 kStrategy strat)
190{
191#ifdef KDEBUG
192  red_count++;
193#ifdef TEST_OPT_DEBUG_RED
194//  if (TEST_OPT_DEBUG)
195//  {
196//    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
197//    PW->wrp();
198//    //printf("\necart(PR)-ecart(PW): %i\n",PR->ecart-PW->ecart);
199//    //pWrite(PR->p);
200//  }
201#endif
202#endif
203  int ret = 0;
204  ring tailRing = PR->tailRing;
205  kTest_L(PR,tailRing);
206  kTest_T(PW);
207
208  poly p1 = PR->GetLmTailRing();   // p2 | p1
209  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
210  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
211  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
212  p_CheckPolyRing(p1, tailRing);
213  p_CheckPolyRing(p2, tailRing);
214
215  pAssume1(p2 != NULL && p1 != NULL &&
216           p_DivisibleBy(p2,  p1, tailRing));
217
218  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
219           (p_GetComp(p2, tailRing) == 0 &&
220            p_MaxComp(pNext(p2),tailRing) == 0));
221
222#ifdef HAVE_PLURAL
223  if (rIsPluralRing(currRing))
224  {
225    // for the time being: we know currRing==strat->tailRing
226    // no exp-bound checking needed
227    // (only needed if exp-bound(tailring)<exp-b(currRing))
228    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
229    else
230    {
231      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
232      assume(_p != NULL);
233      nc_PolyPolyRed(_p, p2,coef, currRing);
234      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
235      PR->pLength=0; // usually not used, GetpLength re-computes it if needed
236    }
237    return 0;
238  }
239#endif
240
241  if (t2==NULL)           // Divisor is just one term, therefore it will
242  {                       // just cancel the leading term
243    PR->LmDeleteAndIter();
244    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
245    return 0;
246  }
247
248  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
249
250  //if (tailRing != currRing)
251  {
252    // check that reduction does not violate exp bound
253    while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
254    {
255      // undo changes of lm
256      p_ExpVectorAdd(lm, p2, tailRing);
257      if (strat == NULL) return 2;
258      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
259      tailRing = strat->tailRing;
260      p1 = PR->GetLmTailRing();
261      p2 = PW->GetLmTailRing();
262      t2 = pNext(p2);
263      lm = p1;
264      p_ExpVectorSub(lm, p2, tailRing);
265      ret = 1;
266    }
267  }
268
269#ifdef HAVE_SHIFTBBA
270  poly lmRight;
271  if (tailRing->isLPring)
272  {
273    assume(PR->shift == 0);
274    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
275    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
276  }
277#endif
278
279  // take care of coef buisness
280  if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
281  {
282    number bn = pGetCoeff(lm);
283    number an = pGetCoeff(p2);
284    int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
285    p_SetCoeff(lm, bn, tailRing);
286    if ((ct == 0) || (ct == 2))
287      PR->Tail_Mult_nn(an);
288    if (coef != NULL) *coef = an;
289    else n_Delete(&an, tailRing->cf);
290  }
291  else
292  {
293    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
294  }
295
296
297  // and finally,
298#ifdef HAVE_SHIFTBBA
299  if (tailRing->isLPring)
300  {
301    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
302  }
303  else
304#endif
305  {
306    PR->Tail_Minus_mm_Mult_qq(lm, t2, pLength(t2) /*PW->GetpLength() - 1*/, spNoether);
307  }
308  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
309  PR->LmDeleteAndIter();
310
311  return ret;
312}
313
314int ksReducePolyGCD(LObject* PR,
315                 TObject* PW,
316                 poly spNoether,
317                 number *coef,
318                 kStrategy strat)
319{
320#ifdef KDEBUG
321  red_count++;
322#ifdef TEST_OPT_DEBUG_RED
323//  if (TEST_OPT_DEBUG)
324//  {
325//    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
326//    PW->wrp();
327//    //printf("\necart(PR)-ecart(PW): %i\n",PR->ecart-PW->ecart);
328//    //pWrite(PR->p);
329//  }
330#endif
331#endif
332  int ret = 0;
333  ring tailRing = PR->tailRing;
334  kTest_L(PR, tailRing);
335  kTest_T(PW);
336
337  poly p1 = PR->GetLmTailRing();
338  poly p2 = PW->GetLmTailRing();
339  poly t2 = pNext(p2), lm = pOne();
340  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
341  p_CheckPolyRing(p1, tailRing);
342  p_CheckPolyRing(p2, tailRing);
343
344  pAssume1(p2 != NULL && p1 != NULL &&
345           p_DivisibleBy(p2,  p1, tailRing));
346
347  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
348           (p_GetComp(p2, tailRing) == 0 &&
349            p_MaxComp(pNext(p2),tailRing) == 0));
350
351#ifdef HAVE_PLURAL
352  if (rIsPluralRing(currRing))
353  {
354    // for the time being: we know currRing==strat->tailRing
355    // no exp-bound checking needed
356    // (only needed if exp-bound(tailring)<exp-b(currRing))
357    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
358    else
359    {
360      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
361      assume(_p != NULL);
362      nc_PolyPolyRed(_p, p2,coef, currRing);
363      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
364      PR->pLength=0; // usually not used, GetpLength re-computes it if needed
365    }
366    return 0;
367  }
368#endif
369  // check that reduction does not violate exp bound
370  while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
371  {
372    // undo changes of lm
373    p_ExpVectorAdd(lm, p2, tailRing);
374    if (strat == NULL) return 2;
375    if (! kStratChangeTailRing(strat, PR, PW)) return -1;
376    tailRing = strat->tailRing;
377    p1 = PR->GetLmTailRing();
378    p2 = PW->GetLmTailRing();
379    t2 = pNext(p2);
380    lm = p1;
381    p_ExpVectorSub(lm, p2, tailRing);
382    ret = 1;
383  }
384
385#ifdef HAVE_SHIFTBBA
386  poly lmRight;
387  if (tailRing->isLPring)
388  {
389    assume(PR->shift == 0);
390    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
391    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
392  }
393#endif
394
395  number ct, an, bn;
396  // take care of coef buisness
397  if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
398  {
399    ct = n_ExtGcd(pGetCoeff(p1), pGetCoeff(p2), &an, &bn, tailRing->cf);    // Calculate GCD
400#ifdef HAVE_SHIFTBBA
401    if (n_IsZero(an, tailRing->cf) || n_IsZero(bn, tailRing->cf))
402    {
403      // NOTE: not sure why this is not checked in the commutative case, this *does* happen and then zero coeff errors are reported
404
405      // NOTE: we are probably leaking memory of lm=pOne(), but we cannot delete it since it could also be lm=p1
406      n_Delete(&an, tailRing->cf);
407      n_Delete(&bn, tailRing->cf);
408      n_Delete(&ct, tailRing->cf);
409      return ret;
410    }
411#endif
412    /* negate bn since we subtract in Tail_Minus_mm_Mult_qq */
413    bn  = n_InpNeg(bn, tailRing->cf);
414    p_SetCoeff(lm, bn, tailRing);
415    PR->Tail_Mult_nn(an);
416  }
417  else
418  {
419    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
420  }
421
422
423  // and finally,
424#ifdef HAVE_SHIFTBBA
425  if (tailRing->isLPring)
426  {
427    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
428  }
429  else
430#endif
431  {
432    PR->Tail_Minus_mm_Mult_qq(lm, t2, pLength(t2) /*PW->GetpLength() - 1*/, spNoether);
433  }
434  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
435  pSetCoeff(PR->p, ct);
436
437  return ret;
438}
439
440/* Computes a reduction of the lead coefficient only. We have already tested
441 * that lm(PW) divides lm(PR), but lc(PW) does not divide lc(PR). We have
442 * computed division with remainder on the lead coefficients, parameter
443 * coef is the corresponding multiple for PW we need. The new lead
444 * coefficient, i.e. the remainder of lc division has already been
445 * set before calling this function. We do not drop the lead term at
446 * the end, but keep the adjusted, correct lead term. */
447int ksReducePolyLC(LObject* PR,
448                 TObject* PW,
449                 poly spNoether,
450                 number *coef,
451                 kStrategy strat)
452{
453#ifdef KDEBUG
454  red_count++;
455#ifdef TEST_OPT_DEBUG_RED
456//  if (TEST_OPT_DEBUG)
457//  {
458//    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
459//    PW->wrp();
460//    //printf("\necart(PR)-ecart(PW): %i\n",PR->ecart-PW->ecart);
461//    //pWrite(PR->p);
462//  }
463#endif
464#endif
465  /* printf("PR->P: ");
466   * p_Write(PR->p, currRing, PR->tailRing); */
467  int ret = 0;
468  ring tailRing = PR->tailRing;
469  kTest_L(PR,tailRing);
470  kTest_T(PW);
471
472  poly p1 = PR->GetLmTailRing();   // p2 | p1
473  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
474  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
475  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
476  p_CheckPolyRing(p1, tailRing);
477  p_CheckPolyRing(p2, tailRing);
478
479  pAssume1(p2 != NULL && p1 != NULL &&
480           p_DivisibleBy(p2,  p1, tailRing));
481
482  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
483           (p_GetComp(p2, tailRing) == 0 &&
484            p_MaxComp(pNext(p2),tailRing) == 0));
485
486#ifdef HAVE_PLURAL
487  if (rIsPluralRing(currRing))
488  {
489    // for the time being: we know currRing==strat->tailRing
490    // no exp-bound checking needed
491    // (only needed if exp-bound(tailring)<exp-b(currRing))
492    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
493    else
494    {
495      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
496      assume(_p != NULL);
497      nc_PolyPolyRed(_p, p2,coef, currRing);
498      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
499      PR->pLength=0; // usually not used, GetpLength re-computes it if needed
500    }
501    return 0;
502  }
503#endif
504
505  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
506  p_SetCoeff(lm, n_Init(1, tailRing), tailRing);
507  while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
508  {
509    // undo changes of lm
510    p_ExpVectorAdd(lm, p2, tailRing);
511    if (strat == NULL) return 2;
512    /* if (! kStratChangeTailRing(strat, PR, PW)) return -1; */
513    tailRing = strat->tailRing;
514    p1 = PR->GetLmTailRing();
515    p2 = PW->GetLmTailRing();
516    t2 = pNext(p2);
517    lm = p1;
518    p_ExpVectorSub(lm, p2, tailRing);
519    ret = 1;
520  }
521
522#ifdef HAVE_SHIFTBBA
523  poly lmRight;
524  if (tailRing->isLPring)
525  {
526    assume(PR->shift == 0);
527    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
528    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
529  }
530#endif
531
532  // and finally,
533#ifdef HAVE_SHIFTBBA
534  if (tailRing->isLPring)
535  {
536    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(p2, lmRight, tailRing), pLength(p2), spNoether);
537  }
538  else
539#endif
540  {
541    PR->Tail_Minus_mm_Mult_qq(lm, p2, pLength(p2) /*PW->GetpLength() - 1*/, spNoether);
542  }
543  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
544
545  PR->LmDeleteAndIter();
546  p_SetCoeff(PR->p, *coef, currRing);
547
548#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
549  if (TEST_OPT_DEBUG)
550  {
551    Print(" to: "); PR->wrp(); Print("\n");
552    //printf("\nt^%i ", PR->ecart);pWrite(pHead(PR->p));
553  }
554#endif
555  return ret;
556}
557
558int ksReducePolyBound(LObject* PR,
559                 TObject* PW,
560                 int bound,
561                 poly spNoether,
562                 number *coef,
563                 kStrategy strat)
564{
565#ifdef KDEBUG
566  red_count++;
567#ifdef TEST_OPT_DEBUG_RED
568  if (TEST_OPT_DEBUG)
569  {
570    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
571    PW->wrp();
572    //printf("\necart(PR)-ecart(PW): %i\n",PR->ecart-PW->ecart);
573    //pWrite(PR->p);
574  }
575#endif
576#endif
577  int ret = 0;
578  ring tailRing = PR->tailRing;
579  kTest_L(PR,tailRing);
580  kTest_T(PW);
581
582  poly p1 = PR->GetLmTailRing();   // p2 | p1
583  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
584  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
585  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
586  p_CheckPolyRing(p1, tailRing);
587  p_CheckPolyRing(p2, tailRing);
588
589  pAssume1(p2 != NULL && p1 != NULL &&
590           p_DivisibleBy(p2,  p1, tailRing));
591
592  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
593           (p_GetComp(p2, tailRing) == 0 &&
594            p_MaxComp(pNext(p2),tailRing) == 0));
595
596#ifdef HAVE_PLURAL
597  if (rIsPluralRing(currRing))
598  {
599    // for the time being: we know currRing==strat->tailRing
600    // no exp-bound checking needed
601    // (only needed if exp-bound(tailring)<exp-b(currRing))
602    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
603    else
604    {
605      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
606      assume(_p != NULL);
607      nc_PolyPolyRed(_p, p2,coef, currRing);
608      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
609      PR->pLength=0; // usually not used, GetpLength re-computes it if needed
610    }
611    return 0;
612  }
613#endif
614
615  if (t2==NULL)           // Divisor is just one term, therefore it will
616  {                       // just cancel the leading term
617    PR->LmDeleteAndIter();
618    if (coef != NULL) *coef = n_Init(1, tailRing);
619    return 0;
620  }
621
622  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
623
624  if (tailRing != currRing)
625  {
626    // check that reduction does not violate exp bound
627    while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
628    {
629      // undo changes of lm
630      p_ExpVectorAdd(lm, p2, tailRing);
631      if (strat == NULL) return 2;
632      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
633      tailRing = strat->tailRing;
634      p1 = PR->GetLmTailRing();
635      p2 = PW->GetLmTailRing();
636      t2 = pNext(p2);
637      lm = p1;
638      p_ExpVectorSub(lm, p2, tailRing);
639      ret = 1;
640    }
641  }
642
643#ifdef HAVE_SHIFTBBA
644  poly lmRight;
645  if (tailRing->isLPring)
646  {
647    assume(PR->shift == 0);
648    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
649    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
650  }
651#endif
652
653  // take care of coef buisness
654  if (! n_IsOne(pGetCoeff(p2), tailRing))
655  {
656    number bn = pGetCoeff(lm);
657    number an = pGetCoeff(p2);
658    int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
659    p_SetCoeff(lm, bn, tailRing);
660    if ((ct == 0) || (ct == 2))
661      PR->Tail_Mult_nn(an);
662    if (coef != NULL) *coef = an;
663    else n_Delete(&an, tailRing);
664  }
665  else
666  {
667    if (coef != NULL) *coef = n_Init(1, tailRing);
668  }
669
670
671  // and finally,
672#ifdef HAVE_SHIFTBBA
673  if (tailRing->isLPring)
674  {
675    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
676  }
677  else
678#endif
679  {
680    PR->Tail_Minus_mm_Mult_qq(lm, t2, pLength(t2) /*PW->GetpLength() - 1*/, spNoether);
681  }
682  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
683  PR->LmDeleteAndIter();
684
685#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
686  if (TEST_OPT_DEBUG)
687  {
688    Print(" to: "); PR->wrp(); Print("\n");
689    //printf("\nt^%i ", PR->ecart);pWrite(pHead(PR->p));
690  }
691#endif
692  return ret;
693}
694
695/***************************************************************
696 *
697 * Reduces PR with PW
698 * Assumes PR != NULL, PW != NULL, Lm(PW) divides Lm(PR)
699 *
700 ***************************************************************/
701
702int ksReducePolySig(LObject* PR,
703                 TObject* PW,
704                 long /*idx*/,
705                 poly spNoether,
706                 number *coef,
707                 kStrategy strat)
708{
709#ifdef KDEBUG
710  red_count++;
711#ifdef TEST_OPT_DEBUG_RED
712  if (TEST_OPT_DEBUG)
713  {
714    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
715    PW->wrp();
716  }
717#endif
718#endif
719  int ret = 0;
720  ring tailRing = PR->tailRing;
721  kTest_L(PR,tailRing);
722  kTest_T(PW);
723
724  // signature-based stuff:
725  // checking for sig-safeness first
726  // NOTE: This has to be done in the current ring
727  //
728  /**********************************************
729   *
730   * TODO:
731   * --------------------------------------------
732   * if strat->sbaOrder == 1
733   * Since we are subdividing lower index and
734   * current index reductions it is enough to
735   * look at the polynomial part of the signature
736   * for a check. This should speed-up checking
737   * a lot!
738   * if !strat->sbaOrder == 0
739   * We are not subdividing lower and current index
740   * due to the fact that we are using the induced
741   * Schreyer order
742   *
743   * nevertheless, this different behaviour is
744   * taken care of by is_sigsafe
745   * => one reduction procedure can be used for
746   * both, the incremental and the non-incremental
747   * attempt!
748   * --------------------------------------------
749   *
750   *********************************************/
751  //printf("COMPARE IDX: %ld -- %ld\n",idx,strat->currIdx);
752  if (!PW->is_sigsafe)
753  {
754    poly sigMult = pCopy(PW->sig);   // copy signature of reducer
755//#if 1
756#ifdef DEBUGF5
757    printf("IN KSREDUCEPOLYSIG: \n");
758    pWrite(pHead(f1));
759    pWrite(pHead(f2));
760    pWrite(sigMult);
761    printf("--------------\n");
762#endif
763    p_ExpVectorAddSub(sigMult,PR->GetLmCurrRing(),PW->GetLmCurrRing(),currRing);
764//#if 1
765#ifdef DEBUGF5
766    printf("------------------- IN KSREDUCEPOLYSIG: --------------------\n");
767    pWrite(pHead(f1));
768    pWrite(pHead(f2));
769    pWrite(sigMult);
770    pWrite(PR->sig);
771    printf("--------------\n");
772#endif
773    int sigSafe = p_LmCmp(PR->sig,sigMult,currRing);
774    // now we can delete the copied polynomial data used for checking for
775    // sig-safeness of the reduction step
776//#if 1
777#ifdef DEBUGF5
778    printf("%d -- %d sig\n",sigSafe,PW->is_sigsafe);
779
780#endif
781    //pDelete(&f1);
782    pDelete(&sigMult);
783    // go on with the computations only if the signature of p2 is greater than the
784    // signature of fm*p1
785    if(sigSafe != 1)
786    {
787      PR->is_redundant = TRUE;
788      return 3;
789    }
790    //PW->is_sigsafe  = TRUE;
791  }
792  PR->is_redundant = FALSE;
793  poly p1 = PR->GetLmTailRing();   // p2 | p1
794  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
795  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
796  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
797  p_CheckPolyRing(p1, tailRing);
798  p_CheckPolyRing(p2, tailRing);
799
800  pAssume1(p2 != NULL && p1 != NULL &&
801      p_DivisibleBy(p2,  p1, tailRing));
802
803  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
804      (p_GetComp(p2, tailRing) == 0 &&
805       p_MaxComp(pNext(p2),tailRing) == 0));
806
807#ifdef HAVE_PLURAL
808  if (rIsPluralRing(currRing))
809  {
810    // for the time being: we know currRing==strat->tailRing
811    // no exp-bound checking needed
812    // (only needed if exp-bound(tailring)<exp-b(currRing))
813    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
814    else
815    {
816      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
817      assume(_p != NULL);
818      nc_PolyPolyRed(_p, p2, coef, currRing);
819      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
820      PR->pLength=0; // usaully not used, GetpLength re-comoutes it if needed
821    }
822    return 0;
823  }
824#endif
825
826  if (t2==NULL)           // Divisor is just one term, therefore it will
827  {                       // just cancel the leading term
828    PR->LmDeleteAndIter();
829    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
830    return 0;
831  }
832
833  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
834
835  if (tailRing != currRing)
836  {
837    // check that reduction does not violate exp bound
838    while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
839    {
840      // undo changes of lm
841      p_ExpVectorAdd(lm, p2, tailRing);
842      if (strat == NULL) return 2;
843      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
844      tailRing = strat->tailRing;
845      p1 = PR->GetLmTailRing();
846      p2 = PW->GetLmTailRing();
847      t2 = pNext(p2);
848      lm = p1;
849      p_ExpVectorSub(lm, p2, tailRing);
850      ret = 1;
851    }
852  }
853
854#ifdef HAVE_SHIFTBBA
855  poly lmRight;
856  if (tailRing->isLPring)
857  {
858    assume(PR->shift == 0);
859    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
860    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
861  }
862#endif
863
864  // take care of coef buisness
865  if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
866  {
867    number bn = pGetCoeff(lm);
868    number an = pGetCoeff(p2);
869    int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
870    p_SetCoeff(lm, bn, tailRing);
871    if ((ct == 0) || (ct == 2))
872      PR->Tail_Mult_nn(an);
873    if (coef != NULL) *coef = an;
874    else n_Delete(&an, tailRing->cf);
875  }
876  else
877  {
878    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
879  }
880
881
882  // and finally,
883#ifdef HAVE_SHIFTBBA
884  if (tailRing->isLPring)
885  {
886    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
887  }
888  else
889#endif
890  {
891    PR->Tail_Minus_mm_Mult_qq(lm, t2, PW->GetpLength() - 1, spNoether);
892  }
893  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
894  PR->LmDeleteAndIter();
895
896#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
897  if (TEST_OPT_DEBUG)
898  {
899    Print(" to: "); PR->wrp(); Print("\n");
900  }
901#endif
902  return ret;
903}
904
905int ksReducePolySigRing(LObject* PR,
906                 TObject* PW,
907                 long /*idx*/,
908                 poly spNoether,
909                 number *coef,
910                 kStrategy strat)
911{
912#ifdef KDEBUG
913  red_count++;
914#ifdef TEST_OPT_DEBUG_RED
915  if (TEST_OPT_DEBUG)
916  {
917    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
918    PW->wrp();
919  }
920#endif
921#endif
922  int ret = 0;
923  ring tailRing = PR->tailRing;
924  kTest_L(PR,tailRing);
925  kTest_T(PW);
926
927  // signature-based stuff:
928  // checking for sig-safeness first
929  // NOTE: This has to be done in the current ring
930  //
931  /**********************************************
932   *
933   * TODO:
934   * --------------------------------------------
935   * if strat->sbaOrder == 1
936   * Since we are subdividing lower index and
937   * current index reductions it is enough to
938   * look at the polynomial part of the signature
939   * for a check. This should speed-up checking
940   * a lot!
941   * if !strat->sbaOrder == 0
942   * We are not subdividing lower and current index
943   * due to the fact that we are using the induced
944   * Schreyer order
945   *
946   * nevertheless, this different behaviour is
947   * taken care of by is_sigsafe
948   * => one reduction procedure can be used for
949   * both, the incremental and the non-incremental
950   * attempt!
951   * --------------------------------------------
952   *
953   *********************************************/
954  //printf("COMPARE IDX: %ld -- %ld\n",idx,strat->currIdx);
955  if (!PW->is_sigsafe)
956  {
957    poly sigMult = pCopy(PW->sig);   // copy signature of reducer
958//#if 1
959#ifdef DEBUGF5
960    printf("IN KSREDUCEPOLYSIG: \n");
961    pWrite(pHead(f1));
962    pWrite(pHead(f2));
963    pWrite(sigMult);
964    printf("--------------\n");
965#endif
966    p_ExpVectorAddSub(sigMult,PR->GetLmCurrRing(),PW->GetLmCurrRing(),currRing);
967    //I have also to set the leading coeficient for sigMult (in the case of rings)
968    if(rField_is_Ring(currRing))
969    {
970      pSetCoeff(sigMult,nMult(nDiv(pGetCoeff(PR->p),pGetCoeff(PW->p)), pGetCoeff(sigMult)));
971      if(nIsZero(pGetCoeff(sigMult)))
972      {
973        sigMult = NULL;
974      }
975    }
976//#if 1
977#ifdef DEBUGF5
978    printf("------------------- IN KSREDUCEPOLYSIG: --------------------\n");
979    pWrite(pHead(f1));
980    pWrite(pHead(f2));
981    pWrite(sigMult);
982    pWrite(PR->sig);
983    printf("--------------\n");
984#endif
985    int sigSafe;
986    if(!rField_is_Ring(currRing))
987      sigSafe = p_LmCmp(PR->sig,sigMult,currRing);
988    // now we can delete the copied polynomial data used for checking for
989    // sig-safeness of the reduction step
990//#if 1
991#ifdef DEBUGF5
992    printf("%d -- %d sig\n",sigSafe,PW->is_sigsafe);
993
994#endif
995    if(rField_is_Ring(currRing))
996    {
997      // Set the sig
998      poly origsig = pCopy(PR->sig);
999      if(sigMult != NULL)
1000        PR->sig = pHead(pSub(PR->sig, sigMult));
1001      //The sigs have the same lm, have to substract
1002      //It may happen that now the signature is 0 (drop)
1003      if(PR->sig == NULL)
1004      {
1005        strat->sigdrop=TRUE;
1006      }
1007      else
1008      {
1009        if(pLtCmp(PR->sig,origsig) == 1)
1010        {
1011          // do not allow this reduction - it will increase it's signature
1012          // and the partially standard basis is just till the old sig, not the new one
1013          PR->is_redundant = TRUE;
1014          pDelete(&PR->sig);
1015          PR->sig = origsig;
1016          strat->blockred++;
1017          return 3;
1018        }
1019        if(pLtCmp(PR->sig,origsig) == -1)
1020        {
1021          strat->sigdrop=TRUE;
1022        }
1023      }
1024      pDelete(&origsig);
1025    }
1026    //pDelete(&f1);
1027    // go on with the computations only if the signature of p2 is greater than the
1028    // signature of fm*p1
1029    if(sigSafe != 1 && !rField_is_Ring(currRing))
1030    {
1031      PR->is_redundant = TRUE;
1032      return 3;
1033    }
1034    //PW->is_sigsafe  = TRUE;
1035  }
1036  PR->is_redundant = FALSE;
1037  poly p1 = PR->GetLmTailRing();   // p2 | p1
1038  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
1039  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
1040  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
1041  p_CheckPolyRing(p1, tailRing);
1042  p_CheckPolyRing(p2, tailRing);
1043
1044  pAssume1(p2 != NULL && p1 != NULL &&
1045      p_DivisibleBy(p2,  p1, tailRing));
1046
1047  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
1048      (p_GetComp(p2, tailRing) == 0 &&
1049       p_MaxComp(pNext(p2),tailRing) == 0));
1050
1051#ifdef HAVE_PLURAL
1052  if (rIsPluralRing(currRing))
1053  {
1054    // for the time being: we know currRing==strat->tailRing
1055    // no exp-bound checking needed
1056    // (only needed if exp-bound(tailring)<exp-b(currRing))
1057    if (PR->bucket!=NULL)  nc_kBucketPolyRed_Z(PR->bucket, p2,coef);
1058    else
1059    {
1060      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
1061      assume(_p != NULL);
1062      nc_PolyPolyRed(_p, p2, coef, currRing);
1063      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
1064      PR->pLength=0; // usaully not used, GetpLength re-comoutes it if needed
1065    }
1066    return 0;
1067  }
1068#endif
1069
1070  if (t2==NULL)           // Divisor is just one term, therefore it will
1071  {                       // just cancel the leading term
1072    PR->LmDeleteAndIter();
1073    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
1074    return 0;
1075  }
1076
1077  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
1078
1079  if (tailRing != currRing)
1080  {
1081    // check that reduction does not violate exp bound
1082    while (PW->max_exp != NULL && !p_LmExpVectorAddIsOk(lm, PW->max_exp, tailRing))
1083    {
1084      // undo changes of lm
1085      p_ExpVectorAdd(lm, p2, tailRing);
1086      if (strat == NULL) return 2;
1087      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
1088      tailRing = strat->tailRing;
1089      p1 = PR->GetLmTailRing();
1090      p2 = PW->GetLmTailRing();
1091      t2 = pNext(p2);
1092      lm = p1;
1093      p_ExpVectorSub(lm, p2, tailRing);
1094      ret = 1;
1095    }
1096  }
1097
1098#ifdef HAVE_SHIFTBBA
1099  poly lmRight;
1100  if (tailRing->isLPring)
1101  {
1102    assume(PR->shift == 0);
1103    assume(PW->shift == si_max(p_mFirstVblock(PW->p, tailRing) - 1, 0));
1104    k_SplitFrame(lm, lmRight, PW->shift + 1, tailRing);
1105  }
1106#endif
1107
1108  // take care of coef buisness
1109  if(rField_is_Ring(currRing))
1110  {
1111    p_SetCoeff(lm, nDiv(pGetCoeff(lm),pGetCoeff(p2)), tailRing);
1112    if (coef != NULL) *coef = n_Init(1, tailRing->cf);
1113  }
1114  else
1115  {
1116    if (! n_IsOne(pGetCoeff(p2), tailRing->cf))
1117    {
1118      number bn = pGetCoeff(lm);
1119      number an = pGetCoeff(p2);
1120      int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
1121      p_SetCoeff(lm, bn, tailRing);
1122      if (((ct == 0) || (ct == 2)))
1123        PR->Tail_Mult_nn(an);
1124      if (coef != NULL) *coef = an;
1125      else n_Delete(&an, tailRing->cf);
1126    }
1127    else
1128    {
1129      if (coef != NULL) *coef = n_Init(1, tailRing->cf);
1130    }
1131  }
1132
1133  // and finally,
1134#ifdef HAVE_SHIFTBBA
1135  if (tailRing->isLPring)
1136  {
1137    PR->Tail_Minus_mm_Mult_qq(lm, tailRing->p_Procs->pp_Mult_mm(t2, lmRight, tailRing), pLength(t2), spNoether);
1138  }
1139  else
1140#endif
1141  {
1142    PR->Tail_Minus_mm_Mult_qq(lm, t2, PW->GetpLength() - 1, spNoether);
1143  }
1144  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
1145  PR->LmDeleteAndIter();
1146
1147#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
1148  if (TEST_OPT_DEBUG)
1149  {
1150    Print(" to: "); PR->wrp(); Print("\n");
1151  }
1152#endif
1153  return ret;
1154}
1155
1156/***************************************************************
1157 *
1158 * Creates S-Poly of p1 and p2
1159 *
1160 *
1161 ***************************************************************/
1162void ksCreateSpoly(LObject* Pair,   poly spNoether,
1163                   int use_buckets, ring tailRing,
1164                   poly m1, poly m2, TObject** R)
1165{
1166#ifdef KDEBUG
1167  create_count++;
1168#endif
1169  kTest_L(Pair,tailRing);
1170  poly p1 = Pair->p1;
1171  poly p2 = Pair->p2;
1172  Pair->tailRing = tailRing;
1173
1174  assume(p1 != NULL);
1175  assume(p2 != NULL);
1176  assume(tailRing != NULL);
1177
1178  poly a1 = pNext(p1), a2 = pNext(p2);
1179  number lc1 = pGetCoeff(p1), lc2 = pGetCoeff(p2);
1180  int co=0/*, ct = ksCheckCoeff(&lc1, &lc2, currRing->cf)*/; // gcd and zero divisors
1181  (void) ksCheckCoeff(&lc1, &lc2, currRing->cf);
1182
1183  int l1=0, l2=0;
1184
1185  if (currRing->pCompIndex >= 0)
1186  {
1187    if (__p_GetComp(p1, currRing)!=__p_GetComp(p2, currRing))
1188    {
1189      if (__p_GetComp(p1, currRing)==0)
1190      {
1191        co=1;
1192        p_SetCompP(p1,__p_GetComp(p2, currRing), currRing, tailRing);
1193      }
1194      else
1195      {
1196        co=2;
1197        p_SetCompP(p2, __p_GetComp(p1, currRing), currRing, tailRing);
1198      }
1199    }
1200  }
1201
1202  // get m1 = LCM(LM(p1), LM(p2))/LM(p1)
1203  //     m2 = LCM(LM(p1), LM(p2))/LM(p2)
1204  if (m1 == NULL)
1205    k_GetLeadTerms(p1, p2, currRing, m1, m2, tailRing);
1206
1207#ifdef HAVE_SHIFTBBA
1208  poly m12, m22;
1209  if (tailRing->isLPring)
1210  {
1211    assume(p_mFirstVblock(p1, tailRing) <= 1 || p_mFirstVblock(p2, tailRing) <= 1);
1212    k_SplitFrame(m1, m12, si_max(p_mFirstVblock(p1, tailRing), 1), tailRing);
1213    k_SplitFrame(m2, m22, si_max(p_mFirstVblock(p2, tailRing), 1), tailRing);
1214    // manually free the coeffs, because pSetCoeff0 is used in the next step
1215    n_Delete(&(m1->coef), tailRing->cf);
1216    n_Delete(&(m2->coef), tailRing->cf);
1217  }
1218#endif
1219
1220  pSetCoeff0(m1, lc2);
1221  pSetCoeff0(m2, lc1);  // and now, m1 * LT(p1) == m2 * LT(p2)
1222
1223  if (R != NULL)
1224  {
1225    if (Pair->i_r1 == -1)
1226    {
1227      l1 = pLength(p1) - 1;
1228    }
1229    else
1230    {
1231      l1 = (R[Pair->i_r1])->GetpLength() - 1;
1232    }
1233    if ((Pair->i_r2 == -1)||(R[Pair->i_r2]==NULL))
1234    {
1235      l2 = pLength(p2) - 1;
1236    }
1237    else
1238    {
1239      l2 = (R[Pair->i_r2])->GetpLength() - 1;
1240    }
1241  }
1242
1243  // get m2 * a2
1244  if (spNoether != NULL)
1245  {
1246    l2 = -1;
1247    a2 = tailRing->p_Procs->pp_Mult_mm_Noether(a2, m2, spNoether, l2, tailRing);
1248    assume(l2 == pLength(a2));
1249  }
1250  else
1251#ifdef HAVE_SHIFTBBA
1252    if (tailRing->isLPring)
1253    {
1254      // m2*a2*m22
1255      a2 = tailRing->p_Procs->pp_Mult_mm(tailRing->p_Procs->pp_mm_Mult(a2, m2, tailRing), m22, tailRing);
1256    }
1257    else
1258#endif
1259    {
1260      a2 = tailRing->p_Procs->pp_Mult_mm(a2, m2, tailRing);
1261    }
1262#ifdef HAVE_RINGS
1263  if (!(rField_is_Domain(currRing))) l2 = pLength(a2);
1264#endif
1265
1266  Pair->SetLmTail(m2, a2, l2, use_buckets, tailRing);
1267
1268#ifdef HAVE_SHIFTBBA
1269  if (tailRing->isLPring)
1270  {
1271    // get m2*a2*m22 - m1*a1*m12
1272    Pair->Tail_Minus_mm_Mult_qq(m1, tailRing->p_Procs->pp_Mult_mm(a1, m12, tailRing), l1, spNoether);
1273  }
1274  else
1275#endif
1276  {
1277    // get m2*a2 - m1*a1
1278    Pair->Tail_Minus_mm_Mult_qq(m1, a1, l1, spNoether);
1279  }
1280
1281  // Clean-up time
1282  Pair->LmDeleteAndIter();
1283  p_LmDelete(m1, tailRing);
1284#ifdef HAVE_SHIFTBBA
1285  if (tailRing->isLPring)
1286  {
1287    // just to be sure, check that the shift is correct
1288    assume(Pair->shift == 0);
1289    assume(si_max(p_mFirstVblock(Pair->p, tailRing) - 1, 0) == Pair->shift); // == 0
1290
1291    p_LmDelete(m12, tailRing);
1292    p_LmDelete(m22, tailRing);
1293    // m2 is already deleted
1294  }
1295#endif
1296
1297  if (co != 0)
1298  {
1299    if (co==1)
1300    {
1301      p_SetCompP(p1,0, currRing, tailRing);
1302    }
1303    else
1304    {
1305      p_SetCompP(p2,0, currRing, tailRing);
1306    }
1307  }
1308}
1309
1310int ksReducePolyTail(LObject* PR, TObject* PW, poly Current, poly spNoether)
1311{
1312  BOOLEAN ret;
1313  number coef;
1314  poly Lp =     PR->GetLmCurrRing();
1315  poly Save =   PW->GetLmCurrRing();
1316
1317  kTest_L(PR,PR->tailRing);
1318  kTest_T(PW);
1319  pAssume(pIsMonomOf(Lp, Current));
1320
1321  assume(Lp != NULL && Current != NULL && pNext(Current) != NULL);
1322  assume(PR->bucket == NULL);
1323
1324  LObject Red(pNext(Current), PR->tailRing);
1325  TObject With(PW, Lp == Save);
1326
1327  pAssume(!pHaveCommonMonoms(Red.p, With.p));
1328  ret = ksReducePoly(&Red, &With, spNoether, &coef);
1329
1330  if (!ret)
1331  {
1332    if (! n_IsOne(coef, currRing->cf))
1333    {
1334      pNext(Current) = NULL;
1335      if (Current == PR->p && PR->t_p != NULL)
1336        pNext(PR->t_p) = NULL;
1337      PR->Mult_nn(coef);
1338    }
1339
1340    n_Delete(&coef, currRing->cf);
1341    pNext(Current) = Red.GetLmTailRing();
1342    if (Current == PR->p && PR->t_p != NULL)
1343      pNext(PR->t_p) = pNext(Current);
1344  }
1345
1346  if (Lp == Save)
1347    With.Delete();
1348
1349  return ret;
1350}
1351
1352int ksReducePolyTailBound(LObject* PR, TObject* PW, int bound, poly Current, poly spNoether)
1353{
1354  BOOLEAN ret;
1355  number coef;
1356  poly Lp =     PR->GetLmCurrRing();
1357  poly Save =   PW->GetLmCurrRing();
1358
1359  kTest_L(PR,PR->tailRing);
1360  kTest_T(PW);
1361  pAssume(pIsMonomOf(Lp, Current));
1362
1363  assume(Lp != NULL && Current != NULL && pNext(Current) != NULL);
1364  assume(PR->bucket == NULL);
1365
1366  LObject Red(pNext(Current), PR->tailRing);
1367  TObject With(PW, Lp == Save);
1368
1369  pAssume(!pHaveCommonMonoms(Red.p, With.p));
1370  ret = ksReducePolyBound(&Red, &With,bound, spNoether, &coef);
1371
1372  if (!ret)
1373  {
1374    if (! n_IsOne(coef, currRing))
1375    {
1376      pNext(Current) = NULL;
1377      if (Current == PR->p && PR->t_p != NULL)
1378        pNext(PR->t_p) = NULL;
1379      PR->Mult_nn(coef);
1380    }
1381
1382    n_Delete(&coef, currRing);
1383    pNext(Current) = Red.GetLmTailRing();
1384    if (Current == PR->p && PR->t_p != NULL)
1385      pNext(PR->t_p) = pNext(Current);
1386  }
1387
1388  if (Lp == Save)
1389    With.Delete();
1390
1391  return ret;
1392}
1393
1394/***************************************************************
1395 *
1396 * Auxillary Routines
1397 *
1398 *
1399 ***************************************************************/
1400
1401/*2
1402* creates the leading term of the S-polynomial of p1 and p2
1403* do not destroy p1 and p2
1404* remarks:
1405*   1. the coefficient is 0 (p_Init)
1406*   1. a) in the case of coefficient ring, the coefficient is calculated
1407*   2. pNext is undefined
1408*/
1409//static void bbb() { int i=0; }
1410poly ksCreateShortSpoly(poly p1, poly p2, ring tailRing)
1411{
1412  poly a1 = pNext(p1), a2 = pNext(p2);
1413#ifdef HAVE_SHIFTBBA
1414  int shift1, shift2;
1415  if (tailRing->isLPring) {
1416    // assume: LM is shifted, tail unshifted
1417    assume(p_FirstVblock(a1, tailRing) <= 1);
1418    assume(p_FirstVblock(a2, tailRing) <= 1);
1419    // save the shift of the LM so we can shift the other monomials on demand
1420    shift1 = p_mFirstVblock(p1, tailRing) - 1;
1421    shift2 = p_mFirstVblock(p2, tailRing) - 1;
1422  }
1423#endif
1424  long c1=p_GetComp(p1, currRing),c2=p_GetComp(p2, currRing);
1425  long c;
1426  poly m1,m2;
1427  number t1 = NULL,t2 = NULL;
1428  int cm,i;
1429  BOOLEAN equal;
1430
1431#ifdef HAVE_RINGS
1432  BOOLEAN is_Ring=rField_is_Ring(currRing);
1433  number lc1 = pGetCoeff(p1), lc2 = pGetCoeff(p2);
1434  if (is_Ring)
1435  {
1436    ksCheckCoeff(&lc1, &lc2, currRing->cf); // gcd and zero divisors
1437    if (a1 != NULL) t2 = nMult(pGetCoeff(a1),lc2);
1438    if (a2 != NULL) t1 = nMult(pGetCoeff(a2),lc1);
1439    while (a1 != NULL && nIsZero(t2))
1440    {
1441      pIter(a1);
1442      nDelete(&t2);
1443      if (a1 != NULL) t2 = nMult(pGetCoeff(a1),lc2);
1444    }
1445    while (a2 != NULL && nIsZero(t1))
1446    {
1447      pIter(a2);
1448      nDelete(&t1);
1449      if (a2 != NULL) t1 = nMult(pGetCoeff(a2),lc1);
1450    }
1451  }
1452#endif
1453
1454#ifdef HAVE_SHIFTBBA
1455  // shift the next monomial on demand
1456  if (tailRing->isLPring)
1457  {
1458    a1 = p_LPCopyAndShiftLM(a1, shift1, tailRing);
1459    a2 = p_LPCopyAndShiftLM(a2, shift2, tailRing);
1460  }
1461#endif
1462  if (a1==NULL)
1463  {
1464    if(a2!=NULL)
1465    {
1466      m2=p_Init(currRing);
1467x2:
1468      for (i = (currRing->N); i; i--)
1469      {
1470        c = p_GetExpDiff(p1, p2,i, currRing);
1471        if (c>0)
1472        {
1473          p_SetExp(m2,i,(c+p_GetExp(a2,i,tailRing)),currRing);
1474        }
1475        else
1476        {
1477          p_SetExp(m2,i,p_GetExp(a2,i,tailRing),currRing);
1478        }
1479      }
1480      if ((c1==c2)||(c2!=0))
1481      {
1482        p_SetComp(m2,p_GetComp(a2,tailRing), currRing);
1483      }
1484      else
1485      {
1486        p_SetComp(m2,c1,currRing);
1487      }
1488      p_Setm(m2, currRing);
1489#ifdef HAVE_RINGS
1490      if (is_Ring)
1491      {
1492          nDelete(&lc1);
1493          nDelete(&lc2);
1494          nDelete(&t2);
1495          pSetCoeff0(m2, t1);
1496      }
1497#endif
1498      return m2;
1499    }
1500    else
1501    {
1502#ifdef HAVE_RINGS
1503      if (is_Ring)
1504      {
1505        nDelete(&lc1);
1506        nDelete(&lc2);
1507        nDelete(&t1);
1508        nDelete(&t2);
1509      }
1510#endif
1511      return NULL;
1512    }
1513  }
1514  if (a2==NULL)
1515  {
1516    m1=p_Init(currRing);
1517x1:
1518    for (i = (currRing->N); i; i--)
1519    {
1520      c = p_GetExpDiff(p2, p1,i,currRing);
1521      if (c>0)
1522      {
1523        p_SetExp(m1,i,(c+p_GetExp(a1,i, tailRing)),currRing);
1524      }
1525      else
1526      {
1527        p_SetExp(m1,i,p_GetExp(a1,i, tailRing), currRing);
1528      }
1529    }
1530    if ((c1==c2)||(c1!=0))
1531    {
1532      p_SetComp(m1,p_GetComp(a1,tailRing),currRing);
1533    }
1534    else
1535    {
1536      p_SetComp(m1,c2,currRing);
1537    }
1538    p_Setm(m1, currRing);
1539#ifdef HAVE_RINGS
1540    if (is_Ring)
1541    {
1542      pSetCoeff0(m1, t2);
1543      nDelete(&lc1);
1544      nDelete(&lc2);
1545      nDelete(&t1);
1546    }
1547#endif
1548    return m1;
1549  }
1550  m1 = p_Init(currRing);
1551  m2 = p_Init(currRing);
1552  loop
1553  {
1554    for (i = (currRing->N); i; i--)
1555    {
1556      c = p_GetExpDiff(p1, p2,i,currRing);
1557      if (c > 0)
1558      {
1559        p_SetExp(m2,i,(c+p_GetExp(a2,i,tailRing)), currRing);
1560        p_SetExp(m1,i,p_GetExp(a1,i, tailRing), currRing);
1561      }
1562      else
1563      {
1564        p_SetExp(m1,i,(p_GetExp(a1,i,tailRing)-c), currRing);
1565        p_SetExp(m2,i,p_GetExp(a2,i, tailRing), currRing);
1566      }
1567    }
1568    if(c1==c2)
1569    {
1570      p_SetComp(m1,p_GetComp(a1, tailRing), currRing);
1571      p_SetComp(m2,p_GetComp(a2, tailRing), currRing);
1572    }
1573    else
1574    {
1575      if(c1!=0)
1576      {
1577        p_SetComp(m1,p_GetComp(a1, tailRing), currRing);
1578        p_SetComp(m2,c1, currRing);
1579      }
1580      else
1581      {
1582        p_SetComp(m2,p_GetComp(a2, tailRing), currRing);
1583        p_SetComp(m1,c2, currRing);
1584      }
1585    }
1586    p_Setm(m1,currRing);
1587    p_Setm(m2,currRing);
1588    cm = p_LmCmp(m1, m2,currRing);
1589    if (cm!=0)
1590    {
1591      if(cm==1)
1592      {
1593        p_LmFree(m2,currRing);
1594#ifdef HAVE_RINGS
1595        if (is_Ring)
1596        {
1597          pSetCoeff0(m1, t2);
1598          nDelete(&lc1);
1599          nDelete(&lc2);
1600          nDelete(&t1);
1601        }
1602#endif
1603        return m1;
1604      }
1605      else
1606      {
1607        p_LmFree(m1,currRing);
1608#ifdef HAVE_RINGS
1609        if (is_Ring)
1610        {
1611          pSetCoeff0(m2, t1);
1612          nDelete(&lc1);
1613          nDelete(&lc2);
1614          nDelete(&t2);
1615        }
1616#endif
1617        return m2;
1618      }
1619    }
1620#ifdef HAVE_RINGS
1621    if (is_Ring)
1622    {
1623      equal = nEqual(t1,t2);
1624    }
1625    else
1626#endif
1627    {
1628      t1 = nMult(pGetCoeff(a2),pGetCoeff(p1));
1629      t2 = nMult(pGetCoeff(a1),pGetCoeff(p2));
1630      equal = nEqual(t1,t2);
1631      nDelete(&t2);
1632      nDelete(&t1);
1633    }
1634    if (!equal)
1635    {
1636      p_LmFree(m2,currRing);
1637#ifdef HAVE_RINGS
1638      if (is_Ring)
1639      {
1640          pSetCoeff0(m1, nSub(t1, t2));
1641          nDelete(&lc1);
1642          nDelete(&lc2);
1643          nDelete(&t1);
1644          nDelete(&t2);
1645      }
1646#endif
1647      return m1;
1648    }
1649    pIter(a1);
1650    pIter(a2);
1651#ifdef HAVE_RINGS
1652    if (is_Ring)
1653    {
1654      if (a2 != NULL)
1655      {
1656        nDelete(&t1);
1657        t1 = nMult(pGetCoeff(a2),lc1);
1658      }
1659      if (a1 != NULL)
1660      {
1661        nDelete(&t2);
1662        t2 = nMult(pGetCoeff(a1),lc2);
1663      }
1664      while ((a1 != NULL) && nIsZero(t2))
1665      {
1666        pIter(a1);
1667        if (a1 != NULL)
1668        {
1669          nDelete(&t2);
1670          t2 = nMult(pGetCoeff(a1),lc2);
1671        }
1672      }
1673      while ((a2 != NULL) && nIsZero(t1))
1674      {
1675        pIter(a2);
1676        if (a2 != NULL)
1677        {
1678          nDelete(&t1);
1679          t1 = nMult(pGetCoeff(a2),lc1);
1680        }
1681      }
1682    }
1683#endif
1684#ifdef HAVE_SHIFTBBA
1685    if (tailRing->isLPring)
1686    {
1687      a1 = p_LPCopyAndShiftLM(a1, shift1, tailRing);
1688      a2 = p_LPCopyAndShiftLM(a2, shift2, tailRing);
1689    }
1690#endif
1691    if (a2==NULL)
1692    {
1693      p_LmFree(m2,currRing);
1694      if (a1==NULL)
1695      {
1696#ifdef HAVE_RINGS
1697        if (is_Ring)
1698        {
1699          nDelete(&lc1);
1700          nDelete(&lc2);
1701          nDelete(&t1);
1702          nDelete(&t2);
1703        }
1704#endif
1705        p_LmFree(m1,currRing);
1706        return NULL;
1707      }
1708      goto x1;
1709    }
1710    if (a1==NULL)
1711    {
1712      p_LmFree(m1,currRing);
1713      goto x2;
1714    }
1715  }
1716}
Note: See TracBrowser for help on using the repository browser.