p_polys.cc
Go to the documentation of this file.
1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /***************************************************************
5  * File: p_polys.cc
6  * Purpose: implementation of ring independent poly procedures?
7  * Author: obachman (Olaf Bachmann)
8  * Created: 8/00
9  *******************************************************************/
10 
11 #include <ctype.h>
12 
13 #include <omalloc/omalloc.h>
14 
15 #include <misc/auxiliary.h>
16 
17 #include <misc/options.h>
18 #include <misc/intvec.h>
19 
20 
21 #include <coeffs/longrat.h> // snumber is needed...
22 #include <coeffs/numbers.h> // ndCopyMap
23 
24 #include <polys/PolyEnumerator.h>
25 
26 #define TRANSEXT_PRIVATES
27 
30 
31 #include <polys/weight.h>
32 #include <polys/simpleideals.h>
33 
34 #include "ring.h"
35 #include "p_polys.h"
36 
40 
41 
42 // #include <???/ideals.h>
43 // #include <???/int64vec.h>
44 
45 #ifndef SING_NDEBUG
46 // #include <???/febase.h>
47 #endif
48 
49 #ifdef HAVE_PLURAL
50 #include "nc/nc.h"
51 #include "nc/sca.h"
52 #endif
53 
54 #include "coeffrings.h"
55 #include "clapsing.h"
56 
57 #define ADIDEBUG 0
58 
59 /*
60  * lift ideal with coeffs over Z (mod N) to Q via Farey
61  */
62 poly p_Farey(poly p, number N, const ring r)
63 {
64  poly h=p_Copy(p,r);
65  poly hh=h;
66  while(h!=NULL)
67  {
68  number c=pGetCoeff(h);
69  pSetCoeff0(h,n_Farey(c,N,r->cf));
70  n_Delete(&c,r->cf);
71  pIter(h);
72  }
73  while((hh!=NULL)&&(n_IsZero(pGetCoeff(hh),r->cf)))
74  {
75  p_LmDelete(&hh,r);
76  }
77  h=hh;
78  while((h!=NULL) && (pNext(h)!=NULL))
79  {
80  if(n_IsZero(pGetCoeff(pNext(h)),r->cf))
81  {
82  p_LmDelete(&pNext(h),r);
83  }
84  else pIter(h);
85  }
86  return hh;
87 }
88 /*2
89 * xx,q: arrays of length 0..rl-1
90 * xx[i]: SB mod q[i]
91 * assume: char=0
92 * assume: q[i]!=0
93 * destroys xx
94 */
95 poly p_ChineseRemainder(poly *xx, number *x,number *q, int rl, CFArray &inv_cache,const ring R)
96 {
97  poly r,h,hh;
98  int j;
99  poly res_p=NULL;
100  loop
101  {
102  /* search the lead term */
103  r=NULL;
104  for(j=rl-1;j>=0;j--)
105  {
106  h=xx[j];
107  if ((h!=NULL)
108  &&((r==NULL)||(p_LmCmp(r,h,R)==-1)))
109  r=h;
110  }
111  /* nothing found -> return */
112  if (r==NULL) break;
113  /* create the monomial in h */
114  h=p_Head(r,R);
115  /* collect the coeffs in x[..]*/
116  for(j=rl-1;j>=0;j--)
117  {
118  hh=xx[j];
119  if ((hh!=NULL) && (p_LmCmp(h,hh,R)==0))
120  {
121  x[j]=pGetCoeff(hh);
122  hh=p_LmFreeAndNext(hh,R);
123  xx[j]=hh;
124  }
125  else
126  x[j]=n_Init(0, R);
127  }
128  number n=n_ChineseRemainderSym(x,q,rl,TRUE,inv_cache,R->cf);
129  for(j=rl-1;j>=0;j--)
130  {
131  x[j]=NULL; // n_Init(0...) takes no memory
132  }
133  if (n_IsZero(n,R)) p_Delete(&h,R);
134  else
135  {
136  //Print("new mon:");pWrite(h);
137  p_SetCoeff(h,n,R);
138  pNext(h)=res_p;
139  res_p=h; // building res_p in reverse order!
140  }
141  }
142  res_p=pReverse(res_p);
143  p_Test(res_p, R);
144  return res_p;
145 }
146 /***************************************************************
147  *
148  * Completing what needs to be set for the monomial
149  *
150  ***************************************************************/
151 // this is special for the syz stuff
152 static int* _components = NULL;
153 static long* _componentsShifted = NULL;
154 static int _componentsExternal = 0;
155 
157 
158 #ifndef SING_NDEBUG
159 # define MYTEST 0
160 #else /* ifndef SING_NDEBUG */
161 # define MYTEST 0
162 #endif /* ifndef SING_NDEBUG */
163 
164 void p_Setm_General(poly p, const ring r)
165 {
166  p_LmCheckPolyRing(p, r);
167  int pos=0;
168  if (r->typ!=NULL)
169  {
170  loop
171  {
172  unsigned long ord=0;
173  sro_ord* o=&(r->typ[pos]);
174  switch(o->ord_typ)
175  {
176  case ro_dp:
177  {
178  int a,e;
179  a=o->data.dp.start;
180  e=o->data.dp.end;
181  for(int i=a;i<=e;i++) ord+=p_GetExp(p,i,r);
182  p->exp[o->data.dp.place]=ord;
183  break;
184  }
185  case ro_wp_neg:
187  // no break;
188  case ro_wp:
189  {
190  int a,e;
191  a=o->data.wp.start;
192  e=o->data.wp.end;
193  int *w=o->data.wp.weights;
194 #if 1
195  for(int i=a;i<=e;i++) ord+=((unsigned long)p_GetExp(p,i,r))*((unsigned long)w[i-a]);
196 #else
197  long ai;
198  int ei,wi;
199  for(int i=a;i<=e;i++)
200  {
201  ei=p_GetExp(p,i,r);
202  wi=w[i-a];
203  ai=ei*wi;
204  if (ai/ei!=wi) pSetm_error=TRUE;
205  ord+=ai;
206  if (ord<ai) pSetm_error=TRUE;
207  }
208 #endif
209  p->exp[o->data.wp.place]=ord;
210  break;
211  }
212  case ro_am:
213  {
214  ord = POLY_NEGWEIGHT_OFFSET;
215  const short a=o->data.am.start;
216  const short e=o->data.am.end;
217  const int * w=o->data.am.weights;
218 #if 1
219  for(short i=a; i<=e; i++, w++)
220  ord += ((*w) * p_GetExp(p,i,r));
221 #else
222  long ai;
223  int ei,wi;
224  for(short i=a;i<=e;i++)
225  {
226  ei=p_GetExp(p,i,r);
227  wi=w[i-a];
228  ai=ei*wi;
229  if (ai/ei!=wi) pSetm_error=TRUE;
230  ord += ai;
231  if (ord<ai) pSetm_error=TRUE;
232  }
233 #endif
234  const int c = p_GetComp(p,r);
235 
236  const short len_gen= o->data.am.len_gen;
237 
238  if ((c > 0) && (c <= len_gen))
239  {
240  assume( w == o->data.am.weights_m );
241  assume( w[0] == len_gen );
242  ord += w[c];
243  }
244 
245  p->exp[o->data.am.place] = ord;
246  break;
247  }
248  case ro_wp64:
249  {
250  int64 ord=0;
251  int a,e;
252  a=o->data.wp64.start;
253  e=o->data.wp64.end;
254  int64 *w=o->data.wp64.weights64;
255  int64 ei,wi,ai;
256  for(int i=a;i<=e;i++)
257  {
258  //Print("exp %d w %d \n",p_GetExp(p,i,r),(int)w[i-a]);
259  //ord+=((int64)p_GetExp(p,i,r))*w[i-a];
260  ei=(int64)p_GetExp(p,i,r);
261  wi=w[i-a];
262  ai=ei*wi;
263  if(ei!=0 && ai/ei!=wi)
264  {
266  #if SIZEOF_LONG == 4
267  Print("ai %lld, wi %lld\n",ai,wi);
268  #else
269  Print("ai %ld, wi %ld\n",ai,wi);
270  #endif
271  }
272  ord+=ai;
273  if (ord<ai)
274  {
276  #if SIZEOF_LONG == 4
277  Print("ai %lld, ord %lld\n",ai,ord);
278  #else
279  Print("ai %ld, ord %ld\n",ai,ord);
280  #endif
281  }
282  }
283  int64 mask=(int64)0x7fffffff;
284  long a_0=(long)(ord&mask); //2^31
285  long a_1=(long)(ord >>31 ); /*(ord/(mask+1));*/
286 
287  //Print("mask: %x, ord: %d, a_0: %d, a_1: %d\n"
288  //,(int)mask,(int)ord,(int)a_0,(int)a_1);
289  //Print("mask: %d",mask);
290 
291  p->exp[o->data.wp64.place]=a_1;
292  p->exp[o->data.wp64.place+1]=a_0;
293 // if(p_Setm_error) PrintS("***************************\n"
294 // "***************************\n"
295 // "**WARNING: overflow error**\n"
296 // "***************************\n"
297 // "***************************\n");
298  break;
299  }
300  case ro_cp:
301  {
302  int a,e;
303  a=o->data.cp.start;
304  e=o->data.cp.end;
305  int pl=o->data.cp.place;
306  for(int i=a;i<=e;i++) { p->exp[pl]=p_GetExp(p,i,r); pl++; }
307  break;
308  }
309  case ro_syzcomp:
310  {
311  long c=p_GetComp(p,r);
312  long sc = c;
313  int* Components = (_componentsExternal ? _components :
314  o->data.syzcomp.Components);
315  long* ShiftedComponents = (_componentsExternal ? _componentsShifted:
316  o->data.syzcomp.ShiftedComponents);
317  if (ShiftedComponents != NULL)
318  {
319  assume(Components != NULL);
320  assume(c == 0 || Components[c] != 0);
321  sc = ShiftedComponents[Components[c]];
322  assume(c == 0 || sc != 0);
323  }
324  p->exp[o->data.syzcomp.place]=sc;
325  break;
326  }
327  case ro_syz:
328  {
329  const unsigned long c = p_GetComp(p, r);
330  const short place = o->data.syz.place;
331  const int limit = o->data.syz.limit;
332 
333  if (c > (unsigned long)limit)
334  p->exp[place] = o->data.syz.curr_index;
335  else if (c > 0)
336  {
337  assume( (1 <= c) && (c <= (unsigned long)limit) );
338  p->exp[place]= o->data.syz.syz_index[c];
339  }
340  else
341  {
342  assume(c == 0);
343  p->exp[place]= 0;
344  }
345  break;
346  }
347  // Prefix for Induced Schreyer ordering
348  case ro_isTemp: // Do nothing?? (to be removed into suffix later on...?)
349  {
350  assume(p != NULL);
351 
352 #ifndef SING_NDEBUG
353 #if MYTEST
354  Print("p_Setm_General: ro_isTemp ord: pos: %d, p: ", pos); p_wrp(p, r);
355 #endif
356 #endif
357  int c = p_GetComp(p, r);
358 
359  assume( c >= 0 );
360 
361  // Let's simulate case ro_syz above....
362  // Should accumulate (by Suffix) and be a level indicator
363  const int* const pVarOffset = o->data.isTemp.pVarOffset;
364 
365  assume( pVarOffset != NULL );
366 
367  // TODO: Can this be done in the suffix???
368  for( int i = 1; i <= r->N; i++ ) // No v[0] here!!!
369  {
370  const int vo = pVarOffset[i];
371  if( vo != -1) // TODO: optimize: can be done once!
372  {
373  // Hans! Please don't break it again! p_SetExp(p, ..., r, vo) is correct:
374  p_SetExp(p, p_GetExp(p, i, r), r, vo); // copy put them verbatim
375  // Hans! Please don't break it again! p_GetExp(p, r, vo) is correct:
376  assume( p_GetExp(p, r, vo) == p_GetExp(p, i, r) ); // copy put them verbatim
377  }
378  }
379 #ifndef SING_NDEBUG
380  for( int i = 1; i <= r->N; i++ ) // No v[0] here!!!
381  {
382  const int vo = pVarOffset[i];
383  if( vo != -1) // TODO: optimize: can be done once!
384  {
385  // Hans! Please don't break it again! p_GetExp(p, r, vo) is correct:
386  assume( p_GetExp(p, r, vo) == p_GetExp(p, i, r) ); // copy put them verbatim
387  }
388  }
389 #if MYTEST
390 // if( p->exp[o->data.isTemp.start] > 0 )
391  PrintS("after Values: "); p_wrp(p, r);
392 #endif
393 #endif
394  break;
395  }
396 
397  // Suffix for Induced Schreyer ordering
398  case ro_is:
399  {
400 #ifndef SING_NDEBUG
401 #if MYTEST
402  Print("p_Setm_General: ro_is ord: pos: %d, p: ", pos); p_wrp(p, r);
403 #endif
404 #endif
405 
406  assume(p != NULL);
407 
408  int c = p_GetComp(p, r);
409 
410  assume( c >= 0 );
411  const ideal F = o->data.is.F;
412  const int limit = o->data.is.limit;
413  assume( limit >= 0 );
414  const int start = o->data.is.start;
415 
416  if( F != NULL && c > limit )
417  {
418 #ifndef SING_NDEBUG
419 #if MYTEST
420  Print("p_Setm_General: ro_is : in rSetm: pos: %d, c: %d > limit: %d\n", c, pos, limit);
421  PrintS("preComputed Values: ");
422  p_wrp(p, r);
423 #endif
424 #endif
425 // if( c > limit ) // BUG???
426  p->exp[start] = 1;
427 // else
428 // p->exp[start] = 0;
429 
430 
431  c -= limit;
432  assume( c > 0 );
433  c--;
434 
435  if( c >= IDELEMS(F) )
436  break;
437 
438  assume( c < IDELEMS(F) ); // What about others???
439 
440  const poly pp = F->m[c]; // get reference monomial!!!
441 
442  if(pp == NULL)
443  break;
444 
445  assume(pp != NULL);
446 
447 #ifndef SING_NDEBUG
448 #if MYTEST
449  Print("Respective F[c - %d: %d] pp: ", limit, c);
450  p_wrp(pp, r);
451 #endif
452 #endif
453 
454  const int end = o->data.is.end;
455  assume(start <= end);
456 
457 
458 // const int st = o->data.isTemp.start;
459 
460 #ifndef SING_NDEBUG
461  Print("p_Setm_General: is(-Temp-) :: c: %d, limit: %d, [st:%d] ===>>> %ld\n", c, limit, start, p->exp[start]);
462 #endif
463 
464  // p_ExpVectorAdd(p, pp, r);
465 
466  for( int i = start; i <= end; i++) // v[0] may be here...
467  p->exp[i] += pp->exp[i]; // !!!!!!!! ADD corresponding LT(F)
468 
469  // p_MemAddAdjust(p, ri);
470  if (r->NegWeightL_Offset != NULL)
471  {
472  for (int i=r->NegWeightL_Size-1; i>=0; i--)
473  {
474  const int _i = r->NegWeightL_Offset[i];
475  if( start <= _i && _i <= end )
476  p->exp[_i] -= POLY_NEGWEIGHT_OFFSET;
477  }
478  }
479 
480 
481 #ifndef SING_NDEBUG
482  const int* const pVarOffset = o->data.is.pVarOffset;
483 
484  assume( pVarOffset != NULL );
485 
486  for( int i = 1; i <= r->N; i++ ) // No v[0] here!!!
487  {
488  const int vo = pVarOffset[i];
489  if( vo != -1) // TODO: optimize: can be done once!
490  // Hans! Please don't break it again! p_GetExp(p/pp, r, vo) is correct:
491  assume( p_GetExp(p, r, vo) == (p_GetExp(p, i, r) + p_GetExp(pp, r, vo)) );
492  }
493  // TODO: how to check this for computed values???
494 #if MYTEST
495  PrintS("Computed Values: "); p_wrp(p, r);
496 #endif
497 #endif
498  } else
499  {
500  p->exp[start] = 0; //!!!!????? where?????
501 
502  const int* const pVarOffset = o->data.is.pVarOffset;
503 
504  // What about v[0] - component: it will be added later by
505  // suffix!!!
506  // TODO: Test it!
507  const int vo = pVarOffset[0];
508  if( vo != -1 )
509  p->exp[vo] = c; // initial component v[0]!
510 
511 #ifndef SING_NDEBUG
512 #if MYTEST
513  Print("ELSE p_Setm_General: ro_is :: c: %d <= limit: %d, vo: %d, exp: %d\n", c, limit, vo, p->exp[vo]);
514  p_wrp(p, r);
515 #endif
516 #endif
517  }
518 
519  break;
520  }
521  default:
522  dReportError("wrong ord in rSetm:%d\n",o->ord_typ);
523  return;
524  }
525  pos++;
526  if (pos == r->OrdSize) return;
527  }
528  }
529 }
530 
531 void p_Setm_Syz(poly p, ring r, int* Components, long* ShiftedComponents)
532 {
533  _components = Components;
534  _componentsShifted = ShiftedComponents;
536  p_Setm_General(p, r);
538 }
539 
540 // dummy for lp, ls, etc
541 void p_Setm_Dummy(poly p, const ring r)
542 {
543  p_LmCheckPolyRing(p, r);
544 }
545 
546 // for dp, Dp, ds, etc
547 void p_Setm_TotalDegree(poly p, const ring r)
548 {
549  p_LmCheckPolyRing(p, r);
550  p->exp[r->pOrdIndex] = p_Totaldegree(p, r);
551 }
552 
553 // for wp, Wp, ws, etc
554 void p_Setm_WFirstTotalDegree(poly p, const ring r)
555 {
556  p_LmCheckPolyRing(p, r);
557  p->exp[r->pOrdIndex] = p_WFirstTotalDegree(p, r);
558 }
559 
561 {
562  // covers lp, rp, ls,
563  if (r->typ == NULL) return p_Setm_Dummy;
564 
565  if (r->OrdSize == 1)
566  {
567  if (r->typ[0].ord_typ == ro_dp &&
568  r->typ[0].data.dp.start == 1 &&
569  r->typ[0].data.dp.end == r->N &&
570  r->typ[0].data.dp.place == r->pOrdIndex)
571  return p_Setm_TotalDegree;
572  if (r->typ[0].ord_typ == ro_wp &&
573  r->typ[0].data.wp.start == 1 &&
574  r->typ[0].data.wp.end == r->N &&
575  r->typ[0].data.wp.place == r->pOrdIndex &&
576  r->typ[0].data.wp.weights == r->firstwv)
578  }
579  return p_Setm_General;
580 }
581 
582 
583 /* -------------------------------------------------------------------*/
584 /* several possibilities for pFDeg: the degree of the head term */
585 
586 /* comptible with ordering */
587 long p_Deg(poly a, const ring r)
588 {
589  p_LmCheckPolyRing(a, r);
590 // assume(p_GetOrder(a, r) == p_WTotaldegree(a, r)); // WRONG assume!
591  return p_GetOrder(a, r);
592 }
593 
594 // p_WTotalDegree for weighted orderings
595 // whose first block covers all variables
596 long p_WFirstTotalDegree(poly p, const ring r)
597 {
598  int i;
599  long sum = 0;
600 
601  for (i=1; i<= r->firstBlockEnds; i++)
602  {
603  sum += p_GetExp(p, i, r)*r->firstwv[i-1];
604  }
605  return sum;
606 }
607 
608 /*2
609 * compute the degree of the leading monomial of p
610 * with respect to weigths from the ordering
611 * the ordering is not compatible with degree so do not use p->Order
612 */
613 long p_WTotaldegree(poly p, const ring r)
614 {
615  p_LmCheckPolyRing(p, r);
616  int i, k;
617  long j =0;
618 
619  // iterate through each block:
620  for (i=0;r->order[i]!=0;i++)
621  {
622  int b0=r->block0[i];
623  int b1=r->block1[i];
624  switch(r->order[i])
625  {
626  case ringorder_M:
627  for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
628  { // in jedem block:
629  j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/]*r->OrdSgn;
630  }
631  break;
632  case ringorder_wp:
633  case ringorder_ws:
634  case ringorder_Wp:
635  case ringorder_Ws:
636  for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
637  { // in jedem block:
638  j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/];
639  }
640  break;
641  case ringorder_lp:
642  case ringorder_ls:
643  case ringorder_rs:
644  case ringorder_dp:
645  case ringorder_ds:
646  case ringorder_Dp:
647  case ringorder_Ds:
648  case ringorder_rp:
649  for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
650  {
651  j+= p_GetExp(p,k,r);
652  }
653  break;
654  case ringorder_a64:
655  {
656  int64* w=(int64*)r->wvhdl[i];
657  for (k=0;k<=(b1 /*r->block1[i]*/ - b0 /*r->block0[i]*/);k++)
658  {
659  //there should be added a line which checks if w[k]>2^31
660  j+= p_GetExp(p,k+1, r)*(long)w[k];
661  }
662  //break;
663  return j;
664  }
665  case ringorder_c:
666  case ringorder_C:
667  case ringorder_S:
668  case ringorder_s:
669  case ringorder_aa:
670  case ringorder_IS:
671  break;
672  case ringorder_a:
673  for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
674  { // only one line
675  j+= p_GetExp(p,k, r)*r->wvhdl[i][ k- b0 /*r->block0[i]*/];
676  }
677  //break;
678  return j;
679 
680 #ifndef SING_NDEBUG
681  default:
682  Print("missing order %d in p_WTotaldegree\n",r->order[i]);
683  break;
684 #endif
685  }
686  }
687  return j;
688 }
689 
690 long p_DegW(poly p, const short *w, const ring R)
691 {
692  p_Test(p, R);
693  assume( w != NULL );
694  long r=-LONG_MAX;
695 
696  while (p!=NULL)
697  {
698  long t=totaldegreeWecart_IV(p,R,w);
699  if (t>r) r=t;
700  pIter(p);
701  }
702  return r;
703 }
704 
705 int p_Weight(int i, const ring r)
706 {
707  if ((r->firstwv==NULL) || (i>r->firstBlockEnds))
708  {
709  return 1;
710  }
711  return r->firstwv[i-1];
712 }
713 
714 long p_WDegree(poly p, const ring r)
715 {
716  if (r->firstwv==NULL) return p_Totaldegree(p, r);
717  p_LmCheckPolyRing(p, r);
718  int i;
719  long j =0;
720 
721  for(i=1;i<=r->firstBlockEnds;i++)
722  j+=p_GetExp(p, i, r)*r->firstwv[i-1];
723 
724  for (;i<=rVar(r);i++)
725  j+=p_GetExp(p,i, r)*p_Weight(i, r);
726 
727  return j;
728 }
729 
730 
731 /* ---------------------------------------------------------------------*/
732 /* several possibilities for pLDeg: the maximal degree of a monomial in p*/
733 /* compute in l also the pLength of p */
734 
735 /*2
736 * compute the length of a polynomial (in l)
737 * and the degree of the monomial with maximal degree: the last one
738 */
739 long pLDeg0(poly p,int *l, const ring r)
740 {
741  p_CheckPolyRing(p, r);
742  long k= p_GetComp(p, r);
743  int ll=1;
744 
745  if (k > 0)
746  {
747  while ((pNext(p)!=NULL) && (p_GetComp(pNext(p), r)==k))
748  {
749  pIter(p);
750  ll++;
751  }
752  }
753  else
754  {
755  while (pNext(p)!=NULL)
756  {
757  pIter(p);
758  ll++;
759  }
760  }
761  *l=ll;
762  return r->pFDeg(p, r);
763 }
764 
765 /*2
766 * compute the length of a polynomial (in l)
767 * and the degree of the monomial with maximal degree: the last one
768 * but search in all components before syzcomp
769 */
770 long pLDeg0c(poly p,int *l, const ring r)
771 {
772  assume(p!=NULL);
773  p_Test(p,r);
774  p_CheckPolyRing(p, r);
775  long o;
776  int ll=1;
777 
778  if (! rIsSyzIndexRing(r))
779  {
780  while (pNext(p) != NULL)
781  {
782  pIter(p);
783  ll++;
784  }
785  o = r->pFDeg(p, r);
786  }
787  else
788  {
789  int curr_limit = rGetCurrSyzLimit(r);
790  poly pp = p;
791  while ((p=pNext(p))!=NULL)
792  {
793  if (p_GetComp(p, r)<=curr_limit/*syzComp*/)
794  ll++;
795  else break;
796  pp = p;
797  }
798  p_Test(pp,r);
799  o = r->pFDeg(pp, r);
800  }
801  *l=ll;
802  return o;
803 }
804 
805 /*2
806 * compute the length of a polynomial (in l)
807 * and the degree of the monomial with maximal degree: the first one
808 * this works for the polynomial case with degree orderings
809 * (both c,dp and dp,c)
810 */
811 long pLDegb(poly p,int *l, const ring r)
812 {
813  p_CheckPolyRing(p, r);
814  long k= p_GetComp(p, r);
815  long o = r->pFDeg(p, r);
816  int ll=1;
817 
818  if (k != 0)
819  {
820  while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
821  {
822  ll++;
823  }
824  }
825  else
826  {
827  while ((p=pNext(p)) !=NULL)
828  {
829  ll++;
830  }
831  }
832  *l=ll;
833  return o;
834 }
835 
836 /*2
837 * compute the length of a polynomial (in l)
838 * and the degree of the monomial with maximal degree:
839 * this is NOT the last one, we have to look for it
840 */
841 long pLDeg1(poly p,int *l, const ring r)
842 {
843  p_CheckPolyRing(p, r);
844  long k= p_GetComp(p, r);
845  int ll=1;
846  long t,max;
847 
848  max=r->pFDeg(p, r);
849  if (k > 0)
850  {
851  while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
852  {
853  t=r->pFDeg(p, r);
854  if (t>max) max=t;
855  ll++;
856  }
857  }
858  else
859  {
860  while ((p=pNext(p))!=NULL)
861  {
862  t=r->pFDeg(p, r);
863  if (t>max) max=t;
864  ll++;
865  }
866  }
867  *l=ll;
868  return max;
869 }
870 
871 /*2
872 * compute the length of a polynomial (in l)
873 * and the degree of the monomial with maximal degree:
874 * this is NOT the last one, we have to look for it
875 * in all components
876 */
877 long pLDeg1c(poly p,int *l, const ring r)
878 {
879  p_CheckPolyRing(p, r);
880  int ll=1;
881  long t,max;
882 
883  max=r->pFDeg(p, r);
884  if (rIsSyzIndexRing(r))
885  {
886  long limit = rGetCurrSyzLimit(r);
887  while ((p=pNext(p))!=NULL)
888  {
889  if (p_GetComp(p, r)<=limit)
890  {
891  if ((t=r->pFDeg(p, r))>max) max=t;
892  ll++;
893  }
894  else break;
895  }
896  }
897  else
898  {
899  while ((p=pNext(p))!=NULL)
900  {
901  if ((t=r->pFDeg(p, r))>max) max=t;
902  ll++;
903  }
904  }
905  *l=ll;
906  return max;
907 }
908 
909 // like pLDeg1, only pFDeg == pDeg
910 long pLDeg1_Deg(poly p,int *l, const ring r)
911 {
912  assume(r->pFDeg == p_Deg);
913  p_CheckPolyRing(p, r);
914  long k= p_GetComp(p, r);
915  int ll=1;
916  long t,max;
917 
918  max=p_GetOrder(p, r);
919  if (k > 0)
920  {
921  while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
922  {
923  t=p_GetOrder(p, r);
924  if (t>max) max=t;
925  ll++;
926  }
927  }
928  else
929  {
930  while ((p=pNext(p))!=NULL)
931  {
932  t=p_GetOrder(p, r);
933  if (t>max) max=t;
934  ll++;
935  }
936  }
937  *l=ll;
938  return max;
939 }
940 
941 long pLDeg1c_Deg(poly p,int *l, const ring r)
942 {
943  assume(r->pFDeg == p_Deg);
944  p_CheckPolyRing(p, r);
945  int ll=1;
946  long t,max;
947 
948  max=p_GetOrder(p, r);
949  if (rIsSyzIndexRing(r))
950  {
951  long limit = rGetCurrSyzLimit(r);
952  while ((p=pNext(p))!=NULL)
953  {
954  if (p_GetComp(p, r)<=limit)
955  {
956  if ((t=p_GetOrder(p, r))>max) max=t;
957  ll++;
958  }
959  else break;
960  }
961  }
962  else
963  {
964  while ((p=pNext(p))!=NULL)
965  {
966  if ((t=p_GetOrder(p, r))>max) max=t;
967  ll++;
968  }
969  }
970  *l=ll;
971  return max;
972 }
973 
974 // like pLDeg1, only pFDeg == pTotoalDegree
975 long pLDeg1_Totaldegree(poly p,int *l, const ring r)
976 {
977  p_CheckPolyRing(p, r);
978  long k= p_GetComp(p, r);
979  int ll=1;
980  long t,max;
981 
982  max=p_Totaldegree(p, r);
983  if (k > 0)
984  {
985  while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
986  {
987  t=p_Totaldegree(p, r);
988  if (t>max) max=t;
989  ll++;
990  }
991  }
992  else
993  {
994  while ((p=pNext(p))!=NULL)
995  {
996  t=p_Totaldegree(p, r);
997  if (t>max) max=t;
998  ll++;
999  }
1000  }
1001  *l=ll;
1002  return max;
1003 }
1004 
1005 long pLDeg1c_Totaldegree(poly p,int *l, const ring r)
1006 {
1007  p_CheckPolyRing(p, r);
1008  int ll=1;
1009  long t,max;
1010 
1011  max=p_Totaldegree(p, r);
1012  if (rIsSyzIndexRing(r))
1013  {
1014  long limit = rGetCurrSyzLimit(r);
1015  while ((p=pNext(p))!=NULL)
1016  {
1017  if (p_GetComp(p, r)<=limit)
1018  {
1019  if ((t=p_Totaldegree(p, r))>max) max=t;
1020  ll++;
1021  }
1022  else break;
1023  }
1024  }
1025  else
1026  {
1027  while ((p=pNext(p))!=NULL)
1028  {
1029  if ((t=p_Totaldegree(p, r))>max) max=t;
1030  ll++;
1031  }
1032  }
1033  *l=ll;
1034  return max;
1035 }
1036 
1037 // like pLDeg1, only pFDeg == p_WFirstTotalDegree
1038 long pLDeg1_WFirstTotalDegree(poly p,int *l, const ring r)
1039 {
1040  p_CheckPolyRing(p, r);
1041  long k= p_GetComp(p, r);
1042  int ll=1;
1043  long t,max;
1044 
1045  max=p_WFirstTotalDegree(p, r);
1046  if (k > 0)
1047  {
1048  while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
1049  {
1050  t=p_WFirstTotalDegree(p, r);
1051  if (t>max) max=t;
1052  ll++;
1053  }
1054  }
1055  else
1056  {
1057  while ((p=pNext(p))!=NULL)
1058  {
1059  t=p_WFirstTotalDegree(p, r);
1060  if (t>max) max=t;
1061  ll++;
1062  }
1063  }
1064  *l=ll;
1065  return max;
1066 }
1067 
1068 long pLDeg1c_WFirstTotalDegree(poly p,int *l, const ring r)
1069 {
1070  p_CheckPolyRing(p, r);
1071  int ll=1;
1072  long t,max;
1073 
1074  max=p_WFirstTotalDegree(p, r);
1075  if (rIsSyzIndexRing(r))
1076  {
1077  long limit = rGetCurrSyzLimit(r);
1078  while ((p=pNext(p))!=NULL)
1079  {
1080  if (p_GetComp(p, r)<=limit)
1081  {
1082  if ((t=p_Totaldegree(p, r))>max) max=t;
1083  ll++;
1084  }
1085  else break;
1086  }
1087  }
1088  else
1089  {
1090  while ((p=pNext(p))!=NULL)
1091  {
1092  if ((t=p_Totaldegree(p, r))>max) max=t;
1093  ll++;
1094  }
1095  }
1096  *l=ll;
1097  return max;
1098 }
1099 
1100 /***************************************************************
1101  *
1102  * Maximal Exponent business
1103  *
1104  ***************************************************************/
1105 
1106 static inline unsigned long
1107 p_GetMaxExpL2(unsigned long l1, unsigned long l2, const ring r,
1108  unsigned long number_of_exp)
1109 {
1110  const unsigned long bitmask = r->bitmask;
1111  unsigned long ml1 = l1 & bitmask;
1112  unsigned long ml2 = l2 & bitmask;
1113  unsigned long max = (ml1 > ml2 ? ml1 : ml2);
1114  unsigned long j = number_of_exp - 1;
1115 
1116  if (j > 0)
1117  {
1118  unsigned long mask = bitmask << r->BitsPerExp;
1119  while (1)
1120  {
1121  ml1 = l1 & mask;
1122  ml2 = l2 & mask;
1123  max |= ((ml1 > ml2 ? ml1 : ml2) & mask);
1124  j--;
1125  if (j == 0) break;
1126  mask = mask << r->BitsPerExp;
1127  }
1128  }
1129  return max;
1130 }
1131 
1132 static inline unsigned long
1133 p_GetMaxExpL2(unsigned long l1, unsigned long l2, const ring r)
1134 {
1135  return p_GetMaxExpL2(l1, l2, r, r->ExpPerLong);
1136 }
1137 
1138 poly p_GetMaxExpP(poly p, const ring r)
1139 {
1140  p_CheckPolyRing(p, r);
1141  if (p == NULL) return p_Init(r);
1142  poly max = p_LmInit(p, r);
1143  pIter(p);
1144  if (p == NULL) return max;
1145  int i, offset;
1146  unsigned long l_p, l_max;
1147  unsigned long divmask = r->divmask;
1148 
1149  do
1150  {
1151  offset = r->VarL_Offset[0];
1152  l_p = p->exp[offset];
1153  l_max = max->exp[offset];
1154  // do the divisibility trick to find out whether l has an exponent
1155  if (l_p > l_max ||
1156  (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
1157  max->exp[offset] = p_GetMaxExpL2(l_max, l_p, r);
1158 
1159  for (i=1; i<r->VarL_Size; i++)
1160  {
1161  offset = r->VarL_Offset[i];
1162  l_p = p->exp[offset];
1163  l_max = max->exp[offset];
1164  // do the divisibility trick to find out whether l has an exponent
1165  if (l_p > l_max ||
1166  (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
1167  max->exp[offset] = p_GetMaxExpL2(l_max, l_p, r);
1168  }
1169  pIter(p);
1170  }
1171  while (p != NULL);
1172  return max;
1173 }
1174 
1175 unsigned long p_GetMaxExpL(poly p, const ring r, unsigned long l_max)
1176 {
1177  unsigned long l_p, divmask = r->divmask;
1178  int i;
1179 
1180  while (p != NULL)
1181  {
1182  l_p = p->exp[r->VarL_Offset[0]];
1183  if (l_p > l_max ||
1184  (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
1185  l_max = p_GetMaxExpL2(l_max, l_p, r);
1186  for (i=1; i<r->VarL_Size; i++)
1187  {
1188  l_p = p->exp[r->VarL_Offset[i]];
1189  // do the divisibility trick to find out whether l has an exponent
1190  if (l_p > l_max ||
1191  (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
1192  l_max = p_GetMaxExpL2(l_max, l_p, r);
1193  }
1194  pIter(p);
1195  }
1196  return l_max;
1197 }
1198 
1199 
1200 
1201 
1202 /***************************************************************
1203  *
1204  * Misc things
1205  *
1206  ***************************************************************/
1207 // returns TRUE, if all monoms have the same component
1208 BOOLEAN p_OneComp(poly p, const ring r)
1209 {
1210  if(p!=NULL)
1211  {
1212  long i = p_GetComp(p, r);
1213  while (pNext(p)!=NULL)
1214  {
1215  pIter(p);
1216  if(i != p_GetComp(p, r)) return FALSE;
1217  }
1218  }
1219  return TRUE;
1220 }
1221 
1222 /*2
1223 *test if a monomial /head term is a pure power,
1224 * i.e. depends on only one variable
1225 */
1226 int p_IsPurePower(const poly p, const ring r)
1227 {
1228  int i,k=0;
1229 
1230  for (i=r->N;i;i--)
1231  {
1232  if (p_GetExp(p,i, r)!=0)
1233  {
1234  if(k!=0) return 0;
1235  k=i;
1236  }
1237  }
1238  return k;
1239 }
1240 
1241 /*2
1242 *test if a polynomial is univariate
1243 * return -1 for constant,
1244 * 0 for not univariate,s
1245 * i if dep. on var(i)
1246 */
1247 int p_IsUnivariate(poly p, const ring r)
1248 {
1249  int i,k=-1;
1250 
1251  while (p!=NULL)
1252  {
1253  for (i=r->N;i;i--)
1254  {
1255  if (p_GetExp(p,i, r)!=0)
1256  {
1257  if((k!=-1)&&(k!=i)) return 0;
1258  k=i;
1259  }
1260  }
1261  pIter(p);
1262  }
1263  return k;
1264 }
1265 
1266 // set entry e[i] to 1 if var(i) occurs in p, ignore var(j) if e[j]>0
1267 int p_GetVariables(poly p, int * e, const ring r)
1268 {
1269  int i;
1270  int n=0;
1271  while(p!=NULL)
1272  {
1273  n=0;
1274  for(i=r->N; i>0; i--)
1275  {
1276  if(e[i]==0)
1277  {
1278  if (p_GetExp(p,i,r)>0)
1279  {
1280  e[i]=1;
1281  n++;
1282  }
1283  }
1284  else
1285  n++;
1286  }
1287  if (n==r->N) break;
1288  pIter(p);
1289  }
1290  return n;
1291 }
1292 
1293 
1294 /*2
1295 * returns a polynomial representing the integer i
1296 */
1297 poly p_ISet(long i, const ring r)
1298 {
1299  poly rc = NULL;
1300  if (i!=0)
1301  {
1302  rc = p_Init(r);
1303  pSetCoeff0(rc,n_Init(i,r->cf));
1304  if (n_IsZero(pGetCoeff(rc),r->cf))
1305  p_LmDelete(&rc,r);
1306  }
1307  return rc;
1308 }
1309 
1310 /*2
1311 * an optimized version of p_ISet for the special case 1
1312 */
1313 poly p_One(const ring r)
1314 {
1315  poly rc = p_Init(r);
1316  pSetCoeff0(rc,n_Init(1,r->cf));
1317  return rc;
1318 }
1319 
1321 {
1322  *h=pNext(p);
1323  pNext(p)=NULL;
1324 }
1325 
1326 /*2
1327 * pair has no common factor ? or is no polynomial
1328 */
1329 BOOLEAN p_HasNotCF(poly p1, poly p2, const ring r)
1330 {
1331 
1332  if (p_GetComp(p1,r) > 0 || p_GetComp(p2,r) > 0)
1333  return FALSE;
1334  int i = rVar(r);
1335  loop
1336  {
1337  if ((p_GetExp(p1, i, r) > 0) && (p_GetExp(p2, i, r) > 0))
1338  return FALSE;
1339  i--;
1340  if (i == 0)
1341  return TRUE;
1342  }
1343 }
1344 
1345 /*2
1346 * convert monomial given as string to poly, e.g. 1x3y5z
1347 */
1348 const char * p_Read(const char *st, poly &rc, const ring r)
1349 {
1350  if (r==NULL) { rc=NULL;return st;}
1351  int i,j;
1352  rc = p_Init(r);
1353  const char *s = n_Read(st,&(p_GetCoeff(rc, r)),r->cf);
1354  if (s==st)
1355  /* i.e. it does not start with a coeff: test if it is a ringvar*/
1356  {
1357  j = r_IsRingVar(s,r->names,r->N);
1358  if (j >= 0)
1359  {
1360  p_IncrExp(rc,1+j,r);
1361  while (*s!='\0') s++;
1362  goto done;
1363  }
1364  }
1365  while (*s!='\0')
1366  {
1367  char ss[2];
1368  ss[0] = *s++;
1369  ss[1] = '\0';
1370  j = r_IsRingVar(ss,r->names,r->N);
1371  if (j >= 0)
1372  {
1373  const char *s_save=s;
1374  s = eati(s,&i);
1375  if (((unsigned long)i) > r->bitmask/2)
1376  {
1377  // exponent to large: it is not a monomial
1378  p_LmDelete(&rc,r);
1379  return s_save;
1380  }
1381  p_AddExp(rc,1+j, (long)i, r);
1382  }
1383  else
1384  {
1385  // 1st char of is not a varname
1386  // We return the parsed polynomial nevertheless. This is needed when
1387  // we are parsing coefficients in a rational function field.
1388  s--;
1389  break;
1390  }
1391  }
1392 done:
1393  if (n_IsZero(pGetCoeff(rc),r->cf)) p_LmDelete(&rc,r);
1394  else
1395  {
1396 #ifdef HAVE_PLURAL
1397  // in super-commutative ring
1398  // squares of anti-commutative variables are zeroes!
1399  if(rIsSCA(r))
1400  {
1401  const unsigned int iFirstAltVar = scaFirstAltVar(r);
1402  const unsigned int iLastAltVar = scaLastAltVar(r);
1403 
1404  assume(rc != NULL);
1405 
1406  for(unsigned int k = iFirstAltVar; k <= iLastAltVar; k++)
1407  if( p_GetExp(rc, k, r) > 1 )
1408  {
1409  p_LmDelete(&rc, r);
1410  goto finish;
1411  }
1412  }
1413 #endif
1414 
1415  p_Setm(rc,r);
1416  }
1417 finish:
1418  return s;
1419 }
1420 poly p_mInit(const char *st, BOOLEAN &ok, const ring r)
1421 {
1422  poly p;
1423  const char *s=p_Read(st,p,r);
1424  if (*s!='\0')
1425  {
1426  if ((s!=st)&&isdigit(st[0]))
1427  {
1429  }
1430  ok=FALSE;
1431  p_Delete(&p,r);
1432  return NULL;
1433  }
1434  p_Test(p,r);
1435  ok=!errorreported;
1436  return p;
1437 }
1438 
1439 /*2
1440 * returns a polynomial representing the number n
1441 * destroys n
1442 */
1443 poly p_NSet(number n, const ring r)
1444 {
1445  if (n_IsZero(n,r->cf))
1446  {
1447  n_Delete(&n, r->cf);
1448  return NULL;
1449  }
1450  else
1451  {
1452  poly rc = p_Init(r);
1453  pSetCoeff0(rc,n);
1454  return rc;
1455  }
1456 }
1457 /*2
1458 * assumes that LM(a) = LM(b)*m, for some monomial m,
1459 * returns the multiplicant m,
1460 * leaves a and b unmodified
1461 */
1462 poly p_Divide(poly a, poly b, const ring r)
1463 {
1464  assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0));
1465  int i;
1466  poly result = p_Init(r);
1467 
1468  for(i=(int)r->N; i; i--)
1469  p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r);
1470  p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r);
1471  p_Setm(result,r);
1472  return result;
1473 }
1474 
1475 poly p_Div_nn(poly p, const number n, const ring r)
1476 {
1477  pAssume(!n_IsZero(n,r->cf));
1478  p_Test(p, r);
1479  poly result = p;
1480  poly prev = NULL;
1481  while (p!=NULL)
1482  {
1483  number nc = n_Div(pGetCoeff(p),n,r->cf);
1484  if (!n_IsZero(nc,r->cf))
1485  {
1486  p_SetCoeff(p,nc,r);
1487  prev=p;
1488  pIter(p);
1489  }
1490  else
1491  {
1492  if (prev==NULL)
1493  {
1494  p_LmDelete(&result,r);
1495  p=result;
1496  }
1497  else
1498  {
1499  p_LmDelete(&pNext(prev),r);
1500  p=pNext(prev);
1501  }
1502  }
1503  }
1504  p_Test(result,r);
1505  return(result);
1506 }
1507 
1508 /*2
1509 * divides a by the monomial b, ignores monomials which are not divisible
1510 * assumes that b is not NULL, destroyes b
1511 */
1512 poly p_DivideM(poly a, poly b, const ring r)
1513 {
1514  if (a==NULL) { p_Delete(&b,r); return NULL; }
1515  poly result=a;
1516  poly prev=NULL;
1517  int i;
1518 #ifdef HAVE_RINGS
1519  number inv=pGetCoeff(b);
1520 #else
1521  number inv=n_Invers(pGetCoeff(b),r->cf);
1522 #endif
1523 
1524  while (a!=NULL)
1525  {
1526  if (p_DivisibleBy(b,a,r))
1527  {
1528  for(i=(int)r->N; i; i--)
1529  p_SubExp(a,i, p_GetExp(b,i,r),r);
1530  p_SubComp(a, p_GetComp(b,r),r);
1531  p_Setm(a,r);
1532  prev=a;
1533  pIter(a);
1534  }
1535  else
1536  {
1537  if (prev==NULL)
1538  {
1539  p_LmDelete(&result,r);
1540  a=result;
1541  }
1542  else
1543  {
1544  p_LmDelete(&pNext(prev),r);
1545  a=pNext(prev);
1546  }
1547  }
1548  }
1549 #ifdef HAVE_RINGS
1550  if (n_IsUnit(inv,r->cf))
1551  {
1552  inv = n_Invers(inv,r->cf);
1553  p_Mult_nn(result,inv,r);
1554  n_Delete(&inv, r->cf);
1555  }
1556  else
1557  {
1558  result = p_Div_nn(result,inv,r);
1559  }
1560 #else
1561  result = p_Mult_nn(result,inv,r);
1562  n_Delete(&inv, r->cf);
1563 #endif
1564  p_Delete(&b, r);
1565  return result;
1566 }
1567 
1568 #ifdef HAVE_RINGS
1569 /* TRUE iff LT(f) | LT(g) */
1571 {
1572  int exponent;
1573  for(int i = (int)rVar(r); i>0; i--)
1574  {
1575  exponent = p_GetExp(g, i, r) - p_GetExp(f, i, r);
1576  if (exponent < 0) return FALSE;
1577  }
1578  return n_DivBy(pGetCoeff(g), pGetCoeff(f), r->cf);
1579 }
1580 #endif
1581 
1582 // returns the LCM of the head terms of a and b in *m
1583 void p_Lcm(const poly a, const poly b, poly m, const ring r)
1584 {
1585  for (int i=rVar(r); i; --i)
1586  p_SetExp(m,i, si_max( p_GetExp(a,i,r), p_GetExp(b,i,r)),r);
1587 
1588  p_SetComp(m, si_max(p_GetComp(a,r), p_GetComp(b,r)),r);
1589  /* Don't do a pSetm here, otherwise hres/lres chockes */
1590 }
1591 
1592 
1593 
1594 #ifdef HAVE_RATGRING
1595 /*2
1596 * returns the rational LCM of the head terms of a and b
1597 * without coefficient!!!
1598 */
1599 poly p_LcmRat(const poly a, const poly b, const long lCompM, const ring r)
1600 {
1601  poly m = // p_One( r);
1602  p_Init(r);
1603 
1604 // const int (currRing->N) = r->N;
1605 
1606  // for (int i = (currRing->N); i>=r->real_var_start; i--)
1607  for (int i = r->real_var_end; i>=r->real_var_start; i--)
1608  {
1609  const int lExpA = p_GetExp (a, i, r);
1610  const int lExpB = p_GetExp (b, i, r);
1611 
1612  p_SetExp (m, i, si_max(lExpA, lExpB), r);
1613  }
1614 
1615  p_SetComp (m, lCompM, r);
1616  p_Setm(m,r);
1617  n_New(&(p_GetCoeff(m, r)), r);
1618 
1619  return(m);
1620 };
1621 
1622 void p_LmDeleteAndNextRat(poly *p, int ishift, ring r)
1623 {
1624  /* modifies p*/
1625  // Print("start: "); Print(" "); p_wrp(*p,r);
1626  p_LmCheckPolyRing2(*p, r);
1627  poly q = p_Head(*p,r);
1628  const long cmp = p_GetComp(*p, r);
1629  while ( ( (*p)!=NULL ) && ( p_Comp_k_n(*p, q, ishift+1, r) ) && (p_GetComp(*p, r) == cmp) )
1630  {
1631  p_LmDelete(p,r);
1632  // Print("while: ");p_wrp(*p,r);Print(" ");
1633  }
1634  // p_wrp(*p,r);Print(" ");
1635  // PrintS("end\n");
1636  p_LmDelete(&q,r);
1637 }
1638 
1639 
1640 /* returns x-coeff of p, i.e. a poly in x, s.t. corresponding xd-monomials
1641 have the same D-part and the component 0
1642 does not destroy p
1643 */
1644 poly p_GetCoeffRat(poly p, int ishift, ring r)
1645 {
1646  poly q = pNext(p);
1647  poly res; // = p_Head(p,r);
1648  res = p_GetExp_k_n(p, ishift+1, r->N, r); // does pSetm internally
1649  p_SetCoeff(res,n_Copy(p_GetCoeff(p,r),r),r);
1650  poly s;
1651  long cmp = p_GetComp(p, r);
1652  while ( (q!= NULL) && (p_Comp_k_n(p, q, ishift+1, r)) && (p_GetComp(q, r) == cmp) )
1653  {
1654  s = p_GetExp_k_n(q, ishift+1, r->N, r);
1655  p_SetCoeff(s,n_Copy(p_GetCoeff(q,r),r),r);
1656  res = p_Add_q(res,s,r);
1657  q = pNext(q);
1658  }
1659  cmp = 0;
1660  p_SetCompP(res,cmp,r);
1661  return res;
1662 }
1663 
1664 
1665 
1666 void p_ContentRat(poly &ph, const ring r)
1667 // changes ph
1668 // for rat coefficients in K(x1,..xN)
1669 {
1670  // init array of RatLeadCoeffs
1671  // poly p_GetCoeffRat(poly p, int ishift, ring r);
1672 
1673  int len=pLength(ph);
1674  poly *C = (poly *)omAlloc0((len+1)*sizeof(poly)); //rat coeffs
1675  poly *LM = (poly *)omAlloc0((len+1)*sizeof(poly)); // rat lead terms
1676  int *D = (int *)omAlloc0((len+1)*sizeof(int)); //degrees of coeffs
1677  int *L = (int *)omAlloc0((len+1)*sizeof(int)); //lengths of coeffs
1678  int k = 0;
1679  poly p = p_Copy(ph, r); // ph will be needed below
1680  int mintdeg = p_Totaldegree(p, r);
1681  int minlen = len;
1682  int dd = 0; int i;
1683  int HasConstantCoef = 0;
1684  int is = r->real_var_start - 1;
1685  while (p!=NULL)
1686  {
1687  LM[k] = p_GetExp_k_n(p,1,is, r); // need LmRat istead of p_HeadRat(p, is, currRing); !
1688  C[k] = p_GetCoeffRat(p, is, r);
1689  D[k] = p_Totaldegree(C[k], r);
1690  mintdeg = si_min(mintdeg,D[k]);
1691  L[k] = pLength(C[k]);
1692  minlen = si_min(minlen,L[k]);
1693  if (p_IsConstant(C[k], r))
1694  {
1695  // C[k] = const, so the content will be numerical
1696  HasConstantCoef = 1;
1697  // smth like goto cleanup and return(pContent(p));
1698  }
1699  p_LmDeleteAndNextRat(&p, is, r);
1700  k++;
1701  }
1702 
1703  // look for 1 element of minimal degree and of minimal length
1704  k--;
1705  poly d;
1706  int mindeglen = len;
1707  if (k<=0) // this poly is not a ratgring poly -> pContent
1708  {
1709  p_Delete(&C[0], r);
1710  p_Delete(&LM[0], r);
1711  p_Content(ph, r);
1712  goto cleanup;
1713  }
1714 
1715  int pmindeglen;
1716  for(i=0; i<=k; i++)
1717  {
1718  if (D[i] == mintdeg)
1719  {
1720  if (L[i] < mindeglen)
1721  {
1722  mindeglen=L[i];
1723  pmindeglen = i;
1724  }
1725  }
1726  }
1727  d = p_Copy(C[pmindeglen], r);
1728  // there are dd>=1 mindeg elements
1729  // and pmideglen is the coordinate of one of the smallest among them
1730 
1731  // poly g = singclap_gcd(p_Copy(p,r),p_Copy(q,r));
1732  // return naGcd(d,d2,currRing);
1733 
1734  // adjoin pContentRat here?
1735  for(i=0; i<=k; i++)
1736  {
1737  d=singclap_gcd(d,p_Copy(C[i], r), r);
1738  if (p_Totaldegree(d, r)==0)
1739  {
1740  // cleanup, pContent, return
1741  p_Delete(&d, r);
1742  for(;k>=0;k--)
1743  {
1744  p_Delete(&C[k], r);
1745  p_Delete(&LM[k], r);
1746  }
1747  p_Content(ph, r);
1748  goto cleanup;
1749  }
1750  }
1751  for(i=0; i<=k; i++)
1752  {
1753  poly h=singclap_pdivide(C[i],d, r);
1754  p_Delete(&C[i], r);
1755  C[i]=h;
1756  }
1757 
1758  // zusammensetzen,
1759  p=NULL; // just to be sure
1760  for(i=0; i<=k; i++)
1761  {
1762  p = p_Add_q(p, p_Mult_q(C[i],LM[i], r), r);
1763  C[i]=NULL; LM[i]=NULL;
1764  }
1765  p_Delete(&ph, r); // do not need it anymore
1766  ph = p;
1767  // aufraeumen, return
1768 cleanup:
1769  omFree(C);
1770  omFree(LM);
1771  omFree(D);
1772  omFree(L);
1773 }
1774 
1775 
1776 #endif
1777 
1778 
1779 /* assumes that p and divisor are univariate polynomials in r,
1780  mentioning the same variable;
1781  assumes divisor != NULL;
1782  p may be NULL;
1783  assumes a global monomial ordering in r;
1784  performs polynomial division of p by divisor:
1785  - afterwards p contains the remainder of the division, i.e.,
1786  p_before = result * divisor + p_afterwards;
1787  - if needResult == TRUE, then the method computes and returns 'result',
1788  otherwise NULL is returned (This parametrization can be used when
1789  one is only interested in the remainder of the division. In this
1790  case, the method will be slightly faster.)
1791  leaves divisor unmodified */
1792 poly p_PolyDiv(poly &p, const poly divisor, const BOOLEAN needResult, const ring r)
1793 {
1794  assume(divisor != NULL);
1795  if (p == NULL) return NULL;
1796 
1797  poly result = NULL;
1798  number divisorLC = p_GetCoeff(divisor, r);
1799  int divisorLE = p_GetExp(divisor, 1, r);
1800  while ((p != NULL) && (p_Deg(p, r) >= p_Deg(divisor, r)))
1801  {
1802  /* determine t = LT(p) / LT(divisor) */
1803  poly t = p_ISet(1, r);
1804  number c = n_Div(p_GetCoeff(p, r), divisorLC, r->cf);
1805  n_Normalize(c,r->cf);
1806  p_SetCoeff(t, c, r);
1807  int e = p_GetExp(p, 1, r) - divisorLE;
1808  p_SetExp(t, 1, e, r);
1809  p_Setm(t, r);
1810  if (needResult) result = p_Add_q(result, p_Copy(t, r), r);
1811  p = p_Add_q(p, p_Neg(p_Mult_q(t, p_Copy(divisor, r), r), r), r);
1812  }
1813  return result;
1814 }
1815 
1816 /*2
1817 * returns the partial differentiate of a by the k-th variable
1818 * does not destroy the input
1819 */
1820 poly p_Diff(poly a, int k, const ring r)
1821 {
1822  poly res, f, last;
1823  number t;
1824 
1825  last = res = NULL;
1826  while (a!=NULL)
1827  {
1828  if (p_GetExp(a,k,r)!=0)
1829  {
1830  f = p_LmInit(a,r);
1831  t = n_Init(p_GetExp(a,k,r),r->cf);
1832  pSetCoeff0(f,n_Mult(t,pGetCoeff(a),r->cf));
1833  n_Delete(&t,r->cf);
1834  if (n_IsZero(pGetCoeff(f),r->cf))
1835  p_LmDelete(&f,r);
1836  else
1837  {
1838  p_DecrExp(f,k,r);
1839  p_Setm(f,r);
1840  if (res==NULL)
1841  {
1842  res=last=f;
1843  }
1844  else
1845  {
1846  pNext(last)=f;
1847  last=f;
1848  }
1849  }
1850  }
1851  pIter(a);
1852  }
1853  return res;
1854 }
1855 
1856 static poly p_DiffOpM(poly a, poly b,BOOLEAN multiply, const ring r)
1857 {
1858  int i,j,s;
1859  number n,h,hh;
1860  poly p=p_One(r);
1861  n=n_Mult(pGetCoeff(a),pGetCoeff(b),r->cf);
1862  for(i=rVar(r);i>0;i--)
1863  {
1864  s=p_GetExp(b,i,r);
1865  if (s<p_GetExp(a,i,r))
1866  {
1867  n_Delete(&n,r->cf);
1868  p_LmDelete(&p,r);
1869  return NULL;
1870  }
1871  if (multiply)
1872  {
1873  for(j=p_GetExp(a,i,r); j>0;j--)
1874  {
1875  h = n_Init(s,r->cf);
1876  hh=n_Mult(n,h,r->cf);
1877  n_Delete(&h,r->cf);
1878  n_Delete(&n,r->cf);
1879  n=hh;
1880  s--;
1881  }
1882  p_SetExp(p,i,s,r);
1883  }
1884  else
1885  {
1886  p_SetExp(p,i,s-p_GetExp(a,i,r),r);
1887  }
1888  }
1889  p_Setm(p,r);
1890  /*if (multiply)*/ p_SetCoeff(p,n,r);
1891  if (n_IsZero(n,r->cf)) p=p_LmDeleteAndNext(p,r); // return NULL as p is a monomial
1892  return p;
1893 }
1894 
1895 poly p_DiffOp(poly a, poly b,BOOLEAN multiply, const ring r)
1896 {
1897  poly result=NULL;
1898  poly h;
1899  for(;a!=NULL;pIter(a))
1900  {
1901  for(h=b;h!=NULL;pIter(h))
1902  {
1903  result=p_Add_q(result,p_DiffOpM(a,h,multiply,r),r);
1904  }
1905  }
1906  return result;
1907 }
1908 /*2
1909 * subtract p2 from p1, p1 and p2 are destroyed
1910 * do not put attention on speed: the procedure is only used in the interpreter
1911 */
1912 poly p_Sub(poly p1, poly p2, const ring r)
1913 {
1914  return p_Add_q(p1, p_Neg(p2,r),r);
1915 }
1916 
1917 /*3
1918 * compute for a monomial m
1919 * the power m^exp, exp > 1
1920 * destroys p
1921 */
1922 static poly p_MonPower(poly p, int exp, const ring r)
1923 {
1924  int i;
1925 
1926  if(!n_IsOne(pGetCoeff(p),r->cf))
1927  {
1928  number x, y;
1929  y = pGetCoeff(p);
1930  n_Power(y,exp,&x,r->cf);
1931  n_Delete(&y,r->cf);
1932  pSetCoeff0(p,x);
1933  }
1934  for (i=rVar(r); i!=0; i--)
1935  {
1936  p_MultExp(p,i, exp,r);
1937  }
1938  p_Setm(p,r);
1939  return p;
1940 }
1941 
1942 /*3
1943 * compute for monomials p*q
1944 * destroys p, keeps q
1945 */
1946 static void p_MonMult(poly p, poly q, const ring r)
1947 {
1948  number x, y;
1949 
1950  y = pGetCoeff(p);
1951  x = n_Mult(y,pGetCoeff(q),r->cf);
1952  n_Delete(&y,r->cf);
1953  pSetCoeff0(p,x);
1954  //for (int i=pVariables; i!=0; i--)
1955  //{
1956  // pAddExp(p,i, pGetExp(q,i));
1957  //}
1958  //p->Order += q->Order;
1959  p_ExpVectorAdd(p,q,r);
1960 }
1961 
1962 /*3
1963 * compute for monomials p*q
1964 * keeps p, q
1965 */
1966 static poly p_MonMultC(poly p, poly q, const ring rr)
1967 {
1968  number x;
1969  poly r = p_Init(rr);
1970 
1971  x = n_Mult(pGetCoeff(p),pGetCoeff(q),rr->cf);
1972  pSetCoeff0(r,x);
1973  p_ExpVectorSum(r,p, q, rr);
1974  return r;
1975 }
1976 
1977 /*3
1978 * create binomial coef.
1979 */
1980 static number* pnBin(int exp, const ring r)
1981 {
1982  int e, i, h;
1983  number x, y, *bin=NULL;
1984 
1985  x = n_Init(exp,r->cf);
1986  if (n_IsZero(x,r->cf))
1987  {
1988  n_Delete(&x,r->cf);
1989  return bin;
1990  }
1991  h = (exp >> 1) + 1;
1992  bin = (number *)omAlloc0(h*sizeof(number));
1993  bin[1] = x;
1994  if (exp < 4)
1995  return bin;
1996  i = exp - 1;
1997  for (e=2; e<h; e++)
1998  {
1999  x = n_Init(i,r->cf);
2000  i--;
2001  y = n_Mult(x,bin[e-1],r->cf);
2002  n_Delete(&x,r->cf);
2003  x = n_Init(e,r->cf);
2004  bin[e] = n_ExactDiv(y,x,r->cf);
2005  n_Delete(&x,r->cf);
2006  n_Delete(&y,r->cf);
2007  }
2008  return bin;
2009 }
2010 
2011 static void pnFreeBin(number *bin, int exp,const coeffs r)
2012 {
2013  int e, h = (exp >> 1) + 1;
2014 
2015  if (bin[1] != NULL)
2016  {
2017  for (e=1; e<h; e++)
2018  n_Delete(&(bin[e]),r);
2019  }
2020  omFreeSize((ADDRESS)bin, h*sizeof(number));
2021 }
2022 
2023 /*
2024 * compute for a poly p = head+tail, tail is monomial
2025 * (head + tail)^exp, exp > 1
2026 * with binomial coef.
2027 */
2028 static poly p_TwoMonPower(poly p, int exp, const ring r)
2029 {
2030  int eh, e;
2031  long al;
2032  poly *a;
2033  poly tail, b, res, h;
2034  number x;
2035  number *bin = pnBin(exp,r);
2036 
2037  tail = pNext(p);
2038  if (bin == NULL)
2039  {
2040  p_MonPower(p,exp,r);
2041  p_MonPower(tail,exp,r);
2042  p_Test(p,r);
2043  return p;
2044  }
2045  eh = exp >> 1;
2046  al = (exp + 1) * sizeof(poly);
2047  a = (poly *)omAlloc(al);
2048  a[1] = p;
2049  for (e=1; e<exp; e++)
2050  {
2051  a[e+1] = p_MonMultC(a[e],p,r);
2052  }
2053  res = a[exp];
2054  b = p_Head(tail,r);
2055  for (e=exp-1; e>eh; e--)
2056  {
2057  h = a[e];
2058  x = n_Mult(bin[exp-e],pGetCoeff(h),r->cf);
2059  p_SetCoeff(h,x,r);
2060  p_MonMult(h,b,r);
2061  res = pNext(res) = h;
2062  p_MonMult(b,tail,r);
2063  }
2064  for (e=eh; e!=0; e--)
2065  {
2066  h = a[e];
2067  x = n_Mult(bin[e],pGetCoeff(h),r->cf);
2068  p_SetCoeff(h,x,r);
2069  p_MonMult(h,b,r);
2070  res = pNext(res) = h;
2071  p_MonMult(b,tail,r);
2072  }
2073  p_LmDelete(&tail,r);
2074  pNext(res) = b;
2075  pNext(b) = NULL;
2076  res = a[exp];
2077  omFreeSize((ADDRESS)a, al);
2078  pnFreeBin(bin, exp, r->cf);
2079 // tail=res;
2080 // while((tail!=NULL)&&(pNext(tail)!=NULL))
2081 // {
2082 // if(nIsZero(pGetCoeff(pNext(tail))))
2083 // {
2084 // pLmDelete(&pNext(tail));
2085 // }
2086 // else
2087 // pIter(tail);
2088 // }
2089  p_Test(res,r);
2090  return res;
2091 }
2092 
2093 static poly p_Pow(poly p, int i, const ring r)
2094 {
2095  poly rc = p_Copy(p,r);
2096  i -= 2;
2097  do
2098  {
2099  rc = p_Mult_q(rc,p_Copy(p,r),r);
2100  p_Normalize(rc,r);
2101  i--;
2102  }
2103  while (i != 0);
2104  return p_Mult_q(rc,p,r);
2105 }
2106 
2107 static poly p_Pow_charp(poly p, int i, const ring r)
2108 {
2109  //assume char_p == i
2110  poly h=p;
2111  while(h!=NULL) { p_MonPower(h,i,r);pIter(h);}
2112  return p;
2113 }
2114 
2115 /*2
2116 * returns the i-th power of p
2117 * p will be destroyed
2118 */
2119 poly p_Power(poly p, int i, const ring r)
2120 {
2121  poly rc=NULL;
2122 
2123  if (i==0)
2124  {
2125  p_Delete(&p,r);
2126  return p_One(r);
2127  }
2128 
2129  if(p!=NULL)
2130  {
2131  if ( (i > 0) && ((unsigned long ) i > (r->bitmask)))
2132  {
2133  Werror("exponent %d is too large, max. is %ld",i,r->bitmask);
2134  return NULL;
2135  }
2136  switch (i)
2137  {
2138 // cannot happen, see above
2139 // case 0:
2140 // {
2141 // rc=pOne();
2142 // pDelete(&p);
2143 // break;
2144 // }
2145  case 1:
2146  rc=p;
2147  break;
2148  case 2:
2149  rc=p_Mult_q(p_Copy(p,r),p,r);
2150  break;
2151  default:
2152  if (i < 0)
2153  {
2154  p_Delete(&p,r);
2155  return NULL;
2156  }
2157  else
2158  {
2159 #ifdef HAVE_PLURAL
2160  if (rIsPluralRing(r)) /* in the NC case nothing helps :-( */
2161  {
2162  int j=i;
2163  rc = p_Copy(p,r);
2164  while (j>1)
2165  {
2166  rc = p_Mult_q(p_Copy(p,r),rc,r);
2167  j--;
2168  }
2169  p_Delete(&p,r);
2170  return rc;
2171  }
2172 #endif
2173  rc = pNext(p);
2174  if (rc == NULL)
2175  return p_MonPower(p,i,r);
2176  /* else: binom ?*/
2177  int char_p=rChar(r);
2178  if ((char_p>0) && (i>char_p)
2179  && ((rField_is_Zp(r,char_p)
2180  || (rField_is_Zp_a(r,char_p)))))
2181  {
2182  poly h=p_Pow_charp(p_Copy(p,r),char_p,r);
2183  int rest=i-char_p;
2184  while (rest>=char_p)
2185  {
2186  rest-=char_p;
2187  h=p_Mult_q(h,p_Pow_charp(p_Copy(p,r),char_p,r),r);
2188  }
2189  poly res=h;
2190  if (rest>0)
2191  res=p_Mult_q(p_Power(p_Copy(p,r),rest,r),h,r);
2192  p_Delete(&p,r);
2193  return res;
2194  }
2195  if ((pNext(rc) != NULL)
2196  || rField_is_Ring(r)
2197  )
2198  return p_Pow(p,i,r);
2199  if ((char_p==0) || (i<=char_p))
2200  return p_TwoMonPower(p,i,r);
2201  return p_Pow(p,i,r);
2202  }
2203  /*end default:*/
2204  }
2205  }
2206  return rc;
2207 }
2208 
2209 /* --------------------------------------------------------------------------------*/
2210 /* content suff */
2211 
2212 static number p_InitContent(poly ph, const ring r);
2213 
2214 #define CLEARENUMERATORS 1
2215 
2216 void p_Content(poly ph, const ring r)
2217 {
2218  assume( ph != NULL );
2219 
2220  assume( r != NULL ); assume( r->cf != NULL );
2221 
2222 
2223 #if CLEARENUMERATORS
2224  if( 0 )
2225  {
2226  const coeffs C = r->cf;
2227  // experimentall (recursive enumerator treatment) of alg. Ext!
2228  CPolyCoeffsEnumerator itr(ph);
2229  n_ClearContent(itr, r->cf);
2230 
2231  p_Test(ph, r); n_Test(pGetCoeff(ph), C);
2232  assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
2233 
2234  // if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2235  return;
2236  }
2237 #endif
2238 
2239 
2240 #ifdef HAVE_RINGS
2241  if (rField_is_Ring(r))
2242  {
2243  if (rField_has_Units(r))
2244  {
2245  number k = n_GetUnit(pGetCoeff(ph),r->cf);
2246  if (!n_IsOne(k,r->cf))
2247  {
2248  number tmpGMP = k;
2249  k = n_Invers(k,r->cf);
2250  n_Delete(&tmpGMP,r->cf);
2251  poly h = pNext(ph);
2252  p_SetCoeff(ph, n_Mult(pGetCoeff(ph), k,r->cf),r);
2253  while (h != NULL)
2254  {
2255  p_SetCoeff(h, n_Mult(pGetCoeff(h), k,r->cf),r);
2256  pIter(h);
2257  }
2258 // assume( n_GreaterZero(pGetCoeff(ph),r->cf) );
2259 // if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2260  }
2261  n_Delete(&k,r->cf);
2262  }
2263  return;
2264  }
2265 #endif
2266  number h,d;
2267  poly p;
2268 
2269  if(TEST_OPT_CONTENTSB) return;
2270  if(pNext(ph)==NULL)
2271  {
2272  p_SetCoeff(ph,n_Init(1,r->cf),r);
2273  }
2274  else
2275  {
2276  assume( pNext(ph) != NULL );
2277 #if CLEARENUMERATORS
2278  if( nCoeff_is_Q(r->cf) )
2279  {
2280  // experimentall (recursive enumerator treatment) of alg. Ext!
2281  CPolyCoeffsEnumerator itr(ph);
2282  n_ClearContent(itr, r->cf);
2283 
2284  p_Test(ph, r); n_Test(pGetCoeff(ph), r->cf);
2285  assume(n_GreaterZero(pGetCoeff(ph), r->cf)); // ??
2286 
2287  // if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2288  return;
2289  }
2290 #endif
2291 
2292  n_Normalize(pGetCoeff(ph),r->cf);
2293  if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2294  if (rField_is_Q(r)||(getCoeffType(r->cf)==n_transExt)) // should not be used anymore if CLEARENUMERATORS is 1
2295  {
2296  h=p_InitContent(ph,r);
2297  p=ph;
2298  }
2299  else
2300  {
2301  h=n_Copy(pGetCoeff(ph),r->cf);
2302  p = pNext(ph);
2303  }
2304  while (p!=NULL)
2305  {
2306  n_Normalize(pGetCoeff(p),r->cf);
2307  d=n_SubringGcd(h,pGetCoeff(p),r->cf);
2308  n_Delete(&h,r->cf);
2309  h = d;
2310  if(n_IsOne(h,r->cf))
2311  {
2312  break;
2313  }
2314  pIter(p);
2315  }
2316  p = ph;
2317  //number tmp;
2318  if(!n_IsOne(h,r->cf))
2319  {
2320  while (p!=NULL)
2321  {
2322  //d = nDiv(pGetCoeff(p),h);
2323  //tmp = nExactDiv(pGetCoeff(p),h);
2324  //if (!nEqual(d,tmp))
2325  //{
2326  // StringSetS("** div0:");nWrite(pGetCoeff(p));StringAppendS("/");
2327  // nWrite(h);StringAppendS("=");nWrite(d);StringAppendS(" int:");
2328  // nWrite(tmp);Print(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
2329  //}
2330  //nDelete(&tmp);
2331  d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2332  p_SetCoeff(p,d,r);
2333  pIter(p);
2334  }
2335  }
2336  n_Delete(&h,r->cf);
2337  if (rField_is_Q_a(r))
2338  {
2339  // special handling for alg. ext.:
2340  if (getCoeffType(r->cf)==n_algExt)
2341  {
2342  h = n_Init(1, r->cf->extRing->cf);
2343  p=ph;
2344  while (p!=NULL)
2345  { // each monom: coeff in Q_a
2346  poly c_n_n=(poly)pGetCoeff(p);
2347  poly c_n=c_n_n;
2348  while (c_n!=NULL)
2349  { // each monom: coeff in Q
2350  d=n_NormalizeHelper(h,pGetCoeff(c_n),r->cf->extRing->cf);
2351  n_Delete(&h,r->cf->extRing->cf);
2352  h=d;
2353  pIter(c_n);
2354  }
2355  pIter(p);
2356  }
2357  /* h contains the 1/lcm of all denominators in c_n_n*/
2358  //n_Normalize(h,r->cf->extRing->cf);
2359  if(!n_IsOne(h,r->cf->extRing->cf))
2360  {
2361  p=ph;
2362  while (p!=NULL)
2363  { // each monom: coeff in Q_a
2364  poly c_n=(poly)pGetCoeff(p);
2365  while (c_n!=NULL)
2366  { // each monom: coeff in Q
2367  d=n_Mult(h,pGetCoeff(c_n),r->cf->extRing->cf);
2368  n_Normalize(d,r->cf->extRing->cf);
2369  n_Delete(&pGetCoeff(c_n),r->cf->extRing->cf);
2370  pGetCoeff(c_n)=d;
2371  pIter(c_n);
2372  }
2373  pIter(p);
2374  }
2375  }
2376  n_Delete(&h,r->cf->extRing->cf);
2377  }
2378  /*else
2379  {
2380  // special handling for rat. functions.:
2381  number hzz =NULL;
2382  p=ph;
2383  while (p!=NULL)
2384  { // each monom: coeff in Q_a (Z_a)
2385  fraction f=(fraction)pGetCoeff(p);
2386  poly c_n=NUM(f);
2387  if (hzz==NULL)
2388  {
2389  hzz=n_Copy(pGetCoeff(c_n),r->cf->extRing->cf);
2390  pIter(c_n);
2391  }
2392  while ((c_n!=NULL)&&(!n_IsOne(hzz,r->cf->extRing->cf)))
2393  { // each monom: coeff in Q (Z)
2394  d=n_Gcd(hzz,pGetCoeff(c_n),r->cf->extRing->cf);
2395  n_Delete(&hzz,r->cf->extRing->cf);
2396  hzz=d;
2397  pIter(c_n);
2398  }
2399  pIter(p);
2400  }
2401  // hzz contains the gcd of all numerators in f
2402  h=n_Invers(hzz,r->cf->extRing->cf);
2403  n_Delete(&hzz,r->cf->extRing->cf);
2404  n_Normalize(h,r->cf->extRing->cf);
2405  if(!n_IsOne(h,r->cf->extRing->cf))
2406  {
2407  p=ph;
2408  while (p!=NULL)
2409  { // each monom: coeff in Q_a (Z_a)
2410  fraction f=(fraction)pGetCoeff(p);
2411  NUM(f)=p_Mult_nn(NUM(f),h,r->cf->extRing);
2412  p_Normalize(NUM(f),r->cf->extRing);
2413  pIter(p);
2414  }
2415  }
2416  n_Delete(&h,r->cf->extRing->cf);
2417  }*/
2418  }
2419  }
2420  if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2421 }
2422 
2423 // Not yet?
2424 #if 1 // currently only used by Singular/janet
2425 void p_SimpleContent(poly ph, int smax, const ring r)
2426 {
2427  if(TEST_OPT_CONTENTSB) return;
2428  if (ph==NULL) return;
2429  if (pNext(ph)==NULL)
2430  {
2431  p_SetCoeff(ph,n_Init(1,r->cf),r);
2432  return;
2433  }
2434  if ((pNext(pNext(ph))==NULL)||(!rField_is_Q(r)))
2435  {
2436  return;
2437  }
2438  number d=p_InitContent(ph,r);
2439  if (n_Size(d,r->cf)<=smax)
2440  {
2441  //if (TEST_OPT_PROT) PrintS("G");
2442  return;
2443  }
2444 
2445 
2446  poly p=ph;
2447  number h=d;
2448  if (smax==1) smax=2;
2449  while (p!=NULL)
2450  {
2451 #if 0
2452  d=n_Gcd(h,pGetCoeff(p),r->cf);
2453  n_Delete(&h,r->cf);
2454  h = d;
2455 #else
2456  STATISTIC(n_Gcd); nlInpGcd(h,pGetCoeff(p),r->cf); // FIXME? TODO? // extern void nlInpGcd(number &a, number b, const coeffs r);
2457 #endif
2458  if(n_Size(h,r->cf)<smax)
2459  {
2460  //if (TEST_OPT_PROT) PrintS("g");
2461  return;
2462  }
2463  pIter(p);
2464  }
2465  p = ph;
2466  if (!n_GreaterZero(pGetCoeff(p),r->cf)) h=n_InpNeg(h,r->cf);
2467  if(n_IsOne(h,r->cf)) return;
2468  //if (TEST_OPT_PROT) PrintS("c");
2469  while (p!=NULL)
2470  {
2471 #if 1
2472  d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2473  p_SetCoeff(p,d,r);
2474 #else
2475  STATISTIC(n_ExactDiv); nlInpExactDiv(pGetCoeff(p),h,r->cf); // no such function... ?
2476 #endif
2477  pIter(p);
2478  }
2479  n_Delete(&h,r->cf);
2480 }
2481 #endif
2482 
2483 static number p_InitContent(poly ph, const ring r)
2484 // only for coefficients in Q and rational functions
2485 #if 0
2486 {
2488  assume(ph!=NULL);
2489  assume(pNext(ph)!=NULL);
2490  assume(rField_is_Q(r));
2491  if (pNext(pNext(ph))==NULL)
2492  {
2493  return n_GetNumerator(pGetCoeff(pNext(ph)),r->cf);
2494  }
2495  poly p=ph;
2496  number n1=n_GetNumerator(pGetCoeff(p),r->cf);
2497  pIter(p);
2498  number n2=n_GetNumerator(pGetCoeff(p),r->cf);
2499  pIter(p);
2500  number d;
2501  number t;
2502  loop
2503  {
2504  nlNormalize(pGetCoeff(p),r->cf);
2505  t=n_GetNumerator(pGetCoeff(p),r->cf);
2506  if (nlGreaterZero(t,r->cf))
2507  d=nlAdd(n1,t,r->cf);
2508  else
2509  d=nlSub(n1,t,r->cf);
2510  nlDelete(&t,r->cf);
2511  nlDelete(&n1,r->cf);
2512  n1=d;
2513  pIter(p);
2514  if (p==NULL) break;
2515  nlNormalize(pGetCoeff(p),r->cf);
2516  t=n_GetNumerator(pGetCoeff(p),r->cf);
2517  if (nlGreaterZero(t,r->cf))
2518  d=nlAdd(n2,t,r->cf);
2519  else
2520  d=nlSub(n2,t,r->cf);
2521  nlDelete(&t,r->cf);
2522  nlDelete(&n2,r->cf);
2523  n2=d;
2524  pIter(p);
2525  if (p==NULL) break;
2526  }
2527  d=nlGcd(n1,n2,r->cf);
2528  nlDelete(&n1,r->cf);
2529  nlDelete(&n2,r->cf);
2530  return d;
2531 }
2532 #else
2533 {
2534  number d=pGetCoeff(ph);
2535  int s;
2536  int s2=-1;
2537  if(rField_is_Q(r))
2538  {
2539  if (SR_HDL(d)&SR_INT) return d;
2540  s=mpz_size1(d->z);
2541  }
2542  else
2543  s=n_Size(d,r);
2544  number d2=d;
2545  loop
2546  {
2547  pIter(ph);
2548  if(ph==NULL)
2549  {
2550  if (s2==-1) return n_Copy(d,r->cf);
2551  break;
2552  }
2553  if (rField_is_Q(r))
2554  {
2555  if (SR_HDL(pGetCoeff(ph))&SR_INT)
2556  {
2557  s2=s;
2558  d2=d;
2559  s=0;
2560  d=pGetCoeff(ph);
2561  if (s2==0) break;
2562  }
2563  else if (mpz_size1((pGetCoeff(ph)->z))<=s)
2564  {
2565  s2=s;
2566  d2=d;
2567  d=pGetCoeff(ph);
2568  s=mpz_size1(d->z);
2569  }
2570  }
2571  else
2572  {
2573  int ns=n_Size(pGetCoeff(ph),r);
2574  if (ns<=3)
2575  {
2576  s2=s;
2577  d2=d;
2578  d=pGetCoeff(ph);
2579  s=ns;
2580  if (s2<=3) break;
2581  }
2582  else if (ns<s)
2583  {
2584  s2=s;
2585  d2=d;
2586  d=pGetCoeff(ph);
2587  s=ns;
2588  }
2589  }
2590  }
2591  return n_SubringGcd(d,d2,r->cf);
2592 }
2593 #endif
2594 
2595 //void pContent(poly ph)
2596 //{
2597 // number h,d;
2598 // poly p;
2599 //
2600 // p = ph;
2601 // if(pNext(p)==NULL)
2602 // {
2603 // pSetCoeff(p,nInit(1));
2604 // }
2605 // else
2606 // {
2607 //#ifdef PDEBUG
2608 // if (!pTest(p)) return;
2609 //#endif
2610 // nNormalize(pGetCoeff(p));
2611 // if(!nGreaterZero(pGetCoeff(ph)))
2612 // {
2613 // ph = pNeg(ph);
2614 // nNormalize(pGetCoeff(p));
2615 // }
2616 // h=pGetCoeff(p);
2617 // pIter(p);
2618 // while (p!=NULL)
2619 // {
2620 // nNormalize(pGetCoeff(p));
2621 // if (nGreater(h,pGetCoeff(p))) h=pGetCoeff(p);
2622 // pIter(p);
2623 // }
2624 // h=nCopy(h);
2625 // p=ph;
2626 // while (p!=NULL)
2627 // {
2628 // d=n_Gcd(h,pGetCoeff(p));
2629 // nDelete(&h);
2630 // h = d;
2631 // if(nIsOne(h))
2632 // {
2633 // break;
2634 // }
2635 // pIter(p);
2636 // }
2637 // p = ph;
2638 // //number tmp;
2639 // if(!nIsOne(h))
2640 // {
2641 // while (p!=NULL)
2642 // {
2643 // d = nExactDiv(pGetCoeff(p),h);
2644 // pSetCoeff(p,d);
2645 // pIter(p);
2646 // }
2647 // }
2648 // nDelete(&h);
2649 // if ( (nGetChar() == 1) || (nGetChar() < 0) ) /* Q[a],Q(a),Zp[a],Z/p(a) */
2650 // {
2651 // pTest(ph);
2652 // singclap_divide_content(ph);
2653 // pTest(ph);
2654 // }
2655 // }
2656 //}
2657 #if 0
2658 void p_Content(poly ph, const ring r)
2659 {
2660  number h,d;
2661  poly p;
2662 
2663  if(pNext(ph)==NULL)
2664  {
2665  p_SetCoeff(ph,n_Init(1,r->cf),r);
2666  }
2667  else
2668  {
2669  n_Normalize(pGetCoeff(ph),r->cf);
2670  if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2671  h=n_Copy(pGetCoeff(ph),r->cf);
2672  p = pNext(ph);
2673  while (p!=NULL)
2674  {
2675  n_Normalize(pGetCoeff(p),r->cf);
2676  d=n_Gcd(h,pGetCoeff(p),r->cf);
2677  n_Delete(&h,r->cf);
2678  h = d;
2679  if(n_IsOne(h,r->cf))
2680  {
2681  break;
2682  }
2683  pIter(p);
2684  }
2685  p = ph;
2686  //number tmp;
2687  if(!n_IsOne(h,r->cf))
2688  {
2689  while (p!=NULL)
2690  {
2691  //d = nDiv(pGetCoeff(p),h);
2692  //tmp = nExactDiv(pGetCoeff(p),h);
2693  //if (!nEqual(d,tmp))
2694  //{
2695  // StringSetS("** div0:");nWrite(pGetCoeff(p));StringAppendS("/");
2696  // nWrite(h);StringAppendS("=");nWrite(d);StringAppendS(" int:");
2697  // nWrite(tmp);Print(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
2698  //}
2699  //nDelete(&tmp);
2700  d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2701  p_SetCoeff(p,d,r->cf);
2702  pIter(p);
2703  }
2704  }
2705  n_Delete(&h,r->cf);
2706  //if ( (n_GetChar(r) == 1) || (n_GetChar(r) < 0) ) /* Q[a],Q(a),Zp[a],Z/p(a) */
2707  //{
2708  // singclap_divide_content(ph);
2709  // if(!n_GreaterZero(pGetCoeff(ph),r)) ph = p_Neg(ph,r);
2710  //}
2711  }
2712 }
2713 #endif
2714 /* ---------------------------------------------------------------------------*/
2715 /* cleardenom suff */
2716 poly p_Cleardenom(poly p, const ring r)
2717 {
2718  if( p == NULL )
2719  return NULL;
2720 
2721  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
2722 
2723 #if CLEARENUMERATORS
2724  if( 0 )
2725  {
2726  CPolyCoeffsEnumerator itr(p);
2727 
2728  n_ClearDenominators(itr, C);
2729 
2730  n_ClearContent(itr, C); // divide out the content
2731 
2732  p_Test(p, r); n_Test(pGetCoeff(p), C);
2733  assume(n_GreaterZero(pGetCoeff(p), C)); // ??
2734 // if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2735 
2736  return p;
2737  }
2738 #endif
2739 
2740 
2741  number d, h;
2742 
2743  if (rField_is_Ring(r))
2744  {
2745  p_Content(p,r);
2746  if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2747  return p;
2748  }
2749 
2751  {
2752  if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2753  return p;
2754  }
2755 
2756  assume(p != NULL);
2757 
2758  if(pNext(p)==NULL)
2759  {
2760  if (!TEST_OPT_CONTENTSB
2761  && !rField_is_Ring(r))
2762  p_SetCoeff(p,n_Init(1,r->cf),r);
2763  else if(!n_GreaterZero(pGetCoeff(p),C))
2764  p = p_Neg(p,r);
2765  return p;
2766  }
2767 
2768  assume(pNext(p)!=NULL);
2769  poly start=p;
2770 
2771 #if 0 && CLEARENUMERATORS
2772 //CF: does not seem to work that well..
2773 
2774  if( nCoeff_is_Q(C) || nCoeff_is_Q_a(C) )
2775  {
2776  CPolyCoeffsEnumerator itr(p);
2777 
2778  n_ClearDenominators(itr, C);
2779 
2780  n_ClearContent(itr, C); // divide out the content
2781 
2782  p_Test(p, r); n_Test(pGetCoeff(p), C);
2783  assume(n_GreaterZero(pGetCoeff(p), C)); // ??
2784 // if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2785 
2786  return start;
2787  }
2788 #endif
2789 
2790  if(1)
2791  {
2792  h = n_Init(1,r->cf);
2793  while (p!=NULL)
2794  {
2795  n_Normalize(pGetCoeff(p),r->cf);
2796  d=n_NormalizeHelper(h,pGetCoeff(p),r->cf);
2797  n_Delete(&h,r->cf);
2798  h=d;
2799  pIter(p);
2800  }
2801  /* contains the 1/lcm of all denominators */
2802  if(!n_IsOne(h,r->cf))
2803  {
2804  p = start;
2805  while (p!=NULL)
2806  {
2807  /* should be: // NOTE: don't use ->coef!!!!
2808  * number hh;
2809  * nGetDenom(p->coef,&hh);
2810  * nMult(&h,&hh,&d);
2811  * nNormalize(d);
2812  * nDelete(&hh);
2813  * nMult(d,p->coef,&hh);
2814  * nDelete(&d);
2815  * nDelete(&(p->coef));
2816  * p->coef =hh;
2817  */
2818  d=n_Mult(h,pGetCoeff(p),r->cf);
2819  n_Normalize(d,r->cf);
2820  p_SetCoeff(p,d,r);
2821  pIter(p);
2822  }
2823  n_Delete(&h,r->cf);
2824  }
2825  n_Delete(&h,r->cf);
2826  p=start;
2827 
2828  p_Content(p,r);
2829 #ifdef HAVE_RATGRING
2830  if (rIsRatGRing(r))
2831  {
2832  /* quick unit detection in the rational case is done in gr_nc_bba */
2833  p_ContentRat(p, r);
2834  start=p;
2835  }
2836 #endif
2837  }
2838 
2839  if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2840 
2841  return start;
2842 }
2843 
2844 void p_Cleardenom_n(poly ph,const ring r,number &c)
2845 {
2846  const coeffs C = r->cf;
2847  number d, h;
2848 
2849  assume( ph != NULL );
2850 
2851  poly p = ph;
2852 
2853 #if CLEARENUMERATORS
2854  if( 0 )
2855  {
2856  CPolyCoeffsEnumerator itr(ph);
2857 
2858  n_ClearDenominators(itr, d, C); // multiply with common denom. d
2859  n_ClearContent(itr, h, C); // divide by the content h
2860 
2861  c = n_Div(d, h, C); // d/h
2862 
2863  n_Delete(&d, C);
2864  n_Delete(&h, C);
2865 
2866  n_Test(c, C);
2867 
2868  p_Test(ph, r); n_Test(pGetCoeff(ph), C);
2869  assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
2870 /*
2871  if(!n_GreaterZero(pGetCoeff(ph),C))
2872  {
2873  ph = p_Neg(ph,r);
2874  c = n_InpNeg(c, C);
2875  }
2876 */
2877  return;
2878  }
2879 #endif
2880 
2881 
2882  if( pNext(p) == NULL )
2883  {
2884  c=n_Invers(pGetCoeff(p), C);
2885  p_SetCoeff(p, n_Init(1, C), r);
2886 
2887  assume( n_GreaterZero(pGetCoeff(ph),C) );
2888  if(!n_GreaterZero(pGetCoeff(ph),C))
2889  {
2890  ph = p_Neg(ph,r);
2891  c = n_InpNeg(c, C);
2892  }
2893 
2894  return;
2895  }
2896 
2897  assume( pNext(p) != NULL );
2898 
2899 #if CLEARENUMERATORS
2900  if( nCoeff_is_Q(C) || nCoeff_is_Q_a(C) )
2901  {
2902  CPolyCoeffsEnumerator itr(ph);
2903 
2904  n_ClearDenominators(itr, d, C); // multiply with common denom. d
2905  n_ClearContent(itr, h, C); // divide by the content h
2906 
2907  c = n_Div(d, h, C); // d/h
2908 
2909  n_Delete(&d, C);
2910  n_Delete(&h, C);
2911 
2912  n_Test(c, C);
2913 
2914  p_Test(ph, r); n_Test(pGetCoeff(ph), C);
2915  assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
2916 /*
2917  if(!n_GreaterZero(pGetCoeff(ph),C))
2918  {
2919  ph = p_Neg(ph,r);
2920  c = n_InpNeg(c, C);
2921  }
2922 */
2923  return;
2924  }
2925 #endif
2926 
2927 
2928 
2929 
2930  if(1)
2931  {
2932  h = n_Init(1,r->cf);
2933  while (p!=NULL)
2934  {
2935  n_Normalize(pGetCoeff(p),r->cf);
2936  d=n_NormalizeHelper(h,pGetCoeff(p),r->cf);
2937  n_Delete(&h,r->cf);
2938  h=d;
2939  pIter(p);
2940  }
2941  c=h;
2942  /* contains the 1/lcm of all denominators */
2943  if(!n_IsOne(h,r->cf))
2944  {
2945  p = ph;
2946  while (p!=NULL)
2947  {
2948  /* should be: // NOTE: don't use ->coef!!!!
2949  * number hh;
2950  * nGetDenom(p->coef,&hh);
2951  * nMult(&h,&hh,&d);
2952  * nNormalize(d);
2953  * nDelete(&hh);
2954  * nMult(d,p->coef,&hh);
2955  * nDelete(&d);
2956  * nDelete(&(p->coef));
2957  * p->coef =hh;
2958  */
2959  d=n_Mult(h,pGetCoeff(p),r->cf);
2960  n_Normalize(d,r->cf);
2961  p_SetCoeff(p,d,r);
2962  pIter(p);
2963  }
2964  if (rField_is_Q_a(r))
2965  {
2966  loop
2967  {
2968  h = n_Init(1,r->cf);
2969  p=ph;
2970  while (p!=NULL)
2971  {
2972  d=n_NormalizeHelper(h,pGetCoeff(p),r->cf);
2973  n_Delete(&h,r->cf);
2974  h=d;
2975  pIter(p);
2976  }
2977  /* contains the 1/lcm of all denominators */
2978  if(!n_IsOne(h,r->cf))
2979  {
2980  p = ph;
2981  while (p!=NULL)
2982  {
2983  /* should be: // NOTE: don't use ->coef!!!!
2984  * number hh;
2985  * nGetDenom(p->coef,&hh);
2986  * nMult(&h,&hh,&d);
2987  * nNormalize(d);
2988  * nDelete(&hh);
2989  * nMult(d,p->coef,&hh);
2990  * nDelete(&d);
2991  * nDelete(&(p->coef));
2992  * p->coef =hh;
2993  */
2994  d=n_Mult(h,pGetCoeff(p),r->cf);
2995  n_Normalize(d,r->cf);
2996  p_SetCoeff(p,d,r);
2997  pIter(p);
2998  }
2999  number t=n_Mult(c,h,r->cf);
3000  n_Delete(&c,r->cf);
3001  c=t;
3002  }
3003  else
3004  {
3005  break;
3006  }
3007  n_Delete(&h,r->cf);
3008  }
3009  }
3010  }
3011  }
3012 
3013  if(!n_GreaterZero(pGetCoeff(ph),C))
3014  {
3015  ph = p_Neg(ph,r);
3016  c = n_InpNeg(c, C);
3017  }
3018 
3019 }
3020 
3021  // normalization: for poly over Q: make poly primitive, integral
3022  // Qa make poly integral with leading
3023  // coefficient minimal in N
3024  // Q(t) make poly primitive, integral
3025 
3026 void p_ProjectiveUnique(poly ph, const ring r)
3027 {
3028  if( ph == NULL )
3029  return;
3030 
3031  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
3032 
3033  number h;
3034  poly p;
3035 
3036  if (rField_is_Ring(r))
3037  {
3038  p_Content(ph,r);
3039  if(!n_GreaterZero(pGetCoeff(ph),C)) ph = p_Neg(ph,r);
3040  assume( n_GreaterZero(pGetCoeff(ph),C) );
3041  return;
3042  }
3043 
3045  {
3046  assume( n_GreaterZero(pGetCoeff(ph),C) );
3047  if(!n_GreaterZero(pGetCoeff(ph),C)) ph = p_Neg(ph,r);
3048  return;
3049  }
3050  p = ph;
3051 
3052  assume(p != NULL);
3053 
3054  if(pNext(p)==NULL) // a monomial
3055  {
3056  p_SetCoeff(p, n_Init(1, C), r);
3057  return;
3058  }
3059 
3060  assume(pNext(p)!=NULL);
3061 
3062  if(!rField_is_Q(r) && !nCoeff_is_transExt(C))
3063  {
3064  h = p_GetCoeff(p, C);
3065  number hInv = n_Invers(h, C);
3066  pIter(p);
3067  while (p!=NULL)
3068  {
3069  p_SetCoeff(p, n_Mult(p_GetCoeff(p, C), hInv, C), r);
3070  pIter(p);
3071  }
3072  n_Delete(&hInv, C);
3073  p = ph;
3074  p_SetCoeff(p, n_Init(1, C), r);
3075  }
3076 
3077  p_Cleardenom(ph, r); //performs also a p_Content
3078 
3079 
3080  /* normalize ph over a transcendental extension s.t.
3081  lead (ph) is > 0 if extRing->cf == Q
3082  or lead (ph) is monic if extRing->cf == Zp*/
3083  if (nCoeff_is_transExt(C))
3084  {
3085  p= ph;
3086  h= p_GetCoeff (p, C);
3087  fraction f = (fraction) h;
3088  number n=p_GetCoeff (NUM (f),C->extRing->cf);
3089  if (rField_is_Q (C->extRing))
3090  {
3091  if (!n_GreaterZero(n,C->extRing->cf))
3092  {
3093  p=p_Neg (p,r);
3094  }
3095  }
3096  else if (rField_is_Zp(C->extRing))
3097  {
3098  if (!n_IsOne (n, C->extRing->cf))
3099  {
3100  n=n_Invers (n,C->extRing->cf);
3101  nMapFunc nMap;
3102  nMap= n_SetMap (C->extRing->cf, C);
3103  number ninv= nMap (n,C->extRing->cf, C);
3104  p=p_Mult_nn (p, ninv, r);
3105  n_Delete (&ninv, C);
3106  n_Delete (&n, C->extRing->cf);
3107  }
3108  }
3109  p= ph;
3110  }
3111 
3112  return;
3113 }
3114 
3115 #if 0 /*unused*/
3116 number p_GetAllDenom(poly ph, const ring r)
3117 {
3118  number d=n_Init(1,r->cf);
3119  poly p = ph;
3120 
3121  while (p!=NULL)
3122  {
3123  number h=n_GetDenom(pGetCoeff(p),r->cf);
3124  if (!n_IsOne(h,r->cf))
3125  {
3126  number dd=n_Mult(d,h,r->cf);
3127  n_Delete(&d,r->cf);
3128  d=dd;
3129  }
3130  n_Delete(&h,r->cf);
3131  pIter(p);
3132  }
3133  return d;
3134 }
3135 #endif
3136 
3137 int p_Size(poly p, const ring r)
3138 {
3139  int count = 0;
3140  if (r->cf->has_simple_Alloc)
3141  return pLength(p);
3142  while ( p != NULL )
3143  {
3144  count+= n_Size( pGetCoeff( p ), r->cf );
3145  pIter( p );
3146  }
3147  return count;
3148 }
3149 
3150 /*2
3151 *make p homogeneous by multiplying the monomials by powers of x_varnum
3152 *assume: deg(var(varnum))==1
3153 */
3154 poly p_Homogen (poly p, int varnum, const ring r)
3155 {
3156  pFDegProc deg;
3157  if (r->pLexOrder && (r->order[0]==ringorder_lp))
3158  deg=p_Totaldegree;
3159  else
3160  deg=r->pFDeg;
3161 
3162  poly q=NULL, qn;
3163  int o,ii;
3164  sBucket_pt bp;
3165 
3166  if (p!=NULL)
3167  {
3168  if ((varnum < 1) || (varnum > rVar(r)))
3169  {
3170  return NULL;
3171  }
3172  o=deg(p,r);
3173  q=pNext(p);
3174  while (q != NULL)
3175  {
3176  ii=deg(q,r);
3177  if (ii>o) o=ii;
3178  pIter(q);
3179  }
3180  q = p_Copy(p,r);
3181  bp = sBucketCreate(r);
3182  while (q != NULL)
3183  {
3184  ii = o-deg(q,r);
3185  if (ii!=0)
3186  {
3187  p_AddExp(q,varnum, (long)ii,r);
3188  p_Setm(q,r);
3189  }
3190  qn = pNext(q);
3191  pNext(q) = NULL;
3192  sBucket_Add_p(bp, q, 1);
3193  q = qn;
3194  }
3195  sBucketDestroyAdd(bp, &q, &ii);
3196  }
3197  return q;
3198 }
3199 
3200 /*2
3201 *tests if p is homogeneous with respect to the actual weigths
3202 */
3204 {
3205  poly qp=p;
3206  int o;
3207 
3208  if ((p == NULL) || (pNext(p) == NULL)) return TRUE;
3209  pFDegProc d;
3210  if (r->pLexOrder && (r->order[0]==ringorder_lp))
3211  d=p_Totaldegree;
3212  else
3213  d=r->pFDeg;
3214  o = d(p,r);
3215  do
3216  {
3217  if (d(qp,r) != o) return FALSE;
3218  pIter(qp);
3219  }
3220  while (qp != NULL);
3221  return TRUE;
3222 }
3223 
3224 /*----------utilities for syzygies--------------*/
3225 BOOLEAN p_VectorHasUnitB(poly p, int * k, const ring r)
3226 {
3227  poly q=p,qq;
3228  int i;
3229 
3230  while (q!=NULL)
3231  {
3232  if (p_LmIsConstantComp(q,r))
3233  {
3234  i = p_GetComp(q,r);
3235  qq = p;
3236  while ((qq != q) && (p_GetComp(qq,r) != i)) pIter(qq);
3237  if (qq == q)
3238  {
3239  *k = i;
3240  return TRUE;
3241  }
3242  else
3243  pIter(q);
3244  }
3245  else pIter(q);
3246  }
3247  return FALSE;
3248 }
3249 
3250 void p_VectorHasUnit(poly p, int * k, int * len, const ring r)
3251 {
3252  poly q=p,qq;
3253  int i,j=0;
3254 
3255  *len = 0;
3256  while (q!=NULL)
3257  {
3258  if (p_LmIsConstantComp(q,r))
3259  {
3260  i = p_GetComp(q,r);
3261  qq = p;
3262  while ((qq != q) && (p_GetComp(qq,r) != i)) pIter(qq);
3263  if (qq == q)
3264  {
3265  j = 0;
3266  while (qq!=NULL)
3267  {
3268  if (p_GetComp(qq,r)==i) j++;
3269  pIter(qq);
3270  }
3271  if ((*len == 0) || (j<*len))
3272  {
3273  *len = j;
3274  *k = i;
3275  }
3276  }
3277  }
3278  pIter(q);
3279  }
3280 }
3281 
3282 poly p_TakeOutComp1(poly * p, int k, const ring r)
3283 {
3284  poly q = *p;
3285 
3286  if (q==NULL) return NULL;
3287 
3288  poly qq=NULL,result = NULL;
3289 
3290  if (p_GetComp(q,r)==k)
3291  {
3292  result = q; /* *p */
3293  while ((q!=NULL) && (p_GetComp(q,r)==k))
3294  {
3295  p_SetComp(q,0,r);
3296  p_SetmComp(q,r);
3297  qq = q;
3298  pIter(q);
3299  }
3300  *p = q;
3301  pNext(qq) = NULL;
3302  }
3303  if (q==NULL) return result;
3304 // if (pGetComp(q) > k) pGetComp(q)--;
3305  while (pNext(q)!=NULL)
3306  {
3307  if (p_GetComp(pNext(q),r)==k)
3308  {
3309  if (result==NULL)
3310  {
3311  result = pNext(q);
3312  qq = result;
3313  }
3314  else
3315  {
3316  pNext(qq) = pNext(q);
3317  pIter(qq);
3318  }
3319  pNext(q) = pNext(pNext(q));
3320  pNext(qq) =NULL;
3321  p_SetComp(qq,0,r);
3322  p_SetmComp(qq,r);
3323  }
3324  else
3325  {
3326  pIter(q);
3327 // if (pGetComp(q) > k) pGetComp(q)--;
3328  }
3329  }
3330  return result;
3331 }
3332 
3333 poly p_TakeOutComp(poly * p, int k, const ring r)
3334 {
3335  poly q = *p,qq=NULL,result = NULL;
3336 
3337  if (q==NULL) return NULL;
3338  BOOLEAN use_setmcomp=rOrd_SetCompRequiresSetm(r);
3339  if (p_GetComp(q,r)==k)
3340  {
3341  result = q;
3342  do
3343  {
3344  p_SetComp(q,0,r);
3345  if (use_setmcomp) p_SetmComp(q,r);
3346  qq = q;
3347  pIter(q);
3348  }
3349  while ((q!=NULL) && (p_GetComp(q,r)==k));
3350  *p = q;
3351  pNext(qq) = NULL;
3352  }
3353  if (q==NULL) return result;
3354  if (p_GetComp(q,r) > k)
3355  {
3356  p_SubComp(q,1,r);
3357  if (use_setmcomp) p_SetmComp(q,r);
3358  }
3359  poly pNext_q;
3360  while ((pNext_q=pNext(q))!=NULL)
3361  {
3362  if (p_GetComp(pNext_q,r)==k)
3363  {
3364  if (result==NULL)
3365  {
3366  result = pNext_q;
3367  qq = result;
3368  }
3369  else
3370  {
3371  pNext(qq) = pNext_q;
3372  pIter(qq);
3373  }
3374  pNext(q) = pNext(pNext_q);
3375  pNext(qq) =NULL;
3376  p_SetComp(qq,0,r);
3377  if (use_setmcomp) p_SetmComp(qq,r);
3378  }
3379  else
3380  {
3381  /*pIter(q);*/ q=pNext_q;
3382  if (p_GetComp(q,r) > k)
3383  {
3384  p_SubComp(q,1,r);
3385  if (use_setmcomp) p_SetmComp(q,r);
3386  }
3387  }
3388  }
3389  return result;
3390 }
3391 
3392 // Splits *p into two polys: *q which consists of all monoms with
3393 // component == comp and *p of all other monoms *lq == pLength(*q)
3394 void p_TakeOutComp(poly *r_p, long comp, poly *r_q, int *lq, const ring r)
3395 {
3396  spolyrec pp, qq;
3397  poly p, q, p_prev;
3398  int l = 0;
3399 
3400 #ifdef HAVE_ASSUME
3401  int lp = pLength(*r_p);
3402 #endif
3403 
3404  pNext(&pp) = *r_p;
3405  p = *r_p;
3406  p_prev = &pp;
3407  q = &qq;
3408 
3409  while(p != NULL)
3410  {
3411  while (p_GetComp(p,r) == comp)
3412  {
3413  pNext(q) = p;
3414  pIter(q);
3415  p_SetComp(p, 0,r);
3416  p_SetmComp(p,r);
3417  pIter(p);
3418  l++;
3419  if (p == NULL)
3420  {
3421  pNext(p_prev) = NULL;
3422  goto Finish;
3423  }
3424  }
3425  pNext(p_prev) = p;
3426  p_prev = p;
3427  pIter(p);
3428  }
3429 
3430  Finish:
3431  pNext(q) = NULL;
3432  *r_p = pNext(&pp);
3433  *r_q = pNext(&qq);
3434  *lq = l;
3435 #ifdef HAVE_ASSUME
3436  assume(pLength(*r_p) + pLength(*r_q) == lp);
3437 #endif
3438  p_Test(*r_p,r);
3439  p_Test(*r_q,r);
3440 }
3441 
3442 void p_DeleteComp(poly * p,int k, const ring r)
3443 {
3444  poly q;
3445 
3446  while ((*p!=NULL) && (p_GetComp(*p,r)==k)) p_LmDelete(p,r);
3447  if (*p==NULL) return;
3448  q = *p;
3449  if (p_GetComp(q,r)>k)
3450  {
3451  p_SubComp(q,1,r);
3452  p_SetmComp(q,r);
3453  }
3454  while (pNext(q)!=NULL)
3455  {
3456  if (p_GetComp(pNext(q),r)==k)
3457  p_LmDelete(&(pNext(q)),r);
3458  else
3459  {
3460  pIter(q);
3461  if (p_GetComp(q,r)>k)
3462  {
3463  p_SubComp(q,1,r);
3464  p_SetmComp(q,r);
3465  }
3466  }
3467  }
3468 }
3469 
3470 /*2
3471 * convert a vector to a set of polys,
3472 * allocates the polyset, (entries 0..(*len)-1)
3473 * the vector will not be changed
3474 */
3475 void p_Vec2Polys(poly v, poly* *p, int *len, const ring r)
3476 {
3477  poly h;
3478  int k;
3479 
3480  *len=p_MaxComp(v,r);
3481  if (*len==0) *len=1;
3482  *p=(poly*)omAlloc0((*len)*sizeof(poly));
3483  while (v!=NULL)
3484  {
3485  h=p_Head(v,r);
3486  k=p_GetComp(h,r);
3487  p_SetComp(h,0,r);
3488  (*p)[k-1]=p_Add_q((*p)[k-1],h,r);
3489  pIter(v);
3490  }
3491 }
3492 
3493 //
3494 // resets the pFDeg and pLDeg: if pLDeg is not given, it is
3495 // set to currRing->pLDegOrig, i.e. to the respective LDegProc which
3496 // only uses pFDeg (and not p_Deg, or pTotalDegree, etc)
3497 void pSetDegProcs(ring r, pFDegProc new_FDeg, pLDegProc new_lDeg)
3498 {
3499  assume(new_FDeg != NULL);
3500  r->pFDeg = new_FDeg;
3501 
3502  if (new_lDeg == NULL)
3503  new_lDeg = r->pLDegOrig;
3504 
3505  r->pLDeg = new_lDeg;
3506 }
3507 
3508 // restores pFDeg and pLDeg:
3509 void pRestoreDegProcs(ring r, pFDegProc old_FDeg, pLDegProc old_lDeg)
3510 {
3511  assume(old_FDeg != NULL && old_lDeg != NULL);
3512  r->pFDeg = old_FDeg;
3513  r->pLDeg = old_lDeg;
3514 }
3515 
3516 /*-------- several access procedures to monomials -------------------- */
3517 /*
3518 * the module weights for std
3519 */
3523 
3524 static long pModDeg(poly p, ring r)
3525 {
3526  long d=pOldFDeg(p, r);
3527  int c=p_GetComp(p, r);
3528  if ((c>0) && ((r->pModW)->range(c-1))) d+= (*(r->pModW))[c-1];
3529  return d;
3530  //return pOldFDeg(p, r)+(*pModW)[p_GetComp(p, r)-1];
3531 }
3532 
3533 void p_SetModDeg(intvec *w, ring r)
3534 {
3535  if (w!=NULL)
3536  {
3537  r->pModW = w;
3538  pOldFDeg = r->pFDeg;
3539  pOldLDeg = r->pLDeg;
3540  pOldLexOrder = r->pLexOrder;
3541  pSetDegProcs(r,pModDeg);
3542  r->pLexOrder = TRUE;
3543  }
3544  else
3545  {
3546  r->pModW = NULL;
3547  pRestoreDegProcs(r,pOldFDeg, pOldLDeg);
3548  r->pLexOrder = pOldLexOrder;
3549  }
3550 }
3551 
3552 /*2
3553 * handle memory request for sets of polynomials (ideals)
3554 * l is the length of *p, increment is the difference (may be negative)
3555 */
3556 void pEnlargeSet(poly* *p, int l, int increment)
3557 {
3558  poly* h;
3559 
3560  if (*p==NULL)
3561  {
3562  if (increment==0) return;
3563  h=(poly*)omAlloc0(increment*sizeof(poly));
3564  }
3565  else
3566  {
3567  h=(poly*)omReallocSize((poly*)*p,l*sizeof(poly),(l+increment)*sizeof(poly));
3568  if (increment>0)
3569  {
3570  //for (i=l; i<l+increment; i++)
3571  // h[i]=NULL;
3572  memset(&(h[l]),0,increment*sizeof(poly));
3573  }
3574  }
3575  *p=h;
3576 }
3577 
3578 /*2
3579 *divides p1 by its leading coefficient
3580 */
3581 void p_Norm(poly p1, const ring r)
3582 {
3583  if (rField_is_Ring(r))
3584  {
3585  if (!n_IsUnit(pGetCoeff(p1), r->cf)) return;
3586  // Werror("p_Norm not possible in the case of coefficient rings.");
3587  }
3588  else if (p1!=NULL)
3589  {
3590  if (pNext(p1)==NULL)
3591  {
3592  p_SetCoeff(p1,n_Init(1,r->cf),r);
3593  return;
3594  }
3595  poly h;
3596  if (!n_IsOne(pGetCoeff(p1),r->cf))
3597  {
3598  number k, c;
3599  n_Normalize(pGetCoeff(p1),r->cf);
3600  k = pGetCoeff(p1);
3601  c = n_Init(1,r->cf);
3602  pSetCoeff0(p1,c);
3603  h = pNext(p1);
3604  while (h!=NULL)
3605  {
3606  c=n_Div(pGetCoeff(h),k,r->cf);
3607  // no need to normalize: Z/p, R
3608  // normalize already in nDiv: Q_a, Z/p_a
3609  // remains: Q
3610  if (rField_is_Q(r) && (!n_IsOne(c,r->cf))) n_Normalize(c,r->cf);
3611  p_SetCoeff(h,c,r);
3612  pIter(h);
3613  }
3614  n_Delete(&k,r->cf);
3615  }
3616  else
3617  {
3618  //if (r->cf->cfNormalize != nDummy2) //TODO: OPTIMIZE
3619  {
3620  h = pNext(p1);
3621  while (h!=NULL)
3622  {
3623  n_Normalize(pGetCoeff(h),r->cf);
3624  pIter(h);
3625  }
3626  }
3627  }
3628  }
3629 }
3630 
3631 /*2
3632 *normalize all coefficients
3633 */
3634 void p_Normalize(poly p,const ring r)
3635 {
3636  if (rField_has_simple_inverse(r)) return; /* Z/p, GF(p,n), R, long R/C */
3637  while (p!=NULL)
3638  {
3639  // no test befor n_Normalize: n_Normalize should fix problems
3640  n_Normalize(pGetCoeff(p),r->cf);
3641  pIter(p);
3642  }
3643 }
3644 
3645 // splits p into polys with Exp(n) == 0 and Exp(n) != 0
3646 // Poly with Exp(n) != 0 is reversed
3647 static void p_SplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero, const ring r)
3648 {
3649  if (p == NULL)
3650  {
3651  *non_zero = NULL;
3652  *zero = NULL;
3653  return;
3654  }
3655  spolyrec sz;
3656  poly z, n_z, next;
3657  z = &sz;
3658  n_z = NULL;
3659 
3660  while(p != NULL)
3661  {
3662  next = pNext(p);
3663  if (p_GetExp(p, n,r) == 0)
3664  {
3665  pNext(z) = p;
3666  pIter(z);
3667  }
3668  else
3669  {
3670  pNext(p) = n_z;
3671  n_z = p;
3672  }
3673  p = next;
3674  }
3675  pNext(z) = NULL;
3676  *zero = pNext(&sz);
3677  *non_zero = n_z;
3678 }
3679 /*3
3680 * substitute the n-th variable by 1 in p
3681 * destroy p
3682 */
3683 static poly p_Subst1 (poly p,int n, const ring r)
3684 {
3685  poly qq=NULL, result = NULL;
3686  poly zero=NULL, non_zero=NULL;
3687 
3688  // reverse, so that add is likely to be linear
3689  p_SplitAndReversePoly(p, n, &non_zero, &zero,r);
3690 
3691  while (non_zero != NULL)
3692  {
3693  assume(p_GetExp(non_zero, n,r) != 0);
3694  qq = non_zero;
3695  pIter(non_zero);
3696  qq->next = NULL;
3697  p_SetExp(qq,n,0,r);
3698  p_Setm(qq,r);
3699  result = p_Add_q(result,qq,r);
3700  }
3701  p = p_Add_q(result, zero,r);
3702  p_Test(p,r);
3703  return p;
3704 }
3705 
3706 /*3
3707 * substitute the n-th variable by number e in p
3708 * destroy p
3709 */
3710 static poly p_Subst2 (poly p,int n, number e, const ring r)
3711 {
3712  assume( ! n_IsZero(e,r->cf) );
3713  poly qq,result = NULL;
3714  number nn, nm;
3715  poly zero, non_zero;
3716 
3717  // reverse, so that add is likely to be linear
3718  p_SplitAndReversePoly(p, n, &non_zero, &zero,r);
3719 
3720  while (non_zero != NULL)
3721  {
3722  assume(p_GetExp(non_zero, n, r) != 0);
3723  qq = non_zero;
3724  pIter(non_zero);
3725  qq->next = NULL;
3726  n_Power(e, p_GetExp(qq, n, r), &nn,r->cf);
3727  nm = n_Mult(nn, pGetCoeff(qq),r->cf);
3728 #ifdef HAVE_RINGS
3729  if (n_IsZero(nm,r->cf))
3730  {
3731  p_LmFree(&qq,r);
3732  n_Delete(&nm,r->cf);
3733  }
3734  else
3735 #endif
3736  {
3737  p_SetCoeff(qq, nm,r);
3738  p_SetExp(qq, n, 0,r);
3739  p_Setm(qq,r);
3740  result = p_Add_q(result,qq,r);
3741  }
3742  n_Delete(&nn,r->cf);
3743  }
3744  p = p_Add_q(result, zero,r);
3745  p_Test(p,r);
3746  return p;
3747 }
3748 
3749 
3750 /* delete monoms whose n-th exponent is different from zero */
3751 static poly p_Subst0(poly p, int n, const ring r)
3752 {
3753  spolyrec res;
3754  poly h = &res;
3755  pNext(h) = p;
3756 
3757  while (pNext(h)!=NULL)
3758  {
3759  if (p_GetExp(pNext(h),n,r)!=0)
3760  {
3761  p_LmDelete(&pNext(h),r);
3762  }
3763  else
3764  {
3765  pIter(h);
3766  }
3767  }
3768  p_Test(pNext(&res),r);
3769  return pNext(&res);
3770 }
3771 
3772 /*2
3773 * substitute the n-th variable by e in p
3774 * destroy p
3775 */
3776 poly p_Subst(poly p, int n, poly e, const ring r)
3777 {
3778  if (e == NULL) return p_Subst0(p, n,r);
3779 
3780  if (p_IsConstant(e,r))
3781  {
3782  if (n_IsOne(pGetCoeff(e),r->cf)) return p_Subst1(p,n,r);
3783  else return p_Subst2(p, n, pGetCoeff(e),r);
3784  }
3785 
3786 #ifdef HAVE_PLURAL
3787  if (rIsPluralRing(r))
3788  {
3789  return nc_pSubst(p,n,e,r);
3790  }
3791 #endif
3792 
3793  int exponent,i;
3794  poly h, res, m;
3795  int *me,*ee;
3796  number nu,nu1;
3797 
3798  me=(int *)omAlloc((rVar(r)+1)*sizeof(int));
3799  ee=(int *)omAlloc((rVar(r)+1)*sizeof(int));
3800  if (e!=NULL) p_GetExpV(e,ee,r);
3801  res=NULL;
3802  h=p;
3803  while (h!=NULL)
3804  {
3805  if ((e!=NULL) || (p_GetExp(h,n,r)==0))
3806  {
3807  m=p_Head(h,r);
3808  p_GetExpV(m,me,r);
3809  exponent=me[n];
3810  me[n]=0;
3811  for(i=rVar(r);i>0;i--)
3812  me[i]+=exponent*ee[i];
3813  p_SetExpV(m,me,r);
3814  if (e!=NULL)
3815  {
3816  n_Power(pGetCoeff(e),exponent,&nu,r->cf);
3817  nu1=n_Mult(pGetCoeff(m),nu,r->cf);
3818  n_Delete(&nu,r->cf);
3819  p_SetCoeff(m,nu1,r);
3820  }
3821  res=p_Add_q(res,m,r);
3822  }
3823  p_LmDelete(&h,r);
3824  }
3825  omFreeSize((ADDRESS)me,(rVar(r)+1)*sizeof(int));
3826  omFreeSize((ADDRESS)ee,(rVar(r)+1)*sizeof(int));
3827  return res;
3828 }
3829 
3830 /*2
3831  * returns a re-ordered convertion of a number as a polynomial,
3832  * with permutation of parameters
3833  * NOTE: this only works for Frank's alg. & trans. fields
3834  */
3835 poly n_PermNumber(const number z, const int *par_perm, const int , const ring src, const ring dst)
3836 {
3837 #if 0
3838  PrintS("\nSource Ring: \n");
3839  rWrite(src);
3840 
3841  if(0)
3842  {
3843  number zz = n_Copy(z, src->cf);
3844  PrintS("z: "); n_Write(zz, src);
3845  n_Delete(&zz, src->cf);
3846  }
3847 
3848  PrintS("\nDestination Ring: \n");
3849  rWrite(dst);
3850 
3851  /*Print("\nOldPar: %d\n", OldPar);
3852  for( int i = 1; i <= OldPar; i++ )
3853  {
3854  Print("par(%d) -> par/var (%d)\n", i, par_perm[i-1]);
3855  }*/
3856 #endif
3857  if( z == NULL )
3858  return NULL;
3859 
3860  const coeffs srcCf = src->cf;
3861  assume( srcCf != NULL );
3862 
3863  assume( !nCoeff_is_GF(srcCf) );
3864  assume( src->cf->extRing!=NULL );
3865 
3866  poly zz = NULL;
3867 
3868  const ring srcExtRing = srcCf->extRing;
3869  assume( srcExtRing != NULL );
3870 
3871  const coeffs dstCf = dst->cf;
3872  assume( dstCf != NULL );
3873 
3874  if( nCoeff_is_algExt(srcCf) ) // nCoeff_is_GF(srcCf)?
3875  {
3876  zz = (poly) z;
3877  if( zz == NULL ) return NULL;
3878  }
3879  else if (nCoeff_is_transExt(srcCf))
3880  {
3881  assume( !IS0(z) );
3882 
3883  zz = NUM((fraction)z);
3884  p_Test (zz, srcExtRing);
3885 
3886  if( zz == NULL ) return NULL;
3887  if( !DENIS1((fraction)z) )
3888  {
3889  if (!p_IsConstant(DEN((fraction)z),srcExtRing))
3890  WarnS("Not defined: Cannot map a rational fraction and make a polynomial out of it! Ignoring the denumerator.");
3891  }
3892  }
3893  else
3894  {
3895  assume (FALSE);
3896  WerrorS("Number permutation is not implemented for this data yet!");
3897  return NULL;
3898  }
3899 
3900  assume( zz != NULL );
3901  p_Test (zz, srcExtRing);
3902 
3903  nMapFunc nMap = n_SetMap(srcExtRing->cf, dstCf);
3904 
3905  assume( nMap != NULL );
3906 
3907  poly qq;
3908  if ((par_perm == NULL) && (rPar(dst) != 0 && rVar (srcExtRing) > 0))
3909  {
3910  int* perm;
3911  perm=(int *)omAlloc0((rVar(srcExtRing)+1)*sizeof(int));
3912  perm[0]= 0;
3913  for(int i=si_min(rVar(srcExtRing),rPar(dst));i>0;i--)
3914  perm[i]=-i;
3915  qq = p_PermPoly(zz, perm, srcExtRing, dst, nMap, NULL, rVar(srcExtRing)-1);
3916  omFreeSize ((ADDRESS)perm, (rVar(srcExtRing)+1)*sizeof(int));
3917  }
3918  else
3919  qq = p_PermPoly(zz, par_perm-1, srcExtRing, dst, nMap, NULL, rVar (srcExtRing)-1);
3920 
3921  if(nCoeff_is_transExt(srcCf)
3922  && (!DENIS1((fraction)z))
3923  && p_IsConstant(DEN((fraction)z),srcExtRing))
3924  {
3925  number n=nMap(pGetCoeff(DEN((fraction)z)),srcExtRing->cf, dstCf);
3926  qq=p_Div_nn(qq,n,dst);
3927  n_Delete(&n,dstCf);
3928  p_Normalize(qq,dst);
3929  }
3930  p_Test (qq, dst);
3931 
3932  return qq;
3933 }
3934 
3935 
3936 /*2
3937 *returns a re-ordered copy of a polynomial, with permutation of the variables
3938 */
3939 poly p_PermPoly (poly p, const int * perm, const ring oldRing, const ring dst,
3940  nMapFunc nMap, const int *par_perm, int OldPar, BOOLEAN use_mult)
3941 {
3942 #if 0
3943  p_Test(p, oldRing);
3944  PrintS("p_PermPoly::p: "); p_Write(p, oldRing, oldRing);
3945 #endif
3946  const int OldpVariables = rVar(oldRing);
3947  poly result = NULL;
3948  poly result_last = NULL;
3949  poly aq = NULL; /* the map coefficient */
3950  poly qq; /* the mapped monomial */
3951  assume(dst != NULL);
3952  assume(dst->cf != NULL);
3953  #ifdef HAVE_PLURAL
3954  poly tmp_mm=p_One(dst);
3955  #endif
3956  while (p != NULL)
3957  {
3958  // map the coefficient
3959  if ( ((OldPar == 0) || (par_perm == NULL) || rField_is_GF(oldRing) || (nMap==ndCopyMap))
3960  && (nMap != NULL) )
3961  {
3962  qq = p_Init(dst);
3963  assume( nMap != NULL );
3964  number n = nMap(p_GetCoeff(p, oldRing), oldRing->cf, dst->cf);
3965  n_Test (n,dst->cf);
3966  if ( nCoeff_is_algExt(dst->cf) )
3967  n_Normalize(n, dst->cf);
3968  p_GetCoeff(qq, dst) = n;// Note: n can be a ZERO!!!
3969  }
3970  else
3971  {
3972  qq = p_One(dst);
3973 // aq = naPermNumber(p_GetCoeff(p, oldRing), par_perm, OldPar, oldRing); // no dst???
3974 // poly n_PermNumber(const number z, const int *par_perm, const int P, const ring src, const ring dst)
3975  aq = n_PermNumber(p_GetCoeff(p, oldRing), par_perm, OldPar, oldRing, dst);
3976  p_Test(aq, dst);
3977  if ( nCoeff_is_algExt(dst->cf) )
3978  p_Normalize(aq,dst);
3979  if (aq == NULL)
3980  p_SetCoeff(qq, n_Init(0, dst->cf),dst); // Very dirty trick!!!
3981  p_Test(aq, dst);
3982  }
3983  if (rRing_has_Comp(dst))
3984  p_SetComp(qq, p_GetComp(p, oldRing), dst);
3985  if ( n_IsZero(pGetCoeff(qq), dst->cf) )
3986  {
3987  p_LmDelete(&qq,dst);
3988  qq = NULL;
3989  }
3990  else
3991  {
3992  // map pars:
3993  int mapped_to_par = 0;
3994  for(int i = 1; i <= OldpVariables; i++)
3995  {
3996  int e = p_GetExp(p, i, oldRing);
3997  if (e != 0)
3998  {
3999  if (perm==NULL)
4000  p_SetExp(qq, i, e, dst);
4001  else if (perm[i]>0)
4002  {
4003  #ifdef HAVE_PLURAL
4004  if(use_mult)
4005  {
4006  p_SetExp(tmp_mm,perm[i],e,dst);
4007  p_Setm(tmp_mm,dst);
4008  qq=p_Mult_mm(qq,tmp_mm,dst);
4009  p_SetExp(tmp_mm,perm[i],0,dst);
4010 
4011  }
4012  else
4013  #endif
4014  p_AddExp(qq,perm[i], e/*p_GetExp( p,i,oldRing)*/, dst);
4015  }
4016  else if (perm[i]<0)
4017  {
4018  number c = p_GetCoeff(qq, dst);
4019  if (rField_is_GF(dst))
4020  {
4021  assume( dst->cf->extRing == NULL );
4022  number ee = n_Param(1, dst);
4023  number eee;
4024  n_Power(ee, e, &eee, dst->cf); //nfDelete(ee,dst);
4025  ee = n_Mult(c, eee, dst->cf);
4026  //nfDelete(c,dst);nfDelete(eee,dst);
4027  pSetCoeff0(qq,ee);
4028  }
4029  else if (nCoeff_is_Extension(dst->cf))
4030  {
4031  const int par = -perm[i];
4032  assume( par > 0 );
4033 // WarnS("longalg missing 3");
4034 #if 1
4035  const coeffs C = dst->cf;
4036  assume( C != NULL );
4037  const ring R = C->extRing;
4038  assume( R != NULL );
4039  assume( par <= rVar(R) );
4040  poly pcn; // = (number)c
4041  assume( !n_IsZero(c, C) );
4042  if( nCoeff_is_algExt(C) )
4043  pcn = (poly) c;
4044  else // nCoeff_is_transExt(C)
4045  pcn = NUM((fraction)c);
4046  if (pNext(pcn) == NULL) // c->z
4047  p_AddExp(pcn, -perm[i], e, R);
4048  else /* more difficult: we have really to multiply: */
4049  {
4050  poly mmc = p_ISet(1, R);
4051  p_SetExp(mmc, -perm[i], e, R);
4052  p_Setm(mmc, R);
4053  number nnc;
4054  // convert back to a number: number nnc = mmc;
4055  if( nCoeff_is_algExt(C) )
4056  nnc = (number) mmc;
4057  else // nCoeff_is_transExt(C)
4058  nnc = ntInit(mmc, C);
4059  p_GetCoeff(qq, dst) = n_Mult((number)c, nnc, C);
4060  n_Delete((number *)&c, C);
4061  n_Delete((number *)&nnc, C);
4062  }
4063  mapped_to_par=1;
4064 #endif
4065  }
4066  }
4067  else
4068  {
4069  /* this variable maps to 0 !*/
4070  p_LmDelete(&qq, dst);
4071  break;
4072  }
4073  }
4074  }
4075  if ( mapped_to_par && (qq!= NULL) && nCoeff_is_algExt(dst->cf) )
4076  {
4077  number n = p_GetCoeff(qq, dst);
4078  n_Normalize(n, dst->cf);
4079  p_GetCoeff(qq, dst) = n;
4080  }
4081  }
4082  pIter(p);
4083 
4084 #if 0
4085  p_Test(aq,dst);
4086  PrintS("aq: "); p_Write(aq, dst, dst);
4087 #endif
4088 
4089 
4090 #if 1
4091  if (qq!=NULL)
4092  {
4093  p_Setm(qq,dst);
4094 
4095  p_Test(aq,dst);
4096  p_Test(qq,dst);
4097 
4098 #if 0
4099  PrintS("qq: "); p_Write(qq, dst, dst);
4100 #endif
4101 
4102  if (aq!=NULL)
4103  qq=p_Mult_q(aq,qq,dst);
4104  aq = qq;
4105  while (pNext(aq) != NULL) pIter(aq);
4106  if (result_last==NULL)
4107  {
4108  result=qq;
4109  }
4110  else
4111  {
4112  pNext(result_last)=qq;
4113  }
4114  result_last=aq;
4115  aq = NULL;
4116  }
4117  else if (aq!=NULL)
4118  {
4119  p_Delete(&aq,dst);
4120  }
4121  }
4122  result=p_SortAdd(result,dst);
4123 #else
4124  // if (qq!=NULL)
4125  // {
4126  // pSetm(qq);
4127  // pTest(qq);
4128  // pTest(aq);
4129  // if (aq!=NULL) qq=pMult(aq,qq);
4130  // aq = qq;
4131  // while (pNext(aq) != NULL) pIter(aq);
4132  // pNext(aq) = result;
4133  // aq = NULL;
4134  // result = qq;
4135  // }
4136  // else if (aq!=NULL)
4137  // {
4138  // pDelete(&aq);
4139  // }
4140  //}
4141  //p = result;
4142  //result = NULL;
4143  //while (p != NULL)
4144  //{
4145  // qq = p;
4146  // pIter(p);
4147  // qq->next = NULL;
4148  // result = pAdd(result, qq);
4149  //}
4150 #endif
4151  p_Test(result,dst);
4152 #if 0
4153  p_Test(result,dst);
4154  PrintS("result: "); p_Write(result,dst,dst);
4155 #endif
4156  #ifdef HAVE_PLURAL
4157  p_LmDelete(&tmp_mm,dst);
4158  #endif
4159  return result;
4160 }
4161 /**************************************************************
4162  *
4163  * Jet
4164  *
4165  **************************************************************/
4166 
4167 poly pp_Jet(poly p, int m, const ring R)
4168 {
4169  poly r=NULL;
4170  poly t=NULL;
4171 
4172  while (p!=NULL)
4173  {
4174  if (p_Totaldegree(p,R)<=m)
4175  {
4176  if (r==NULL)
4177  r=p_Head(p,R);
4178  else
4179  if (t==NULL)
4180  {
4181  pNext(r)=p_Head(p,R);
4182  t=pNext(r);
4183  }
4184  else
4185  {
4186  pNext(t)=p_Head(p,R);
4187  pIter(t);
4188  }
4189  }
4190  pIter(p);
4191  }
4192  return r;
4193 }
4194 
4195 poly p_Jet(poly p, int m,const ring R)
4196 {
4197  while((p!=NULL) && (p_Totaldegree(p,R)>m)) p_LmDelete(&p,R);
4198  if (p==NULL) return NULL;
4199  poly r=p;
4200  while (pNext(p)!=NULL)
4201  {
4202  if (p_Totaldegree(pNext(p),R)>m)
4203  {
4204  p_LmDelete(&pNext(p),R);
4205  }
4206  else
4207  pIter(p);
4208  }
4209  return r;
4210 }
4211 
4212 poly pp_JetW(poly p, int m, short *w, const ring R)
4213 {
4214  poly r=NULL;
4215  poly t=NULL;
4216  while (p!=NULL)
4217  {
4218  if (totaldegreeWecart_IV(p,R,w)<=m)
4219  {
4220  if (r==NULL)
4221  r=p_Head(p,R);
4222  else
4223  if (t==NULL)
4224  {
4225  pNext(r)=p_Head(p,R);
4226  t=pNext(r);
4227  }
4228  else
4229  {
4230  pNext(t)=p_Head(p,R);
4231  pIter(t);
4232  }
4233  }
4234  pIter(p);
4235  }
4236  return r;
4237 }
4238 
4239 poly p_JetW(poly p, int m, short *w, const ring R)
4240 {
4241  while((p!=NULL) && (totaldegreeWecart_IV(p,R,w)>m)) p_LmDelete(&p,R);
4242  if (p==NULL) return NULL;
4243  poly r=p;
4244  while (pNext(p)!=NULL)
4245  {
4246  if (totaldegreeWecart_IV(pNext(p),R,w)>m)
4247  {
4248  p_LmDelete(&pNext(p),R);
4249  }
4250  else
4251  pIter(p);
4252  }
4253  return r;
4254 }
4255 
4256 /*************************************************************/
4257 int p_MinDeg(poly p,intvec *w, const ring R)
4258 {
4259  if(p==NULL)
4260  return -1;
4261  int d=-1;
4262  while(p!=NULL)
4263  {
4264  int d0=0;
4265  for(int j=0;j<rVar(R);j++)
4266  if(w==NULL||j>=w->length())
4267  d0+=p_GetExp(p,j+1,R);
4268  else
4269  d0+=(*w)[j]*p_GetExp(p,j+1,R);
4270  if(d0<d||d==-1)
4271  d=d0;
4272  pIter(p);
4273  }
4274  return d;
4275 }
4276 
4277 /***************************************************************/
4278 
4279 poly p_Series(int n,poly p,poly u, intvec *w, const ring R)
4280 {
4281  short *ww=iv2array(w,R);
4282  if(p!=NULL)
4283  {
4284  if(u==NULL)
4285  p=p_JetW(p,n,ww,R);
4286  else
4287  p=p_JetW(p_Mult_q(p,p_Invers(n-p_MinDeg(p,w,R),u,w,R),R),n,ww,R);
4288  }
4289  omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(short));
4290  return p;
4291 }
4292 
4293 poly p_Invers(int n,poly u,intvec *w, const ring R)
4294 {
4295  if(n<0)
4296  return NULL;
4297  number u0=n_Invers(pGetCoeff(u),R->cf);
4298  poly v=p_NSet(u0,R);
4299  if(n==0)
4300  return v;
4301  short *ww=iv2array(w,R);
4302  poly u1=p_JetW(p_Sub(p_One(R),p_Mult_nn(u,u0,R),R),n,ww,R);
4303  if(u1==NULL)
4304  {
4305  omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(short));
4306  return v;
4307  }
4308  poly v1=p_Mult_nn(p_Copy(u1,R),u0,R);
4309  v=p_Add_q(v,p_Copy(v1,R),R);
4310  for(int i=n/p_MinDeg(u1,w,R);i>1;i--)
4311  {
4312  v1=p_JetW(p_Mult_q(v1,p_Copy(u1,R),R),n,ww,R);
4313  v=p_Add_q(v,p_Copy(v1,R),R);
4314  }
4315  p_Delete(&u1,R);
4316  p_Delete(&v1,R);
4317  omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(short));
4318  return v;
4319 }
4320 
4321 BOOLEAN p_EqualPolys(poly p1,poly p2, const ring r)
4322 {
4323  while ((p1 != NULL) && (p2 != NULL))
4324  {
4325  if (! p_LmEqual(p1, p2,r))
4326  return FALSE;
4327  if (! n_Equal(p_GetCoeff(p1,r), p_GetCoeff(p2,r),r->cf ))
4328  return FALSE;
4329  pIter(p1);
4330  pIter(p2);
4331  }
4332  return (p1==p2);
4333 }
4334 
4335 static inline BOOLEAN p_ExpVectorEqual(poly p1, poly p2, const ring r1, const ring r2)
4336 {
4337  assume( r1 == r2 || rSamePolyRep(r1, r2) );
4338 
4339  p_LmCheckPolyRing1(p1, r1);
4340  p_LmCheckPolyRing1(p2, r2);
4341 
4342  int i = r1->ExpL_Size;
4343 
4344  assume( r1->ExpL_Size == r2->ExpL_Size );
4345 
4346  unsigned long *ep = p1->exp;
4347  unsigned long *eq = p2->exp;
4348 
4349  do
4350  {
4351  i--;
4352  if (ep[i] != eq[i]) return FALSE;
4353  }
4354  while (i);
4355 
4356  return TRUE;
4357 }
4358 
4359 BOOLEAN p_EqualPolys(poly p1,poly p2, const ring r1, const ring r2)
4360 {
4361  assume( r1 == r2 || rSamePolyRep(r1, r2) ); // will be used in rEqual!
4362  assume( r1->cf == r2->cf );
4363 
4364  while ((p1 != NULL) && (p2 != NULL))
4365  {
4366  // returns 1 if ExpVector(p)==ExpVector(q): does not compare numbers !!
4367  // #define p_LmEqual(p1, p2, r) p_ExpVectorEqual(p1, p2, r)
4368 
4369  if (! p_ExpVectorEqual(p1, p2, r1, r2))
4370  return FALSE;
4371 
4372  if (! n_Equal(p_GetCoeff(p1,r1), p_GetCoeff(p2,r2), r1->cf ))
4373  return FALSE;
4374 
4375  pIter(p1);
4376  pIter(p2);
4377  }
4378  return (p1==p2);
4379 }
4380 
4381 /*2
4382 *returns TRUE if p1 is a skalar multiple of p2
4383 *assume p1 != NULL and p2 != NULL
4384 */
4385 BOOLEAN p_ComparePolys(poly p1,poly p2, const ring r)
4386 {
4387  number n,nn;
4388  pAssume(p1 != NULL && p2 != NULL);
4389 
4390  if (!p_LmEqual(p1,p2,r)) //compare leading mons
4391  return FALSE;
4392  if ((pNext(p1)==NULL) && (pNext(p2)!=NULL))
4393  return FALSE;
4394  if ((pNext(p2)==NULL) && (pNext(p1)!=NULL))
4395  return FALSE;
4396  if (pLength(p1) != pLength(p2))
4397  return FALSE;
4398  #ifdef HAVE_RINGS
4399  if (rField_is_Ring(r))
4400  {
4401  if (!n_DivBy(p_GetCoeff(p1, r), p_GetCoeff(p2, r), r)) return FALSE;
4402  }
4403  #endif
4404  n=n_Div(p_GetCoeff(p1,r),p_GetCoeff(p2,r),r);
4405  while ((p1 != NULL) /*&& (p2 != NULL)*/)
4406  {
4407  if ( ! p_LmEqual(p1, p2,r))
4408  {
4409  n_Delete(&n, r);
4410  return FALSE;
4411  }
4412  if (!n_Equal(p_GetCoeff(p1, r), nn = n_Mult(p_GetCoeff(p2, r),n, r->cf), r->cf))
4413  {
4414  n_Delete(&n, r);
4415  n_Delete(&nn, r);
4416  return FALSE;
4417  }
4418  n_Delete(&nn, r);
4419  pIter(p1);
4420  pIter(p2);
4421  }
4422  n_Delete(&n, r);
4423  return TRUE;
4424 }
4425 
4426 /*2
4427 * returns the length of a (numbers of monomials)
4428 * respect syzComp
4429 */
4430 poly p_Last(const poly p, int &l, const ring r)
4431 {
4432  if (p == NULL)
4433  {
4434  l = 0;
4435  return NULL;
4436  }
4437  l = 1;
4438  poly a = p;
4439  if (! rIsSyzIndexRing(r))
4440  {
4441  poly next = pNext(a);
4442  while (next!=NULL)
4443  {
4444  a = next;
4445  next = pNext(a);
4446  l++;
4447  }
4448  }
4449  else
4450  {
4451  int curr_limit = rGetCurrSyzLimit(r);
4452  poly pp = a;
4453  while ((a=pNext(a))!=NULL)
4454  {
4455  if (p_GetComp(a,r)<=curr_limit/*syzComp*/)
4456  l++;
4457  else break;
4458  pp = a;
4459  }
4460  a=pp;
4461  }
4462  return a;
4463 }
4464 
4465 int p_Var(poly m,const ring r)
4466 {
4467  if (m==NULL) return 0;
4468  if (pNext(m)!=NULL) return 0;
4469  int i,e=0;
4470  for (i=rVar(r); i>0; i--)
4471  {
4472  int exp=p_GetExp(m,i,r);
4473  if (exp==1)
4474  {
4475  if (e==0) e=i;
4476  else return 0;
4477  }
4478  else if (exp!=0)
4479  {
4480  return 0;
4481  }
4482  }
4483  return e;
4484 }
4485 
4486 /*2
4487 *the minimal index of used variables - 1
4488 */
4489 int p_LowVar (poly p, const ring r)
4490 {
4491  int k,l,lex;
4492 
4493  if (p == NULL) return -1;
4494 
4495  k = 32000;/*a very large dummy value*/
4496  while (p != NULL)
4497  {
4498  l = 1;
4499  lex = p_GetExp(p,l,r);
4500  while ((l < (rVar(r))) && (lex == 0))
4501  {
4502  l++;
4503  lex = p_GetExp(p,l,r);
4504  }
4505  l--;
4506  if (l < k) k = l;
4507  pIter(p);
4508  }
4509  return k;
4510 }
4511 
4512 /*2
4513 * verschiebt die Indizees der Modulerzeugenden um i
4514 */
4515 void p_Shift (poly * p,int i, const ring r)
4516 {
4517  poly qp1 = *p,qp2 = *p;/*working pointers*/
4518  int j = p_MaxComp(*p,r),k = p_MinComp(*p,r);
4519 
4520  if (j+i < 0) return ;
4521  while (qp1 != NULL)
4522  {
4523  if ((p_GetComp(qp1,r)+i > 0) || ((j == -i) && (j == k)))
4524  {
4525  p_AddComp(qp1,i,r);
4526  p_SetmComp(qp1,r);
4527  qp2 = qp1;
4528  pIter(qp1);
4529  }
4530  else
4531  {
4532  if (qp2 == *p)
4533  {
4534  pIter(*p);
4535  p_LmDelete(&qp2,r);
4536  qp2 = *p;
4537  qp1 = *p;
4538  }
4539  else
4540  {
4541  qp2->next = qp1->next;
4542  if (qp1!=NULL) p_LmDelete(&qp1,r);
4543  qp1 = qp2->next;
4544  }
4545  }
4546  }
4547 }
4548 
4549 /***************************************************************
4550  *
4551  * Storage Managament Routines
4552  *
4553  ***************************************************************/
4554 
4555 
4556 static inline unsigned long GetBitFields(const long e,
4557  const unsigned int s, const unsigned int n)
4558 {
4559 #define Sy_bit_L(x) (((unsigned long)1L)<<(x))
4560  unsigned int i = 0;
4561  unsigned long ev = 0L;
4562  assume(n > 0 && s < BIT_SIZEOF_LONG);
4563  do
4564  {
4565  assume(s+i < BIT_SIZEOF_LONG);
4566  if (e > (long) i) ev |= Sy_bit_L(s+i);
4567  else break;
4568  i++;
4569  }
4570  while (i < n);
4571  return ev;
4572 }
4573 
4574 // Short Exponent Vectors are used for fast divisibility tests
4575 // ShortExpVectors "squeeze" an exponent vector into one word as follows:
4576 // Let n = BIT_SIZEOF_LONG / pVariables.
4577 // If n == 0 (i.e. pVariables > BIT_SIZE_OF_LONG), let m == the number
4578 // of non-zero exponents. If (m>BIT_SIZEOF_LONG), then sev = ~0, else
4579 // first m bits of sev are set to 1.
4580 // Otherwise (i.e. pVariables <= BIT_SIZE_OF_LONG)
4581 // represented by a bit-field of length n (resp. n+1 for some
4582 // exponents). If the value of an exponent is greater or equal to n, then
4583 // all of its respective n bits are set to 1. If the value of an exponent
4584 // is smaller than n, say m, then only the first m bits of the respective
4585 // n bits are set to 1, the others are set to 0.
4586 // This way, we have:
4587 // exp1 / exp2 ==> (ev1 & ~ev2) == 0, i.e.,
4588 // if (ev1 & ~ev2) then exp1 does not divide exp2
4589 unsigned long p_GetShortExpVector(const poly p, const ring r)
4590 {
4591  assume(p != NULL);
4592  unsigned long ev = 0; // short exponent vector
4593  unsigned int n = BIT_SIZEOF_LONG / r->N; // number of bits per exp
4594  unsigned int m1; // highest bit which is filled with (n+1)
4595  int i=0,j=1;
4596 
4597  if (n == 0)
4598  {
4599  if (r->N <2*BIT_SIZEOF_LONG)
4600  {
4601  n=1;
4602  m1=0;
4603  }
4604  else
4605  {
4606  for (; j<=r->N; j++)
4607  {
4608  if (p_GetExp(p,j,r) > 0) i++;
4609  if (i == BIT_SIZEOF_LONG) break;
4610  }
4611  if (i>0)
4612  ev = ~(0UL) >> (BIT_SIZEOF_LONG - i);
4613  return ev;
4614  }
4615  }
4616  else
4617  {
4618  m1 = (n+1)*(BIT_SIZEOF_LONG - n*r->N);
4619  }
4620 
4621  n++;
4622  while (i<m1)
4623  {
4624  ev |= GetBitFields(p_GetExp(p, j,r), i, n);
4625  i += n;
4626  j++;
4627  }
4628 
4629  n--;
4630  while (i<BIT_SIZEOF_LONG)
4631  {
4632  ev |= GetBitFields(p_GetExp(p, j,r), i, n);
4633  i += n;
4634  j++;
4635  }
4636  return ev;
4637 }
4638 
4639 
4640 /// p_GetShortExpVector of p * pp
4641 unsigned long p_GetShortExpVector(const poly p, const poly pp, const ring r)
4642 {
4643  assume(p != NULL);
4644  assume(pp != NULL);
4645 
4646  unsigned long ev = 0; // short exponent vector
4647  unsigned int n = BIT_SIZEOF_LONG / r->N; // number of bits per exp
4648  unsigned int m1; // highest bit which is filled with (n+1)
4649  int j=1;
4650  unsigned long i = 0L;
4651 
4652  if (n == 0)
4653  {
4654  if (r->N <2*BIT_SIZEOF_LONG)
4655  {
4656  n=1;
4657  m1=0;
4658  }
4659  else
4660  {
4661  for (; j<=r->N; j++)
4662  {
4663  if (p_GetExp(p,j,r) > 0 || p_GetExp(pp,j,r) > 0) i++;
4664  if (i == BIT_SIZEOF_LONG) break;
4665  }
4666  if (i>0)
4667  ev = ~(0UL) >> (BIT_SIZEOF_LONG - i);
4668  return ev;
4669  }
4670  }
4671  else
4672  {
4673  m1 = (n+1)*(BIT_SIZEOF_LONG - n*r->N);
4674  }
4675 
4676  n++;
4677  while (i<m1)
4678  {
4679  ev |= GetBitFields(p_GetExp(p, j,r) + p_GetExp(pp, j,r), i, n);
4680  i += n;
4681  j++;
4682  }
4683 
4684  n--;
4685  while (i<BIT_SIZEOF_LONG)
4686  {
4687  ev |= GetBitFields(p_GetExp(p, j,r) + p_GetExp(pp, j,r), i, n);
4688  i += n;
4689  j++;
4690  }
4691  return ev;
4692 }
4693 
4694 
4695 
4696 /***************************************************************
4697  *
4698  * p_ShallowDelete
4699  *
4700  ***************************************************************/
4701 #undef LINKAGE
4702 #define LINKAGE
4703 #undef p_Delete__T
4704 #define p_Delete__T p_ShallowDelete
4705 #undef n_Delete__T
4706 #define n_Delete__T(n, r) do {} while (0)
4707 
4709 
#define p_LmCheckPolyRing2(p, r)
Definition: monomials.h:207
#define n_New(n, r)
Definition: coeffs.h:444
int status int void size_t count
Definition: si_signals.h:59
BOOLEAN p_VectorHasUnitB(poly p, int *k, const ring r)
Definition: p_polys.cc:3225
void p_SetModDeg(intvec *w, ring r)
Definition: p_polys.cc:3533
for idElimination, like a, except pFDeg, pWeigths ignore it
Definition: ring.h:99
#define STATISTIC(f)
Definition: numstats.h:16
unsigned long p_GetMaxExpL(poly p, const ring r, unsigned long l_max)
return the maximal exponent of p in form of the maximal long var
Definition: p_polys.cc:1175
static FORCE_INLINE number n_GetNumerator(number &n, const coeffs r)
return the numerator of n (if elements of r are by nature not fractional, result is n) ...
Definition: coeffs.h:612
static FORCE_INLINE number n_GetUnit(number n, const coeffs r)
in Z: 1 in Z/kZ (where k is not a prime): largest divisor of n (taken in Z) that is co-prime with k i...
Definition: coeffs.h:536
LINLINE number nlSub(number la, number li, const coeffs r)
Definition: longrat.cc:2579
static unsigned long p_AddComp(poly p, unsigned long v, ring r)
Definition: p_polys.h:442
static FORCE_INLINE number n_Gcd(number a, number b, const coeffs r)
in Z: return the gcd of &#39;a&#39; and &#39;b&#39; in Z/nZ, Z/2^kZ: computed as in the case Z in Z/pZ...
Definition: coeffs.h:690
void p_Setm_General(poly p, const ring r)
Definition: p_polys.cc:164
const CanonicalForm int s
Definition: facAbsFact.cc:55
poly p_Diff(poly a, int k, const ring r)
Definition: p_polys.cc:1820
static poly p_MonPower(poly p, int exp, const ring r)
Definition: p_polys.cc:1922
const char * eati(const char *s, int *i)
Definition: reporter.cc:373
const CanonicalForm int const CFList const Variable & y
Definition: facAbsFact.cc:57
int p_GetVariables(poly p, int *e, const ring r)
set entry e[i] to 1 if var(i) occurs in p, ignore var(j) if e[j]>0 return #(e[i]>0) ...
Definition: p_polys.cc:1267
#define D(A)
Definition: gentable.cc:119
for int64 weights
Definition: ring.h:79
static BOOLEAN p_LmIsConstantComp(const poly p, const ring r)
Definition: p_polys.h:932
static FORCE_INLINE BOOLEAN n_IsUnit(number n, const coeffs r)
TRUE iff n has a multiplicative inverse in the given coeff field/ring r.
Definition: coeffs.h:519
Definition: ring.h:68
const poly a
Definition: syzextra.cc:212
static FORCE_INLINE const char * n_Read(const char *s, number *a, const coeffs r)
!!! Recommendation: This method is too cryptic to be part of the user- !!! interface. As defined here, it is merely a helper !!! method for parsing number input strings.
Definition: coeffs.h:602
#define POLY_NEGWEIGHT_OFFSET
Definition: monomials.h:244
#define Print
Definition: emacs.cc:83
long pLDeg1(poly p, int *l, const ring r)
Definition: p_polys.cc:841
static poly p_LmDeleteAndNext(poly p, const ring r)
Definition: p_polys.h:720
return
Definition: syzextra.cc:280
Definition: ring.h:61
static BOOLEAN rField_is_Zp_a(const ring r)
Definition: ring.h:518
long pLDeg1c_Totaldegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1005
BOOLEAN p_ComparePolys(poly p1, poly p2, const ring r)
returns TRUE if p1 is a skalar multiple of p2 assume p1 != NULL and p2 != NULL
Definition: p_polys.cc:4385
p_SetmProc p_GetSetmProc(const ring r)
Definition: p_polys.cc:560
BOOLEAN nlGreaterZero(number za, const coeffs r)
Definition: longrat.cc:1153
used for all transcendental extensions, i.e., the top-most extension in an extension tower is transce...
Definition: coeffs.h:39
loop
Definition: myNF.cc:98
if(0 > strat->sl)
Definition: myNF.cc:73
static int si_min(const int a, const int b)
Definition: auxiliary.h:124
long pLDeg1c(poly p, int *l, const ring r)
Definition: p_polys.cc:877
#define FALSE
Definition: auxiliary.h:97
static BOOLEAN pOldLexOrder
Definition: p_polys.cc:3522
poly p_Homogen(poly p, int varnum, const ring r)
Definition: p_polys.cc:3154
void p_Split(poly p, poly *h)
Definition: p_polys.cc:1320
return P p
Definition: myNF.cc:203
opposite of ls
Definition: ring.h:100
void p_VectorHasUnit(poly p, int *k, int *len, const ring r)
Definition: p_polys.cc:3250
short * iv2array(intvec *iv, const ring R)
Definition: weight.cc:208
static FORCE_INLINE BOOLEAN n_IsOne(number n, const coeffs r)
TRUE iff &#39;n&#39; represents the one element.
Definition: coeffs.h:472
BOOLEAN p_IsHomogeneous(poly p, const ring r)
Definition: p_polys.cc:3203
static int rPar(const ring r)
(r->cf->P)
Definition: ring.h:587
static poly p_Mult_mm(poly p, poly m, const ring r)
Definition: p_polys.h:968
poly p_NSet(number n, const ring r)
returns the poly representing the number n, destroys n
Definition: p_polys.cc:1443
#define pAssume(cond)
Definition: monomials.h:98
static long p_IncrExp(poly p, int v, ring r)
Definition: p_polys.h:586
void p_Setm_WFirstTotalDegree(poly p, const ring r)
Definition: p_polys.cc:554
int p_MinDeg(poly p, intvec *w, const ring R)
Definition: p_polys.cc:4257
number ndCopyMap(number a, const coeffs aRing, const coeffs r)
Definition: numbers.cc:239
static unsigned long p_SetComp(poly p, unsigned long c, ring r)
Definition: p_polys.h:242
static BOOLEAN rIsSyzIndexRing(const ring r)
Definition: ring.h:708
static int _componentsExternal
Definition: p_polys.cc:154
static poly p_DiffOpM(poly a, poly b, BOOLEAN multiply, const ring r)
Definition: p_polys.cc:1856
#define p_GetComp(p, r)
Definition: monomials.h:72
static int rGetCurrSyzLimit(const ring r)
Definition: ring.h:711
static poly last
Definition: hdegree.cc:1077
static BOOLEAN rIsRatGRing(const ring r)
Definition: ring.h:415
#define TEST_OPT_CONTENTSB
Definition: options.h:121
long p_WDegree(poly p, const ring r)
Definition: p_polys.cc:714
static FORCE_INLINE number n_Init(long i, const coeffs r)
a number representing i in the given coeff field/ring r
Definition: coeffs.h:542
static long pModDeg(poly p, ring r)
Definition: p_polys.cc:3524
static void p_GetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1443
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
LINLINE number nlAdd(number la, number li, const coeffs r)
Definition: longrat.cc:2513
static poly p_Subst1(poly p, int n, const ring r)
Definition: p_polys.cc:3683
int rChar(ring r)
Definition: ring.cc:684
void nlInpGcd(number &a, number b, const coeffs r)
Definition: longrat.cc:2757
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:580
poly singclap_gcd(poly f, poly g, const ring r)
destroys f and g
Definition: clapsing.cc:287
long pLDeg0c(poly p, int *l, const ring r)
Definition: p_polys.cc:770
void sBucketDestroyAdd(sBucket_pt bucket, poly *p, int *length)
Definition: sbuckets.h:72
poly p_Div_nn(poly p, const number n, const ring r)
Definition: p_polys.cc:1475
long int64
Definition: auxiliary.h:69
void p_Lcm(const poly a, const poly b, poly m, const ring r)
Definition: p_polys.cc:1583
static BOOLEAN rField_is_Q_a(const ring r)
Definition: ring.h:528
number nlGcd(number a, number b, const coeffs r)
Definition: longrat.cc:1190
#define TRUE
Definition: auxiliary.h:101
static void p_MonMult(poly p, poly q, const ring r)
Definition: p_polys.cc:1946
static long p_Totaldegree(poly p, const ring r)
Definition: p_polys.h:1430
static FORCE_INLINE void n_Normalize(number &n, const coeffs r)
inplace-normalization of n; produces some canonical representation of n;
Definition: coeffs.h:582
void sBucket_Add_p(sBucket_pt bucket, poly p, int length)
adds poly p to bucket destroys p!
Definition: sbuckets.cc:201
void * ADDRESS
Definition: auxiliary.h:118
poly p_Subst(poly p, int n, poly e, const ring r)
Definition: p_polys.cc:3776
g
Definition: cfModGcd.cc:4031
void WerrorS(const char *s)
Definition: feFopen.cc:24
int k
Definition: cfEzgcd.cc:93
void p_Setm_TotalDegree(poly p, const ring r)
Definition: p_polys.cc:547
static FORCE_INLINE number n_NormalizeHelper(number a, number b, const coeffs r)
assume that r is a quotient field (otherwise, return 1) for arguments (a1/a2,b1/b2) return (lcm(a1...
Definition: coeffs.h:721
poly p_TakeOutComp1(poly *p, int k, const ring r)
Definition: p_polys.cc:3282
static BOOLEAN rField_is_GF(const ring r)
Definition: ring.h:510
void p_Norm(poly p1, const ring r)
Definition: p_polys.cc:3581
static FORCE_INLINE BOOLEAN nCoeff_is_Q(const coeffs r)
Definition: coeffs.h:840
void p_SimpleContent(poly ph, int smax, const ring r)
Definition: p_polys.cc:2425
static long p_MultExp(poly p, int v, long ee, ring r)
Definition: p_polys.h:616
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
Definition: monomials.h:51
long pLDeg1_Deg(poly p, int *l, const ring r)
Definition: p_polys.cc:910
#define WarnS
Definition: emacs.cc:81
Definition: ring.h:66
poly singclap_pdivide(poly f, poly g, const ring r)
Definition: clapsing.cc:547
#define omAlloc(size)
Definition: omAllocDecl.h:210
static poly p_TwoMonPower(poly p, int exp, const ring r)
Definition: p_polys.cc:2028
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:407
poly p_Sub(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1912
long(* pLDegProc)(poly p, int *length, ring r)
Definition: ring.h:45
static void p_LmFree(poly p, ring)
Definition: p_polys.h:678
Definition: ring.h:64
static int pLength(poly a)
Definition: p_polys.h:189
static void p_SetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1451
union sro_ord::@0 data
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition: p_polys.h:804
poly pp
Definition: myNF.cc:296
static FORCE_INLINE BOOLEAN nCoeff_is_Q_a(const coeffs r)
Definition: coeffs.h:899
static long p_SubExp(poly p, int v, long ee, ring r)
Definition: p_polys.h:608
long totaldegreeWecart_IV(poly p, ring r, const short *w)
Definition: weight.cc:239
long pLDeg1c_Deg(poly p, int *l, const ring r)
Definition: p_polys.cc:941
static BOOLEAN rField_has_simple_inverse(const ring r)
Definition: ring.h:537
int comp(const CanonicalForm &A, const CanonicalForm &B)
compare polynomials
static FORCE_INLINE number n_Param(const int iParameter, const coeffs r)
return the (iParameter^th) parameter as a NEW number NOTE: parameter numbering: 1..n_NumberOfParameters(...)
Definition: coeffs.h:817
void p_Cleardenom_n(poly ph, const ring r, number &c)
Definition: p_polys.cc:2844
int p_IsUnivariate(poly p, const ring r)
return i, if poly depends only on var(i)
Definition: p_polys.cc:1247
#define pIter(p)
Definition: monomials.h:44
poly pp_JetW(poly p, int m, short *w, const ring R)
Definition: p_polys.cc:4212
poly res
Definition: myNF.cc:322
static poly p_SortAdd(poly p, const ring r, BOOLEAN revert=FALSE)
Definition: p_polys.h:1142
static FORCE_INLINE number n_Mult(number a, number b, const coeffs r)
return the product of &#39;a&#39; and &#39;b&#39;, i.e., a*b
Definition: coeffs.h:640
BOOLEAN p_CheckPolyRing(poly p, ring r)
Definition: pDebug.cc:111
int p_Weight(int i, const ring r)
Definition: p_polys.cc:705
#define omReallocSize(addr, o_size, size)
Definition: omAllocDecl.h:220
const char * p_Read(const char *st, poly &rc, const ring r)
Definition: p_polys.cc:1348
ro_typ ord_typ
Definition: ring.h:228
void p_ContentRat(poly &ph, const ring r)
Definition: p_polys.cc:1666
static FORCE_INLINE void n_ClearContent(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs r)
Computes the content and (inplace) divides it out on a collection of numbers number c is the content ...
Definition: coeffs.h:942
static poly p_Head(poly p, const ring r)
Definition: p_polys.h:812
long p_DegW(poly p, const short *w, const ring R)
Definition: p_polys.cc:690
void p_Setm_Syz(poly p, ring r, int *Components, long *ShiftedComponents)
Definition: p_polys.cc:531
int r_IsRingVar(const char *n, char **names, int N)
Definition: ring.cc:222
long p_Deg(poly a, const ring r)
Definition: p_polys.cc:587
const ring r
Definition: syzextra.cc:208
poly p_LcmRat(const poly a, const poly b, const long lCompM, const ring r)
Definition: p_polys.cc:1599
poly p_PermPoly(poly p, const int *perm, const ring oldRing, const ring dst, nMapFunc nMap, const int *par_perm, int OldPar, BOOLEAN use_mult)
Definition: p_polys.cc:3939
int p_Size(poly p, const ring r)
Definition: p_polys.cc:3137
poly p_Invers(int n, poly u, intvec *w, const ring R)
Definition: p_polys.cc:4293
static int p_Comp_k_n(poly a, poly b, int k, ring r)
Definition: p_polys.h:635
poly p_Farey(poly p, number N, const ring r)
Definition: p_polys.cc:62
#define TEST_OPT_INTSTRATEGY
Definition: options.h:105
static FORCE_INLINE BOOLEAN nCoeff_is_algExt(const coeffs r)
TRUE iff r represents an algebraic extension field.
Definition: coeffs.h:924
static void p_SetCompP(poly p, int i, ring r)
Definition: p_polys.h:249
Definition: intvec.h:14
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:49
poly p_One(const ring r)
Definition: p_polys.cc:1313
Concrete implementation of enumerators over polynomials.
for(int i=0;i< R->ExpL_Size;i++) Print("%09lx "
Definition: cfEzgcd.cc:66
static long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset)
get a single variable exponent : the integer VarOffset encodes:
Definition: p_polys.h:464
int j
Definition: myNF.cc:70
static int max(int a, int b)
Definition: fast_mult.cc:264
This is a polynomial enumerator for simple iteration over coefficients of polynomials.
void pSetDegProcs(ring r, pFDegProc new_FDeg, pLDegProc new_lDeg)
Definition: p_polys.cc:3497
#define omFree(addr)
Definition: omAllocDecl.h:261
number ntInit(long i, const coeffs cf)
Definition: transext.cc:613
#define assume(x)
Definition: mod2.h:403
static BOOLEAN rIsPluralRing(const ring r)
we must always have this test!
Definition: ring.h:404
static BOOLEAN p_IsConstant(const poly p, const ring r)
Definition: p_polys.h:1876
The main handler for Singular numbers which are suitable for Singular polynomials.
void p_LmDeleteAndNextRat(poly *p, int ishift, ring r)
Definition: p_polys.cc:1622
static poly p_MonMultC(poly p, poly q, const ring rr)
Definition: p_polys.cc:1966
long p_WFirstTotalDegree(poly p, const ring r)
Definition: p_polys.cc:596
static poly p_Subst0(poly p, int n, const ring r)
Definition: p_polys.cc:3751
void p_Vec2Polys(poly v, poly **p, int *len, const ring r)
Definition: p_polys.cc:3475
number(* nMapFunc)(number a, const coeffs src, const coeffs dst)
maps "a", which lives in src, into dst
Definition: coeffs.h:73
static FORCE_INLINE BOOLEAN n_DivBy(number a, number b, const coeffs r)
test whether &#39;a&#39; is divisible &#39;b&#39;; for r encoding a field: TRUE iff &#39;b&#39; does not represent zero in Z:...
Definition: coeffs.h:787
long pLDeg0(poly p, int *l, const ring r)
Definition: p_polys.cc:739
static FORCE_INLINE number n_ChineseRemainderSym(number *a, number *b, int rl, BOOLEAN sym, CFArray &inv_cache, const coeffs r)
Definition: coeffs.h:798
sBucket_pt sBucketCreate(const ring r)
Definition: sbuckets.cc:120
const ring R
Definition: DebugPrint.cc:36
static FORCE_INLINE void n_Write(number &n, const coeffs r, const BOOLEAN bShortOut=TRUE)
Definition: coeffs.h:595
#define n_Test(a, r)
BOOLEAN n_Test(number a, const coeffs r)
Definition: coeffs.h:742
poly p_Jet(poly p, int m, const ring R)
Definition: p_polys.cc:4195
static BOOLEAN p_DivisibleBy(poly a, poly b, const ring r)
Definition: p_polys.h:1777
void p_Setm_Dummy(poly p, const ring r)
Definition: p_polys.cc:541
Definition: ring.h:226
All the auxiliary stuff.
static int p_LmCmp(poly p, poly q, const ring r)
Definition: p_polys.h:1467
static FORCE_INLINE number n_Invers(number a, const coeffs r)
return the multiplicative inverse of &#39;a&#39;; raise an error if &#39;a&#39; is not invertible ...
Definition: coeffs.h:568
int m
Definition: cfEzgcd.cc:119
BOOLEAN rSamePolyRep(ring r1, ring r2)
returns TRUE, if r1 and r2 represents the monomials in the same way FALSE, otherwise this is an analo...
Definition: ring.cc:1682
static FORCE_INLINE number n_InpNeg(number n, const coeffs r)
in-place negation of n MUST BE USED: n = n_InpNeg(n) (no copy is returned)
Definition: coeffs.h:561
static int si_max(const int a, const int b)
Definition: auxiliary.h:123
static FORCE_INLINE BOOLEAN nCoeff_is_transExt(const coeffs r)
TRUE iff r represents a transcendental extension field.
Definition: coeffs.h:932
static number * pnBin(int exp, const ring r)
Definition: p_polys.cc:1980
static number p_InitContent(poly ph, const ring r)
Definition: p_polys.cc:2483
poly p_GetCoeffRat(poly p, int ishift, ring r)
Definition: p_polys.cc:1644
FILE * f
Definition: checklibs.c:7
int i
Definition: cfEzgcd.cc:123
Induced (Schreyer) ordering.
Definition: ring.h:101
void PrintS(const char *s)
Definition: reporter.cc:284
static poly p_Mult_nn(poly p, number n, const ring r)
Definition: p_polys.h:895
void(* p_SetmProc)(poly p, const ring r)
Definition: ring.h:47
static BOOLEAN rField_is_Q(const ring r)
Definition: ring.h:501
int p_LowVar(poly p, const ring r)
the minimal index of used variables - 1
Definition: p_polys.cc:4489
#define p_LmCheckPolyRing1(p, r)
Definition: monomials.h:185
static long p_MinComp(poly p, ring lmRing, ring tailRing)
Definition: p_polys.h:308
poly p_Divide(poly a, poly b, const ring r)
Definition: p_polys.cc:1462
static void p_ExpVectorAdd(poly p1, poly p2, const ring r)
Definition: p_polys.h:1334
S?
Definition: ring.h:83
BOOLEAN rOrd_SetCompRequiresSetm(const ring r)
return TRUE if p_SetComp requires p_Setm
Definition: ring.cc:1876
static poly p_LmFreeAndNext(poly p, ring)
Definition: p_polys.h:698
BOOLEAN pSetm_error
Definition: p_polys.cc:156
void rWrite(ring r, BOOLEAN details)
Definition: ring.cc:236
void p_Content(poly ph, const ring r)
Definition: p_polys.cc:2216
#define IDELEMS(i)
Definition: simpleideals.h:24
static FORCE_INLINE BOOLEAN n_IsZero(number n, const coeffs r)
TRUE iff &#39;n&#39; represents the zero element.
Definition: coeffs.h:468
static FORCE_INLINE BOOLEAN nCoeff_is_GF(const coeffs r)
Definition: coeffs.h:853
static FORCE_INLINE nMapFunc n_SetMap(const coeffs src, const coeffs dst)
set the mapping function pointers for translating numbers from src to dst
Definition: coeffs.h:725
static poly p_Pow(poly p, int i, const ring r)
Definition: p_polys.cc:2093
poly p_JetW(poly p, int m, short *w, const ring R)
Definition: p_polys.cc:4239
poly p_PolyDiv(poly &p, const poly divisor, const BOOLEAN needResult, const ring r)
assumes that p and divisor are univariate polynomials in r, mentioning the same variable; assumes div...
Definition: p_polys.cc:1792
static short scaFirstAltVar(ring r)
Definition: sca.h:18
static poly pReverse(poly p)
Definition: p_polys.h:330
static FORCE_INLINE n_coeffType getCoeffType(const coeffs r)
Returns the type of coeffs domain.
Definition: coeffs.h:425
Definition: ring.h:69
poly p_Series(int n, poly p, poly u, intvec *w, const ring R)
Definition: p_polys.cc:4279
BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r)
Definition: p_polys.cc:4321
#define p_Test(p, r)
Definition: p_polys.h:160
short errorreported
Definition: feFopen.cc:23
static unsigned long p_SubComp(poly p, unsigned long v, ring r)
Definition: p_polys.h:448
void pRestoreDegProcs(ring r, pFDegProc old_FDeg, pLDegProc old_lDeg)
Definition: p_polys.cc:3509
Definition: ring.h:69
static long p_GetOrder(poly p, ring r)
Definition: p_polys.h:416
static BOOLEAN rField_is_Zp(const ring r)
Definition: ring.h:495
static FORCE_INLINE number n_Farey(number a, number b, const coeffs r)
Definition: coeffs.h:801
BOOLEAN p_HasNotCF(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1329
void p_Shift(poly *p, int i, const ring r)
shifts components of the vector p by i
Definition: p_polys.cc:4515
void p_Normalize(poly p, const ring r)
Definition: p_polys.cc:3634
#define rRing_has_Comp(r)
Definition: monomials.h:274
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:843
#define p_SetmComp
Definition: p_polys.h:239
void nlNormalize(number &x, const coeffs r)
Definition: longrat.cc:1332
poly p_mInit(const char *st, BOOLEAN &ok, const ring r)
Definition: p_polys.cc:1420
unsigned long p_GetShortExpVector(const poly p, const ring r)
Definition: p_polys.cc:4589
#define p_LmEqual(p1, p2, r)
Definition: p_polys.h:1611
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:37
static unsigned long GetBitFields(const long e, const unsigned int s, const unsigned int n)
Definition: p_polys.cc:4556
#define mpz_size1(A)
Definition: si_gmp.h:12
static unsigned long p_SetExp(poly p, const unsigned long e, const unsigned long iBitmask, const int VarOffset)
set a single variable exponent : VarOffset encodes the position in p->exp
Definition: p_polys.h:483
poly p_DivideM(poly a, poly b, const ring r)
Definition: p_polys.cc:1512
BOOLEAN p_DivisibleByRingCase(poly f, poly g, const ring r)
divisibility check over ground ring (which may contain zero divisors); TRUE iff LT(f) divides LT(g)...
Definition: p_polys.cc:1570
int p_IsPurePower(const poly p, const ring r)
return i, if head depends only on var(i)
Definition: p_polys.cc:1226
static FORCE_INLINE void n_Power(number a, int b, number *res, const coeffs r)
fill res with the power a^b
Definition: coeffs.h:636
static pLDegProc pOldLDeg
Definition: p_polys.cc:3521
long pLDegb(poly p, int *l, const ring r)
Definition: p_polys.cc:811
static BOOLEAN rField_is_Ring(const ring r)
Definition: ring.h:477
#define NULL
Definition: omList.c:10
void pEnlargeSet(poly **p, int l, int increment)
Definition: p_polys.cc:3556
long pLDeg1_Totaldegree(poly p, int *l, const ring r)
Definition: p_polys.cc:975
int length() const
Definition: intvec.h:86
static FORCE_INLINE number n_Copy(number n, const coeffs r)
return a copy of &#39;n&#39;
Definition: coeffs.h:455
BOOLEAN p_LmCheckPolyRing(poly p, ring r)
Definition: pDebug.cc:119
poly n_PermNumber(const number z, const int *par_perm, const int, const ring src, const ring dst)
Definition: p_polys.cc:3835
static int * _components
Definition: p_polys.cc:152
used for all algebraic extensions, i.e., the top-most extension in an extension tower is algebraic ...
Definition: coeffs.h:36
LINLINE void nlDelete(number *a, const coeffs r)
Definition: longrat.cc:2478
poly p_Last(const poly p, int &l, const ring r)
Definition: p_polys.cc:4430
static FORCE_INLINE number n_Div(number a, number b, const coeffs r)
return the quotient of &#39;a&#39; and &#39;b&#39;, i.e., a/b; raises an error if &#39;b&#39; is not invertible in r exceptio...
Definition: coeffs.h:619
long(* pFDegProc)(poly p, ring r)
Definition: ring.h:46
static poly p_Pow_charp(poly p, int i, const ring r)
Definition: p_polys.cc:2107
static pFDegProc pOldFDeg
Definition: p_polys.cc:3520
#define SR_INT
Definition: longrat.h:68
const CanonicalForm & w
Definition: facAbsFact.cc:55
static short scaLastAltVar(ring r)
Definition: sca.h:25
poly p_ChineseRemainder(poly *xx, number *x, number *q, int rl, CFArray &inv_cache, const ring R)
Definition: p_polys.cc:95
Variable x
Definition: cfModGcd.cc:4023
Definition: ring.h:63
static FORCE_INLINE number n_GetDenom(number &n, const coeffs r)
return the denominator of n (if elements of r are by nature not fractional, result is 1) ...
Definition: coeffs.h:607
static bool rIsSCA(const ring r)
Definition: nc.h:206
Definition: ring.h:60
#define pNext(p)
Definition: monomials.h:43
static FORCE_INLINE BOOLEAN n_Equal(number a, number b, const coeffs r)
TRUE iff &#39;a&#39; and &#39;b&#39; represent the same number; they may have different representations.
Definition: coeffs.h:464
static poly p_LmInit(poly p, const ring r)
Definition: p_polys.h:1258
#define Sy_bit_L(x)
static FORCE_INLINE number n_ExactDiv(number a, number b, const coeffs r)
assume that there is a canonical subring in cf and we know that division is possible for these a and ...
Definition: coeffs.h:626
static void p_Setm(poly p, const ring r)
Definition: p_polys.h:228
#define pSetCoeff0(p, n)
Definition: monomials.h:67
poly p_DiffOp(poly a, poly b, BOOLEAN multiply, const ring r)
Definition: p_polys.cc:1895
static long p_AddExp(poly p, int v, long ee, ring r)
Definition: p_polys.h:601
#define p_GetCoeff(p, r)
Definition: monomials.h:57
long pLDeg1_WFirstTotalDegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1038
static poly p_GetExp_k_n(poly p, int l, int k, const ring r)
Definition: p_polys.h:1295
static FORCE_INLINE number n_SubringGcd(number a, number b, const coeffs r)
Definition: coeffs.h:692
long pLDeg1c_WFirstTotalDegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1068
Definition: ring.h:62
int dReportError(const char *fmt,...)
Definition: dError.cc:45
static FORCE_INLINE BOOLEAN nCoeff_is_Extension(const coeffs r)
Definition: coeffs.h:860
static long * _componentsShifted
Definition: p_polys.cc:153
static unsigned long p_GetMaxExpL2(unsigned long l1, unsigned long l2, const ring r, unsigned long number_of_exp)
Definition: p_polys.cc:1107
p exp[i]
Definition: DebugPrint.cc:39
static poly p_Neg(poly p, const ring r)
Definition: p_polys.h:1013
static void p_LmDelete(poly p, const ring r)
Definition: p_polys.h:706
int exponent(const CanonicalForm &f, int q)
int exponent ( const CanonicalForm & f, int q )
static poly p_Subst2(poly p, int n, number e, const ring r)
Definition: p_polys.cc:3710
#define BIT_SIZEOF_LONG
Definition: auxiliary.h:81
long p_WTotaldegree(poly p, const ring r)
Definition: p_polys.cc:613
#define SR_HDL(A)
Definition: tgb.cc:35
END_NAMESPACE const void * p2
Definition: syzextra.cc:202
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete &#39;p&#39;
Definition: coeffs.h:459
void p_wrp(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:237
poly pp_Jet(poly p, int m, const ring R)
Definition: p_polys.cc:4167
poly nc_pSubst(poly p, int n, poly e, const ring r)
substitute the n-th variable by e in p destroy p e is not a constant
Definition: old.gring.cc:3288
void p_Write(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:206
static FORCE_INLINE BOOLEAN n_GreaterZero(number n, const coeffs r)
ordered fields: TRUE iff &#39;n&#39; is positive; in Z/pZ: TRUE iff 0 < m <= roundedBelow(p/2), where m is the long representing n in C: TRUE iff (Im(n) != 0 and Im(n) >= 0) or (Im(n) == 0 and Re(n) >= 0) in K(a)/<p(a)>: TRUE iff (n != 0 and (LC(n) > 0 or deg(n) > 0)) in K(t_1, ..., t_n): TRUE iff (LC(numerator(n) is a constant and > 0) or (LC(numerator(n) is not a constant) in Z/2^kZ: TRUE iff 0 < n <= 2^(k-1) in Z/mZ: TRUE iff the internal mpz is greater than zero in Z: TRUE iff n > 0
Definition: coeffs.h:498
static void p_ExpVectorSum(poly pr, poly p1, poly p2, const ring r)
Definition: p_polys.h:1348
static void p_SplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero, const ring r)
Definition: p_polys.cc:3647
static BOOLEAN p_ExpVectorEqual(poly p1, poly p2, const ring r1, const ring r2)
Definition: p_polys.cc:4335
static long p_DecrExp(poly p, int v, ring r)
Definition: p_polys.h:593
polyrec * poly
Definition: hilb.h:10
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:877
static BOOLEAN rField_has_Units(const ring r)
Definition: ring.h:483
poly p_TakeOutComp(poly *p, int k, const ring r)
Definition: p_polys.cc:3333
int offset
Definition: libparse.cc:1091
int perm[100]
static Poly * h
Definition: janet.cc:978
s?
Definition: ring.h:84
int BOOLEAN
Definition: auxiliary.h:88
static poly p_Init(const ring r, omBin bin)
Definition: p_polys.h:1243
const poly b
Definition: syzextra.cc:213
poly p_Cleardenom(poly p, const ring r)
Definition: p_polys.cc:2716
int p_Var(poly m, const ring r)
Definition: p_polys.cc:4465
static FORCE_INLINE int n_Size(number n, const coeffs r)
return a non-negative measure for the complexity of n; return 0 only when n represents zero; (used fo...
Definition: coeffs.h:574
void p_ProjectiveUnique(poly ph, const ring r)
Definition: p_polys.cc:3026
poly p_ISet(long i, const ring r)
returns the poly representing the integer i
Definition: p_polys.cc:1297
static poly p_Mult_q(poly p, poly q, const ring r)
Definition: p_polys.h:1020
poly p_Power(poly p, int i, const ring r)
Definition: p_polys.cc:2119
void Werror(const char *fmt,...)
Definition: reporter.cc:189
void p_DeleteComp(poly *p, int k, const ring r)
Definition: p_polys.cc:3442
#define omAlloc0(size)
Definition: omAllocDecl.h:211
return result
Definition: facAbsBiFact.cc:76
int l
Definition: cfEzgcd.cc:94
static FORCE_INLINE void n_ClearDenominators(ICoeffsEnumerator &numberCollectionEnumerator, number &d, const coeffs r)
(inplace) Clears denominators on a collection of numbers number d is the LCM of all the coefficient d...
Definition: coeffs.h:949
static long p_MaxComp(poly p, ring lmRing, ring tailRing)
Definition: p_polys.h:287
BOOLEAN p_OneComp(poly p, const ring r)
return TRUE if all monoms have the same component
Definition: p_polys.cc:1208
poly p_GetMaxExpP(poly p, const ring r)
return monomial r such that GetExp(r,i) is maximum of all monomials in p; coeff == 0...
Definition: p_polys.cc:1138
ListNode * next
Definition: janet.h:31
static void pnFreeBin(number *bin, int exp, const coeffs r)
Definition: p_polys.cc:2011