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

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