source: git/kernel/GBEngine/kspoly.cc @ a3f0fea

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