source: git/kernel/GBEngine/kspoly.cc @ 036a5e

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