My Project
Loading...
Searching...
No Matches
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 "misc/auxiliary.h"
14
15#include "misc/options.h"
16#include "misc/intvec.h"
17
18
19#include "coeffs/longrat.h" // snumber is needed...
20#include "coeffs/numbers.h" // ndCopyMap
21
23
24#define TRANSEXT_PRIVATES
25
28
29#include "polys/weight.h"
30#include "polys/simpleideals.h"
31
32#include "ring.h"
33#include "p_polys.h"
34
38
39
40#ifdef HAVE_PLURAL
41#include "nc/nc.h"
42#include "nc/sca.h"
43#endif
44
45#ifdef HAVE_SHIFTBBA
46#include "polys/shiftop.h"
47#endif
48
49#include "clapsing.h"
50
51/*
52 * lift ideal with coeffs over Z (mod N) to Q via Farey
53 */
54poly p_Farey(poly p, number N, const ring r)
55{
56 poly h=p_Copy(p,r);
57 poly hh=h;
58 while(h!=NULL)
59 {
60 number c=pGetCoeff(h);
61 pSetCoeff0(h,n_Farey(c,N,r->cf));
62 n_Delete(&c,r->cf);
63 pIter(h);
64 }
65 while((hh!=NULL)&&(n_IsZero(pGetCoeff(hh),r->cf)))
66 {
67 p_LmDelete(&hh,r);
68 }
69 h=hh;
70 while((h!=NULL) && (pNext(h)!=NULL))
71 {
72 if(n_IsZero(pGetCoeff(pNext(h)),r->cf))
73 {
74 p_LmDelete(&pNext(h),r);
75 }
76 else pIter(h);
77 }
78 return hh;
79}
80/*2
81* xx,q: arrays of length 0..rl-1
82* xx[i]: SB mod q[i]
83* assume: char=0
84* assume: q[i]!=0
85* x: work space
86* destroys xx
87*/
88poly p_ChineseRemainder(poly *xx, number *x,number *q, int rl, CFArray &inv_cache,const ring R)
89{
90 poly r,h,hh;
91 int j;
92 poly res_p=NULL;
93 loop
94 {
95 /* search the lead term */
96 r=NULL;
97 for(j=rl-1;j>=0;j--)
98 {
99 h=xx[j];
100 if ((h!=NULL)
101 &&((r==NULL)||(p_LmCmp(r,h,R)==-1)))
102 r=h;
103 }
104 /* nothing found -> return */
105 if (r==NULL) break;
106 /* create the monomial in h */
107 h=p_Head(r,R);
108 /* collect the coeffs in x[..]*/
109 for(j=rl-1;j>=0;j--)
110 {
111 hh=xx[j];
112 if ((hh!=NULL) && (p_LmCmp(h,hh,R)==0))
113 {
114 x[j]=pGetCoeff(hh);
115 hh=p_LmFreeAndNext(hh,R);
116 xx[j]=hh;
117 }
118 else
119 x[j]=n_Init(0, R->cf);
120 }
121 number n=n_ChineseRemainderSym(x,q,rl,TRUE,inv_cache,R->cf);
122 for(j=rl-1;j>=0;j--)
123 {
124 x[j]=NULL; // n_Init(0...) takes no memory
125 }
126 if (n_IsZero(n,R->cf)) p_Delete(&h,R);
127 else
128 {
129 //Print("new mon:");pWrite(h);
130 p_SetCoeff(h,n,R);
131 pNext(h)=res_p;
132 res_p=h; // building res_p in reverse order!
133 }
134 }
135 res_p=pReverse(res_p);
136 p_Test(res_p, R);
137 return res_p;
138}
139
140/***************************************************************
141 *
142 * Completing what needs to be set for the monomial
143 *
144 ***************************************************************/
145// this is special for the syz stuff
149
151
152#ifndef SING_NDEBUG
153# define MYTEST 0
154#else /* ifndef SING_NDEBUG */
155# define MYTEST 0
156#endif /* ifndef SING_NDEBUG */
157
158void p_Setm_General(poly p, const ring r)
159{
161 int pos=0;
162 if (r->typ!=NULL)
163 {
164 loop
165 {
166 unsigned long ord=0;
167 sro_ord* o=&(r->typ[pos]);
168 switch(o->ord_typ)
169 {
170 case ro_dp:
171 {
172 int a,e;
173 a=o->data.dp.start;
174 e=o->data.dp.end;
175 for(int i=a;i<=e;i++) ord+=p_GetExp(p,i,r);
176 p->exp[o->data.dp.place]=ord;
177 break;
178 }
179 case ro_wp_neg:
181 // no break;
182 case ro_wp:
183 {
184 int a,e;
185 a=o->data.wp.start;
186 e=o->data.wp.end;
187 int *w=o->data.wp.weights;
188#if 1
189 for(int i=a;i<=e;i++) ord+=((unsigned long)p_GetExp(p,i,r))*((unsigned long)w[i-a]);
190#else
191 long ai;
192 int ei,wi;
193 for(int i=a;i<=e;i++)
194 {
195 ei=p_GetExp(p,i,r);
196 wi=w[i-a];
197 ai=ei*wi;
198 if (ai/ei!=wi) pSetm_error=TRUE;
199 ord+=ai;
200 if (ord<ai) pSetm_error=TRUE;
201 }
202#endif
203 p->exp[o->data.wp.place]=ord;
204 break;
205 }
206 case ro_am:
207 {
209 const short a=o->data.am.start;
210 const short e=o->data.am.end;
211 const int * w=o->data.am.weights;
212#if 1
213 for(short i=a; i<=e; i++, w++)
214 ord += ((*w) * p_GetExp(p,i,r));
215#else
216 long ai;
217 int ei,wi;
218 for(short i=a;i<=e;i++)
219 {
220 ei=p_GetExp(p,i,r);
221 wi=w[i-a];
222 ai=ei*wi;
223 if (ai/ei!=wi) pSetm_error=TRUE;
224 ord += ai;
225 if (ord<ai) pSetm_error=TRUE;
226 }
227#endif
228 const int c = p_GetComp(p,r);
229
230 const short len_gen= o->data.am.len_gen;
231
232 if ((c > 0) && (c <= len_gen))
233 {
234 assume( w == o->data.am.weights_m );
235 assume( w[0] == len_gen );
236 ord += w[c];
237 }
238
239 p->exp[o->data.am.place] = ord;
240 break;
241 }
242 case ro_wp64:
243 {
244 int64 ord=0;
245 int a,e;
246 a=o->data.wp64.start;
247 e=o->data.wp64.end;
248 int64 *w=o->data.wp64.weights64;
249 int64 ei,wi,ai;
250 for(int i=a;i<=e;i++)
251 {
252 //Print("exp %d w %d \n",p_GetExp(p,i,r),(int)w[i-a]);
253 //ord+=((int64)p_GetExp(p,i,r))*w[i-a];
254 ei=(int64)p_GetExp(p,i,r);
255 wi=w[i-a];
256 ai=ei*wi;
257 if(ei!=0 && ai/ei!=wi)
258 {
260 #if SIZEOF_LONG == 4
261 Print("ai %lld, wi %lld\n",ai,wi);
262 #else
263 Print("ai %ld, wi %ld\n",ai,wi);
264 #endif
265 }
266 ord+=ai;
267 if (ord<ai)
268 {
270 #if SIZEOF_LONG == 4
271 Print("ai %lld, ord %lld\n",ai,ord);
272 #else
273 Print("ai %ld, ord %ld\n",ai,ord);
274 #endif
275 }
276 }
277 #if SIZEOF_LONG == 4
278 int64 mask=(int64)0x7fffffff;
279 long a_0=(long)(ord&mask); //2^31
280 long a_1=(long)(ord >>31 ); /*(ord/(mask+1));*/
281
282 //Print("mask: %x, ord: %d, a_0: %d, a_1: %d\n"
283 //,(int)mask,(int)ord,(int)a_0,(int)a_1);
284 //Print("mask: %d",mask);
285
286 p->exp[o->data.wp64.place]=a_1;
287 p->exp[o->data.wp64.place+1]=a_0;
288 #elif SIZEOF_LONG == 8
289 p->exp[o->data.wp64.place]=ord;
290 #endif
291// if(p_Setm_error) PrintS("***************************\n"
292// "***************************\n"
293// "**WARNING: overflow error**\n"
294// "***************************\n"
295// "***************************\n");
296 break;
297 }
298 case ro_cp:
299 {
300 int a,e;
301 a=o->data.cp.start;
302 e=o->data.cp.end;
303 int pl=o->data.cp.place;
304 for(int i=a;i<=e;i++) { p->exp[pl]=p_GetExp(p,i,r); pl++; }
305 break;
306 }
307 case ro_syzcomp:
308 {
309 long c=__p_GetComp(p,r);
310 long sc = c;
311 int* Components = (_componentsExternal ? _components :
312 o->data.syzcomp.Components);
313 long* ShiftedComponents = (_componentsExternal ? _componentsShifted:
314 o->data.syzcomp.ShiftedComponents);
315 if (ShiftedComponents != NULL)
316 {
317 assume(Components != NULL);
318 assume(c == 0 || Components[c] != 0);
319 sc = ShiftedComponents[Components[c]];
320 assume(c == 0 || sc != 0);
321 }
322 p->exp[o->data.syzcomp.place]=sc;
323 break;
324 }
325 case ro_syz:
326 {
327 const unsigned long c = __p_GetComp(p, r);
328 const short place = o->data.syz.place;
329 const int limit = o->data.syz.limit;
330
331 if (c > (unsigned long)limit)
332 p->exp[place] = o->data.syz.curr_index;
333 else if (c > 0)
334 {
335 assume( (1 <= c) && (c <= (unsigned long)limit) );
336 p->exp[place]= o->data.syz.syz_index[c];
337 }
338 else
339 {
340 assume(c == 0);
341 p->exp[place]= 0;
342 }
343 break;
344 }
345 // Prefix for Induced Schreyer ordering
346 case ro_isTemp: // Do nothing?? (to be removed into suffix later on...?)
347 {
348 assume(p != NULL);
349
350#ifndef SING_NDEBUG
351#if MYTEST
352 Print("p_Setm_General: ro_isTemp ord: pos: %d, p: ", pos); p_wrp(p, r);
353#endif
354#endif
355 int c = p_GetComp(p, r);
356
357 assume( c >= 0 );
358
359 // Let's simulate case ro_syz above....
360 // Should accumulate (by Suffix) and be a level indicator
361 const int* const pVarOffset = o->data.isTemp.pVarOffset;
362
363 assume( pVarOffset != NULL );
364
365 // TODO: Can this be done in the suffix???
366 for( int i = 1; i <= r->N; i++ ) // No v[0] here!!!
367 {
368 const int vo = pVarOffset[i];
369 if( vo != -1) // TODO: optimize: can be done once!
370 {
371 // Hans! Please don't break it again! p_SetExp(p, ..., r, vo) is correct:
372 p_SetExp(p, p_GetExp(p, i, r), r, vo); // copy put them verbatim
373 // Hans! Please don't break it again! p_GetExp(p, r, vo) is correct:
374 assume( p_GetExp(p, r, vo) == p_GetExp(p, i, r) ); // copy put them verbatim
375 }
376 }
377#ifndef SING_NDEBUG
378 for( int i = 1; i <= r->N; i++ ) // No v[0] here!!!
379 {
380 const int vo = pVarOffset[i];
381 if( vo != -1) // TODO: optimize: can be done once!
382 {
383 // Hans! Please don't break it again! p_GetExp(p, r, vo) is correct:
384 assume( p_GetExp(p, r, vo) == p_GetExp(p, i, r) ); // copy put them verbatim
385 }
386 }
387#if MYTEST
388// if( p->exp[o->data.isTemp.start] > 0 )
389 PrintS("after Values: "); p_wrp(p, r);
390#endif
391#endif
392 break;
393 }
394
395 // Suffix for Induced Schreyer ordering
396 case ro_is:
397 {
398#ifndef SING_NDEBUG
399#if MYTEST
400 Print("p_Setm_General: ro_is ord: pos: %d, p: ", pos); p_wrp(p, r);
401#endif
402#endif
403
404 assume(p != NULL);
405
406 int c = p_GetComp(p, r);
407
408 assume( c >= 0 );
409 const ideal F = o->data.is.F;
410 const int limit = o->data.is.limit;
411 assume( limit >= 0 );
412 const int start = o->data.is.start;
413
414 if( F != NULL && c > limit )
415 {
416#ifndef SING_NDEBUG
417#if MYTEST
418 Print("p_Setm_General: ro_is : in rSetm: pos: %d, c: %d > limit: %d\n", c, pos, limit);
419 PrintS("preComputed Values: ");
420 p_wrp(p, r);
421#endif
422#endif
423// if( c > limit ) // BUG???
424 p->exp[start] = 1;
425// else
426// p->exp[start] = 0;
427
428
429 c -= limit;
430 assume( c > 0 );
431 c--;
432
433 if( c >= IDELEMS(F) )
434 break;
435
436 assume( c < IDELEMS(F) ); // What about others???
437
438 const poly pp = F->m[c]; // get reference monomial!!!
439
440 if(pp == NULL)
441 break;
442
443 assume(pp != NULL);
444
445#ifndef SING_NDEBUG
446#if MYTEST
447 Print("Respective F[c - %d: %d] pp: ", limit, c);
448 p_wrp(pp, r);
449#endif
450#endif
451
452 const int end = o->data.is.end;
453 assume(start <= end);
454
455
456// const int st = o->data.isTemp.start;
457
458#ifndef SING_NDEBUG
459#if MYTEST
460 Print("p_Setm_General: is(-Temp-) :: c: %d, limit: %d, [st:%d] ===>>> %ld\n", c, limit, start, p->exp[start]);
461#endif
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
531void 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
541void p_Setm_Dummy(poly p, const ring r)
542{
544}
545
546// for dp, Dp, ds, etc
547void p_Setm_TotalDegree(poly p, const ring r)
548{
550 p->exp[r->pOrdIndex] = p_Totaldegree(p, r);
551}
552
553// for wp, Wp, ws, etc
554void p_Setm_WFirstTotalDegree(poly p, const ring r)
555{
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 */
587long 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
596long 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*/
613long p_WTotaldegree(poly p, const ring r)
614{
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_am:
633 b1=si_min(b1,r->N);
634 /* no break, continue as ringorder_a*/
635 case ringorder_a:
636 for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
637 { // only one line
638 j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/];
639 }
640 return j*r->OrdSgn;
641 case ringorder_wp:
642 case ringorder_ws:
643 case ringorder_Wp:
644 case ringorder_Ws:
645 for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
646 { // in jedem block:
647 j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/];
648 }
649 break;
650 case ringorder_lp:
651 case ringorder_ls:
652 case ringorder_rs:
653 case ringorder_dp:
654 case ringorder_ds:
655 case ringorder_Dp:
656 case ringorder_Ds:
657 case ringorder_rp:
658 for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
659 {
660 j+= p_GetExp(p,k,r);
661 }
662 break;
663 case ringorder_a64:
664 {
665 int64* w=(int64*)r->wvhdl[i];
666 for (k=0;k<=(b1 /*r->block1[i]*/ - b0 /*r->block0[i]*/);k++)
667 {
668 //there should be added a line which checks if w[k]>2^31
669 j+= p_GetExp(p,k+1, r)*(long)w[k];
670 }
671 //break;
672 return j;
673 }
674 case ringorder_c: /* nothing to do*/
675 case ringorder_C: /* nothing to do*/
676 case ringorder_S: /* nothing to do*/
677 case ringorder_s: /* nothing to do*/
678 case ringorder_IS: /* nothing to do */
679 case ringorder_unspec: /* to make clang happy, does not occur*/
680 case ringorder_no: /* to make clang happy, does not occur*/
681 case ringorder_L: /* to make clang happy, does not occur*/
682 case ringorder_aa: /* ignored by p_WTotaldegree*/
683 break;
684 /* no default: all orderings covered */
685 }
686 }
687 return j;
688}
689
690long p_DegW(poly p, const int *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
705int 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
714long p_WDegree(poly p, const ring r)
715{
716 if (r->firstwv==NULL) return p_Totaldegree(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*/
739long pLDeg0(poly p,int *l, const ring r)
740{
741 p_CheckPolyRing(p, r);
742 long unsigned 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*/
770long 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 long unsigned 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*/
811long pLDegb(poly p,int *l, const ring r)
812{
813 p_CheckPolyRing(p, r);
814 long unsigned 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*/
841long pLDeg1(poly p,int *l, const ring r)
842{
843 p_CheckPolyRing(p, r);
844 long unsigned 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*/
877long 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 unsigned 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
910long pLDeg1_Deg(poly p,int *l, const ring r)
911{
912 assume(r->pFDeg == p_Deg);
913 p_CheckPolyRing(p, r);
914 long unsigned 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
941long 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 unsigned 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
975long pLDeg1_Totaldegree(poly p,int *l, const ring r)
976{
977 p_CheckPolyRing(p, r);
978 long unsigned 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
1005long 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 unsigned 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
1038long pLDeg1_WFirstTotalDegree(poly p,int *l, const ring r)
1039{
1040 p_CheckPolyRing(p, r);
1041 long unsigned k= p_GetComp(p, r);
1042 int ll=1;
1043 long t,max;
1044
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
1068long pLDeg1c_WFirstTotalDegree(poly p,int *l, const ring r)
1069{
1070 p_CheckPolyRing(p, r);
1071 int ll=1;
1072 long t,max;
1073
1075 if (rIsSyzIndexRing(r))
1076 {
1077 long unsigned 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
1106static inline unsigned long
1107p_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
1132static inline unsigned long
1133p_GetMaxExpL2(unsigned long l1, unsigned long l2, const ring r)
1134{
1135 return p_GetMaxExpL2(l1, l2, r, r->ExpPerLong);
1136}
1137
1138poly 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
1175unsigned 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
1208BOOLEAN 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*/
1226int 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*/
1247int 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
1267int 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*/
1297poly 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*/
1313poly 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
1320void p_Split(poly p, poly *h)
1321{
1322 *h=pNext(p);
1323 pNext(p)=NULL;
1324}
1325
1326/*2
1327* pair has no common factor ? or is no polynomial
1328*/
1329BOOLEAN 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
1345BOOLEAN p_HasNotCFRing(poly p1, poly p2, const ring r)
1346{
1347
1348 if (p_GetComp(p1,r) > 0 || p_GetComp(p2,r) > 0)
1349 return FALSE;
1350 int i = rVar(r);
1351 loop
1352 {
1353 if ((p_GetExp(p1, i, r) > 0) && (p_GetExp(p2, i, r) > 0))
1354 return FALSE;
1355 i--;
1356 if (i == 0) {
1357 if (n_DivBy(pGetCoeff(p1), pGetCoeff(p2), r->cf) ||
1358 n_DivBy(pGetCoeff(p2), pGetCoeff(p1), r->cf)) {
1359 return FALSE;
1360 } else {
1361 return TRUE;
1362 }
1363 }
1364 }
1365}
1366
1367/*2
1368* convert monomial given as string to poly, e.g. 1x3y5z
1369*/
1370const char * p_Read(const char *st, poly &rc, const ring r)
1371{
1372 if (r==NULL) { rc=NULL;return st;}
1373 int i,j;
1374 rc = p_Init(r);
1375 const char *s = n_Read(st,&(p_GetCoeff(rc, r)),r->cf);
1376 if (s==st)
1377 /* i.e. it does not start with a coeff: test if it is a ringvar*/
1378 {
1379 j = r_IsRingVar(s,r->names,r->N);
1380 if (j >= 0)
1381 {
1382 p_IncrExp(rc,1+j,r);
1383 while (*s!='\0') s++;
1384 goto done;
1385 }
1386 }
1387 while (*s!='\0')
1388 {
1389 char ss[2];
1390 ss[0] = *s++;
1391 ss[1] = '\0';
1392 j = r_IsRingVar(ss,r->names,r->N);
1393 if (j >= 0)
1394 {
1395 const char *s_save=s;
1396 s = eati(s,&i);
1397 if (((unsigned long)i) > r->bitmask/2)
1398 {
1399 // exponent to large: it is not a monomial
1400 p_LmDelete(&rc,r);
1401 return s_save;
1402 }
1403 p_AddExp(rc,1+j, (long)i, r);
1404 }
1405 else
1406 {
1407 // 1st char of is not a varname
1408 // We return the parsed polynomial nevertheless. This is needed when
1409 // we are parsing coefficients in a rational function field.
1410 s--;
1411 break;
1412 }
1413 }
1414done:
1415 if (n_IsZero(pGetCoeff(rc),r->cf)) p_LmDelete(&rc,r);
1416 else
1417 {
1418#ifdef HAVE_PLURAL
1419 // in super-commutative ring
1420 // squares of anti-commutative variables are zeroes!
1421 if(rIsSCA(r))
1422 {
1423 const unsigned int iFirstAltVar = scaFirstAltVar(r);
1424 const unsigned int iLastAltVar = scaLastAltVar(r);
1425
1426 assume(rc != NULL);
1427
1428 for(unsigned int k = iFirstAltVar; k <= iLastAltVar; k++)
1429 if( p_GetExp(rc, k, r) > 1 )
1430 {
1431 p_LmDelete(&rc, r);
1432 goto finish;
1433 }
1434 }
1435#endif
1436
1437 p_Setm(rc,r);
1438 }
1439finish:
1440 return s;
1441}
1442poly p_mInit(const char *st, BOOLEAN &ok, const ring r)
1443{
1444 poly p;
1445 const char *s=p_Read(st,p,r);
1446 if (*s!='\0')
1447 {
1448 if ((s!=st)&&isdigit(st[0]))
1449 {
1451 }
1452 ok=FALSE;
1453 if (p!=NULL)
1454 {
1455 if (pGetCoeff(p)==NULL) p_LmFree(p,r);
1456 else p_LmDelete(p,r);
1457 }
1458 return NULL;
1459 }
1460 p_Test(p,r);
1461 ok=!errorreported;
1462 return p;
1463}
1464
1465/*2
1466* returns a polynomial representing the number n
1467* destroys n
1468*/
1469poly p_NSet(number n, const ring r)
1470{
1471 if (n_IsZero(n,r->cf))
1472 {
1473 n_Delete(&n, r->cf);
1474 return NULL;
1475 }
1476 else
1477 {
1478 poly rc = p_Init(r);
1479 pSetCoeff0(rc,n);
1480 return rc;
1481 }
1482}
1483/*2
1484* assumes that LM(a) = LM(b)*m, for some monomial m,
1485* returns the multiplicant m,
1486* leaves a and b unmodified
1487*/
1488poly p_MDivide(poly a, poly b, const ring r)
1489{
1490 assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0));
1491 int i;
1492 poly result = p_Init(r);
1493
1494 for(i=(int)r->N; i; i--)
1495 p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r);
1496 p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r);
1497 p_Setm(result,r);
1498 return result;
1499}
1500
1501poly p_Div_nn(poly p, const number n, const ring r)
1502{
1503 pAssume(!n_IsZero(n,r->cf));
1504 p_Test(p, r);
1505 poly result = p;
1506 poly prev = NULL;
1507 while (p!=NULL)
1508 {
1509 number nc = n_Div(pGetCoeff(p),n,r->cf);
1510 if (!n_IsZero(nc,r->cf))
1511 {
1512 p_SetCoeff(p,nc,r);
1513 prev=p;
1514 pIter(p);
1515 }
1516 else
1517 {
1518 if (prev==NULL)
1519 {
1520 p_LmDelete(&result,r);
1521 p=result;
1522 }
1523 else
1524 {
1525 p_LmDelete(&pNext(prev),r);
1526 p=pNext(prev);
1527 }
1528 }
1529 }
1530 p_Test(result,r);
1531 return(result);
1532}
1533
1534poly p_Div_mm(poly p, const poly m, const ring r)
1535{
1536 p_Test(p, r);
1537 p_Test(m, r);
1538 poly result = p;
1539 poly prev = NULL;
1540 number n=pGetCoeff(m);
1541 while (p!=NULL)
1542 {
1543 number nc = n_Div(pGetCoeff(p),n,r->cf);
1544 n_Normalize(nc,r->cf);
1545 if (!n_IsZero(nc,r->cf))
1546 {
1547 p_SetCoeff(p,nc,r);
1548 prev=p;
1549 p_ExpVectorSub(p,m,r);
1550 pIter(p);
1551 }
1552 else
1553 {
1554 if (prev==NULL)
1555 {
1556 p_LmDelete(&result,r);
1557 p=result;
1558 }
1559 else
1560 {
1561 p_LmDelete(&pNext(prev),r);
1562 p=pNext(prev);
1563 }
1564 }
1565 }
1566 p_Test(result,r);
1567 return(result);
1568}
1569
1570/*2
1571* divides a by the monomial b, ignores monomials which are not divisible
1572* assumes that b is not NULL, destroyes b
1573*/
1574poly p_DivideM(poly a, poly b, const ring r)
1575{
1576 if (a==NULL) { p_Delete(&b,r); return NULL; }
1577 poly result=a;
1578
1579 if(!p_IsConstant(b,r))
1580 {
1581 if (rIsNCRing(r))
1582 {
1583 WerrorS("p_DivideM not implemented for non-commuative rings");
1584 return NULL;
1585 }
1586 poly prev=NULL;
1587 while (a!=NULL)
1588 {
1589 if (p_DivisibleBy(b,a,r))
1590 {
1591 p_ExpVectorSub(a,b,r);
1592 prev=a;
1593 pIter(a);
1594 }
1595 else
1596 {
1597 if (prev==NULL)
1598 {
1599 p_LmDelete(&result,r);
1600 a=result;
1601 }
1602 else
1603 {
1604 p_LmDelete(&pNext(prev),r);
1605 a=pNext(prev);
1606 }
1607 }
1608 }
1609 }
1610 if (result!=NULL)
1611 {
1612 number inv=pGetCoeff(b);
1613 //if ((!rField_is_Ring(r)) || n_IsUnit(inv,r->cf))
1614 if (rField_is_Zp(r))
1615 {
1616 inv = n_Invers(inv,r->cf);
1617 __p_Mult_nn(result,inv,r);
1618 n_Delete(&inv, r->cf);
1619 }
1620 else
1621 {
1622 result = p_Div_nn(result,inv,r);
1623 }
1624 }
1625 p_Delete(&b, r);
1626 return result;
1627}
1628
1629poly pp_DivideM(poly a, poly b, const ring r)
1630{
1631 if (a==NULL) { return NULL; }
1632 // TODO: better implementation without copying a,b
1633 return p_DivideM(p_Copy(a,r),p_Head(b,r),r);
1634}
1635
1636#ifdef HAVE_RINGS
1637/* TRUE iff LT(f) | LT(g) */
1638BOOLEAN p_DivisibleByRingCase(poly f, poly g, const ring r)
1639{
1640 int exponent;
1641 for(int i = (int)rVar(r); i>0; i--)
1642 {
1643 exponent = p_GetExp(g, i, r) - p_GetExp(f, i, r);
1644 if (exponent < 0) return FALSE;
1645 }
1646 return n_DivBy(pGetCoeff(g), pGetCoeff(f), r->cf);
1647}
1648#endif
1649
1650// returns the LCM of the head terms of a and b in *m
1651void p_Lcm(const poly a, const poly b, poly m, const ring r)
1652{
1653 for (int i=r->N; i; --i)
1654 p_SetExp(m,i, si_max( p_GetExp(a,i,r), p_GetExp(b,i,r)),r);
1655
1656 p_SetComp(m, si_max(p_GetComp(a,r), p_GetComp(b,r)),r);
1657 /* Don't do a pSetm here, otherwise hres/lres chockes */
1658}
1659
1660poly p_Lcm(const poly a, const poly b, const ring r)
1661{
1662 poly m=p_Init(r);
1663 p_Lcm(a, b, m, r);
1664 p_Setm(m,r);
1665 return(m);
1666}
1667
1668#ifdef HAVE_RATGRING
1669/*2
1670* returns the rational LCM of the head terms of a and b
1671* without coefficient!!!
1672*/
1673poly p_LcmRat(const poly a, const poly b, const long lCompM, const ring r)
1674{
1675 poly m = // p_One( r);
1676 p_Init(r);
1677
1678// const int (currRing->N) = r->N;
1679
1680 // for (int i = (currRing->N); i>=r->real_var_start; i--)
1681 for (int i = r->real_var_end; i>=r->real_var_start; i--)
1682 {
1683 const int lExpA = p_GetExp (a, i, r);
1684 const int lExpB = p_GetExp (b, i, r);
1685
1686 p_SetExp (m, i, si_max(lExpA, lExpB), r);
1687 }
1688
1689 p_SetComp (m, lCompM, r);
1690 p_Setm(m,r);
1691 n_New(&(p_GetCoeff(m, r)), r);
1692
1693 return(m);
1694};
1695
1696void p_LmDeleteAndNextRat(poly *p, int ishift, ring r)
1697{
1698 /* modifies p*/
1699 // Print("start: "); Print(" "); p_wrp(*p,r);
1700 p_LmCheckPolyRing2(*p, r);
1701 poly q = p_Head(*p,r);
1702 const long cmp = p_GetComp(*p, r);
1703 while ( ( (*p)!=NULL ) && ( p_Comp_k_n(*p, q, ishift+1, r) ) && (p_GetComp(*p, r) == cmp) )
1704 {
1705 p_LmDelete(p,r);
1706 // Print("while: ");p_wrp(*p,r);Print(" ");
1707 }
1708 // p_wrp(*p,r);Print(" ");
1709 // PrintS("end\n");
1710 p_LmDelete(&q,r);
1711}
1712
1713
1714/* returns x-coeff of p, i.e. a poly in x, s.t. corresponding xd-monomials
1715have the same D-part and the component 0
1716does not destroy p
1717*/
1718poly p_GetCoeffRat(poly p, int ishift, ring r)
1719{
1720 poly q = pNext(p);
1721 poly res; // = p_Head(p,r);
1722 res = p_GetExp_k_n(p, ishift+1, r->N, r); // does pSetm internally
1723 p_SetCoeff(res,n_Copy(p_GetCoeff(p,r),r),r);
1724 poly s;
1725 long cmp = p_GetComp(p, r);
1726 while ( (q!= NULL) && (p_Comp_k_n(p, q, ishift+1, r)) && (p_GetComp(q, r) == cmp) )
1727 {
1728 s = p_GetExp_k_n(q, ishift+1, r->N, r);
1729 p_SetCoeff(s,n_Copy(p_GetCoeff(q,r),r),r);
1730 res = p_Add_q(res,s,r);
1731 q = pNext(q);
1732 }
1733 cmp = 0;
1734 p_SetCompP(res,cmp,r);
1735 return res;
1736}
1737
1738
1739
1740void p_ContentRat(poly &ph, const ring r)
1741// changes ph
1742// for rat coefficients in K(x1,..xN)
1743{
1744 // init array of RatLeadCoeffs
1745 // poly p_GetCoeffRat(poly p, int ishift, ring r);
1746
1747 int len=pLength(ph);
1748 poly *C = (poly *)omAlloc0((len+1)*sizeof(poly)); //rat coeffs
1749 poly *LM = (poly *)omAlloc0((len+1)*sizeof(poly)); // rat lead terms
1750 int *D = (int *)omAlloc0((len+1)*sizeof(int)); //degrees of coeffs
1751 int *L = (int *)omAlloc0((len+1)*sizeof(int)); //lengths of coeffs
1752 int k = 0;
1753 poly p = p_Copy(ph, r); // ph will be needed below
1754 int mintdeg = p_Totaldegree(p, r);
1755 int minlen = len;
1756 int dd = 0; int i;
1757 int HasConstantCoef = 0;
1758 int is = r->real_var_start - 1;
1759 while (p!=NULL)
1760 {
1761 LM[k] = p_GetExp_k_n(p,1,is, r); // need LmRat istead of p_HeadRat(p, is, currRing); !
1762 C[k] = p_GetCoeffRat(p, is, r);
1763 D[k] = p_Totaldegree(C[k], r);
1764 mintdeg = si_min(mintdeg,D[k]);
1765 L[k] = pLength(C[k]);
1766 minlen = si_min(minlen,L[k]);
1767 if (p_IsConstant(C[k], r))
1768 {
1769 // C[k] = const, so the content will be numerical
1770 HasConstantCoef = 1;
1771 // smth like goto cleanup and return(pContent(p));
1772 }
1773 p_LmDeleteAndNextRat(&p, is, r);
1774 k++;
1775 }
1776
1777 // look for 1 element of minimal degree and of minimal length
1778 k--;
1779 poly d;
1780 int mindeglen = len;
1781 if (k<=0) // this poly is not a ratgring poly -> pContent
1782 {
1783 p_Delete(&C[0], r);
1784 p_Delete(&LM[0], r);
1785 p_ContentForGB(ph, r);
1786 goto cleanup;
1787 }
1788
1789 int pmindeglen;
1790 for(i=0; i<=k; i++)
1791 {
1792 if (D[i] == mintdeg)
1793 {
1794 if (L[i] < mindeglen)
1795 {
1796 mindeglen=L[i];
1797 pmindeglen = i;
1798 }
1799 }
1800 }
1801 d = p_Copy(C[pmindeglen], r);
1802 // there are dd>=1 mindeg elements
1803 // and pmideglen is the coordinate of one of the smallest among them
1804
1805 // poly g = singclap_gcd(p_Copy(p,r),p_Copy(q,r));
1806 // return naGcd(d,d2,currRing);
1807
1808 // adjoin pContentRat here?
1809 for(i=0; i<=k; i++)
1810 {
1811 d=singclap_gcd(d,p_Copy(C[i], r), r);
1812 if (p_Totaldegree(d, r)==0)
1813 {
1814 // cleanup, pContent, return
1815 p_Delete(&d, r);
1816 for(;k>=0;k--)
1817 {
1818 p_Delete(&C[k], r);
1819 p_Delete(&LM[k], r);
1820 }
1821 p_ContentForGB(ph, r);
1822 goto cleanup;
1823 }
1824 }
1825 for(i=0; i<=k; i++)
1826 {
1827 poly h=singclap_pdivide(C[i],d, r);
1828 p_Delete(&C[i], r);
1829 C[i]=h;
1830 }
1831
1832 // zusammensetzen,
1833 p=NULL; // just to be sure
1834 for(i=0; i<=k; i++)
1835 {
1836 p = p_Add_q(p, p_Mult_q(C[i],LM[i], r), r);
1837 C[i]=NULL; LM[i]=NULL;
1838 }
1839 p_Delete(&ph, r); // do not need it anymore
1840 ph = p;
1841 // aufraeumen, return
1842cleanup:
1843 omFree(C);
1844 omFree(LM);
1845 omFree(D);
1846 omFree(L);
1847}
1848
1849
1850#endif
1851
1852
1853/* assumes that p and divisor are univariate polynomials in r,
1854 mentioning the same variable;
1855 assumes divisor != NULL;
1856 p may be NULL;
1857 assumes a global monomial ordering in r;
1858 performs polynomial division of p by divisor:
1859 - afterwards p contains the remainder of the division, i.e.,
1860 p_before = result * divisor + p_afterwards;
1861 - if needResult == TRUE, then the method computes and returns 'result',
1862 otherwise NULL is returned (This parametrization can be used when
1863 one is only interested in the remainder of the division. In this
1864 case, the method will be slightly faster.)
1865 leaves divisor unmodified */
1866poly p_PolyDiv(poly &p, const poly divisor, const BOOLEAN needResult, const ring r)
1867{
1868 assume(divisor != NULL);
1869 if (p == NULL) return NULL;
1870
1871 poly result = NULL;
1872 number divisorLC = p_GetCoeff(divisor, r);
1873 int divisorLE = p_GetExp(divisor, 1, r);
1874 while ((p != NULL) && (p_Deg(p, r) >= p_Deg(divisor, r)))
1875 {
1876 /* determine t = LT(p) / LT(divisor) */
1877 poly t = p_ISet(1, r);
1878 number c = n_Div(p_GetCoeff(p, r), divisorLC, r->cf);
1879 n_Normalize(c,r->cf);
1880 p_SetCoeff(t, c, r);
1881 int e = p_GetExp(p, 1, r) - divisorLE;
1882 p_SetExp(t, 1, e, r);
1883 p_Setm(t, r);
1884 if (needResult) result = p_Add_q(result, p_Copy(t, r), r);
1885 p = p_Add_q(p, p_Neg(p_Mult_q(t, p_Copy(divisor, r), r), r), r);
1886 }
1887 return result;
1888}
1889
1890/*2
1891* returns the partial differentiate of a by the k-th variable
1892* does not destroy the input
1893*/
1894poly p_Diff(poly a, int k, const ring r)
1895{
1896 poly res, f, last;
1897 number t;
1898
1899 last = res = NULL;
1900 while (a!=NULL)
1901 {
1902 if (p_GetExp(a,k,r)!=0)
1903 {
1904 f = p_LmInit(a,r);
1905 t = n_Init(p_GetExp(a,k,r),r->cf);
1906 pSetCoeff0(f,n_Mult(t,pGetCoeff(a),r->cf));
1907 n_Delete(&t,r->cf);
1908 if (n_IsZero(pGetCoeff(f),r->cf))
1909 p_LmDelete(&f,r);
1910 else
1911 {
1912 p_DecrExp(f,k,r);
1913 p_Setm(f,r);
1914 if (res==NULL)
1915 {
1916 res=last=f;
1917 }
1918 else
1919 {
1920 pNext(last)=f;
1921 last=f;
1922 }
1923 }
1924 }
1925 pIter(a);
1926 }
1927 return res;
1928}
1929
1930static poly p_DiffOpM(poly a, poly b,BOOLEAN multiply, const ring r)
1931{
1932 int i,j,s;
1933 number n,h,hh;
1934 poly p=p_One(r);
1935 n=n_Mult(pGetCoeff(a),pGetCoeff(b),r->cf);
1936 for(i=rVar(r);i>0;i--)
1937 {
1938 s=p_GetExp(b,i,r);
1939 if (s<p_GetExp(a,i,r))
1940 {
1941 n_Delete(&n,r->cf);
1942 p_LmDelete(&p,r);
1943 return NULL;
1944 }
1945 if (multiply)
1946 {
1947 for(j=p_GetExp(a,i,r); j>0;j--)
1948 {
1949 h = n_Init(s,r->cf);
1950 hh=n_Mult(n,h,r->cf);
1951 n_Delete(&h,r->cf);
1952 n_Delete(&n,r->cf);
1953 n=hh;
1954 s--;
1955 }
1956 p_SetExp(p,i,s,r);
1957 }
1958 else
1959 {
1960 p_SetExp(p,i,s-p_GetExp(a,i,r),r);
1961 }
1962 }
1963 p_Setm(p,r);
1964 /*if (multiply)*/ p_SetCoeff(p,n,r);
1965 if (n_IsZero(n,r->cf)) p=p_LmDeleteAndNext(p,r); // return NULL as p is a monomial
1966 return p;
1967}
1968
1969poly p_DiffOp(poly a, poly b,BOOLEAN multiply, const ring r)
1970{
1971 poly result=NULL;
1972 poly h;
1973 for(;a!=NULL;pIter(a))
1974 {
1975 for(h=b;h!=NULL;pIter(h))
1976 {
1977 result=p_Add_q(result,p_DiffOpM(a,h,multiply,r),r);
1978 }
1979 }
1980 return result;
1981}
1982/*2
1983* subtract p2 from p1, p1 and p2 are destroyed
1984* do not put attention on speed: the procedure is only used in the interpreter
1985*/
1986poly p_Sub(poly p1, poly p2, const ring r)
1987{
1988 return p_Add_q(p1, p_Neg(p2,r),r);
1989}
1990
1991/*3
1992* compute for a monomial m
1993* the power m^exp, exp > 1
1994* destroys p
1995*/
1996static poly p_MonPower(poly p, int exp, const ring r)
1997{
1998 int i;
1999
2000 if(!n_IsOne(pGetCoeff(p),r->cf))
2001 {
2002 number x, y;
2003 y = pGetCoeff(p);
2004 n_Power(y,exp,&x,r->cf);
2005 n_Delete(&y,r->cf);
2006 pSetCoeff0(p,x);
2007 }
2008 for (i=rVar(r); i!=0; i--)
2009 {
2010 p_MultExp(p,i, exp,r);
2011 }
2012 p_Setm(p,r);
2013 return p;
2014}
2015
2016/*3
2017* compute for monomials p*q
2018* destroys p, keeps q
2019*/
2020static void p_MonMult(poly p, poly q, const ring r)
2021{
2022 number x, y;
2023
2024 y = pGetCoeff(p);
2025 x = n_Mult(y,pGetCoeff(q),r->cf);
2026 n_Delete(&y,r->cf);
2027 pSetCoeff0(p,x);
2028 //for (int i=pVariables; i!=0; i--)
2029 //{
2030 // pAddExp(p,i, pGetExp(q,i));
2031 //}
2032 //p->Order += q->Order;
2033 p_ExpVectorAdd(p,q,r);
2034}
2035
2036/*3
2037* compute for monomials p*q
2038* keeps p, q
2039*/
2040static poly p_MonMultC(poly p, poly q, const ring rr)
2041{
2042 number x;
2043 poly r = p_Init(rr);
2044
2045 x = n_Mult(pGetCoeff(p),pGetCoeff(q),rr->cf);
2046 pSetCoeff0(r,x);
2047 p_ExpVectorSum(r,p, q, rr);
2048 return r;
2049}
2050
2051/*3
2052* create binomial coef.
2053*/
2054static number* pnBin(int exp, const ring r)
2055{
2056 int e, i, h;
2057 number x, y, *bin=NULL;
2058
2059 x = n_Init(exp,r->cf);
2060 if (n_IsZero(x,r->cf))
2061 {
2062 n_Delete(&x,r->cf);
2063 return bin;
2064 }
2065 h = (exp >> 1) + 1;
2066 bin = (number *)omAlloc0(h*sizeof(number));
2067 bin[1] = x;
2068 if (exp < 4)
2069 return bin;
2070 i = exp - 1;
2071 for (e=2; e<h; e++)
2072 {
2073 x = n_Init(i,r->cf);
2074 i--;
2075 y = n_Mult(x,bin[e-1],r->cf);
2076 n_Delete(&x,r->cf);
2077 x = n_Init(e,r->cf);
2078 bin[e] = n_ExactDiv(y,x,r->cf);
2079 n_Delete(&x,r->cf);
2080 n_Delete(&y,r->cf);
2081 }
2082 return bin;
2083}
2084
2085static void pnFreeBin(number *bin, int exp,const coeffs r)
2086{
2087 int e, h = (exp >> 1) + 1;
2088
2089 if (bin[1] != NULL)
2090 {
2091 for (e=1; e<h; e++)
2092 n_Delete(&(bin[e]),r);
2093 }
2094 omFreeSize((ADDRESS)bin, h*sizeof(number));
2095}
2096
2097/*
2098* compute for a poly p = head+tail, tail is monomial
2099* (head + tail)^exp, exp > 1
2100* with binomial coef.
2101*/
2102static poly p_TwoMonPower(poly p, int exp, const ring r)
2103{
2104 int eh, e;
2105 long al;
2106 poly *a;
2107 poly tail, b, res, h;
2108 number x;
2109 number *bin = pnBin(exp,r);
2110
2111 tail = pNext(p);
2112 if (bin == NULL)
2113 {
2114 p_MonPower(p,exp,r);
2115 p_MonPower(tail,exp,r);
2116 p_Test(p,r);
2117 return p;
2118 }
2119 eh = exp >> 1;
2120 al = (exp + 1) * sizeof(poly);
2121 a = (poly *)omAlloc(al);
2122 a[1] = p;
2123 for (e=1; e<exp; e++)
2124 {
2125 a[e+1] = p_MonMultC(a[e],p,r);
2126 }
2127 res = a[exp];
2128 b = p_Head(tail,r);
2129 for (e=exp-1; e>eh; e--)
2130 {
2131 h = a[e];
2132 x = n_Mult(bin[exp-e],pGetCoeff(h),r->cf);
2133 p_SetCoeff(h,x,r);
2134 p_MonMult(h,b,r);
2135 res = pNext(res) = h;
2136 p_MonMult(b,tail,r);
2137 }
2138 for (e=eh; e!=0; e--)
2139 {
2140 h = a[e];
2141 x = n_Mult(bin[e],pGetCoeff(h),r->cf);
2142 p_SetCoeff(h,x,r);
2143 p_MonMult(h,b,r);
2144 res = pNext(res) = h;
2145 p_MonMult(b,tail,r);
2146 }
2147 p_LmDelete(&tail,r);
2148 pNext(res) = b;
2149 pNext(b) = NULL;
2150 res = a[exp];
2151 omFreeSize((ADDRESS)a, al);
2152 pnFreeBin(bin, exp, r->cf);
2153// tail=res;
2154// while((tail!=NULL)&&(pNext(tail)!=NULL))
2155// {
2156// if(nIsZero(pGetCoeff(pNext(tail))))
2157// {
2158// pLmDelete(&pNext(tail));
2159// }
2160// else
2161// pIter(tail);
2162// }
2163 p_Test(res,r);
2164 return res;
2165}
2166
2167static poly p_Pow(poly p, int i, const ring r)
2168{
2169 poly rc = p_Copy(p,r);
2170 i -= 2;
2171 do
2172 {
2173 rc = p_Mult_q(rc,p_Copy(p,r),r);
2174 p_Normalize(rc,r);
2175 i--;
2176 }
2177 while (i != 0);
2178 return p_Mult_q(rc,p,r);
2179}
2180
2181static poly p_Pow_charp(poly p, int i, const ring r)
2182{
2183 //assume char_p == i
2184 poly h=p;
2185 while(h!=NULL) { p_MonPower(h,i,r);pIter(h);}
2186 return p;
2187}
2188
2189/*2
2190* returns the i-th power of p
2191* p will be destroyed
2192*/
2193poly p_Power(poly p, int i, const ring r)
2194{
2195 poly rc=NULL;
2196
2197 if (i==0)
2198 {
2199 p_Delete(&p,r);
2200 return p_One(r);
2201 }
2202
2203 if(p!=NULL)
2204 {
2205 if ( (i > 0) && ((unsigned long ) i > (r->bitmask))
2206 #ifdef HAVE_SHIFTBBA
2207 && (!rIsLPRing(r))
2208 #endif
2209 )
2210 {
2211 Werror("exponent %d is too large, max. is %ld",i,r->bitmask);
2212 return NULL;
2213 }
2214 switch (i)
2215 {
2216// cannot happen, see above
2217// case 0:
2218// {
2219// rc=pOne();
2220// pDelete(&p);
2221// break;
2222// }
2223 case 1:
2224 rc=p;
2225 break;
2226 case 2:
2227 rc=p_Mult_q(p_Copy(p,r),p,r);
2228 break;
2229 default:
2230 if (i < 0)
2231 {
2232 p_Delete(&p,r);
2233 return NULL;
2234 }
2235 else
2236 {
2237#ifdef HAVE_PLURAL
2238 if (rIsNCRing(r)) /* in the NC case nothing helps :-( */
2239 {
2240 int j=i;
2241 rc = p_Copy(p,r);
2242 while (j>1)
2243 {
2244 rc = p_Mult_q(p_Copy(p,r),rc,r);
2245 j--;
2246 }
2247 p_Delete(&p,r);
2248 return rc;
2249 }
2250#endif
2251 rc = pNext(p);
2252 if (rc == NULL)
2253 return p_MonPower(p,i,r);
2254 /* else: binom ?*/
2255 int char_p=rInternalChar(r);
2256 if ((char_p>0) && (i>char_p)
2257 && ((rField_is_Zp(r,char_p)
2258 || (rField_is_Zp_a(r,char_p)))))
2259 {
2260 poly h=p_Pow_charp(p_Copy(p,r),char_p,r);
2261 int rest=i-char_p;
2262 while (rest>=char_p)
2263 {
2264 rest-=char_p;
2265 h=p_Mult_q(h,p_Pow_charp(p_Copy(p,r),char_p,r),r);
2266 }
2267 poly res=h;
2268 if (rest>0)
2269 res=p_Mult_q(p_Power(p_Copy(p,r),rest,r),h,r);
2270 p_Delete(&p,r);
2271 return res;
2272 }
2273 if ((pNext(rc) != NULL)
2274 || rField_is_Ring(r)
2275 )
2276 return p_Pow(p,i,r);
2277 if ((char_p==0) || (i<=char_p))
2278 return p_TwoMonPower(p,i,r);
2279 return p_Pow(p,i,r);
2280 }
2281 /*end default:*/
2282 }
2283 }
2284 return rc;
2285}
2286
2287/* --------------------------------------------------------------------------------*/
2288/* content suff */
2289//number p_InitContent(poly ph, const ring r);
2290
2291void p_Content(poly ph, const ring r)
2292{
2293 if (ph==NULL) return;
2294 const coeffs cf=r->cf;
2295 if (pNext(ph)==NULL)
2296 {
2297 p_SetCoeff(ph,n_Init(1,cf),r);
2298 return;
2299 }
2300 if ((cf->cfSubringGcd==ndGcd)
2301 || (cf->cfGcd==ndGcd)) /* trivial gcd*/
2302 return;
2303 number h;
2304 if ((rField_is_Q(r))
2305 || (rField_is_Q_a(r))
2306 || (rField_is_Zp_a)(r)
2307 || (rField_is_Z(r))
2308 )
2309 {
2310 h=p_InitContent(ph,r); /* first guess of a gcd of all coeffs */
2311 }
2312 else
2313 {
2314 h=n_Copy(pGetCoeff(ph),cf);
2315 }
2316 poly p;
2317 if(n_IsOne(h,cf))
2318 {
2319 goto content_finish;
2320 }
2321 p=ph;
2322 // take the SubringGcd of all coeffs
2323 while (p!=NULL)
2324 {
2326 number d=n_SubringGcd(h,pGetCoeff(p),cf);
2327 n_Delete(&h,cf);
2328 h = d;
2329 if(n_IsOne(h,cf))
2330 {
2331 goto content_finish;
2332 }
2333 pIter(p);
2334 }
2335 // if found<>1, divide by it
2336 p = ph;
2337 while (p!=NULL)
2338 {
2339 number d = n_ExactDiv(pGetCoeff(p),h,cf);
2340 p_SetCoeff(p,d,r);
2341 pIter(p);
2342 }
2343content_finish:
2344 n_Delete(&h,r->cf);
2345 // and last: check leading sign:
2346 if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2347}
2348
2349void p_Content_n(poly ph, number &c,const ring r)
2350{
2351 const coeffs cf=r->cf;
2352 if (ph==NULL)
2353 {
2354 c=n_Init(1,cf);
2355 return;
2356 }
2357 if (pNext(ph)==NULL)
2358 {
2359 c=pGetCoeff(ph);
2360 p_SetCoeff0(ph,n_Init(1,cf),r);
2361 }
2362 if ((cf->cfSubringGcd==ndGcd)
2363 || (cf->cfGcd==ndGcd)) /* trivial gcd*/
2364 {
2365 c=n_Init(1,r->cf);
2366 return;
2367 }
2368 number h;
2369 if ((rField_is_Q(r))
2370 || (rField_is_Q_a(r))
2371 || (rField_is_Zp_a)(r)
2372 || (rField_is_Z(r))
2373 )
2374 {
2375 h=p_InitContent(ph,r); /* first guess of a gcd of all coeffs */
2376 }
2377 else
2378 {
2379 h=n_Copy(pGetCoeff(ph),cf);
2380 }
2381 poly p;
2382 if(n_IsOne(h,cf))
2383 {
2384 goto content_finish;
2385 }
2386 p=ph;
2387 // take the SubringGcd of all coeffs
2388 while (p!=NULL)
2389 {
2391 number d=n_SubringGcd(h,pGetCoeff(p),cf);
2392 n_Delete(&h,cf);
2393 h = d;
2394 if(n_IsOne(h,cf))
2395 {
2396 goto content_finish;
2397 }
2398 pIter(p);
2399 }
2400 // if found<>1, divide by it
2401 p = ph;
2402 while (p!=NULL)
2403 {
2404 number d = n_ExactDiv(pGetCoeff(p),h,cf);
2405 p_SetCoeff(p,d,r);
2406 pIter(p);
2407 }
2408content_finish:
2409 c=h;
2410 // and last: check leading sign:
2411 if(!n_GreaterZero(pGetCoeff(ph),r->cf))
2412 {
2413 c = n_InpNeg(c,r->cf);
2414 ph = p_Neg(ph,r);
2415 }
2416}
2417
2418#define CLEARENUMERATORS 1
2419
2420void p_ContentForGB(poly ph, const ring r)
2421{
2422 if(TEST_OPT_CONTENTSB) return;
2423 assume( ph != NULL );
2424
2425 assume( r != NULL ); assume( r->cf != NULL );
2426
2427
2428#if CLEARENUMERATORS
2429 if( 0 )
2430 {
2431 const coeffs C = r->cf;
2432 // experimentall (recursive enumerator treatment) of alg. Ext!
2433 CPolyCoeffsEnumerator itr(ph);
2434 n_ClearContent(itr, r->cf);
2435
2436 p_Test(ph, r); n_Test(pGetCoeff(ph), C);
2437 assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
2438
2439 // if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2440 return;
2441 }
2442#endif
2443
2444
2445#ifdef HAVE_RINGS
2446 if (rField_is_Ring(r))
2447 {
2448 if (rField_has_Units(r))
2449 {
2450 number k = n_GetUnit(pGetCoeff(ph),r->cf);
2451 if (!n_IsOne(k,r->cf))
2452 {
2453 number tmpGMP = k;
2454 k = n_Invers(k,r->cf);
2455 n_Delete(&tmpGMP,r->cf);
2456 poly h = pNext(ph);
2457 p_SetCoeff(ph, n_Mult(pGetCoeff(ph), k,r->cf),r);
2458 while (h != NULL)
2459 {
2460 p_SetCoeff(h, n_Mult(pGetCoeff(h), k,r->cf),r);
2461 pIter(h);
2462 }
2463// assume( n_GreaterZero(pGetCoeff(ph),r->cf) );
2464// if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2465 }
2466 n_Delete(&k,r->cf);
2467 }
2468 return;
2469 }
2470#endif
2471 number h,d;
2472 poly p;
2473
2474 if(pNext(ph)==NULL)
2475 {
2476 p_SetCoeff(ph,n_Init(1,r->cf),r);
2477 }
2478 else
2479 {
2480 assume( pNext(ph) != NULL );
2481#if CLEARENUMERATORS
2482 if( nCoeff_is_Q(r->cf) )
2483 {
2484 // experimentall (recursive enumerator treatment) of alg. Ext!
2485 CPolyCoeffsEnumerator itr(ph);
2486 n_ClearContent(itr, r->cf);
2487
2488 p_Test(ph, r); n_Test(pGetCoeff(ph), r->cf);
2489 assume(n_GreaterZero(pGetCoeff(ph), r->cf)); // ??
2490
2491 // if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2492 return;
2493 }
2494#endif
2495
2496 n_Normalize(pGetCoeff(ph),r->cf);
2497 if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2498 if (rField_is_Q(r)||(getCoeffType(r->cf)==n_transExt)) // should not be used anymore if CLEARENUMERATORS is 1
2499 {
2500 h=p_InitContent(ph,r);
2501 p=ph;
2502 }
2503 else
2504 {
2505 h=n_Copy(pGetCoeff(ph),r->cf);
2506 p = pNext(ph);
2507 }
2508 while (p!=NULL)
2509 {
2510 n_Normalize(pGetCoeff(p),r->cf);
2511 d=n_SubringGcd(h,pGetCoeff(p),r->cf);
2512 n_Delete(&h,r->cf);
2513 h = d;
2514 if(n_IsOne(h,r->cf))
2515 {
2516 break;
2517 }
2518 pIter(p);
2519 }
2520 //number tmp;
2521 if(!n_IsOne(h,r->cf))
2522 {
2523 p = ph;
2524 while (p!=NULL)
2525 {
2526 //d = nDiv(pGetCoeff(p),h);
2527 //tmp = nExactDiv(pGetCoeff(p),h);
2528 //if (!nEqual(d,tmp))
2529 //{
2530 // StringSetS("** div0:");nWrite(pGetCoeff(p));StringAppendS("/");
2531 // nWrite(h);StringAppendS("=");nWrite(d);StringAppendS(" int:");
2532 // nWrite(tmp);Print(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
2533 //}
2534 //nDelete(&tmp);
2535 d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2536 p_SetCoeff(p,d,r);
2537 pIter(p);
2538 }
2539 }
2540 n_Delete(&h,r->cf);
2541 if (rField_is_Q_a(r))
2542 {
2543 // special handling for alg. ext.:
2544 if (getCoeffType(r->cf)==n_algExt)
2545 {
2546 h = n_Init(1, r->cf->extRing->cf);
2547 p=ph;
2548 while (p!=NULL)
2549 { // each monom: coeff in Q_a
2550 poly c_n_n=(poly)pGetCoeff(p);
2551 poly c_n=c_n_n;
2552 while (c_n!=NULL)
2553 { // each monom: coeff in Q
2554 d=n_NormalizeHelper(h,pGetCoeff(c_n),r->cf->extRing->cf);
2555 n_Delete(&h,r->cf->extRing->cf);
2556 h=d;
2557 pIter(c_n);
2558 }
2559 pIter(p);
2560 }
2561 /* h contains the 1/lcm of all denominators in c_n_n*/
2562 //n_Normalize(h,r->cf->extRing->cf);
2563 if(!n_IsOne(h,r->cf->extRing->cf))
2564 {
2565 p=ph;
2566 while (p!=NULL)
2567 { // each monom: coeff in Q_a
2568 poly c_n=(poly)pGetCoeff(p);
2569 while (c_n!=NULL)
2570 { // each monom: coeff in Q
2571 d=n_Mult(h,pGetCoeff(c_n),r->cf->extRing->cf);
2572 n_Normalize(d,r->cf->extRing->cf);
2573 n_Delete(&pGetCoeff(c_n),r->cf->extRing->cf);
2574 pGetCoeff(c_n)=d;
2575 pIter(c_n);
2576 }
2577 pIter(p);
2578 }
2579 }
2580 n_Delete(&h,r->cf->extRing->cf);
2581 }
2582 /*else
2583 {
2584 // special handling for rat. functions.:
2585 number hzz =NULL;
2586 p=ph;
2587 while (p!=NULL)
2588 { // each monom: coeff in Q_a (Z_a)
2589 fraction f=(fraction)pGetCoeff(p);
2590 poly c_n=NUM(f);
2591 if (hzz==NULL)
2592 {
2593 hzz=n_Copy(pGetCoeff(c_n),r->cf->extRing->cf);
2594 pIter(c_n);
2595 }
2596 while ((c_n!=NULL)&&(!n_IsOne(hzz,r->cf->extRing->cf)))
2597 { // each monom: coeff in Q (Z)
2598 d=n_Gcd(hzz,pGetCoeff(c_n),r->cf->extRing->cf);
2599 n_Delete(&hzz,r->cf->extRing->cf);
2600 hzz=d;
2601 pIter(c_n);
2602 }
2603 pIter(p);
2604 }
2605 // hzz contains the gcd of all numerators in f
2606 h=n_Invers(hzz,r->cf->extRing->cf);
2607 n_Delete(&hzz,r->cf->extRing->cf);
2608 n_Normalize(h,r->cf->extRing->cf);
2609 if(!n_IsOne(h,r->cf->extRing->cf))
2610 {
2611 p=ph;
2612 while (p!=NULL)
2613 { // each monom: coeff in Q_a (Z_a)
2614 fraction f=(fraction)pGetCoeff(p);
2615 NUM(f)=__p_Mult_nn(NUM(f),h,r->cf->extRing);
2616 p_Normalize(NUM(f),r->cf->extRing);
2617 pIter(p);
2618 }
2619 }
2620 n_Delete(&h,r->cf->extRing->cf);
2621 }*/
2622 }
2623 }
2624 if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2625}
2626
2627// Not yet?
2628#if 1 // currently only used by Singular/janet
2629void p_SimpleContent(poly ph, int smax, const ring r)
2630{
2631 if(TEST_OPT_CONTENTSB) return;
2632 if (ph==NULL) return;
2633 if (pNext(ph)==NULL)
2634 {
2635 p_SetCoeff(ph,n_Init(1,r->cf),r);
2636 return;
2637 }
2638 if (pNext(pNext(ph))==NULL)
2639 {
2640 return;
2641 }
2642 if (!(rField_is_Q(r))
2643 && (!rField_is_Q_a(r))
2644 && (!rField_is_Zp_a(r))
2645 && (!rField_is_Z(r))
2646 )
2647 {
2648 return;
2649 }
2650 number d=p_InitContent(ph,r);
2651 number h=d;
2652 if (n_Size(d,r->cf)<=smax)
2653 {
2654 n_Delete(&h,r->cf);
2655 //if (TEST_OPT_PROT) PrintS("G");
2656 return;
2657 }
2658
2659 poly p=ph;
2660 if (smax==1) smax=2;
2661 while (p!=NULL)
2662 {
2663#if 1
2664 d=n_SubringGcd(h,pGetCoeff(p),r->cf);
2665 n_Delete(&h,r->cf);
2666 h = d;
2667#else
2668 n_InpGcd(h,pGetCoeff(p),r->cf);
2669#endif
2670 if(n_Size(h,r->cf)<smax)
2671 {
2672 //if (TEST_OPT_PROT) PrintS("g");
2673 n_Delete(&h,r->cf);
2674 return;
2675 }
2676 pIter(p);
2677 }
2678 p = ph;
2679 if (!n_GreaterZero(pGetCoeff(p),r->cf)) h=n_InpNeg(h,r->cf);
2680 if(n_IsOne(h,r->cf))
2681 {
2682 n_Delete(&h,r->cf);
2683 return;
2684 }
2685 if (TEST_OPT_PROT) PrintS("c");
2686 while (p!=NULL)
2687 {
2688#if 1
2689 d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2690 p_SetCoeff(p,d,r);
2691#else
2692 STATISTIC(n_ExactDiv); nlInpExactDiv(pGetCoeff(p),h,r->cf); // no such function... ?
2693#endif
2694 pIter(p);
2695 }
2696 n_Delete(&h,r->cf);
2697}
2698#endif
2699
2700number p_InitContent(poly ph, const ring r)
2701// only for coefficients in Q and rational functions
2702#if 0
2703{
2705 assume(ph!=NULL);
2706 assume(pNext(ph)!=NULL);
2707 assume(rField_is_Q(r));
2708 if (pNext(pNext(ph))==NULL)
2709 {
2710 return n_GetNumerator(pGetCoeff(pNext(ph)),r->cf);
2711 }
2712 poly p=ph;
2713 number n1=n_GetNumerator(pGetCoeff(p),r->cf);
2714 pIter(p);
2715 number n2=n_GetNumerator(pGetCoeff(p),r->cf);
2716 pIter(p);
2717 number d;
2718 number t;
2719 loop
2720 {
2721 nlNormalize(pGetCoeff(p),r->cf);
2722 t=n_GetNumerator(pGetCoeff(p),r->cf);
2723 if (nlGreaterZero(t,r->cf))
2724 d=nlAdd(n1,t,r->cf);
2725 else
2726 d=nlSub(n1,t,r->cf);
2727 nlDelete(&t,r->cf);
2728 nlDelete(&n1,r->cf);
2729 n1=d;
2730 pIter(p);
2731 if (p==NULL) break;
2732 nlNormalize(pGetCoeff(p),r->cf);
2733 t=n_GetNumerator(pGetCoeff(p),r->cf);
2734 if (nlGreaterZero(t,r->cf))
2735 d=nlAdd(n2,t,r->cf);
2736 else
2737 d=nlSub(n2,t,r->cf);
2738 nlDelete(&t,r->cf);
2739 nlDelete(&n2,r->cf);
2740 n2=d;
2741 pIter(p);
2742 if (p==NULL) break;
2743 }
2744 d=nlGcd(n1,n2,r->cf);
2745 nlDelete(&n1,r->cf);
2746 nlDelete(&n2,r->cf);
2747 return d;
2748}
2749#else
2750{
2751 /* ph has al least 2 terms */
2752 number d=pGetCoeff(ph);
2753 int s=n_Size(d,r->cf);
2754 pIter(ph);
2755 number d2=pGetCoeff(ph);
2756 int s2=n_Size(d2,r->cf);
2757 pIter(ph);
2758 if (ph==NULL)
2759 {
2760 if (s<s2) return n_Copy(d,r->cf);
2761 else return n_Copy(d2,r->cf);
2762 }
2763 do
2764 {
2765 number nd=pGetCoeff(ph);
2766 int ns=n_Size(nd,r->cf);
2767 if (ns<=2)
2768 {
2769 s2=s;
2770 d2=d;
2771 d=nd;
2772 s=ns;
2773 break;
2774 }
2775 else if (ns<s)
2776 {
2777 s2=s;
2778 d2=d;
2779 d=nd;
2780 s=ns;
2781 }
2782 pIter(ph);
2783 }
2784 while(ph!=NULL);
2785 return n_SubringGcd(d,d2,r->cf);
2786}
2787#endif
2788
2789//void pContent(poly ph)
2790//{
2791// number h,d;
2792// poly p;
2793//
2794// p = ph;
2795// if(pNext(p)==NULL)
2796// {
2797// pSetCoeff(p,nInit(1));
2798// }
2799// else
2800// {
2801//#ifdef PDEBUG
2802// if (!pTest(p)) return;
2803//#endif
2804// nNormalize(pGetCoeff(p));
2805// if(!nGreaterZero(pGetCoeff(ph)))
2806// {
2807// ph = pNeg(ph);
2808// nNormalize(pGetCoeff(p));
2809// }
2810// h=pGetCoeff(p);
2811// pIter(p);
2812// while (p!=NULL)
2813// {
2814// nNormalize(pGetCoeff(p));
2815// if (nGreater(h,pGetCoeff(p))) h=pGetCoeff(p);
2816// pIter(p);
2817// }
2818// h=nCopy(h);
2819// p=ph;
2820// while (p!=NULL)
2821// {
2822// d=n_Gcd(h,pGetCoeff(p));
2823// nDelete(&h);
2824// h = d;
2825// if(nIsOne(h))
2826// {
2827// break;
2828// }
2829// pIter(p);
2830// }
2831// p = ph;
2832// //number tmp;
2833// if(!nIsOne(h))
2834// {
2835// while (p!=NULL)
2836// {
2837// d = nExactDiv(pGetCoeff(p),h);
2838// pSetCoeff(p,d);
2839// pIter(p);
2840// }
2841// }
2842// nDelete(&h);
2843// if ( (nGetChar() == 1) || (nGetChar() < 0) ) /* Q[a],Q(a),Zp[a],Z/p(a) */
2844// {
2845// pTest(ph);
2846// singclap_divide_content(ph);
2847// pTest(ph);
2848// }
2849// }
2850//}
2851#if 0
2852void p_Content(poly ph, const ring r)
2853{
2854 number h,d;
2855 poly p;
2856
2857 if(pNext(ph)==NULL)
2858 {
2859 p_SetCoeff(ph,n_Init(1,r->cf),r);
2860 }
2861 else
2862 {
2863 n_Normalize(pGetCoeff(ph),r->cf);
2864 if(!n_GreaterZero(pGetCoeff(ph),r->cf)) ph = p_Neg(ph,r);
2865 h=n_Copy(pGetCoeff(ph),r->cf);
2866 p = pNext(ph);
2867 while (p!=NULL)
2868 {
2869 n_Normalize(pGetCoeff(p),r->cf);
2870 d=n_Gcd(h,pGetCoeff(p),r->cf);
2871 n_Delete(&h,r->cf);
2872 h = d;
2873 if(n_IsOne(h,r->cf))
2874 {
2875 break;
2876 }
2877 pIter(p);
2878 }
2879 p = ph;
2880 //number tmp;
2881 if(!n_IsOne(h,r->cf))
2882 {
2883 while (p!=NULL)
2884 {
2885 //d = nDiv(pGetCoeff(p),h);
2886 //tmp = nExactDiv(pGetCoeff(p),h);
2887 //if (!nEqual(d,tmp))
2888 //{
2889 // StringSetS("** div0:");nWrite(pGetCoeff(p));StringAppendS("/");
2890 // nWrite(h);StringAppendS("=");nWrite(d);StringAppendS(" int:");
2891 // nWrite(tmp);Print(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
2892 //}
2893 //nDelete(&tmp);
2894 d = n_ExactDiv(pGetCoeff(p),h,r->cf);
2895 p_SetCoeff(p,d,r->cf);
2896 pIter(p);
2897 }
2898 }
2899 n_Delete(&h,r->cf);
2900 //if ( (n_GetChar(r) == 1) || (n_GetChar(r) < 0) ) /* Q[a],Q(a),Zp[a],Z/p(a) */
2901 //{
2902 // singclap_divide_content(ph);
2903 // if(!n_GreaterZero(pGetCoeff(ph),r)) ph = p_Neg(ph,r);
2904 //}
2905 }
2906}
2907#endif
2908/* ---------------------------------------------------------------------------*/
2909/* cleardenom suff */
2910poly p_Cleardenom(poly p, const ring r)
2911{
2912 if( p == NULL )
2913 return NULL;
2914
2915 assume( r != NULL );
2916 assume( r->cf != NULL );
2917 const coeffs C = r->cf;
2918
2919#if CLEARENUMERATORS
2920 if( 0 )
2921 {
2923 n_ClearDenominators(itr, C);
2924 n_ClearContent(itr, C); // divide out the content
2925 p_Test(p, r); n_Test(pGetCoeff(p), C);
2926 assume(n_GreaterZero(pGetCoeff(p), C)); // ??
2927// if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2928 return p;
2929 }
2930#endif
2931
2932 number d, h;
2933
2934 if (rField_is_Ring(r))
2935 {
2936 if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2937 return p;
2938 }
2939
2941 {
2942 if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2943 return p;
2944 }
2945
2946 assume(p != NULL);
2947
2948 if(pNext(p)==NULL)
2949 {
2950 if (!TEST_OPT_CONTENTSB)
2951 p_SetCoeff(p,n_Init(1,C),r);
2952 else if(!n_GreaterZero(pGetCoeff(p),C))
2953 p = p_Neg(p,r);
2954 return p;
2955 }
2956
2957 assume(pNext(p)!=NULL);
2958 poly start=p;
2959
2960#if 0 && CLEARENUMERATORS
2961//CF: does not seem to work that well..
2962
2963 if( nCoeff_is_Q(C) || nCoeff_is_Q_a(C) )
2964 {
2966 n_ClearDenominators(itr, C);
2967 n_ClearContent(itr, C); // divide out the content
2968 p_Test(p, r); n_Test(pGetCoeff(p), C);
2969 assume(n_GreaterZero(pGetCoeff(p), C)); // ??
2970// if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
2971 return start;
2972 }
2973#endif
2974
2975 if(1)
2976 {
2977 // get lcm of all denominators ----------------------------------
2978 h = n_Init(1,C);
2979 while (p!=NULL)
2980 {
2983 n_Delete(&h,C);
2984 h=d;
2985 pIter(p);
2986 }
2987 /* h now contains the 1/lcm of all denominators */
2988 if(!n_IsOne(h,C))
2989 {
2990 // multiply by the lcm of all denominators
2991 p = start;
2992 while (p!=NULL)
2993 {
2994 d=n_Mult(h,pGetCoeff(p),C);
2995 n_Normalize(d,C);
2996 p_SetCoeff(p,d,r);
2997 pIter(p);
2998 }
2999 }
3000 n_Delete(&h,C);
3001 p=start;
3002
3003 p_ContentForGB(p,r);
3004#ifdef HAVE_RATGRING
3005 if (rIsRatGRing(r))
3006 {
3007 /* quick unit detection in the rational case is done in gr_nc_bba */
3008 p_ContentRat(p, r);
3009 start=p;
3010 }
3011#endif
3012 }
3013
3014 if(!n_GreaterZero(pGetCoeff(p),C)) p = p_Neg(p,r);
3015
3016 return start;
3017}
3018
3019void p_Cleardenom_n(poly ph,const ring r,number &c)
3020{
3021 const coeffs C = r->cf;
3022 number d, h;
3023
3024 assume( ph != NULL );
3025
3026 poly p = ph;
3027
3028#if CLEARENUMERATORS
3029 if( 0 )
3030 {
3031 CPolyCoeffsEnumerator itr(ph);
3032
3033 n_ClearDenominators(itr, d, C); // multiply with common denom. d
3034 n_ClearContent(itr, h, C); // divide by the content h
3035
3036 c = n_Div(d, h, C); // d/h
3037
3038 n_Delete(&d, C);
3039 n_Delete(&h, C);
3040
3041 n_Test(c, C);
3042
3043 p_Test(ph, r); n_Test(pGetCoeff(ph), C);
3044 assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
3045/*
3046 if(!n_GreaterZero(pGetCoeff(ph),C))
3047 {
3048 ph = p_Neg(ph,r);
3049 c = n_InpNeg(c, C);
3050 }
3051*/
3052 return;
3053 }
3054#endif
3055
3056
3057 if( pNext(p) == NULL )
3058 {
3060 {
3061 c=n_Invers(pGetCoeff(p), C);
3062 p_SetCoeff(p, n_Init(1, C), r);
3063 }
3064 else
3065 {
3066 c=n_Init(1,C);
3067 }
3068
3069 if(!n_GreaterZero(pGetCoeff(ph),C))
3070 {
3071 ph = p_Neg(ph,r);
3072 c = n_InpNeg(c, C);
3073 }
3074
3075 return;
3076 }
3077 if (TEST_OPT_CONTENTSB) { c=n_Init(1,C); return; }
3078
3079 assume( pNext(p) != NULL );
3080
3081#if CLEARENUMERATORS
3082 if( nCoeff_is_Q(C) || nCoeff_is_Q_a(C) )
3083 {
3084 CPolyCoeffsEnumerator itr(ph);
3085
3086 n_ClearDenominators(itr, d, C); // multiply with common denom. d
3087 n_ClearContent(itr, h, C); // divide by the content h
3088
3089 c = n_Div(d, h, C); // d/h
3090
3091 n_Delete(&d, C);
3092 n_Delete(&h, C);
3093
3094 n_Test(c, C);
3095
3096 p_Test(ph, r); n_Test(pGetCoeff(ph), C);
3097 assume(n_GreaterZero(pGetCoeff(ph), C)); // ??
3098/*
3099 if(!n_GreaterZero(pGetCoeff(ph),C))
3100 {
3101 ph = p_Neg(ph,r);
3102 c = n_InpNeg(c, C);
3103 }
3104*/
3105 return;
3106 }
3107#endif
3108
3109
3110
3111
3112 if(1)
3113 {
3114 h = n_Init(1,C);
3115 while (p!=NULL)
3116 {
3119 n_Delete(&h,C);
3120 h=d;
3121 pIter(p);
3122 }
3123 c=h;
3124 /* contains the 1/lcm of all denominators */
3125 if(!n_IsOne(h,C))
3126 {
3127 p = ph;
3128 while (p!=NULL)
3129 {
3130 /* should be: // NOTE: don't use ->coef!!!!
3131 * number hh;
3132 * nGetDenom(p->coef,&hh);
3133 * nMult(&h,&hh,&d);
3134 * nNormalize(d);
3135 * nDelete(&hh);
3136 * nMult(d,p->coef,&hh);
3137 * nDelete(&d);
3138 * nDelete(&(p->coef));
3139 * p->coef =hh;
3140 */
3141 d=n_Mult(h,pGetCoeff(p),C);
3142 n_Normalize(d,C);
3143 p_SetCoeff(p,d,r);
3144 pIter(p);
3145 }
3146 if (rField_is_Q_a(r))
3147 {
3148 loop
3149 {
3150 h = n_Init(1,C);
3151 p=ph;
3152 while (p!=NULL)
3153 {
3155 n_Delete(&h,C);
3156 h=d;
3157 pIter(p);
3158 }
3159 /* contains the 1/lcm of all denominators */
3160 if(!n_IsOne(h,C))
3161 {
3162 p = ph;
3163 while (p!=NULL)
3164 {
3165 /* should be: // NOTE: don't use ->coef!!!!
3166 * number hh;
3167 * nGetDenom(p->coef,&hh);
3168 * nMult(&h,&hh,&d);
3169 * nNormalize(d);
3170 * nDelete(&hh);
3171 * nMult(d,p->coef,&hh);
3172 * nDelete(&d);
3173 * nDelete(&(p->coef));
3174 * p->coef =hh;
3175 */
3176 d=n_Mult(h,pGetCoeff(p),C);
3177 n_Normalize(d,C);
3178 p_SetCoeff(p,d,r);
3179 pIter(p);
3180 }
3181 number t=n_Mult(c,h,C);
3182 n_Delete(&c,C);
3183 c=t;
3184 }
3185 else
3186 {
3187 break;
3188 }
3189 n_Delete(&h,C);
3190 }
3191 }
3192 }
3193 }
3194
3195 if(!n_GreaterZero(pGetCoeff(ph),C))
3196 {
3197 ph = p_Neg(ph,r);
3198 c = n_InpNeg(c, C);
3199 }
3200
3201}
3202
3203 // normalization: for poly over Q: make poly primitive, integral
3204 // Qa make poly integral with leading
3205 // coefficient minimal in N
3206 // Q(t) make poly primitive, integral
3207
3208void p_ProjectiveUnique(poly ph, const ring r)
3209{
3210 if( ph == NULL )
3211 return;
3212
3213 const coeffs C = r->cf;
3214
3215 number h;
3216 poly p;
3217
3218 if (nCoeff_is_Ring(C))
3219 {
3220 p_ContentForGB(ph,r);
3221 if(!n_GreaterZero(pGetCoeff(ph),C)) ph = p_Neg(ph,r);
3222 assume( n_GreaterZero(pGetCoeff(ph),C) );
3223 return;
3224 }
3225
3227 {
3228 if(!n_GreaterZero(pGetCoeff(ph),C)) ph = p_Neg(ph,r);
3229 return;
3230 }
3231 p = ph;
3232
3233 assume(p != NULL);
3234
3235 if(pNext(p)==NULL) // a monomial
3236 {
3237 p_SetCoeff(p, n_Init(1, C), r);
3238 return;
3239 }
3240
3241 assume(pNext(p)!=NULL);
3242
3243 if(!nCoeff_is_Q(C) && !nCoeff_is_transExt(C))
3244 {
3245 h = p_GetCoeff(p, C);
3246 number hInv = n_Invers(h, C);
3247 pIter(p);
3248 while (p!=NULL)
3249 {
3250 p_SetCoeff(p, n_Mult(p_GetCoeff(p, C), hInv, C), r);
3251 pIter(p);
3252 }
3253 n_Delete(&hInv, C);
3254 p = ph;
3255 p_SetCoeff(p, n_Init(1, C), r);
3256 }
3257
3258 p_Cleardenom(ph, r); //removes also Content
3259
3260
3261 /* normalize ph over a transcendental extension s.t.
3262 lead (ph) is > 0 if extRing->cf == Q
3263 or lead (ph) is monic if extRing->cf == Zp*/
3264 if (nCoeff_is_transExt(C))
3265 {
3266 p= ph;
3267 h= p_GetCoeff (p, C);
3268 fraction f = (fraction) h;
3269 number n=p_GetCoeff (NUM (f),C->extRing->cf);
3270 if (rField_is_Q (C->extRing))
3271 {
3272 if (!n_GreaterZero(n,C->extRing->cf))
3273 {
3274 p=p_Neg (p,r);
3275 }
3276 }
3277 else if (rField_is_Zp(C->extRing))
3278 {
3279 if (!n_IsOne (n, C->extRing->cf))
3280 {
3281 n=n_Invers (n,C->extRing->cf);
3282 nMapFunc nMap;
3283 nMap= n_SetMap (C->extRing->cf, C);
3284 number ninv= nMap (n,C->extRing->cf, C);
3285 p=__p_Mult_nn (p, ninv, r);
3286 n_Delete (&ninv, C);
3287 n_Delete (&n, C->extRing->cf);
3288 }
3289 }
3290 p= ph;
3291 }
3292
3293 return;
3294}
3295
3296#if 0 /*unused*/
3297number p_GetAllDenom(poly ph, const ring r)
3298{
3299 number d=n_Init(1,r->cf);
3300 poly p = ph;
3301
3302 while (p!=NULL)
3303 {
3304 number h=n_GetDenom(pGetCoeff(p),r->cf);
3305 if (!n_IsOne(h,r->cf))
3306 {
3307 number dd=n_Mult(d,h,r->cf);
3308 n_Delete(&d,r->cf);
3309 d=dd;
3310 }
3311 n_Delete(&h,r->cf);
3312 pIter(p);
3313 }
3314 return d;
3315}
3316#endif
3317
3318int p_Size(poly p, const ring r)
3319{
3320 int count = 0;
3321 if (r->cf->has_simple_Alloc)
3322 return pLength(p);
3323 while ( p != NULL )
3324 {
3325 count+= n_Size( pGetCoeff( p ), r->cf );
3326 pIter( p );
3327 }
3328 return count;
3329}
3330
3331/*2
3332*make p homogeneous by multiplying the monomials by powers of x_varnum
3333*assume: deg(var(varnum))==1
3334*/
3335poly p_Homogen (poly p, int varnum, const ring r)
3336{
3337 pFDegProc deg;
3338 if (r->pLexOrder && (r->order[0]==ringorder_lp))
3339 deg=p_Totaldegree;
3340 else
3341 deg=r->pFDeg;
3342
3343 poly q=NULL, qn;
3344 int o,ii;
3345 sBucket_pt bp;
3346
3347 if (p!=NULL)
3348 {
3349 if ((varnum < 1) || (varnum > rVar(r)))
3350 {
3351 return NULL;
3352 }
3353 o=deg(p,r);
3354 q=pNext(p);
3355 while (q != NULL)
3356 {
3357 ii=deg(q,r);
3358 if (ii>o) o=ii;
3359 pIter(q);
3360 }
3361 q = p_Copy(p,r);
3362 bp = sBucketCreate(r);
3363 while (q != NULL)
3364 {
3365 ii = o-deg(q,r);
3366 if (ii!=0)
3367 {
3368 p_AddExp(q,varnum, (long)ii,r);
3369 p_Setm(q,r);
3370 }
3371 qn = pNext(q);
3372 pNext(q) = NULL;
3373 sBucket_Add_m(bp, q);
3374 q = qn;
3375 }
3376 sBucketDestroyAdd(bp, &q, &ii);
3377 }
3378 return q;
3379}
3380
3381/*2
3382*tests if p is homogeneous with respect to the actual weigths
3383*/
3384BOOLEAN p_IsHomogeneous (poly p, const ring r)
3385{
3386 poly qp=p;
3387 int o;
3388
3389 if ((p == NULL) || (pNext(p) == NULL)) return TRUE;
3390 pFDegProc d;
3391 if (r->pLexOrder && (r->order[0]==ringorder_lp))
3392 d=p_Totaldegree;
3393 else
3394 d=r->pFDeg;
3395 o = d(p,r);
3396 do
3397 {
3398 if (d(qp,r) != o) return FALSE;
3399 pIter(qp);
3400 }
3401 while (qp != NULL);
3402 return TRUE;
3403}
3404
3405/*2
3406*tests if p is homogeneous with respect to the given weigths
3407*/
3408BOOLEAN p_IsHomogeneousW (poly p, const intvec *w, const ring r)
3409{
3410 poly qp=p;
3411 long o;
3412
3413 if ((p == NULL) || (pNext(p) == NULL)) return TRUE;
3414 pIter(qp);
3415 o = totaldegreeWecart_IV(p,r,w->ivGetVec());
3416 do
3417 {
3418 if (totaldegreeWecart_IV(qp,r,w->ivGetVec()) != o) return FALSE;
3419 pIter(qp);
3420 }
3421 while (qp != NULL);
3422 return TRUE;
3423}
3424
3425BOOLEAN p_IsHomogeneousW (poly p, const intvec *w, const intvec *module_w, const ring r)
3426{
3427 poly qp=p;
3428 long o;
3429
3430 if ((p == NULL) || (pNext(p) == NULL)) return TRUE;
3431 pIter(qp);
3432 o = totaldegreeWecart_IV(p,r,w->ivGetVec())+(*module_w)[p_GetComp(p,r)];
3433 do
3434 {
3435 long oo=totaldegreeWecart_IV(qp,r,w->ivGetVec())+(*module_w)[p_GetComp(qp,r)];
3436 if (oo != o) return FALSE;
3437 pIter(qp);
3438 }
3439 while (qp != NULL);
3440 return TRUE;
3441}
3442
3443/*----------utilities for syzygies--------------*/
3444BOOLEAN p_VectorHasUnitB(poly p, int * k, const ring r)
3445{
3446 poly q=p,qq;
3447 long unsigned i;
3448
3449 while (q!=NULL)
3450 {
3451 if (p_LmIsConstantComp(q,r))
3452 {
3453 i = __p_GetComp(q,r);
3454 qq = p;
3455 while ((qq != q) && (__p_GetComp(qq,r) != i)) pIter(qq);
3456 if (qq == q)
3457 {
3458 *k = i;
3459 return TRUE;
3460 }
3461 }
3462 pIter(q);
3463 }
3464 return FALSE;
3465}
3466
3467void p_VectorHasUnit(poly p, int * k, int * len, const ring r)
3468{
3469 poly q=p,qq;
3470 int j=0;
3471 long unsigned i;
3472
3473 *len = 0;
3474 while (q!=NULL)
3475 {
3476 if (p_LmIsConstantComp(q,r))
3477 {
3478 i = __p_GetComp(q,r);
3479 qq = p;
3480 while ((qq != q) && (__p_GetComp(qq,r) != i)) pIter(qq);
3481 if (qq == q)
3482 {
3483 j = 0;
3484 while (qq!=NULL)
3485 {
3486 if (__p_GetComp(qq,r)==i) j++;
3487 pIter(qq);
3488 }
3489 if ((*len == 0) || (j<*len))
3490 {
3491 *len = j;
3492 *k = i;
3493 }
3494 }
3495 }
3496 pIter(q);
3497 }
3498}
3499
3500poly p_TakeOutComp1(poly * p, int k, const ring r)
3501{
3502 poly q = *p;
3503
3504 if (q==NULL) return NULL;
3505
3506 poly qq=NULL,result = NULL;
3507 long unsigned kk=k;
3508 if (__p_GetComp(q,r)==kk)
3509 {
3510 result = q; /* *p */
3511 while ((q!=NULL) && (__p_GetComp(q,r)==kk))
3512 {
3513 p_SetComp(q,0,r);
3514 p_SetmComp(q,r);
3515 qq = q;
3516 pIter(q);
3517 }
3518 *p = q;
3519 pNext(qq) = NULL;
3520 }
3521 if (q==NULL) return result;
3522// if (pGetComp(q) > k) pGetComp(q)--;
3523 while (pNext(q)!=NULL)
3524 {
3525 if (__p_GetComp(pNext(q),r)==kk)
3526 {
3527 if (result==NULL)
3528 {
3529 result = pNext(q);
3530 qq = result;
3531 }
3532 else
3533 {
3534 pNext(qq) = pNext(q);
3535 pIter(qq);
3536 }
3537 pNext(q) = pNext(pNext(q));
3538 pNext(qq) =NULL;
3539 p_SetComp(qq,0,r);
3540 p_SetmComp(qq,r);
3541 }
3542 else
3543 {
3544 pIter(q);
3545// if (pGetComp(q) > k) pGetComp(q)--;
3546 }
3547 }
3548 return result;
3549}
3550
3551poly p_TakeOutComp(poly * p, int k, const ring r)
3552{
3553 poly q = *p,qq=NULL,result = NULL;
3554
3555 if (q==NULL) return NULL;
3556 BOOLEAN use_setmcomp=rOrd_SetCompRequiresSetm(r);
3557 if (__p_GetComp(q,r)==k)
3558 {
3559 result = q;
3560 do
3561 {
3562 p_SetComp(q,0,r);
3563 if (use_setmcomp) p_SetmComp(q,r);
3564 qq = q;
3565 pIter(q);
3566 }
3567 while ((q!=NULL) && (__p_GetComp(q,r)==k));
3568 *p = q;
3569 pNext(qq) = NULL;
3570 }
3571 if (q==NULL) return result;
3572 if (__p_GetComp(q,r) > k)
3573 {
3574 p_SubComp(q,1,r);
3575 if (use_setmcomp) p_SetmComp(q,r);
3576 }
3577 poly pNext_q;
3578 while ((pNext_q=pNext(q))!=NULL)
3579 {
3580 if (__p_GetComp(pNext_q,r)==k)
3581 {
3582 if (result==NULL)
3583 {
3584 result = pNext_q;
3585 qq = result;
3586 }
3587 else
3588 {
3589 pNext(qq) = pNext_q;
3590 pIter(qq);
3591 }
3592 pNext(q) = pNext(pNext_q);
3593 pNext(qq) =NULL;
3594 p_SetComp(qq,0,r);
3595 if (use_setmcomp) p_SetmComp(qq,r);
3596 }
3597 else
3598 {
3599 /*pIter(q);*/ q=pNext_q;
3600 if (__p_GetComp(q,r) > k)
3601 {
3602 p_SubComp(q,1,r);
3603 if (use_setmcomp) p_SetmComp(q,r);
3604 }
3605 }
3606 }
3607 return result;
3608}
3609
3610// Splits *p into two polys: *q which consists of all monoms with
3611// component == comp and *p of all other monoms *lq == pLength(*q)
3612void p_TakeOutComp(poly *r_p, long comp, poly *r_q, int *lq, const ring r)
3613{
3614 spolyrec pp, qq;
3615 poly p, q, p_prev;
3616 int l = 0;
3617
3618#ifndef SING_NDEBUG
3619 int lp = pLength(*r_p);
3620#endif
3621
3622 pNext(&pp) = *r_p;
3623 p = *r_p;
3624 p_prev = &pp;
3625 q = &qq;
3626
3627 while(p != NULL)
3628 {
3629 while (__p_GetComp(p,r) == comp)
3630 {
3631 pNext(q) = p;
3632 pIter(q);
3633 p_SetComp(p, 0,r);
3634 p_SetmComp(p,r);
3635 pIter(p);
3636 l++;
3637 if (p == NULL)
3638 {
3639 pNext(p_prev) = NULL;
3640 goto Finish;
3641 }
3642 }
3643 pNext(p_prev) = p;
3644 p_prev = p;
3645 pIter(p);
3646 }
3647
3648 Finish:
3649 pNext(q) = NULL;
3650 *r_p = pNext(&pp);
3651 *r_q = pNext(&qq);
3652 *lq = l;
3653#ifndef SING_NDEBUG
3654 assume(pLength(*r_p) + pLength(*r_q) == (unsigned)lp);
3655#endif
3656 p_Test(*r_p,r);
3657 p_Test(*r_q,r);
3658}
3659
3660void p_DeleteComp(poly * p,int k, const ring r)
3661{
3662 poly q;
3663 long unsigned kk=k;
3664
3665 while ((*p!=NULL) && (__p_GetComp(*p,r)==kk)) p_LmDelete(p,r);
3666 if (*p==NULL) return;
3667 q = *p;
3668 if (__p_GetComp(q,r)>kk)
3669 {
3670 p_SubComp(q,1,r);
3671 p_SetmComp(q,r);
3672 }
3673 while (pNext(q)!=NULL)
3674 {
3675 if (__p_GetComp(pNext(q),r)==kk)
3676 p_LmDelete(&(pNext(q)),r);
3677 else
3678 {
3679 pIter(q);
3680 if (__p_GetComp(q,r)>kk)
3681 {
3682 p_SubComp(q,1,r);
3683 p_SetmComp(q,r);
3684 }
3685 }
3686 }
3687}
3688
3689poly p_Vec2Poly(poly v, int k, const ring r)
3690{
3691 poly h;
3692 poly res=NULL;
3693 long unsigned kk=k;
3694
3695 while (v!=NULL)
3696 {
3697 if (__p_GetComp(v,r)==kk)
3698 {
3699 h=p_Head(v,r);
3700 p_SetComp(h,0,r);
3701 pNext(h)=res;res=h;
3702 }
3703 pIter(v);
3704 }
3705 if (res!=NULL) res=pReverse(res);
3706 return res;
3707}
3708
3709/// vector to already allocated array (len>=p_MaxComp(v,r))
3710// also used for p_Vec2Polys
3711void p_Vec2Array(poly v, poly *p, int len, const ring r)
3712{
3713 poly h;
3714 int k;
3715
3716 for(int i=len-1;i>=0;i--) p[i]=NULL;
3717 while (v!=NULL)
3718 {
3719 h=p_Head(v,r);
3720 k=__p_GetComp(h,r);
3721 if (k>len) { Werror("wrong rank:%d, should be %d",len,k); }
3722 else
3723 {
3724 p_SetComp(h,0,r);
3725 p_Setm(h,r);
3726 pNext(h)=p[k-1];p[k-1]=h;
3727 }
3728 pIter(v);
3729 }
3730 for(int i=len-1;i>=0;i--)
3731 {
3732 if (p[i]!=NULL) p[i]=pReverse(p[i]);
3733 }
3734}
3735
3736/*2
3737* convert a vector to a set of polys,
3738* allocates the polyset, (entries 0..(*len)-1)
3739* the vector will not be changed
3740*/
3741void p_Vec2Polys(poly v, poly* *p, int *len, const ring r)
3742{
3743 *len=p_MaxComp(v,r);
3744 if (*len==0) *len=1;
3745 *p=(poly*)omAlloc((*len)*sizeof(poly));
3746 p_Vec2Array(v,*p,*len,r);
3747}
3748
3749//
3750// resets the pFDeg and pLDeg: if pLDeg is not given, it is
3751// set to currRing->pLDegOrig, i.e. to the respective LDegProc which
3752// only uses pFDeg (and not p_Deg, or pTotalDegree, etc)
3753void pSetDegProcs(ring r, pFDegProc new_FDeg, pLDegProc new_lDeg)
3754{
3755 assume(new_FDeg != NULL);
3756 r->pFDeg = new_FDeg;
3757
3758 if (new_lDeg == NULL)
3759 new_lDeg = r->pLDegOrig;
3760
3761 r->pLDeg = new_lDeg;
3762}
3763
3764// restores pFDeg and pLDeg:
3765void pRestoreDegProcs(ring r, pFDegProc old_FDeg, pLDegProc old_lDeg)
3766{
3767 assume(old_FDeg != NULL && old_lDeg != NULL);
3768 r->pFDeg = old_FDeg;
3769 r->pLDeg = old_lDeg;
3770}
3771
3772/*-------- several access procedures to monomials -------------------- */
3773/*
3774* the module weights for std
3775*/
3779
3780static long pModDeg(poly p, ring r)
3781{
3782 long d=pOldFDeg(p, r);
3783 int c=__p_GetComp(p, r);
3784 if ((c>0) && ((r->pModW)->range(c-1))) d+= (*(r->pModW))[c-1];
3785 return d;
3786 //return pOldFDeg(p, r)+(*pModW)[p_GetComp(p, r)-1];
3787}
3788
3789void p_SetModDeg(intvec *w, ring r)
3790{
3791 if (w!=NULL)
3792 {
3793 r->pModW = w;
3794 pOldFDeg = r->pFDeg;
3795 pOldLDeg = r->pLDeg;
3796 pOldLexOrder = r->pLexOrder;
3798 r->pLexOrder = TRUE;
3799 }
3800 else
3801 {
3802 r->pModW = NULL;
3804 r->pLexOrder = pOldLexOrder;
3805 }
3806}
3807
3808/*2
3809* handle memory request for sets of polynomials (ideals)
3810* l is the length of *p, increment is the difference (may be negative)
3811*/
3812void pEnlargeSet(poly* *p, int l, int increment)
3813{
3814 poly* h;
3815
3816 if (*p==NULL)
3817 {
3818 if (increment==0) return;
3819 h=(poly*)omAlloc0(increment*sizeof(poly));
3820 }
3821 else
3822 {
3823 h=(poly*)omReallocSize((poly*)*p,l*sizeof(poly),(l+increment)*sizeof(poly));
3824 if (increment>0)
3825 {
3826 memset(&(h[l]),0,increment*sizeof(poly));
3827 }
3828 }
3829 *p=h;
3830}
3831
3832/*2
3833*divides p1 by its leading coefficient
3834*/
3835void p_Norm(poly p1, const ring r)
3836{
3837 if (LIKELY(rField_is_Ring(r)))
3838 {
3839 if(!n_GreaterZero(pGetCoeff(p1),r->cf)) p1 = p_Neg(p1,r);
3840 if (!n_IsUnit(pGetCoeff(p1), r->cf)) return;
3841 // Werror("p_Norm not possible in the case of coefficient rings.");
3842 }
3843 else if (LIKELY(p1!=NULL))
3844 {
3845 if (UNLIKELY(pNext(p1)==NULL))
3846 {
3847 p_SetCoeff(p1,n_Init(1,r->cf),r);
3848 return;
3849 }
3850 if (!n_IsOne(pGetCoeff(p1),r->cf))
3851 {
3852 number k = pGetCoeff(p1);
3853 pSetCoeff0(p1,n_Init(1,r->cf));
3854 poly h = pNext(p1);
3855 if (LIKELY(rField_is_Zp(r)))
3856 {
3857 if (r->cf->ch>32003)
3858 {
3859 number inv=n_Invers(k,r->cf);
3860 while (h!=NULL)
3861 {
3862 number c=n_Mult(pGetCoeff(h),inv,r->cf);
3863 // no need to normalize
3864 p_SetCoeff(h,c,r);
3865 pIter(h);
3866 }
3867 // no need for n_Delete for Zp: n_Delete(&inv,r->cf);
3868 }
3869 else
3870 {
3871 while (h!=NULL)
3872 {
3873 number c=n_Div(pGetCoeff(h),k,r->cf);
3874 // no need to normalize
3875 p_SetCoeff(h,c,r);
3876 pIter(h);
3877 }
3878 }
3879 }
3880 else if(getCoeffType(r->cf)==n_algExt)
3881 {
3882 n_Normalize(k,r->cf);
3883 number inv=n_Invers(k,r->cf);
3884 while (h!=NULL)
3885 {
3886 number c=n_Mult(pGetCoeff(h),inv,r->cf);
3887 // no need to normalize
3888 // normalize already in nMult: Zp_a, Q_a
3889 p_SetCoeff(h,c,r);
3890 pIter(h);
3891 }
3892 n_Delete(&inv,r->cf);
3893 n_Delete(&k,r->cf);
3894 }
3895 else
3896 {
3897 n_Normalize(k,r->cf);
3898 while (h!=NULL)
3899 {
3900 number c=n_Div(pGetCoeff(h),k,r->cf);
3901 // no need to normalize: Z/p, R
3902 // remains: Q
3903 if (rField_is_Q(r)) n_Normalize(c,r->cf);
3904 p_SetCoeff(h,c,r);
3905 pIter(h);
3906 }
3907 n_Delete(&k,r->cf);
3908 }
3909 }
3910 else
3911 {
3912 //if (r->cf->cfNormalize != nDummy2) //TODO: OPTIMIZE
3913 if (rField_is_Q(r))
3914 {
3915 poly h = pNext(p1);
3916 while (h!=NULL)
3917 {
3918 n_Normalize(pGetCoeff(h),r->cf);
3919 pIter(h);
3920 }
3921 }
3922 }
3923 }
3924}
3925
3926/*2
3927*normalize all coefficients
3928*/
3929void p_Normalize(poly p,const ring r)
3930{
3931 const coeffs cf=r->cf;
3932 /* Z/p, GF(p,n), R, long R/C, Nemo rings */
3933 if (cf->cfNormalize==ndNormalize)
3934 return;
3935 while (p!=NULL)
3936 {
3937 // no test befor n_Normalize: n_Normalize should fix problems
3939 pIter(p);
3940 }
3941}
3942
3943// splits p into polys with Exp(n) == 0 and Exp(n) != 0
3944// Poly with Exp(n) != 0 is reversed
3945static void p_SplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero, const ring r)
3946{
3947 if (p == NULL)
3948 {
3949 *non_zero = NULL;
3950 *zero = NULL;
3951 return;
3952 }
3953 spolyrec sz;
3954 poly z, n_z, next;
3955 z = &sz;
3956 n_z = NULL;
3957
3958 while(p != NULL)
3959 {
3960 next = pNext(p);
3961 if (p_GetExp(p, n,r) == 0)
3962 {
3963 pNext(z) = p;
3964 pIter(z);
3965 }
3966 else
3967 {
3968 pNext(p) = n_z;
3969 n_z = p;
3970 }
3971 p = next;
3972 }
3973 pNext(z) = NULL;
3974 *zero = pNext(&sz);
3975 *non_zero = n_z;
3976}
3977/*3
3978* substitute the n-th variable by 1 in p
3979* destroy p
3980*/
3981static poly p_Subst1 (poly p,int n, const ring r)
3982{
3983 poly qq=NULL, result = NULL;
3984 poly zero=NULL, non_zero=NULL;
3985
3986 // reverse, so that add is likely to be linear
3987 p_SplitAndReversePoly(p, n, &non_zero, &zero,r);
3988
3989 while (non_zero != NULL)
3990 {
3991 assume(p_GetExp(non_zero, n,r) != 0);
3992 qq = non_zero;
3993 pIter(non_zero);
3994 qq->next = NULL;
3995 p_SetExp(qq,n,0,r);
3996 p_Setm(qq,r);
3997 result = p_Add_q(result,qq,r);
3998 }
3999 p = p_Add_q(result, zero,r);
4000 p_Test(p,r);
4001 return p;
4002}
4003
4004/*3
4005* substitute the n-th variable by number e in p
4006* destroy p
4007*/
4008static poly p_Subst2 (poly p,int n, number e, const ring r)
4009{
4010 assume( ! n_IsZero(e,r->cf) );
4011 poly qq,result = NULL;
4012 number nn, nm;
4013 poly zero, non_zero;
4014
4015 // reverse, so that add is likely to be linear
4016 p_SplitAndReversePoly(p, n, &non_zero, &zero,r);
4017
4018 while (non_zero != NULL)
4019 {
4020 assume(p_GetExp(non_zero, n, r) != 0);
4021 qq = non_zero;
4022 pIter(non_zero);
4023 qq->next = NULL;
4024 n_Power(e, p_GetExp(qq, n, r), &nn,r->cf);
4025 nm = n_Mult(nn, pGetCoeff(qq),r->cf);
4026#ifdef HAVE_RINGS
4027 if (n_IsZero(nm,r->cf))
4028 {
4029 p_LmFree(&qq,r);
4030 n_Delete(&nm,r->cf);
4031 }
4032 else
4033#endif
4034 {
4035 p_SetCoeff(qq, nm,r);
4036 p_SetExp(qq, n, 0,r);
4037 p_Setm(qq,r);
4038 result = p_Add_q(result,qq,r);
4039 }
4040 n_Delete(&nn,r->cf);
4041 }
4042 p = p_Add_q(result, zero,r);
4043 p_Test(p,r);
4044 return p;
4045}
4046
4047
4048/* delete monoms whose n-th exponent is different from zero */
4049static poly p_Subst0(poly p, int n, const ring r)
4050{
4051 spolyrec res;
4052 poly h = &res;
4053 pNext(h) = p;
4054
4055 while (pNext(h)!=NULL)
4056 {
4057 if (p_GetExp(pNext(h),n,r)!=0)
4058 {
4059 p_LmDelete(&pNext(h),r);
4060 }
4061 else
4062 {
4063 pIter(h);
4064 }
4065 }
4066 p_Test(pNext(&res),r);
4067 return pNext(&res);
4068}
4069
4070/*2
4071* substitute the n-th variable by e in p
4072* destroy p
4073*/
4074poly p_Subst(poly p, int n, poly e, const ring r)
4075{
4076#ifdef HAVE_SHIFTBBA
4077 // also don't even use p_Subst0 for Letterplace
4078 if (rIsLPRing(r))
4079 {
4080 poly subst = p_LPSubst(p, n, e, r);
4081 p_Delete(&p, r);
4082 return subst;
4083 }
4084#endif
4085
4086 if (e == NULL) return p_Subst0(p, n,r);
4087
4088 if (p_IsConstant(e,r))
4089 {
4090 if (n_IsOne(pGetCoeff(e),r->cf)) return p_Subst1(p,n,r);
4091 else return p_Subst2(p, n, pGetCoeff(e),r);
4092 }
4093
4094#ifdef HAVE_PLURAL
4095 if (rIsPluralRing(r))
4096 {
4097 return nc_pSubst(p,n,e,r);
4098 }
4099#endif
4100
4101 int exponent,i;
4102 poly h, res, m;
4103 int *me,*ee;
4104 number nu,nu1;
4105
4106 me=(int *)omAlloc((rVar(r)+1)*sizeof(int));
4107 ee=(int *)omAlloc((rVar(r)+1)*sizeof(int));
4108 if (e!=NULL) p_GetExpV(e,ee,r);
4109 res=NULL;
4110 h=p;
4111 while (h!=NULL)
4112 {
4113 if ((e!=NULL) || (p_GetExp(h,n,r)==0))
4114 {
4115 m=p_Head(h,r);
4116 p_GetExpV(m,me,r);
4117 exponent=me[n];
4118 me[n]=0;
4119 for(i=rVar(r);i>0;i--)
4120 me[i]+=exponent*ee[i];
4121 p_SetExpV(m,me,r);
4122 if (e!=NULL)
4123 {
4124 n_Power(pGetCoeff(e),exponent,&nu,r->cf);
4125 nu1=n_Mult(pGetCoeff(m),nu,r->cf);
4126 n_Delete(&nu,r->cf);
4127 p_SetCoeff(m,nu1,r);
4128 }
4129 res=p_Add_q(res,m,r);
4130 }
4131 p_LmDelete(&h,r);
4132 }
4133 omFreeSize((ADDRESS)me,(rVar(r)+1)*sizeof(int));
4134 omFreeSize((ADDRESS)ee,(rVar(r)+1)*sizeof(int));
4135 return res;
4136}
4137
4138/*2
4139 * returns a re-ordered convertion of a number as a polynomial,
4140 * with permutation of parameters
4141 * NOTE: this only works for Frank's alg. & trans. fields
4142 */
4143poly n_PermNumber(const number z, const int *par_perm, const int , const ring src, const ring dst)
4144{
4145#if 0
4146 PrintS("\nSource Ring: \n");
4147 rWrite(src);
4148
4149 if(0)
4150 {
4151 number zz = n_Copy(z, src->cf);
4152 PrintS("z: "); n_Write(zz, src);
4153 n_Delete(&zz, src->cf);
4154 }
4155
4156 PrintS("\nDestination Ring: \n");
4157 rWrite(dst);
4158
4159 /*Print("\nOldPar: %d\n", OldPar);
4160 for( int i = 1; i <= OldPar; i++ )
4161 {
4162 Print("par(%d) -> par/var (%d)\n", i, par_perm[i-1]);
4163 }*/
4164#endif
4165 if( z == NULL )
4166 return NULL;
4167
4168 const coeffs srcCf = src->cf;
4169 assume( srcCf != NULL );
4170
4171 assume( !nCoeff_is_GF(srcCf) );
4172 assume( src->cf->extRing!=NULL );
4173
4174 poly zz = NULL;
4175
4176 const ring srcExtRing = srcCf->extRing;
4177 assume( srcExtRing != NULL );
4178
4179 const coeffs dstCf = dst->cf;
4180 assume( dstCf != NULL );
4181
4182 if( nCoeff_is_algExt(srcCf) ) // nCoeff_is_GF(srcCf)?
4183 {
4184 zz = (poly) z;
4185 if( zz == NULL ) return NULL;
4186 }
4187 else if (nCoeff_is_transExt(srcCf))
4188 {
4189 assume( !IS0(z) );
4190
4191 zz = NUM((fraction)z);
4192 p_Test (zz, srcExtRing);
4193
4194 if( zz == NULL ) return NULL;
4195 if( !DENIS1((fraction)z) )
4196 {
4197 if (!p_IsConstant(DEN((fraction)z),srcExtRing))
4198 WarnS("Not defined: Cannot map a rational fraction and make a polynomial out of it! Ignoring the denominator.");
4199 }
4200 }
4201 else
4202 {
4203 assume (FALSE);
4204 WerrorS("Number permutation is not implemented for this data yet!");
4205 return NULL;
4206 }
4207
4208 assume( zz != NULL );
4209 p_Test (zz, srcExtRing);
4210
4211 nMapFunc nMap = n_SetMap(srcExtRing->cf, dstCf);
4212
4213 assume( nMap != NULL );
4214
4215 poly qq;
4216 if ((par_perm == NULL) && (rPar(dst) != 0 && rVar (srcExtRing) > 0))
4217 {
4218 int* perm;
4219 perm=(int *)omAlloc0((rVar(srcExtRing)+1)*sizeof(int));
4220 for(int i=si_min(rVar(srcExtRing),rPar(dst));i>0;i--)
4221 perm[i]=-i;
4222 qq = p_PermPoly(zz, perm, srcExtRing, dst, nMap, NULL, rVar(srcExtRing)-1);
4223 omFreeSize ((ADDRESS)perm, (rVar(srcExtRing)+1)*sizeof(int));
4224 }
4225 else
4226 qq = p_PermPoly(zz, par_perm-1, srcExtRing, dst, nMap, NULL, rVar (srcExtRing)-1);
4227
4228 if(nCoeff_is_transExt(srcCf)
4229 && (!DENIS1((fraction)z))
4230 && p_IsConstant(DEN((fraction)z),srcExtRing))
4231 {
4232 number n=nMap(pGetCoeff(DEN((fraction)z)),srcExtRing->cf, dstCf);
4233 qq=p_Div_nn(qq,n,dst);
4234 n_Delete(&n,dstCf);
4235 p_Normalize(qq,dst);
4236 }
4237 p_Test (qq, dst);
4238
4239 return qq;
4240}
4241
4242
4243/*2
4244*returns a re-ordered copy of a polynomial, with permutation of the variables
4245*/
4246poly p_PermPoly (poly p, const int * perm, const ring oldRing, const ring dst,
4247 nMapFunc nMap, const int *par_perm, int OldPar, BOOLEAN use_mult)
4248{
4249#if 0
4250 p_Test(p, oldRing);
4251 PrintS("p_PermPoly::p: "); p_Write(p, oldRing, oldRing);
4252#endif
4253 const int OldpVariables = rVar(oldRing);
4254 poly result = NULL;
4255 poly result_last = NULL;
4256 poly aq = NULL; /* the map coefficient */
4257 poly qq; /* the mapped monomial */
4258 assume(dst != NULL);
4259 assume(dst->cf != NULL);
4260 #ifdef HAVE_PLURAL
4261 poly tmp_mm=p_One(dst);
4262 #endif
4263 while (p != NULL)
4264 {
4265 // map the coefficient
4266 if ( ((OldPar == 0) || (par_perm == NULL) || rField_is_GF(oldRing) || (nMap==ndCopyMap))
4267 && (nMap != NULL) )
4268 {
4269 qq = p_Init(dst);
4270 assume( nMap != NULL );
4271 number n = nMap(p_GetCoeff(p, oldRing), oldRing->cf, dst->cf);
4272 n_Test (n,dst->cf);
4273 if ( nCoeff_is_algExt(dst->cf) )
4274 n_Normalize(n, dst->cf);
4275 p_GetCoeff(qq, dst) = n;// Note: n can be a ZERO!!!
4276 }
4277 else
4278 {
4279 qq = p_One(dst);
4280// aq = naPermNumber(p_GetCoeff(p, oldRing), par_perm, OldPar, oldRing); // no dst???
4281// poly n_PermNumber(const number z, const int *par_perm, const int P, const ring src, const ring dst)
4282 aq = n_PermNumber(p_GetCoeff(p, oldRing), par_perm, OldPar, oldRing, dst);
4283 p_Test(aq, dst);
4284 if ( nCoeff_is_algExt(dst->cf) )
4285 p_Normalize(aq,dst);
4286 if (aq == NULL)
4287 p_SetCoeff(qq, n_Init(0, dst->cf),dst); // Very dirty trick!!!
4288 p_Test(aq, dst);
4289 }
4290 if (rRing_has_Comp(dst))
4291 p_SetComp(qq, p_GetComp(p, oldRing), dst);
4292 if ( n_IsZero(pGetCoeff(qq), dst->cf) )
4293 {
4294 p_LmDelete(&qq,dst);
4295 qq = NULL;
4296 }
4297 else
4298 {
4299 // map pars:
4300 int mapped_to_par = 0;
4301 for(int i = 1; i <= OldpVariables; i++)
4302 {
4303 int e = p_GetExp(p, i, oldRing);
4304 if (e != 0)
4305 {
4306 if (perm==NULL)
4307 p_SetExp(qq, i, e, dst);
4308 else if (perm[i]>0)
4309 {
4310 #ifdef HAVE_PLURAL
4311 if(use_mult)
4312 {
4313 p_SetExp(tmp_mm,perm[i],e,dst);
4314 p_Setm(tmp_mm,dst);
4315 qq=p_Mult_mm(qq,tmp_mm,dst);
4316 p_SetExp(tmp_mm,perm[i],0,dst);
4317
4318 }
4319 else
4320 #endif
4321 p_AddExp(qq,perm[i], e/*p_GetExp( p,i,oldRing)*/, dst);
4322 }
4323 else if (perm[i]<0)
4324 {
4325 number c = p_GetCoeff(qq, dst);
4326 if (rField_is_GF(dst))
4327 {
4328 assume( dst->cf->extRing == NULL );
4329 number ee = n_Param(1, dst);
4330 number eee;
4331 n_Power(ee, e, &eee, dst->cf); //nfDelete(ee,dst);
4332 ee = n_Mult(c, eee, dst->cf);
4333 //nfDelete(c,dst);nfDelete(eee,dst);
4334 pSetCoeff0(qq,ee);
4335 }
4336 else if (nCoeff_is_Extension(dst->cf))
4337 {
4338 const int par = -perm[i];
4339 assume( par > 0 );
4340// WarnS("longalg missing 3");
4341#if 1
4342 const coeffs C = dst->cf;
4343 assume( C != NULL );
4344 const ring R = C->extRing;
4345 assume( R != NULL );
4346 assume( par <= rVar(R) );
4347 poly pcn; // = (number)c
4348 assume( !n_IsZero(c, C) );
4349 if( nCoeff_is_algExt(C) )
4350 pcn = (poly) c;
4351 else // nCoeff_is_transExt(C)
4352 pcn = NUM((fraction)c);
4353 if (pNext(pcn) == NULL) // c->z
4354 p_AddExp(pcn, -perm[i], e, R);
4355 else /* more difficult: we have really to multiply: */
4356 {
4357 poly mmc = p_ISet(1, R);
4358 p_SetExp(mmc, -perm[i], e, R);
4359 p_Setm(mmc, R);
4360 number nnc;
4361 // convert back to a number: number nnc = mmc;
4362 if( nCoeff_is_algExt(C) )
4363 nnc = (number) mmc;
4364 else // nCoeff_is_transExt(C)
4365 nnc = ntInit(mmc, C);
4366 p_GetCoeff(qq, dst) = n_Mult((number)c, nnc, C);
4367 n_Delete((number *)&c, C);
4368 n_Delete((number *)&nnc, C);
4369 }
4370 mapped_to_par=1;
4371#endif
4372 }
4373 }
4374 else
4375 {
4376 /* this variable maps to 0 !*/
4377 p_LmDelete(&qq, dst);
4378 break;
4379 }
4380 }
4381 }
4382 if ( mapped_to_par && (qq!= NULL) && nCoeff_is_algExt(dst->cf) )
4383 {
4384 number n = p_GetCoeff(qq, dst);
4385 n_Normalize(n, dst->cf);
4386 p_GetCoeff(qq, dst) = n;
4387 }
4388 }
4389 pIter(p);
4390
4391#if 0
4392 p_Test(aq,dst);
4393 PrintS("aq: "); p_Write(aq, dst, dst);
4394#endif
4395
4396
4397#if 1
4398 if (qq!=NULL)
4399 {
4400 p_Setm(qq,dst);
4401
4402 p_Test(aq,dst);
4403 p_Test(qq,dst);
4404
4405#if 0
4406 PrintS("qq: "); p_Write(qq, dst, dst);
4407#endif
4408
4409 if (aq!=NULL)
4410 qq=p_Mult_q(aq,qq,dst);
4411 aq = qq;
4412 while (pNext(aq) != NULL) pIter(aq);
4413 if (result_last==NULL)
4414 {
4415 result=qq;
4416 }
4417 else
4418 {
4419 pNext(result_last)=qq;
4420 }
4421 result_last=aq;
4422 aq = NULL;
4423 }
4424 else if (aq!=NULL)
4425 {
4426 p_Delete(&aq,dst);
4427 }
4428 }
4429 result=p_SortAdd(result,dst);
4430#else
4431 // if (qq!=NULL)
4432 // {
4433 // pSetm(qq);
4434 // pTest(qq);
4435 // pTest(aq);
4436 // if (aq!=NULL) qq=pMult(aq,qq);
4437 // aq = qq;
4438 // while (pNext(aq) != NULL) pIter(aq);
4439 // pNext(aq) = result;
4440 // aq = NULL;
4441 // result = qq;
4442 // }
4443 // else if (aq!=NULL)
4444 // {
4445 // pDelete(&aq);
4446 // }
4447 //}
4448 //p = result;
4449 //result = NULL;
4450 //while (p != NULL)
4451 //{
4452 // qq = p;
4453 // pIter(p);
4454 // qq->next = NULL;
4455 // result = pAdd(result, qq);
4456 //}
4457#endif
4458 p_Test(result,dst);
4459#if 0
4460 p_Test(result,dst);
4461 PrintS("result: "); p_Write(result,dst,dst);
4462#endif
4463 #ifdef HAVE_PLURAL
4464 p_LmDelete(&tmp_mm,dst);
4465 #endif
4466 return result;
4467}
4468/**************************************************************
4469 *
4470 * Jet
4471 *
4472 **************************************************************/
4473
4474poly pp_Jet(poly p, int m, const ring R)
4475{
4476 poly r=NULL;
4477 poly t=NULL;
4478
4479 while (p!=NULL)
4480 {
4481 if (p_Totaldegree(p,R)<=m)
4482 {
4483 if (r==NULL)
4484 r=p_Head(p,R);
4485 else
4486 if (t==NULL)
4487 {
4488 pNext(r)=p_Head(p,R);
4489 t=pNext(r);
4490 }
4491 else
4492 {
4493 pNext(t)=p_Head(p,R);
4494 pIter(t);
4495 }
4496 }
4497 pIter(p);
4498 }
4499 return r;
4500}
4501
4502poly p_Jet(poly p, int m,const ring R)
4503{
4504 while((p!=NULL) && (p_Totaldegree(p,R)>m)) p_LmDelete(&p,R);
4505 if (p==NULL) return NULL;
4506 poly r=p;
4507 while (pNext(p)!=NULL)
4508 {
4509 if (p_Totaldegree(pNext(p),R)>m)
4510 {
4511 p_LmDelete(&pNext(p),R);
4512 }
4513 else
4514 pIter(p);
4515 }
4516 return r;
4517}
4518
4519poly pp_JetW(poly p, int m, int *w, const ring R)
4520{
4521 poly r=NULL;
4522 poly t=NULL;
4523 while (p!=NULL)
4524 {
4525 if (totaldegreeWecart_IV(p,R,w)<=m)
4526 {
4527 if (r==NULL)
4528 r=p_Head(p,R);
4529 else
4530 if (t==NULL)
4531 {
4532 pNext(r)=p_Head(p,R);
4533 t=pNext(r);
4534 }
4535 else
4536 {
4537 pNext(t)=p_Head(p,R);
4538 pIter(t);
4539 }
4540 }
4541 pIter(p);
4542 }
4543 return r;
4544}
4545
4546poly p_JetW(poly p, int m, int *w, const ring R)
4547{
4548 while((p!=NULL) && (totaldegreeWecart_IV(p,R,w)>m)) p_LmDelete(&p,R);
4549 if (p==NULL) return NULL;
4550 poly r=p;
4551 while (pNext(p)!=NULL)
4552 {
4554 {
4555 p_LmDelete(&pNext(p),R);
4556 }
4557 else
4558 pIter(p);
4559 }
4560 return r;
4561}
4562
4563/*************************************************************/
4564int p_MinDeg(poly p,intvec *w, const ring R)
4565{
4566 if(p==NULL)
4567 return -1;
4568 int d=-1;
4569 while(p!=NULL)
4570 {
4571 int d0=0;
4572 for(int j=0;j<rVar(R);j++)
4573 if(w==NULL||j>=w->length())
4574 d0+=p_GetExp(p,j+1,R);
4575 else
4576 d0+=(*w)[j]*p_GetExp(p,j+1,R);
4577 if(d0<d||d==-1)
4578 d=d0;
4579 pIter(p);
4580 }
4581 return d;
4582}
4583
4584/***************************************************************/
4585static poly p_Invers(int n,poly u,intvec *w, const ring R)
4586{
4587 if(n<0)
4588 return NULL;
4589 number u0=n_Invers(pGetCoeff(u),R->cf);
4590 poly v=p_NSet(u0,R);
4591 if(n==0)
4592 return v;
4593 int *ww=iv2array(w,R);
4594 poly u1=p_JetW(p_Sub(p_One(R),__p_Mult_nn(u,u0,R),R),n,ww,R);
4595 if(u1==NULL)
4596 {
4597 omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(int));
4598 return v;
4599 }
4600 poly v1=__p_Mult_nn(p_Copy(u1,R),u0,R);
4601 v=p_Add_q(v,p_Copy(v1,R),R);
4602 for(int i=n/p_MinDeg(u1,w,R);i>1;i--)
4603 {
4604 v1=p_JetW(p_Mult_q(v1,p_Copy(u1,R),R),n,ww,R);
4605 v=p_Add_q(v,p_Copy(v1,R),R);
4606 }
4607 p_Delete(&u1,R);
4608 p_Delete(&v1,R);
4609 omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(int));
4610 return v;
4611}
4612
4613
4614poly p_Series(int n,poly p,poly u, intvec *w, const ring R)
4615{
4616 int *ww=iv2array(w,R);
4617 if(p!=NULL)
4618 {
4619 if(u==NULL)
4620 p=p_JetW(p,n,ww,R);
4621 else
4622 p=p_JetW(p_Mult_q(p,p_Invers(n-p_MinDeg(p,w,R),u,w,R),R),n,ww,R);
4623 }
4624 omFreeSize((ADDRESS)ww,(rVar(R)+1)*sizeof(int));
4625 return p;
4626}
4627
4628BOOLEAN p_EqualPolys(poly p1,poly p2, const ring r)
4629{
4630 while ((p1 != NULL) && (p2 != NULL))
4631 {
4632 if (! p_LmEqual(p1, p2,r))
4633 return FALSE;
4634 if (! n_Equal(p_GetCoeff(p1,r), p_GetCoeff(p2,r),r->cf ))
4635 return FALSE;
4636 pIter(p1);
4637 pIter(p2);
4638 }
4639 return (p1==p2);
4640}
4641
4642static inline BOOLEAN p_ExpVectorEqual(poly p1, poly p2, const ring r1, const ring r2)
4643{
4644 assume( r1 == r2 || rSamePolyRep(r1, r2) );
4645
4646 p_LmCheckPolyRing1(p1, r1);
4647 p_LmCheckPolyRing1(p2, r2);
4648
4649 int i = r1->ExpL_Size;
4650
4651 assume( r1->ExpL_Size == r2->ExpL_Size );
4652
4653 unsigned long *ep = p1->exp;
4654 unsigned long *eq = p2->exp;
4655
4656 do
4657 {
4658 i--;
4659 if (ep[i] != eq[i]) return FALSE;
4660 }
4661 while (i);
4662
4663 return TRUE;
4664}
4665
4666BOOLEAN p_EqualPolys(poly p1,poly p2, const ring r1, const ring r2)
4667{
4668 assume( r1 == r2 || rSamePolyRep(r1, r2) ); // will be used in rEqual!
4669 assume( r1->cf == r2->cf );
4670
4671 while ((p1 != NULL) && (p2 != NULL))
4672 {
4673 // returns 1 if ExpVector(p)==ExpVector(q): does not compare numbers !!
4674 // #define p_LmEqual(p1, p2, r) p_ExpVectorEqual(p1, p2, r)
4675
4676 if (! p_ExpVectorEqual(p1, p2, r1, r2))
4677 return FALSE;
4678
4679 if (! n_Equal(p_GetCoeff(p1,r1), p_GetCoeff(p2,r2), r1->cf ))
4680 return FALSE;
4681
4682 pIter(p1);
4683 pIter(p2);
4684 }
4685 return (p1==p2);
4686}
4687
4688/*2
4689*returns TRUE if p1 is a skalar multiple of p2
4690*assume p1 != NULL and p2 != NULL
4691*/
4692BOOLEAN p_ComparePolys(poly p1,poly p2, const ring r)
4693{
4694 number n,nn;
4695 pAssume(p1 != NULL && p2 != NULL);
4696
4697 if (!p_LmEqual(p1,p2,r)) //compare leading mons
4698 return FALSE;
4699 if ((pNext(p1)==NULL) && (pNext(p2)!=NULL))
4700 return FALSE;
4701 if ((pNext(p2)==NULL) && (pNext(p1)!=NULL))
4702 return FALSE;
4703 if (pLength(p1) != pLength(p2))
4704 return FALSE;
4705 #ifdef HAVE_RINGS
4706 if (rField_is_Ring(r))
4707 {
4708 if (!n_DivBy(pGetCoeff(p1), pGetCoeff(p2), r->cf)) return FALSE;
4709 }
4710 #endif
4711 n=n_Div(pGetCoeff(p1),pGetCoeff(p2),r->cf);
4712 while ((p1 != NULL) /*&& (p2 != NULL)*/)
4713 {
4714 if ( ! p_LmEqual(p1, p2,r))
4715 {
4716 n_Delete(&n, r->cf);
4717 return FALSE;
4718 }
4719 if (!n_Equal(pGetCoeff(p1), nn = n_Mult(pGetCoeff(p2),n, r->cf), r->cf))
4720 {
4721 n_Delete(&n, r->cf);
4722 n_Delete(&nn, r->cf);
4723 return FALSE;
4724 }
4725 n_Delete(&nn, r->cf);
4726 pIter(p1);
4727 pIter(p2);
4728 }
4729 n_Delete(&n, r->cf);
4730 return TRUE;
4731}
4732
4733/*2
4734* returns the length of a (numbers of monomials)
4735* respect syzComp
4736*/
4737poly p_Last(const poly p, int &l, const ring r)
4738{
4739 if (p == NULL)
4740 {
4741 l = 0;
4742 return NULL;
4743 }
4744 l = 1;
4745 poly a = p;
4746 if (! rIsSyzIndexRing(r))
4747 {
4748 poly next = pNext(a);
4749 while (next!=NULL)
4750 {
4751 a = next;
4752 next = pNext(a);
4753 l++;
4754 }
4755 }
4756 else
4757 {
4758 long unsigned curr_limit = rGetCurrSyzLimit(r);
4759 poly pp = a;
4760 while ((a=pNext(a))!=NULL)
4761 {
4762 if (__p_GetComp(a,r)<=curr_limit/*syzComp*/)
4763 l++;
4764 else break;
4765 pp = a;
4766 }
4767 a=pp;
4768 }
4769 return a;
4770}
4771
4772int p_Var(poly m,const ring r)
4773{
4774 if (m==NULL) return 0;
4775 if (pNext(m)!=NULL) return 0;
4776 int i,e=0;
4777 for (i=rVar(r); i>0; i--)
4778 {
4779 int exp=p_GetExp(m,i,r);
4780 if (exp==1)
4781 {
4782 if (e==0) e=i;
4783 else return 0;
4784 }
4785 else if (exp!=0)
4786 {
4787 return 0;
4788 }
4789 }
4790 return e;
4791}
4792
4793/*2
4794*the minimal index of used variables - 1
4795*/
4796int p_LowVar (poly p, const ring r)
4797{
4798 int k,l,lex;
4799
4800 if (p == NULL) return -1;
4801
4802 k = 32000;/*a very large dummy value*/
4803 while (p != NULL)
4804 {
4805 l = 1;
4806 lex = p_GetExp(p,l,r);
4807 while ((l < (rVar(r))) && (lex == 0))
4808 {
4809 l++;
4810 lex = p_GetExp(p,l,r);
4811 }
4812 l--;
4813 if (l < k) k = l;
4814 pIter(p);
4815 }
4816 return k;
4817}
4818
4819/*2
4820* verschiebt die Indizees der Modulerzeugenden um i
4821*/
4822void p_Shift (poly * p,int i, const ring r)
4823{
4824 poly qp1 = *p,qp2 = *p;/*working pointers*/
4825 int j = p_MaxComp(*p,r),k = p_MinComp(*p,r);
4826
4827 if (j+i < 0) return ;
4828 BOOLEAN toPoly= ((j == -i) && (j == k));
4829 while (qp1 != NULL)
4830 {
4831 if (toPoly || (__p_GetComp(qp1,r)+i > 0))
4832 {
4833 p_AddComp(qp1,i,r);
4834 p_SetmComp(qp1,r);
4835 qp2 = qp1;
4836 pIter(qp1);
4837 }
4838 else
4839 {
4840 if (qp2 == *p)
4841 {
4842 pIter(*p);
4843 p_LmDelete(&qp2,r);
4844 qp2 = *p;
4845 qp1 = *p;
4846 }
4847 else
4848 {
4849 qp2->next = qp1->next;
4850 if (qp1!=NULL) p_LmDelete(&qp1,r);
4851 qp1 = qp2->next;
4852 }
4853 }
4854 }
4855}
4856
4857/***************************************************************
4858 *
4859 * Storage Managament Routines
4860 *
4861 ***************************************************************/
4862
4863
4864static inline unsigned long GetBitFields(const long e,
4865 const unsigned int s, const unsigned int n)
4866{
4867#define Sy_bit_L(x) (((unsigned long)1L)<<(x))
4868 unsigned int i = 0;
4869 unsigned long ev = 0L;
4870 assume(n > 0 && s < BIT_SIZEOF_LONG);
4871 do
4872 {
4874 if (e > (long) i) ev |= Sy_bit_L(s+i);
4875 else break;
4876 i++;
4877 }
4878 while (i < n);
4879 return ev;
4880}
4881
4882// Short Exponent Vectors are used for fast divisibility tests
4883// ShortExpVectors "squeeze" an exponent vector into one word as follows:
4884// Let n = BIT_SIZEOF_LONG / pVariables.
4885// If n == 0 (i.e. pVariables > BIT_SIZE_OF_LONG), let m == the number
4886// of non-zero exponents. If (m>BIT_SIZEOF_LONG), then sev = ~0, else
4887// first m bits of sev are set to 1.
4888// Otherwise (i.e. pVariables <= BIT_SIZE_OF_LONG)
4889// represented by a bit-field of length n (resp. n+1 for some
4890// exponents). If the value of an exponent is greater or equal to n, then
4891// all of its respective n bits are set to 1. If the value of an exponent
4892// is smaller than n, say m, then only the first m bits of the respective
4893// n bits are set to 1, the others are set to 0.
4894// This way, we have:
4895// exp1 / exp2 ==> (ev1 & ~ev2) == 0, i.e.,
4896// if (ev1 & ~ev2) then exp1 does not divide exp2
4897unsigned long p_GetShortExpVector(const poly p, const ring r)
4898{
4899 assume(p != NULL);
4900 unsigned long ev = 0; // short exponent vector
4901 unsigned int n = BIT_SIZEOF_LONG / r->N; // number of bits per exp
4902 unsigned int m1; // highest bit which is filled with (n+1)
4903 unsigned int i=0;
4904 int j=1;
4905
4906 if (n == 0)
4907 {
4908 if (r->N <2*BIT_SIZEOF_LONG)
4909 {
4910 n=1;
4911 m1=0;
4912 }
4913 else
4914 {
4915 for (; j<=r->N; j++)
4916 {
4917 if (p_GetExp(p,j,r) > 0) i++;
4918 if (i == BIT_SIZEOF_LONG) break;
4919 }
4920 if (i>0)
4921 ev = ~(0UL) >> (BIT_SIZEOF_LONG - i);
4922 return ev;
4923 }
4924 }
4925 else
4926 {
4927 m1 = (n+1)*(BIT_SIZEOF_LONG - n*r->N);
4928 }
4929
4930 n++;
4931 while (i<m1)
4932 {
4933 ev |= GetBitFields(p_GetExp(p, j,r), i, n);
4934 i += n;
4935 j++;
4936 }
4937
4938 n--;
4939 while (i<BIT_SIZEOF_LONG)
4940 {
4941 ev |= GetBitFields(p_GetExp(p, j,r), i, n);
4942 i += n;
4943 j++;
4944 }
4945 return ev;
4946}
4947
4948
4949/// p_GetShortExpVector of p * pp
4950unsigned long p_GetShortExpVector(const poly p, const poly pp, const ring r)
4951{
4952 assume(p != NULL);
4953 assume(pp != NULL);
4954
4955 unsigned long ev = 0; // short exponent vector
4956 unsigned int n = BIT_SIZEOF_LONG / r->N; // number of bits per exp
4957 unsigned int m1; // highest bit which is filled with (n+1)
4958 int j=1;
4959 unsigned long i = 0L;
4960
4961 if (n == 0)
4962 {
4963 if (r->N <2*BIT_SIZEOF_LONG)
4964 {
4965 n=1;
4966 m1=0;
4967 }
4968 else
4969 {
4970 for (; j<=r->N; j++)
4971 {
4972 if (p_GetExp(p,j,r) > 0 || p_GetExp(pp,j,r) > 0) i++;
4973 if (i == BIT_SIZEOF_LONG) break;
4974 }
4975 if (i>0)
4976 ev = ~(0UL) >> (BIT_SIZEOF_LONG - i);
4977 return ev;
4978 }
4979 }
4980 else
4981 {
4982 m1 = (n+1)*(BIT_SIZEOF_LONG - n*r->N);
4983 }
4984
4985 n++;
4986 while (i<m1)
4987 {
4988 ev |= GetBitFields(p_GetExp(p, j,r) + p_GetExp(pp, j,r), i, n);
4989 i += n;
4990 j++;
4991 }
4992
4993 n--;
4994 while (i<BIT_SIZEOF_LONG)
4995 {
4996 ev |= GetBitFields(p_GetExp(p, j,r) + p_GetExp(pp, j,r), i, n);
4997 i += n;
4998 j++;
4999 }
5000 return ev;
5001}
5002
5003
5004
5005/***************************************************************
5006 *
5007 * p_ShallowDelete
5008 *
5009 ***************************************************************/
5010#undef LINKAGE
5011#define LINKAGE
5012#undef p_Delete__T
5013#define p_Delete__T p_ShallowDelete
5014#undef n_Delete__T
5015#define n_Delete__T(n, r) do {} while (0)
5016
5018
5019/***************************************************************/
5020/*
5021* compare a and b
5022*/
5023int p_Compare(const poly a, const poly b, const ring R)
5024{
5025 int r=p_Cmp(a,b,R);
5026 if ((r==0)&&(a!=NULL))
5027 {
5028 number h=n_Sub(pGetCoeff(a),pGetCoeff(b),R->cf);
5029 /* compare lead coeffs */
5030 r = -1+n_IsZero(h,R->cf)+2*n_GreaterZero(h,R->cf); /* -1: <, 0:==, 1: > */
5031 n_Delete(&h,R->cf);
5032 }
5033 else if (a==NULL)
5034 {
5035 if (b==NULL)
5036 {
5037 /* compare 0, 0 */
5038 r=0;
5039 }
5040 else if(p_IsConstant(b,R))
5041 {
5042 /* compare 0, const */
5043 r = 1-2*n_GreaterZero(pGetCoeff(b),R->cf); /* -1: <, 1: > */
5044 }
5045 }
5046 else if (b==NULL)
5047 {
5048 if (p_IsConstant(a,R))
5049 {
5050 /* compare const, 0 */
5051 r = -1+2*n_GreaterZero(pGetCoeff(a),R->cf); /* -1: <, 1: > */
5052 }
5053 }
5054 return(r);
5055}
5056
5057poly p_GcdMon(poly f, poly g, const ring r)
5058{
5059 assume(f!=NULL);
5060 assume(g!=NULL);
5061 assume(pNext(f)==NULL);
5062 poly G=p_Head(f,r);
5063 poly h=g;
5064 int *mf=(int*)omAlloc((r->N+1)*sizeof(int));
5065 p_GetExpV(f,mf,r);
5066 int *mh=(int*)omAlloc((r->N+1)*sizeof(int));
5067 BOOLEAN const_mon;
5068 BOOLEAN one_coeff=n_IsOne(pGetCoeff(G),r->cf);
5069 loop
5070 {
5071 if (h==NULL) break;
5072 if(!one_coeff)
5073 {
5074 number n=n_SubringGcd(pGetCoeff(G),pGetCoeff(h),r->cf);
5075 one_coeff=n_IsOne(n,r->cf);
5076 p_SetCoeff(G,n,r);
5077 }
5078 p_GetExpV(h,mh,r);
5079 const_mon=TRUE;
5080 for(unsigned j=r->N;j!=0;j--)
5081 {
5082 if (mh[j]<mf[j]) mf[j]=mh[j];
5083 if (mf[j]>0) const_mon=FALSE;
5084 }
5085 if (one_coeff && const_mon) break;
5086 pIter(h);
5087 }
5088 mf[0]=0;
5089 p_SetExpV(G,mf,r); // included is p_SetComp, p_Setm
5090 omFreeSize(mf,(r->N+1)*sizeof(int));
5091 omFreeSize(mh,(r->N+1)*sizeof(int));
5092 return G;
5093}
5094
5095poly p_CopyPowerProduct0(const poly p, number n, const ring r)
5096{
5098 poly np;
5099 omTypeAllocBin(poly, np, r->PolyBin);
5100 p_SetRingOfLm(np, r);
5101 memcpy(np->exp, p->exp, r->ExpL_Size*sizeof(long));
5102 pNext(np) = NULL;
5103 pSetCoeff0(np, n);
5104 return np;
5105}
5106
5107poly p_CopyPowerProduct(const poly p, const ring r)
5108{
5109 if (p == NULL) return NULL;
5110 return p_CopyPowerProduct0(p,n_Init(1,r->cf),r);
5111}
5112
5113poly p_Head0(const poly p, const ring r)
5114{
5115 if (p==NULL) return NULL;
5116 if (pGetCoeff(p)==NULL) return p_CopyPowerProduct0(p,NULL,r);
5117 return p_Head(p,r);
5118}
5119int p_MaxExpPerVar(poly p, int i, const ring r)
5120{
5121 int m=0;
5122 while(p!=NULL)
5123 {
5124 int mm=p_GetExp(p,i,r);
5125 if (mm>m) m=mm;
5126 pIter(p);
5127 }
5128 return m;
5129}
5130
Concrete implementation of enumerators over polynomials.
All the auxiliary stuff.
long int64
Definition: auxiliary.h:68
static int si_max(const int a, const int b)
Definition: auxiliary.h:124
#define BIT_SIZEOF_LONG
Definition: auxiliary.h:80
#define UNLIKELY(X)
Definition: auxiliary.h:404
int BOOLEAN
Definition: auxiliary.h:87
#define TRUE
Definition: auxiliary.h:100
#define FALSE
Definition: auxiliary.h:96
#define LIKELY(X)
Definition: auxiliary.h:403
void * ADDRESS
Definition: auxiliary.h:119
static int si_min(const int a, const int b)
Definition: auxiliary.h:125
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f )
Definition: cf_gcd.cc:676
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
int l
Definition: cfEzgcd.cc:100
int m
Definition: cfEzgcd.cc:128
int i
Definition: cfEzgcd.cc:132
int k
Definition: cfEzgcd.cc:99
return
Definition: cfGcdAlgExt.cc:218
Variable x
Definition: cfModGcd.cc:4081
int p
Definition: cfModGcd.cc:4077
g
Definition: cfModGcd.cc:4089
CanonicalForm cf
Definition: cfModGcd.cc:4082
CanonicalForm b
Definition: cfModGcd.cc:4102
FILE * f
Definition: checklibs.c:9
poly singclap_pdivide(poly f, poly g, const ring r)
Definition: clapsing.cc:624
This is a polynomial enumerator for simple iteration over coefficients of polynomials.
Definition: intvec.h:23
static FORCE_INLINE number n_Mult(number a, number b, const coeffs r)
return the product of 'a' and 'b', i.e., a*b
Definition: coeffs.h:636
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....
Definition: coeffs.h:783
static FORCE_INLINE number n_Copy(number n, const coeffs r)
return a copy of 'n'
Definition: coeffs.h:451
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:695
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:603
static FORCE_INLINE BOOLEAN nCoeff_is_GF(const coeffs r)
Definition: coeffs.h:839
static FORCE_INLINE BOOLEAN nCoeff_is_Extension(const coeffs r)
Definition: coeffs.h:846
number ndCopyMap(number a, const coeffs src, const coeffs dst)
Definition: numbers.cc:282
#define n_Test(a, r)
BOOLEAN n_Test(number a, const coeffs r)
Definition: coeffs.h:712
@ n_algExt
used for all algebraic extensions, i.e., the top-most extension in an extension tower is algebraic
Definition: coeffs.h:35
@ n_transExt
used for all transcendental extensions, i.e., the top-most extension in an extension tower is transce...
Definition: coeffs.h:38
static FORCE_INLINE number n_Gcd(number a, number b, const coeffs r)
in Z: return the gcd of 'a' and 'b' in Z/nZ, Z/2^kZ: computed as in the case Z in Z/pZ,...
Definition: coeffs.h:664
#define n_New(n, r)
Definition: coeffs.h:440
static FORCE_INLINE number n_Invers(number a, const coeffs r)
return the multiplicative inverse of 'a'; raise an error if 'a' is not invertible
Definition: coeffs.h:564
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:515
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:622
static FORCE_INLINE BOOLEAN n_GreaterZero(number n, const coeffs r)
ordered fields: TRUE iff 'n' is positive; in Z/pZ: TRUE iff 0 < m <= roundedBelow(p/2),...
Definition: coeffs.h:494
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:700
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:557
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:632
static FORCE_INLINE number n_Farey(number a, number b, const coeffs r)
Definition: coeffs.h:767
static FORCE_INLINE number n_Div(number a, number b, const coeffs r)
return the quotient of 'a' and 'b', i.e., a/b; raises an error if 'b' is not invertible in r exceptio...
Definition: coeffs.h:615
static FORCE_INLINE BOOLEAN nCoeff_is_Q(const coeffs r)
Definition: coeffs.h:806
static FORCE_INLINE BOOLEAN n_IsZero(number n, const coeffs r)
TRUE iff 'n' represents the zero element.
Definition: coeffs.h:464
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:570
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:532
static FORCE_INLINE number n_Sub(number a, number b, const coeffs r)
return the difference of 'a' and 'b', i.e., a-b
Definition: coeffs.h:655
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:935
static FORCE_INLINE BOOLEAN nCoeff_is_Ring(const coeffs r)
Definition: coeffs.h:730
static FORCE_INLINE n_coeffType getCoeffType(const coeffs r)
Returns the type of coeffs domain.
Definition: coeffs.h:421
static FORCE_INLINE number n_ChineseRemainderSym(number *a, number *b, int rl, BOOLEAN sym, CFArray &inv_cache, const coeffs r)
Definition: coeffs.h:764
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete 'p'
Definition: coeffs.h:455
static FORCE_INLINE void n_Write(number n, const coeffs r, const BOOLEAN bShortOut=TRUE)
Definition: coeffs.h:591
static FORCE_INLINE BOOLEAN nCoeff_is_Zp(const coeffs r)
Definition: coeffs.h:800
static FORCE_INLINE BOOLEAN nCoeff_is_Q_a(const coeffs r)
Definition: coeffs.h:885
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:538
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:928
static FORCE_INLINE BOOLEAN n_DivBy(number a, number b, const coeffs r)
test whether 'a' is divisible 'b'; for r encoding a field: TRUE iff 'b' does not represent zero in Z:...
Definition: coeffs.h:753
static FORCE_INLINE BOOLEAN nCoeff_is_algExt(const coeffs r)
TRUE iff r represents an algebraic extension field.
Definition: coeffs.h:910
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....
Definition: coeffs.h:598
static FORCE_INLINE BOOLEAN n_Equal(number a, number b, const coeffs r)
TRUE iff 'a' and 'b' represent the same number; they may have different representations.
Definition: coeffs.h:460
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:608
static FORCE_INLINE number n_SubringGcd(number a, number b, const coeffs r)
Definition: coeffs.h:666
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 void n_Normalize(number &n, const coeffs r)
inplace-normalization of n; produces some canonical representation of n;
Definition: coeffs.h:578
static FORCE_INLINE BOOLEAN n_IsOne(number n, const coeffs r)
TRUE iff 'n' represents the one element.
Definition: coeffs.h:468
static FORCE_INLINE BOOLEAN nCoeff_is_transExt(const coeffs r)
TRUE iff r represents a transcendental extension field.
Definition: coeffs.h:918
#define Print
Definition: emacs.cc:80
#define WarnS
Definition: emacs.cc:78
return result
Definition: facAbsBiFact.cc:75
const CanonicalForm int s
Definition: facAbsFact.cc:51
const CanonicalForm int const CFList const Variable & y
Definition: facAbsFact.cc:53
CanonicalForm res
Definition: facAbsFact.cc:60
const CanonicalForm & w
Definition: facAbsFact.cc:51
CanonicalForm subst(const CanonicalForm &f, const CFList &a, const CFList &b, const CanonicalForm &Rstar, bool isFunctionField)
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
for(j=0;j< factors.length();j++)
Definition: facHensel.cc:129
int j
Definition: facHensel.cc:110
int comp(const CanonicalForm &A, const CanonicalForm &B)
compare polynomials
static int max(int a, int b)
Definition: fast_mult.cc:264
VAR short errorreported
Definition: feFopen.cc:23
void WerrorS(const char *s)
Definition: feFopen.cc:24
if(!FE_OPT_NO_SHELL_FLAG)(void) system(sys)
const char * eati(const char *s, int *i)
Definition: reporter.cc:373
#define D(A)
Definition: gentable.cc:131
#define STATIC_VAR
Definition: globaldefs.h:7
#define VAR
Definition: globaldefs.h:5
STATIC_VAR poly last
Definition: hdegree.cc:1173
#define exponent
STATIC_VAR int offset
Definition: janet.cc:29
STATIC_VAR TreeM * G
Definition: janet.cc:31
STATIC_VAR Poly * h
Definition: janet.cc:971
ListNode * next
Definition: janet.h:31
static bool rIsSCA(const ring r)
Definition: nc.h:190
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:3203
LINLINE number nlAdd(number la, number li, const coeffs r)
Definition: longrat.cc:2702
LINLINE number nlSub(number la, number li, const coeffs r)
Definition: longrat.cc:2768
LINLINE void nlDelete(number *a, const coeffs r)
Definition: longrat.cc:2667
BOOLEAN nlGreaterZero(number za, const coeffs r)
Definition: longrat.cc:1309
number nlGcd(number a, number b, const coeffs r)
Definition: longrat.cc:1346
void nlNormalize(number &x, const coeffs r)
Definition: longrat.cc:1487
#define assume(x)
Definition: mod2.h:389
int dReportError(const char *fmt,...)
Definition: dError.cc:43
#define p_SetCoeff0(p, n, r)
Definition: monomials.h:60
#define p_GetComp(p, r)
Definition: monomials.h:64
#define pIter(p)
Definition: monomials.h:37
#define pNext(p)
Definition: monomials.h:36
#define p_LmCheckPolyRing1(p, r)
Definition: monomials.h:177
#define p_LmCheckPolyRing2(p, r)
Definition: monomials.h:199
#define pSetCoeff0(p, n)
Definition: monomials.h:59
#define p_GetCoeff(p, r)
Definition: monomials.h:50
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
Definition: monomials.h:44
#define POLY_NEGWEIGHT_OFFSET
Definition: monomials.h:236
#define __p_GetComp(p, r)
Definition: monomials.h:63
#define p_SetRingOfLm(p, r)
Definition: monomials.h:144
#define rRing_has_Comp(r)
Definition: monomials.h:266
#define pAssume(cond)
Definition: monomials.h:90
gmp_float exp(const gmp_float &a)
Definition: mpr_complex.cc:357
The main handler for Singular numbers which are suitable for Singular polynomials.
Definition: lq.h:40
number ndGcd(number, number, const coeffs r)
Definition: numbers.cc:192
void ndNormalize(number &, const coeffs)
Definition: numbers.cc:190
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define omReallocSize(addr, o_size, size)
Definition: omAllocDecl.h:220
#define omTypeAllocBin(type, addr, bin)
Definition: omAllocDecl.h:203
#define omFree(addr)
Definition: omAllocDecl.h:261
#define omAlloc0(size)
Definition: omAllocDecl.h:211
#define NULL
Definition: omList.c:12
#define TEST_OPT_INTSTRATEGY
Definition: options.h:110
#define TEST_OPT_PROT
Definition: options.h:103
#define TEST_OPT_CONTENTSB
Definition: options.h:127
poly p_Diff(poly a, int k, const ring r)
Definition: p_polys.cc:1894
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
poly p_DivideM(poly a, poly b, const ring r)
Definition: p_polys.cc:1574
int p_IsPurePower(const poly p, const ring r)
return i, if head depends only on var(i)
Definition: p_polys.cc:1226
void p_Setm_WFirstTotalDegree(poly p, const ring r)
Definition: p_polys.cc:554
poly pp_Jet(poly p, int m, const ring R)
Definition: p_polys.cc:4474
STATIC_VAR pLDegProc pOldLDeg
Definition: p_polys.cc:3777
void p_Cleardenom_n(poly ph, const ring r, number &c)
Definition: p_polys.cc:3019
long pLDegb(poly p, int *l, const ring r)
Definition: p_polys.cc:811
long pLDeg1_Totaldegree(poly p, int *l, const ring r)
Definition: p_polys.cc:975
long p_WFirstTotalDegree(poly p, const ring r)
Definition: p_polys.cc:596
poly p_Farey(poly p, number N, const ring r)
Definition: p_polys.cc:54
void p_Content_n(poly ph, number &c, const ring r)
Definition: p_polys.cc:2349
long pLDeg1_WFirstTotalDegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1038
void pRestoreDegProcs(ring r, pFDegProc old_FDeg, pLDegProc old_lDeg)
Definition: p_polys.cc:3765
long pLDeg1c_WFirstTotalDegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1068
poly n_PermNumber(const number z, const int *par_perm, const int, const ring src, const ring dst)
Definition: p_polys.cc:4143
static poly p_DiffOpM(poly a, poly b, BOOLEAN multiply, const ring r)
Definition: p_polys.cc:1930
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:1866
int p_Size(poly p, const ring r)
Definition: p_polys.cc:3318
void p_Setm_Dummy(poly p, const ring r)
Definition: p_polys.cc:541
static poly p_Invers(int n, poly u, intvec *w, const ring R)
Definition: p_polys.cc:4585
poly p_GcdMon(poly f, poly g, const ring r)
polynomial gcd for f=mon
Definition: p_polys.cc:5057
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:4692
int p_LowVar(poly p, const ring r)
the minimal index of used variables - 1
Definition: p_polys.cc:4796
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:1638
poly p_Homogen(poly p, int varnum, const ring r)
Definition: p_polys.cc:3335
poly p_Subst(poly p, int n, poly e, const ring r)
Definition: p_polys.cc:4074
static BOOLEAN p_ExpVectorEqual(poly p1, poly p2, const ring r1, const ring r2)
Definition: p_polys.cc:4642
BOOLEAN p_HasNotCF(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1329
void p_Content(poly ph, const ring r)
Definition: p_polys.cc:2291
int p_Weight(int i, const ring r)
Definition: p_polys.cc:705
void p_Setm_TotalDegree(poly p, const ring r)
Definition: p_polys.cc:547
poly p_CopyPowerProduct(const poly p, const ring r)
like p_Head, but with coefficient 1
Definition: p_polys.cc:5107
poly pp_DivideM(poly a, poly b, const ring r)
Definition: p_polys.cc:1629
STATIC_VAR int _componentsExternal
Definition: p_polys.cc:148
void p_SimpleContent(poly ph, int smax, const ring r)
Definition: p_polys.cc:2629
poly p_ISet(long i, const ring r)
returns the poly representing the integer i
Definition: p_polys.cc:1297
STATIC_VAR long * _componentsShifted
Definition: p_polys.cc:147
void p_Vec2Polys(poly v, poly **p, int *len, const ring r)
Definition: p_polys.cc:3741
static poly p_Subst0(poly p, int n, const ring r)
Definition: p_polys.cc:4049
poly p_DiffOp(poly a, poly b, BOOLEAN multiply, const ring r)
Definition: p_polys.cc:1969
static unsigned long p_GetMaxExpL2(unsigned long l1, unsigned long l2, const ring r, unsigned long number_of_exp)
Definition: p_polys.cc:1107
poly p_Jet(poly p, int m, const ring R)
Definition: p_polys.cc:4502
poly p_TakeOutComp1(poly *p, int k, const ring r)
Definition: p_polys.cc:3500
poly p_TakeOutComp(poly *p, int k, const ring r)
Definition: p_polys.cc:3551
long pLDeg1c_Deg(poly p, int *l, const ring r)
Definition: p_polys.cc:941
long pLDeg1(poly p, int *l, const ring r)
Definition: p_polys.cc:841
static number * pnBin(int exp, const ring r)
Definition: p_polys.cc:2054
void p_Shift(poly *p, int i, const ring r)
shifts components of the vector p by i
Definition: p_polys.cc:4822
static void pnFreeBin(number *bin, int exp, const coeffs r)
Definition: p_polys.cc:2085
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:4246
poly p_Power(poly p, int i, const ring r)
Definition: p_polys.cc:2193
poly p_Div_nn(poly p, const number n, const ring r)
Definition: p_polys.cc:1501
void p_Normalize(poly p, const ring r)
Definition: p_polys.cc:3929
void p_DeleteComp(poly *p, int k, const ring r)
Definition: p_polys.cc:3660
poly p_mInit(const char *st, BOOLEAN &ok, const ring r)
Definition: p_polys.cc:1442
poly p_MDivide(poly a, poly b, const ring r)
Definition: p_polys.cc:1488
void p_ContentRat(poly &ph, const ring r)
Definition: p_polys.cc:1740
void p_Norm(poly p1, const ring r)
Definition: p_polys.cc:3835
poly p_Div_mm(poly p, const poly m, const ring r)
divide polynomial by monomial
Definition: p_polys.cc:1534
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
int p_MinDeg(poly p, intvec *w, const ring R)
Definition: p_polys.cc:4564
int p_MaxExpPerVar(poly p, int i, const ring r)
max exponent of variable x_i in p
Definition: p_polys.cc:5119
STATIC_VAR BOOLEAN pOldLexOrder
Definition: p_polys.cc:3778
int p_Compare(const poly a, const poly b, const ring R)
Definition: p_polys.cc:5023
void p_Setm_Syz(poly p, ring r, int *Components, long *ShiftedComponents)
Definition: p_polys.cc:531
STATIC_VAR pFDegProc pOldFDeg
Definition: p_polys.cc:3776
void p_LmDeleteAndNextRat(poly *p, int ishift, ring r)
Definition: p_polys.cc:1696
unsigned long p_GetShortExpVector(const poly p, const ring r)
Definition: p_polys.cc:4897
BOOLEAN p_IsHomogeneousW(poly p, const intvec *w, const ring r)
Definition: p_polys.cc:3408
VAR BOOLEAN pSetm_error
Definition: p_polys.cc:150
long pLDeg1_Deg(poly p, int *l, const ring r)
Definition: p_polys.cc:910
poly p_Series(int n, poly p, poly u, intvec *w, const ring R)
Definition: p_polys.cc:4614
void p_ProjectiveUnique(poly ph, const ring r)
Definition: p_polys.cc:3208
long p_WTotaldegree(poly p, const ring r)
Definition: p_polys.cc:613
long p_DegW(poly p, const int *w, const ring R)
Definition: p_polys.cc:690
p_SetmProc p_GetSetmProc(const ring r)
Definition: p_polys.cc:560
void p_Setm_General(poly p, const ring r)
Definition: p_polys.cc:158
BOOLEAN p_OneComp(poly p, const ring r)
return TRUE if all monoms have the same component
Definition: p_polys.cc:1208
poly p_Cleardenom(poly p, const ring r)
Definition: p_polys.cc:2910
long pLDeg1c(poly p, int *l, const ring r)
Definition: p_polys.cc:877
void p_Split(poly p, poly *h)
Definition: p_polys.cc:1320
long pLDeg1c_Totaldegree(poly p, int *l, const ring r)
Definition: p_polys.cc:1005
poly p_GetCoeffRat(poly p, int ishift, ring r)
Definition: p_polys.cc:1718
BOOLEAN p_VectorHasUnitB(poly p, int *k, const ring r)
Definition: p_polys.cc:3444
long pLDeg0c(poly p, int *l, const ring r)
Definition: p_polys.cc:770
poly p_Vec2Poly(poly v, int k, const ring r)
Definition: p_polys.cc:3689
poly p_LcmRat(const poly a, const poly b, const long lCompM, const ring r)
Definition: p_polys.cc:1673
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 poly p_TwoMonPower(poly p, int exp, const ring r)
Definition: p_polys.cc:2102
void p_SetModDeg(intvec *w, ring r)
Definition: p_polys.cc:3789
BOOLEAN p_HasNotCFRing(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1345
long pLDeg0(poly p, int *l, const ring r)
Definition: p_polys.cc:739
int p_Var(poly m, const ring r)
Definition: p_polys.cc:4772
poly p_One(const ring r)
Definition: p_polys.cc:1313
poly p_Sub(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1986
void p_VectorHasUnit(poly p, int *k, int *len, const ring r)
Definition: p_polys.cc:3467
static void p_SplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero, const ring r)
Definition: p_polys.cc:3945
int p_IsUnivariate(poly p, const ring r)
return i, if poly depends only on var(i)
Definition: p_polys.cc:1247
STATIC_VAR int * _components
Definition: p_polys.cc:146
poly p_NSet(number n, const ring r)
returns the poly representing the number n, destroys n
Definition: p_polys.cc:1469
void pSetDegProcs(ring r, pFDegProc new_FDeg, pLDegProc new_lDeg)
Definition: p_polys.cc:3753
void pEnlargeSet(poly **p, int l, int increment)
Definition: p_polys.cc:3812
long p_WDegree(poly p, const ring r)
Definition: p_polys.cc:714
static long pModDeg(poly p, ring r)
Definition: p_polys.cc:3780
BOOLEAN p_IsHomogeneous(poly p, const ring r)
Definition: p_polys.cc:3384
poly p_Head0(const poly p, const ring r)
like p_Head, but allow NULL coeff
Definition: p_polys.cc:5113
static poly p_MonMultC(poly p, poly q, const ring rr)
Definition: p_polys.cc:2040
static poly p_Pow_charp(poly p, int i, const ring r)
Definition: p_polys.cc:2181
poly pp_JetW(poly p, int m, int *w, const ring R)
Definition: p_polys.cc:4519
long p_Deg(poly a, const ring r)
Definition: p_polys.cc:587
static poly p_Subst1(poly p, int n, const ring r)
Definition: p_polys.cc:3981
poly p_Last(const poly p, int &l, const ring r)
Definition: p_polys.cc:4737
poly p_CopyPowerProduct0(const poly p, number n, const ring r)
like p_Head, but with coefficient n
Definition: p_polys.cc:5095
static void p_MonMult(poly p, poly q, const ring r)
Definition: p_polys.cc:2020
number p_InitContent(poly ph, const ring r)
Definition: p_polys.cc:2700
void p_Vec2Array(poly v, poly *p, int len, const ring r)
vector to already allocated array (len>=p_MaxComp(v,r))
Definition: p_polys.cc:3711
static poly p_MonPower(poly p, int exp, const ring r)
Definition: p_polys.cc:1996
void p_ContentForGB(poly ph, const ring r)
Definition: p_polys.cc:2420
static poly p_Subst2(poly p, int n, number e, const ring r)
Definition: p_polys.cc:4008
void p_Lcm(const poly a, const poly b, poly m, const ring r)
Definition: p_polys.cc:1651
static unsigned long GetBitFields(const long e, const unsigned int s, const unsigned int n)
Definition: p_polys.cc:4864
poly p_ChineseRemainder(poly *xx, number *x, number *q, int rl, CFArray &inv_cache, const ring R)
Definition: p_polys.cc:88
const char * p_Read(const char *st, poly &rc, const ring r)
Definition: p_polys.cc:1370
#define Sy_bit_L(x)
poly p_JetW(poly p, int m, int *w, const ring R)
Definition: p_polys.cc:4546
BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r)
Definition: p_polys.cc:4628
static poly p_Pow(poly p, int i, const ring r)
Definition: p_polys.cc:2167
static poly p_Neg(poly p, const ring r)
Definition: p_polys.h:1109
static void p_ExpVectorSum(poly pr, poly p1, poly p2, const ring r)
Definition: p_polys.h:1427
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:938
static void p_LmDelete(poly p, const ring r)
Definition: p_polys.h:725
static poly p_Mult_q(poly p, poly q, const ring r)
Definition: p_polys.h:1116
BOOLEAN p_LmCheckPolyRing(poly p, ring r)
Definition: pDebug.cc:120
static void p_ExpVectorAdd(poly p1, poly p2, const ring r)
Definition: p_polys.h:1413
static unsigned long p_SubComp(poly p, unsigned long v, ring r)
Definition: p_polys.h:455
static long p_AddExp(poly p, int v, long ee, ring r)
Definition: p_polys.h:608
static poly p_LmInit(poly p, const ring r)
Definition: p_polys.h:1337
#define p_LmEqual(p1, p2, r)
Definition: p_polys.h:1725
static int p_Cmp(poly p1, poly p2, ring r)
Definition: p_polys.h:1729
void p_Write(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:342
static void p_SetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1546
static int p_Comp_k_n(poly a, poly b, int k, ring r)
Definition: p_polys.h:642
static void p_SetCompP(poly p, int i, ring r)
Definition: p_polys.h:256
static unsigned long p_SetExp(poly p, const unsigned long e, const unsigned long iBitmask, const int VarOffset)
set a single variable exponent @Note: VarOffset encodes the position in p->exp
Definition: p_polys.h:490
static long p_MinComp(poly p, ring lmRing, ring tailRing)
Definition: p_polys.h:315
static unsigned long p_SetComp(poly p, unsigned long c, ring r)
Definition: p_polys.h:249
static long p_IncrExp(poly p, int v, ring r)
Definition: p_polys.h:593
static void p_ExpVectorSub(poly p1, poly p2, const ring r)
Definition: p_polys.h:1442
static unsigned long p_AddComp(poly p, unsigned long v, ring r)
Definition: p_polys.h:449
static void p_Setm(poly p, const ring r)
Definition: p_polys.h:235
#define p_SetmComp
Definition: p_polys.h:246
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:414
static poly pReverse(poly p)
Definition: p_polys.h:337
static BOOLEAN p_LmIsConstantComp(const poly p, const ring r)
Definition: p_polys.h:1008
static poly p_Head(const poly p, const ring r)
copy the (leading) term of p
Definition: p_polys.h:862
static int p_LmCmp(poly p, poly q, const ring r)
Definition: p_polys.h:1582
static long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset)
get a single variable exponent @Note: the integer VarOffset encodes:
Definition: p_polys.h:471
static long p_MultExp(poly p, int v, long ee, ring r)
Definition: p_polys.h:623
static BOOLEAN p_IsConstant(const poly p, const ring r)
Definition: p_polys.h:2005
static poly p_GetExp_k_n(poly p, int l, int k, const ring r)
Definition: p_polys.h:1374
static BOOLEAN p_DivisibleBy(poly a, poly b, const ring r)
Definition: p_polys.h:1906
static long p_MaxComp(poly p, ring lmRing, ring tailRing)
Definition: p_polys.h:294
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:903
static long p_DecrExp(poly p, int v, ring r)
Definition: p_polys.h:600
static unsigned pLength(poly a)
Definition: p_polys.h:191
static void p_GetExpV(poly p, int *ev, const ring r)
Definition: p_polys.h:1522
BOOLEAN p_CheckPolyRing(poly p, ring r)
Definition: pDebug.cc:112
static long p_GetOrder(poly p, ring r)
Definition: p_polys.h:423
static poly p_LmFreeAndNext(poly p, ring)
Definition: p_polys.h:713
static poly p_Mult_mm(poly p, poly m, const ring r)
Definition: p_polys.h:1053
static void p_LmFree(poly p, ring)
Definition: p_polys.h:685
static poly p_Init(const ring r, omBin bin)
Definition: p_polys.h:1322
static poly p_LmDeleteAndNext(poly p, const ring r)
Definition: p_polys.h:757
static poly p_SortAdd(poly p, const ring r, BOOLEAN revert=FALSE)
Definition: p_polys.h:1221
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition: p_polys.h:848
static long p_Totaldegree(poly p, const ring r)
Definition: p_polys.h:1509
#define p_Test(p, r)
Definition: p_polys.h:162
#define __p_Mult_nn(p, n, r)
Definition: p_polys.h:973
void p_wrp(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:373
poly singclap_gcd(poly f, poly g, const ring r)
polynomial gcd via singclap_gcd_r resp. idSyzygies destroys f and g
Definition: polys.cc:380
#define NUM
Definition: readcf.cc:180
void PrintS(const char *s)
Definition: reporter.cc:284
void Werror(const char *fmt,...)
Definition: reporter.cc:189
BOOLEAN rOrd_SetCompRequiresSetm(const ring r)
return TRUE if p_SetComp requires p_Setm
Definition: ring.cc:1993
void rWrite(ring r, BOOLEAN details)
Definition: ring.cc:226
int r_IsRingVar(const char *n, char **names, int N)
Definition: ring.cc:212
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:1799
static BOOLEAN rField_is_Zp_a(const ring r)
Definition: ring.h:530
static BOOLEAN rField_is_Z(const ring r)
Definition: ring.h:510
static BOOLEAN rField_is_Zp(const ring r)
Definition: ring.h:501
void(* p_SetmProc)(poly p, const ring r)
Definition: ring.h:39
static BOOLEAN rIsPluralRing(const ring r)
we must always have this test!
Definition: ring.h:400
ro_typ ord_typ
Definition: ring.h:220
long(* pFDegProc)(poly p, ring r)
Definition: ring.h:38
static int rGetCurrSyzLimit(const ring r)
Definition: ring.h:724
long(* pLDegProc)(poly p, int *length, ring r)
Definition: ring.h:37
static BOOLEAN rIsRatGRing(const ring r)
Definition: ring.h:427
static int rPar(const ring r)
(r->cf->P)
Definition: ring.h:600
@ ro_wp64
Definition: ring.h:55
@ ro_syz
Definition: ring.h:60
@ ro_cp
Definition: ring.h:58
@ ro_dp
Definition: ring.h:52
@ ro_is
Definition: ring.h:61
@ ro_wp_neg
Definition: ring.h:56
@ ro_wp
Definition: ring.h:53
@ ro_isTemp
Definition: ring.h:61
@ ro_am
Definition: ring.h:54
@ ro_syzcomp
Definition: ring.h:59
static int rInternalChar(const ring r)
Definition: ring.h:690
static BOOLEAN rIsLPRing(const ring r)
Definition: ring.h:411
@ ringorder_lp
Definition: ring.h:77
@ ringorder_a
Definition: ring.h:70
@ ringorder_am
Definition: ring.h:88
@ ringorder_a64
for int64 weights
Definition: ring.h:71
@ ringorder_rs
opposite of ls
Definition: ring.h:92
@ ringorder_C
Definition: ring.h:73
@ ringorder_S
S?
Definition: ring.h:75
@ ringorder_ds
Definition: ring.h:84
@ ringorder_Dp
Definition: ring.h:80
@ ringorder_unspec
Definition: ring.h:94
@ ringorder_L
Definition: ring.h:89
@ ringorder_Ds
Definition: ring.h:85
@ ringorder_dp
Definition: ring.h:78
@ ringorder_c
Definition: ring.h:72
@ ringorder_rp
Definition: ring.h:79
@ ringorder_aa
for idElimination, like a, except pFDeg, pWeigths ignore it
Definition: ring.h:91
@ ringorder_no
Definition: ring.h:69
@ ringorder_Wp
Definition: ring.h:82
@ ringorder_ws
Definition: ring.h:86
@ ringorder_Ws
Definition: ring.h:87
@ ringorder_IS
Induced (Schreyer) ordering.
Definition: ring.h:93
@ ringorder_ls
Definition: ring.h:83
@ ringorder_s
s?
Definition: ring.h:76
@ ringorder_wp
Definition: ring.h:81
@ ringorder_M
Definition: ring.h:74
static BOOLEAN rField_is_Q_a(const ring r)
Definition: ring.h:540
static BOOLEAN rField_is_Q(const ring r)
Definition: ring.h:507
static BOOLEAN rField_has_Units(const ring r)
Definition: ring.h:491
static BOOLEAN rIsNCRing(const ring r)
Definition: ring.h:421
static BOOLEAN rIsSyzIndexRing(const ring r)
Definition: ring.h:721
static BOOLEAN rField_is_GF(const ring r)
Definition: ring.h:522
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:593
union sro_ord::@1 data
#define rField_is_Ring(R)
Definition: ring.h:486
Definition: ring.h:219
void sBucket_Add_m(sBucket_pt bucket, poly p)
Definition: sbuckets.cc:173
sBucket_pt sBucketCreate(const ring r)
Definition: sbuckets.cc:96
void sBucketDestroyAdd(sBucket_pt bucket, poly *p, int *length)
Definition: sbuckets.h:68
static short scaLastAltVar(ring r)
Definition: sca.h:25
static short scaFirstAltVar(ring r)
Definition: sca.h:18
poly p_LPSubst(poly p, int n, poly e, const ring r)
Definition: shiftop.cc:912
int status int void size_t count
Definition: si_signals.h:59
#define IDELEMS(i)
Definition: simpleideals.h:23
#define R
Definition: sirandom.c:27
#define loop
Definition: structs.h:75
number ntInit(long i, const coeffs cf)
Definition: transext.cc:704
int * iv2array(intvec *iv, const ring R)
Definition: weight.cc:200
long totaldegreeWecart_IV(poly p, ring r, const int *w)
Definition: weight.cc:231