source: git/kernel/ratgring.cc @ 4732f8

spielwiese
Last change on this file since 4732f8 was 4732f8, checked in by Viktor Levandovskyy <levandov@…>, 15 years ago
*levandov: ncratCreateSpoly added git-svn-id: file:///usr/local/Singular/svn/trunk@11432 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 16.7 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/***************************************************************
5 *  File:    ratgring.cc
6 *  Purpose: Ore-noncommutative kernel procedures
7 *  Author:  levandov (Viktor Levandovsky)
8 *  Created: 8/00 - 11/00
9 *  Version: $Id: ratgring.cc,v 1.14 2009-02-21 19:30:56 levandov Exp $
10 *******************************************************************/
11#include "mod2.h"
12#include "ratgring.h"
13#ifdef HAVE_RATGRING
14#include "gring.h"
15#include "febase.h"
16#include "ring.h"
17#include "polys.h"
18#include "numbers.h"
19#include "ideals.h"
20#include "matpol.h"
21#include "kbuckets.h"
22#include "kstd1.h"
23#include "sbuckets.h"
24#include "prCopy.h"
25#include "p_Mult_q.h"
26#include "clapsing.h"
27
28void pLcmRat(poly a, poly b, poly m, int rat_shift)
29{
30  /* rat_shift is the last exp one should count with */
31  int i;
32  for (i=pVariables; i>=rat_shift; i--)
33  {
34    pSetExp(m,i, si_max( pGetExp(a,i), pGetExp(b,i)));
35  }
36  pSetComp(m, si_max(pGetComp(a), pGetComp(b)));
37  /* Don't do a pSetm here, otherwise hres/lres chockes */ 
38}
39
40// void pLcmRat(poly a, poly b, poly m, poly pshift)
41// {
42//   /* shift is the exp of rational elements */
43//   int i;
44//   for (i=pVariables; i; i--)
45//   {
46//     if (!pGetExp(pshift,i))
47//     {
48//       pSetExp(m,i, si_max( pGetExp(a,i), pGetExp(b,i)));
49//     }
50//     else
51//     {
52//       /* do we really need it? */
53//       pSetExp(m,i,0);
54//     }
55//   }
56//   pSetComp(m, si_max(pGetComp(a), pGetComp(b)));
57//   /* Don't do a pSetm here, otherwise hres/lres chockes */
58// }
59
60/* returns a subpoly of p, s.t. its monomials have the same D-part */
61
62poly p_HeadRat(poly p, int ishift, ring r)
63{
64  poly q   = pNext(p);
65  if (q == NULL) return p;
66  poly res = p_Head(p,r);
67  while ( (q!=NULL) && (p_Comp_k_n(p, q, ishift+1, r)))
68  {
69    res = p_Add_q(res,p_Head(q,r),r);
70    q   = pNext(q);
71  }
72  return res;
73}
74
75/* TO TEST!!! */
76/* returns x-coeff of p, i.e. a poly in x, s.t. corresponding xd-monomials
77have the same D-part
78does not destroy p
79*/
80
81poly p_GetCoeffRat(poly p, int ishift, ring r)
82{
83  poly q   = pNext(p);
84  poly res; // = p_Head(p,r);
85  res = p_GetExp_k_n(p, ishift+1, r->N, r);
86  p_SetCoeff(res,n_Copy(p_GetCoeff(p,r),r),r);
87  poly s;
88  while ((q!= NULL) && (p_Comp_k_n(p, q, ishift+1, r)))
89  {
90    s   = p_GetExp_k_n(q, ishift+1, r->N, r);
91    p_SetCoeff(s,n_Copy(p_GetCoeff(q,r),r),r);
92    res = p_Add_q(res,s,r);
93    q   = pNext(q);
94  }
95  return res;
96}
97
98void p_LmDeleteAndNextRat(poly *p, int ishift, ring r)
99{
100  /* modifies p*/
101  //  Print("start: "); Print(" "); p_wrp(*p,r);
102  p_LmCheckPolyRing2(*p, r);
103  poly q = p_Head(*p,r);
104  // in the next line ishift is correct
105  while ( ( (*p)!=NULL ) && ( p_Comp_k_n(*p, q, ishift, r) ))
106  {
107    p_LmDelete(p,r);
108    //    Print("while: ");p_wrp(*p,r);Print(" ");
109  }
110  //  p_wrp(*p,r);Print(" ");
111  //  PrintS("end\n");
112  p_LmDelete(&q,r);
113}
114
115/* to test!!! */
116/* ExpVector(pr) = ExpVector(p1) - ExpVector(p2) */
117void p_ExpVectorDiffRat(poly pr, poly p1, poly p2, int ishift, ring r)
118{
119  p_LmCheckPolyRing1(p1, r);
120  p_LmCheckPolyRing1(p2, r);
121  p_LmCheckPolyRing1(pr, r); 
122  int i;
123  poly t=pr;
124  Exponent_t e1,e2;
125  for (i=ishift+1; i<=r->N; i++)
126  {
127    e1 = p_GetExp(p1, i, r);
128    e2 = p_GetExp(p2, i, r);
129    //    pAssume1(p_GetExp(p1, i, r) >= p_GetExp(p2, i, r));
130    if (e1 < e2)
131    {
132#ifdef PDEBUG
133      Print("negative ExpVectorDiff\n");
134#endif   
135      p_Delete(&t,r);
136      break;
137    }
138    else
139    {
140      p_SetExp(t,i, e1-e2,r);
141    }
142  }
143  p_Setm(t,r);
144}
145
146/* returns ideal (u,v) s.t. up + vq = 0 */
147
148ideal ncGCD2(poly p, poly q, const ring r)
149{
150  // todo: must destroy p,q
151  intvec *w = NULL;
152  ideal h = idInit(2,1);
153  h->m[0] = p_Copy(p,r);
154  h->m[1] = p_Copy(q,r);
155#ifdef PDEBUG
156  Print("running syzygy comp. for nc_GCD:\n");
157#endif
158  ideal sh = idSyzygies(h, testHomog, &w);
159#ifdef PDEBUG
160  Print("done syzygy comp. for nc_GCD\n");
161#endif
162  /* in comm case, there is only 1 syzygy */
163  /*   singclap_gcd(); */
164  poly K, K1, K2;
165  K  = sh->m[0]; /* take just the first element - to be enhanced later */
166  K1 = pTakeOutComp(&K, 1); // 1st component is taken out from K
167//  pShift(&K,-2); // 2nd component to 0th comp.
168  K2 = pTakeOutComp(&K, 1);
169//  K2 = K;
170
171  Print("syz1: "); p_wrp(K1,r);
172  Print("syz2: "); p_wrp(K2,r);
173
174  /* checking signs before multiplying */   
175  number ck1 = p_GetCoeff(K1,r);
176  number ck2 = p_GetCoeff(K2,r);
177  BOOLEAN bck1, bck2;
178  bck1 = n_GreaterZero(ck1,r);
179  bck2 = n_GreaterZero(ck2,r);
180  /* K1 <0, K2 <0 (-K1,-K2)    */
181//   if ( !(bck1 && bck2) ) /* - , - */
182//   {
183//     K1 = p_Neg(K1,r);
184//     K2 = p_Neg(K2,r);
185//   }
186  id_Delete(&h,r);
187  h = idInit(2,1);
188  h->m[0] = p_Copy(K1,r);
189  h->m[1] = p_Copy(K2,r);
190  id_Delete(&sh,r);
191  return(h);
192}
193
194/* returns ideal (u,v) s.t. up + vq = 0 */
195
196ideal ncGCD(poly p, poly q, const ring r)
197{
198  // destroys p and q
199  // assume: p,q are in the comm. ring
200  // to be used in the coeff business
201#ifdef PDEBUG
202  PrintS(" GCD_start:");
203#endif
204  poly g = singclap_gcd(p_Copy(p,r),p_Copy(q,r));
205#ifdef PDEBUG
206  p_wrp(g,r);
207  PrintS(" GCD_end;\n");
208#endif
209  poly u = singclap_pdivide(q,g); //q/g
210  poly v = singclap_pdivide(p,g); //p/g
211  v = p_Neg(v,r);
212  p_Delete(&p,r);
213  p_Delete(&q,r);
214  ideal h = idInit(2,1);
215  h->m[0] = u; // p_Copy(u,r);
216  h->m[1] = v; // p_Copy(v,r);
217  return(h);
218}
219
220/* PINLINE1 void p_ExpVectorDiff
221   remains as is -> BUT we can do memory shift on smaller number of exp's */
222
223
224/*4 - follow the numbering of gring.cc
225* creates the S-polynomial of p1 and p2
226* do not destroy p1 and p2
227*/
228// poly nc_rat_CreateSpoly(poly p1, poly p2, poly spNoether, int ishift, const ring r)
229// {
230//   if ((p_GetComp(p1,r)!=p_GetComp(p2,r))
231//   && (p_GetComp(p1,r)!=0)
232//   && (p_GetComp(p2,r)!=0))
233//   {
234// #ifdef PDEBUG
235//     Print("nc_CreateSpoly : different components!");
236// #endif
237//     return(NULL);
238//   }
239//   /* prod. crit does not apply yet */
240// //   if ((r->nc->type==nc_lie) && pHasNotCF(p1,p2)) /* prod crit */
241// //   {
242// //     return(nc_p_Bracket_qq(pCopy(p2),p1));
243// //   }
244//   poly pL=pOne();
245//   poly m1=pOne();
246//   poly m2=pOne();
247//   /* define shift */
248//   int is = ishift; /* TODO */
249//   pLcmRat(p1,p2,pL,is);
250//   p_Setm(pL,r);
251//   poly pr1 = p_GetExp_k_n(p1,1,ishift-1,r); /* rat D-exp of p1 */
252//   poly pr2 = p_GetExp_k_n(p2,1,ishift-1,r); /* rat D-exp of p2 */
253// #ifdef PDEBUG
254//   p_Test(pL,r);
255// #endif
256//   p_ExpVectorDiff(m1,pL,p1,r); /* purely in D part by construction */
257//   //p_SetComp(m1,0,r);
258//   //p_Setm(m1,r);
259// #ifdef PDEBUG
260//   p_Test(m1,r);
261// #endif
262//   p_ExpVectorDiff(m2,pL,p2,r); /* purely in D part by construction */
263//   //p_SetComp(m2,0,r);
264//   //p_Setm(m2,r);
265// #ifdef PDEBUG
266//   p_Test(m2,r);
267// #endif
268//   p_Delete(&pL,r);
269//   /* zero exponents ! */
270
271//   /* EXTRACT LEADCOEF */
272
273//   poly H1  = p_HeadRat(p1,is,r);
274//   poly M1  = r->nc->p_Procs.mm_Mult_p(m1,p_Copy(H1,r),r);
275
276//   /* POLY:  number C1  = n_Copy(p_GetCoeff(M1,r),r); */
277//   /* RAT: */
278
279//   poly C1  = p_GetCoeffRat(M1,ishift,r);
280
281//   poly H2  = p_HeadRat(p2,is,r);
282//   poly M2  = r->nc->p_Procs.mm_Mult_p(m2,p_Copy(H2,r),r);
283
284//   /* POLY:  number C2  = n_Copy(p_GetCoeff(M2,r),r); */
285//   /* RAT: */
286
287//   poly C2  = p_GetCoeffRat(M2,ishift,r);
288
289// /* we do not assume that X's commute */
290// /* we just run NC syzygies */
291
292// /* NEW IDEA: change the ring to K<X>, map things there
293//    and return the result back; seems to be a good optimization */
294// /* to be done later */
295// /* problem: map to subalgebra. contexts, induced (non-unique) orderings etc. */
296
297//   intvec *w = NULL;
298//   ideal h = idInit(2,1);
299//   h->m[0] = p_Copy(C1,r);
300//   h->m[1] = p_Copy(C2,r);
301// #ifdef PDEBUG
302//   Print("running syzygy comp. for coeffs");
303// #endif
304//   ideal sh = idSyzygies(h, testHomog, &w);
305//   /* in comm case, there is only 1 syzygy */
306//   /*   singclap_gcd(); */
307//   poly K,K1,K2;
308//   K  = sh->m[0];
309//   K1 = pTakeOutComp(&K, 1); // 1st component is taken out from K
310//   pShift(&K,-2); // 2nd component to 0th comp.
311//   K2 = K;
312
313//   /* checking signs before multiplying */   
314//   number ck1 = p_GetCoeff(K1,r);
315//   number ck2 = p_GetCoeff(K2,r);
316//   BOOLEAN bck1, bck2;
317//   bck1 = n_GreaterZero(ck1,r);
318//   bck2 = n_GreaterZero(ck2,r);
319//   /* K1 >0, K2 >0 (K1,-K2)    */
320//   /* K1 >0, K2 <0 (K1,-K2)    */
321//   /* K1 <0, K2 >0 (-K1,K2)    */
322//   /* K1 <0, K2 <0 (-K1,K2)    */
323//   if ( (bck1) && (bck2) ) /* +, + */
324//   {
325//     K2 = p_Neg(K2,r);
326//   }
327//   if ( (bck1) && (!bck2) ) /* + , - */
328//   {
329//     K2 = p_Neg(K2,r);
330//   }
331//   if ( (!bck1) && (bck2) ) /* - , + */
332//   {
333//     K1 = p_Neg(K1,r);
334//   }
335//   if ( !(bck1 && bck2) ) /* - , - */
336//   {
337//     K1 = p_Neg(K1,r);
338//   }
339
340//   poly P1,P2;
341
342//   //  p_LmDeleteRat(M1,ishift,r); // get tail(D^(gamma-alpha) * lm(p1)) = h_f
343//   P1 = p_Copy(p1,r);
344//   p_LmDeleteAndNextRat(P1,ishift,r); // get tail(p1) = t_f
345//   P1 = r->nc->p_Procs.mm_Mult_p(m1,P1,r);
346//   P1 = p_Add_q(P1,M1,r);
347
348//   //  p_LmDeleteRat(M2,ishift,r);
349//   P2 = p_Copy(p2,r);
350//   p_LmDeleteAndNextRat(P2,ishift,r);// get tail(p2)=t_g
351//   P2 = r->nc->p_Procs.mm_Mult_p(m2,P2,r);
352//   P2 = p_Add_q(P2,M2,r);
353
354//   /* coeff business */
355
356//   P1 = p_Mult_q(P1,K1,r);
357//   P2 = p_Mult_q(P2,K2,r);
358//   P1 = p_Add_q(P1,P2,r);
359
360//   /* cleaning up */
361
362// #ifdef PDEBUG
363//   p_Test(p1,r);
364// #endif
365//   /* questionable: */
366//   if (P1!=NULL) pCleardenom(P1);
367//   if (P1!=NULL) pContent(P1);
368//   return(P1);
369// }
370
371
372/*4 - follow the numbering of gring.cc
373* creates the S-polynomial of p1 and p2
374* do not destroy p1 and p2
375*/
376poly nc_rat_CreateSpoly(poly pp1, poly pp2, int ishift, const ring r)
377{
378
379  poly p1 = p_Copy(pp1,r);
380  poly p2 = p_Copy(pp2,r);
381
382  const long lCompP1 = p_GetComp(p1,r);
383  const long lCompP2 = p_GetComp(p2,r);
384
385  if ((lCompP1!=lCompP2) && (lCompP1!=0) && (lCompP2!=0))
386  {
387#ifdef PDEBUG
388    Werror("nc_rat_CreateSpoly: different non-zero components!");
389#endif
390    return(NULL);
391  }
392
393/* note: prod. crit does not apply! */
394  poly pL=pOne();
395  poly m1=pOne();
396  poly m2=pOne();
397  int is = ishift; /* TODO */
398  pLcmRat(p1,p2,pL,is);
399  p_Setm(pL,r);
400#ifdef PDEBUG
401  p_Test(pL,r);
402#endif
403  poly pr1 = p_GetExp_k_n(p1,1,ishift,r); /* rat D-exp of p1 */
404  poly pr2 = p_GetExp_k_n(p2,1,ishift,r); /* rat D-exp of p2 */
405  p_ExpVectorDiff(m1,pL,pr1,r); /* purely in D part by construction */
406  p_ExpVectorDiff(m2,pL,pr2,r); /* purely in D part by construction */
407  p_Delete(&pr1,r);
408  p_Delete(&pr2,r);
409  p_Delete(&pL,r);
410#ifdef PDEBUG
411  p_Test(m1,r);
412  PrintS("d^{gamma-alpha} = "); p_wrp(m1,r); PrintLn();
413  p_Test(m2,r);
414  PrintS("d^{gamma-beta} = "); p_wrp(m2,r); PrintLn();
415#endif
416
417  poly HF = NULL;
418  HF = p_HeadRat(p1,is,r); // lm_D(f)
419  HF  = nc_mm_Mult_p(m1, HF, r); // // d^{gamma-alpha} lm_D(f)
420  poly C  = p_GetCoeffRat(HF,  is, r); // c = lc_D(h_f) in the paper
421
422  poly HG = NULL;
423  HG = p_HeadRat(p2,is,r); // lm_D(g)
424  HG  = nc_mm_Mult_p(m2, HG, r); // // d^{gamma-beta} lm_D(g)
425  poly K  = p_GetCoeffRat(HG,  is, r); // k = lc_D(h_g) in the paper
426
427#ifdef PDEBUG
428  PrintS("f: "); p_wrp(p1,r); PrintS("\n");
429  PrintS("c: "); p_wrp(C,r); PrintS("\n");
430  PrintS("g: "); p_wrp(p2,r); PrintS("\n");
431  PrintS("k: "); p_wrp(K,r); PrintS("\n");
432#endif
433 
434  ideal ncsyz = ncGCD(C,K,r);
435  poly KK = ncsyz->m[0]; ncsyz->m[0]=NULL; //p_Copy(ncsyz->m[0],r); // k'
436  poly CC = ncsyz->m[1]; ncsyz->m[1]= NULL; //p_Copy(ncsyz->m[1],r); // c'
437  id_Delete(&ncsyz,r);
438
439  p_LmDeleteAndNextRat(&p1, is+1, r); // t_f
440  p_LmDeleteAndNextRat(&HF, is+1, r); // r_f = h_f - lt_D(h_f)
441
442  p_LmDeleteAndNextRat(&p2, is+1, r); // t_g
443  p_LmDeleteAndNextRat(&HG, is+1, r); // r_g = h_g - lt_D(h_g)
444
445
446#ifdef PDEBUG
447  PrintS(" t_f: "); p_wrp(p1,r); PrintS("\n"); 
448  PrintS(" t_g: "); p_wrp(p2,r); PrintS("\n"); 
449  PrintS(" r_f: "); p_wrp(HF,r); PrintS("\n"); 
450  PrintS(" r_g: "); p_wrp(HG,r); PrintS("\n"); 
451  PrintS(" c': "); p_wrp(CC,r); PrintS("\n"); 
452  PrintS(" k': "); p_wrp(KK,r); PrintS("\n"); 
453
454#endif
455
456  // k'(r_f + d^{gamma-alpha} t_f)
457
458  p1 = p_Mult_q(m1, p1, r); // p1 = d^{gamma-alpha} t_f
459  p1 = p_Add_q(p1,HF,r); // p1 = r_f + d^{gamma-alpha} t_f
460  p1 = p_Mult_q(KK,p1,r); // p1 = k'(r_f + d^{gamma-alpha} t_f)
461
462  // c'(r_f + d^{gamma-beta} t_g)
463
464  p2 = p_Mult_q(m2, p2, r); // p2 = d^{gamma-beta} t_g
465  p2 = p_Add_q(p2,HG,r); // p2 = r_g + d^{gamma-beta} t_g
466  p2 = p_Mult_q(CC,p2,r); // p2 = c'(r_g + d^{gamma-beta} t_g)
467
468#ifdef PDEBUG
469  p_Test(p1,r);
470  p_Test(p2,r);
471  PrintS(" k'(r_f + d^{gamma-alpha} t_f): "); p_wrp(p1,r);
472  PrintS(" c'(r_g + d^{gamma-beta} t_g): "); p_wrp(p2,r);
473#endif
474
475  poly out = p_Add_q(p1,p2,r); // delete p1, p2; // the sum
476
477#ifdef PDEBUG
478  p_Test(out,r);
479#endif
480
481  if ( out!=NULL ) pContent(out);
482  return(out);
483}
484
485
486/*2
487* reduction of p2 with p1
488* do not destroy p1, but p2
489* p1 divides p2 -> for use in NF algorithm
490* works in an integer fashion
491*/
492
493poly nc_rat_ReduceSpolyNew(const poly p1, poly p2, int ishift, const ring r)
494{
495  const long lCompP1 = p_GetComp(p1,r);
496  const long lCompP2 = p_GetComp(p2,r);
497
498  if ((lCompP1!=lCompP2) && (lCompP1!=0) && (lCompP2!=0))
499  {
500#ifdef PDEBUG
501    Werror("nc_rat_ReduceSpolyNew: different non-zero components!");
502#endif
503    return(NULL);
504  }
505
506  int is = ishift; /* TODO */
507
508  poly m = pOne();
509  p_ExpVectorDiffRat(m, p2, p1, ishift, r); // includes X and D parts
510  //p_Setm(m,r);
511  //  m = p_GetExp_k_n(m,1,ishift,r); /* rat D-exp of m */
512#ifdef PDEBUG
513  p_Test(m,r);
514  PrintS("d^alpha = "); p_wrp(m,r); PrintLn();
515#endif
516
517  /* pSetComp(m,r)=0? */
518  poly HH = NULL;
519  poly H  = NULL;
520  HH = p_HeadRat(p1,is,r); //p_Copy(p_HeadRat(p1,is,r),r); // lm_D(g)
521//  H  = r->nc->p_Procs.mm_Mult_p(m, p_Copy(HH, r), r); // d^aplha lm_D(g)
522  H  = nc_mm_Mult_p(m, HH, r); // d^aplha lm_D(g) == h_g in the paper
523
524  poly K  = p_GetCoeffRat(H,  is, r); //p_Copy( p_GetCoeffRat(H,  is, r), r); // k in the paper
525  poly P  = p_GetCoeffRat(p2, is, r); //p_Copy( p_GetCoeffRat(p2, is, r), r); // lc_D(p_2) == lc_D(f)
526
527#ifdef PDEBUG
528  PrintS("k: "); p_wrp(K,r); PrintS("\n");
529  PrintS("p: "); p_wrp(P,r); PrintS("\n");
530  PrintS("f: "); p_wrp(p2,r); PrintS("\n");
531  PrintS("g: "); p_wrp(p1,r); PrintS("\n");
532#endif
533  // alt:
534  poly out = p_Copy(p1,r);
535  p_LmDeleteAndNextRat(&out, is+1, r); // out == t_g
536
537  ideal ncsyz = ncGCD(P,K,r);
538  poly KK = ncsyz->m[0]; ncsyz->m[0]=NULL; //p_Copy(ncsyz->m[0],r); // k'
539  poly PP = ncsyz->m[1]; ncsyz->m[1]= NULL; //p_Copy(ncsyz->m[1],r); // p'
540
541#ifdef PDEBUG
542  PrintS("t_g: "); p_wrp(out,r);
543  PrintS("k': "); p_wrp(KK,r); PrintS("\n"); 
544  PrintS("p': "); p_wrp(PP,r); PrintS("\n"); 
545#endif
546  id_Delete(&ncsyz,r);
547  p_LmDeleteAndNextRat(&p2, is+1, r); // t_f
548  p_LmDeleteAndNextRat(&H, is+1, r); // r_g = h_g - lt_D(h_g)
549
550#ifdef PDEBUG
551  PrintS(" t_f: "); p_wrp(p2,r);
552  PrintS(" r_g: "); p_wrp(H,r);
553#endif
554
555  p2 = p_Mult_q(KK, p2, r); // p2 = k' t_f
556
557#ifdef PDEBUG
558  p_Test(p2,r);
559  PrintS(" k' t_f: "); p_wrp(p2,r);
560#endif
561
562//  out = r->nc->p_Procs.mm_Mult_p(m, out, r); // d^aplha t_g
563  out = nc_mm_Mult_p(m, out, r); // d^aplha t_g 
564  p_Delete(&m,r);
565
566#ifdef PDEBUG
567  PrintS(" d^a t_g: "); p_wrp(out,r);
568  PrintS(" end reduction\n");
569#endif
570
571  out = p_Add_q(H, out, r); // r_g + d^a t_g
572
573#ifdef PDEBUG
574  p_Test(out,r);
575#endif
576  out = p_Mult_q(PP, out, r); // p' (r_g + d^a t_g)
577  out = p_Add_q(p2,out,r); // delete out, p2; // the sum
578
579#ifdef PDEBUG
580  p_Test(out,r);
581#endif
582
583  if ( out!=NULL ) pContent(out);
584  return(out);
585}
586
587// return: FALSE, if there exists i in ishift..r->N,
588//                 such that a->exp[i] > b->exp[i]
589//         TRUE, otherwise
590
591BOOLEAN p_DivisibleByRat(poly a, poly b, int ishift, const ring r)
592{
593#ifdef PDEBUG
594  PrintS("invoke p_DivByRat with a = ");
595  p_wrp(p_Head(a,r),r);
596  PrintS(" and b= ");
597  p_wrp(p_Head(b,r),r); 
598  PrintLn();
599#endif
600  int i;
601  for(i=r->N; i>ishift; i--)
602  {
603#ifdef PDEBUG
604    Print("i=%d,",i);
605#endif
606    if (p_GetExp(a,i,r) > p_GetExp(b,i,r)) return FALSE;
607  }
608  return ((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(a,r)==0));
609}
610/*2
611*reduces h with elements from reducer choosing the best possible
612* element in t with respect to the given red_length
613* arrays reducer and red_length are [0..(rl-1)]
614*/
615int redRat (poly* h, poly *reducer, int *red_length, int rl, int ishift, ring r)
616{
617  if ((*h)==NULL) return 0;
618
619  int j,i,l;
620
621  loop
622  {
623    j=rl;l=MAX_INT_VAL;
624    for(i=rl-1;i>=0;i--)
625    {
626      //      Print("test %d, l=%d (curr=%d, l=%d\n",i,red_length[i],j,l);
627      if ((l>red_length[i]) && (p_DivisibleByRat(reducer[i],*h,ishift,r)))
628      {
629        j=i; l=red_length[i];
630        //        PrintS(" yes\n");
631      }
632      //      else PrintS(" no\n");
633    }
634    if (j >=rl)
635    {
636      return 1; // not reducible
637    }
638
639    if (TEST_OPT_DEBUG)
640    {
641      PrintS("reduce ");
642      p_wrp(*h,r);
643      PrintS(" with ");
644      p_wrp(reducer[j],r);
645    }
646    poly hh=nc_rat_ReduceSpolyNew(reducer[j], *h, ishift, r);
647    //    p_Delete(h,r);
648    *h=hh;
649    if (TEST_OPT_DEBUG)
650    {
651      PrintS(" to ");
652      p_wrp(*h,r);
653      PrintLn();
654    }
655    if ((*h)==NULL)
656    {
657      return 0;
658    }
659  }
660}
661#endif
Note: See TracBrowser for help on using the repository browser.