1 | /**************************************** |
---|
2 | * Computer Algebra System SINGULAR * |
---|
3 | ****************************************/ |
---|
4 | /*************************************************************** |
---|
5 | * File: p_polys.h |
---|
6 | * Purpose: declaration of poly stuf which are independent of |
---|
7 | * currRing |
---|
8 | * Author: obachman (Olaf Bachmann) |
---|
9 | * Created: 9/00 |
---|
10 | *******************************************************************/ |
---|
11 | /*************************************************************** |
---|
12 | * Purpose: implementation of poly procs which iter over ExpVector |
---|
13 | * Author: obachman (Olaf Bachmann) |
---|
14 | * Created: 8/00 |
---|
15 | *******************************************************************/ |
---|
16 | #ifndef P_POLYS_H |
---|
17 | #define P_POLYS_H |
---|
18 | |
---|
19 | #include "misc/mylimits.h" |
---|
20 | #include "misc/intvec.h" |
---|
21 | #include "coeffs/coeffs.h" |
---|
22 | |
---|
23 | #include "polys/monomials/monomials.h" |
---|
24 | #include "polys/monomials/ring.h" |
---|
25 | |
---|
26 | #include "polys/templates/p_MemAdd.h" |
---|
27 | #include "polys/templates/p_MemCmp.h" |
---|
28 | #include "polys/templates/p_Procs.h" |
---|
29 | |
---|
30 | #include "polys/sbuckets.h" |
---|
31 | |
---|
32 | #ifdef HAVE_PLURAL |
---|
33 | #include "polys/nc/nc.h" |
---|
34 | #endif |
---|
35 | |
---|
36 | poly p_Farey(poly p, number N, const ring r); |
---|
37 | /* |
---|
38 | * xx,q: arrays of length 0..rl-1 |
---|
39 | * xx[i]: SB mod q[i] |
---|
40 | * assume: char=0 |
---|
41 | * assume: q[i]!=0 |
---|
42 | * destroys xx |
---|
43 | */ |
---|
44 | poly p_ChineseRemainder(poly *xx, number *x,number *q, int rl, CFArray &inv_cache, const ring R); |
---|
45 | /*************************************************************** |
---|
46 | * |
---|
47 | * Divisiblity tests, args must be != NULL, except for |
---|
48 | * pDivisbleBy |
---|
49 | * |
---|
50 | ***************************************************************/ |
---|
51 | unsigned long p_GetShortExpVector(const poly a, const ring r); |
---|
52 | |
---|
53 | /// p_GetShortExpVector of p * pp |
---|
54 | unsigned long p_GetShortExpVector(const poly p, const poly pp, const ring r); |
---|
55 | |
---|
56 | #ifdef HAVE_RINGS |
---|
57 | /*! divisibility check over ground ring (which may contain zero divisors); |
---|
58 | TRUE iff LT(f) divides LT(g), i.e., LT(f)*c*m = LT(g), for some |
---|
59 | coefficient c and some monomial m; |
---|
60 | does not take components into account |
---|
61 | */ |
---|
62 | BOOLEAN p_DivisibleByRingCase(poly f, poly g, const ring r); |
---|
63 | #endif |
---|
64 | |
---|
65 | /*************************************************************** |
---|
66 | * |
---|
67 | * Misc things on polys |
---|
68 | * |
---|
69 | ***************************************************************/ |
---|
70 | |
---|
71 | poly p_One(const ring r); |
---|
72 | |
---|
73 | int p_MinDeg(poly p,intvec *w, const ring R); |
---|
74 | |
---|
75 | long p_DegW(poly p, const int *w, const ring R); |
---|
76 | |
---|
77 | /// return TRUE if all monoms have the same component |
---|
78 | BOOLEAN p_OneComp(poly p, const ring r); |
---|
79 | |
---|
80 | /// return i, if head depends only on var(i) |
---|
81 | int p_IsPurePower(const poly p, const ring r); |
---|
82 | |
---|
83 | /// return i, if poly depends only on var(i) |
---|
84 | int p_IsUnivariate(poly p, const ring r); |
---|
85 | |
---|
86 | /// set entry e[i] to 1 if var(i) occurs in p, ignore var(j) if e[j]>0 |
---|
87 | /// return #(e[i]>0) |
---|
88 | int p_GetVariables(poly p, int * e, const ring r); |
---|
89 | |
---|
90 | /// returns the poly representing the integer i |
---|
91 | poly p_ISet(long i, const ring r); |
---|
92 | |
---|
93 | /// returns the poly representing the number n, destroys n |
---|
94 | poly p_NSet(number n, const ring r); |
---|
95 | |
---|
96 | void p_Vec2Polys(poly v, poly**p, int *len, const ring r); |
---|
97 | poly p_Vec2Poly(poly v, int k, const ring r); |
---|
98 | |
---|
99 | /// julia: vector to already allocated array (len=p_MaxComp(v,r)) |
---|
100 | void p_Vec2Array(poly v, poly *p, int len, const ring r); |
---|
101 | |
---|
102 | /*************************************************************** |
---|
103 | * |
---|
104 | * Copying/Deletion of polys: args may be NULL |
---|
105 | * |
---|
106 | ***************************************************************/ |
---|
107 | |
---|
108 | // simply deletes monomials, does not free coeffs |
---|
109 | void p_ShallowDelete(poly *p, const ring r); |
---|
110 | |
---|
111 | |
---|
112 | |
---|
113 | /*************************************************************** |
---|
114 | * |
---|
115 | * Copying/Deleteion of polys: args may be NULL |
---|
116 | * - p/q as arg mean a poly |
---|
117 | * - m a monomial |
---|
118 | * - n a number |
---|
119 | * - pp (resp. qq, mm, nn) means arg is constant |
---|
120 | * - p (resp, q, m, n) means arg is destroyed |
---|
121 | * |
---|
122 | ***************************************************************/ |
---|
123 | |
---|
124 | poly p_Sub(poly a, poly b, const ring r); |
---|
125 | |
---|
126 | poly p_Power(poly p, int i, const ring r); |
---|
127 | |
---|
128 | |
---|
129 | /*************************************************************** |
---|
130 | * |
---|
131 | * PDEBUG stuff |
---|
132 | * |
---|
133 | ***************************************************************/ |
---|
134 | #ifdef PDEBUG |
---|
135 | // Returns TRUE if m is monom of p, FALSE otherwise |
---|
136 | BOOLEAN pIsMonomOf(poly p, poly m); |
---|
137 | // Returns TRUE if p and q have common monoms |
---|
138 | BOOLEAN pHaveCommonMonoms(poly p, poly q); |
---|
139 | |
---|
140 | // p_Check* routines return TRUE if everything is ok, |
---|
141 | // else, they report error message and return false |
---|
142 | |
---|
143 | // check if Lm(p) is from ring r |
---|
144 | BOOLEAN p_LmCheckIsFromRing(poly p, ring r); |
---|
145 | // check if Lm(p) != NULL, r != NULL and initialized && Lm(p) is from r |
---|
146 | BOOLEAN p_LmCheckPolyRing(poly p, ring r); |
---|
147 | // check if all monoms of p are from ring r |
---|
148 | BOOLEAN p_CheckIsFromRing(poly p, ring r); |
---|
149 | // check r != NULL and initialized && all monoms of p are from r |
---|
150 | BOOLEAN p_CheckPolyRing(poly p, ring r); |
---|
151 | // check if r != NULL and initialized |
---|
152 | BOOLEAN p_CheckRing(ring r); |
---|
153 | // only do check if cond |
---|
154 | |
---|
155 | |
---|
156 | #define pIfThen(cond, check) do {if (cond) {check;}} while (0) |
---|
157 | |
---|
158 | BOOLEAN _p_Test(poly p, ring r, int level); |
---|
159 | BOOLEAN _p_LmTest(poly p, ring r, int level); |
---|
160 | BOOLEAN _pp_Test(poly p, ring lmRing, ring tailRing, int level); |
---|
161 | |
---|
162 | #define p_Test(p,r) _p_Test(p, r, PDEBUG) |
---|
163 | #define p_LmTest(p,r) _p_LmTest(p, r, PDEBUG) |
---|
164 | #define pp_Test(p, lmRing, tailRing) _pp_Test(p, lmRing, tailRing, PDEBUG) |
---|
165 | |
---|
166 | #else // ! PDEBUG |
---|
167 | |
---|
168 | #define pIsMonomOf(p, q) (TRUE) |
---|
169 | #define pHaveCommonMonoms(p, q) (TRUE) |
---|
170 | #define p_LmCheckIsFromRing(p,r) (TRUE) |
---|
171 | #define p_LmCheckPolyRing(p,r) (TRUE) |
---|
172 | #define p_CheckIsFromRing(p,r) (TRUE) |
---|
173 | #define p_CheckPolyRing(p,r) (TRUE) |
---|
174 | #define p_CheckRing(r) (TRUE) |
---|
175 | #define P_CheckIf(cond, check) (TRUE) |
---|
176 | |
---|
177 | #define p_Test(p,r) (TRUE) |
---|
178 | #define p_LmTest(p,r) (TRUE) |
---|
179 | #define pp_Test(p, lmRing, tailRing) (TRUE) |
---|
180 | |
---|
181 | #endif |
---|
182 | |
---|
183 | /*************************************************************** |
---|
184 | * |
---|
185 | * Misc stuff |
---|
186 | * |
---|
187 | ***************************************************************/ |
---|
188 | /*2 |
---|
189 | * returns the length of a polynomial (numbers of monomials) |
---|
190 | */ |
---|
191 | static inline unsigned pLength(poly a) |
---|
192 | { |
---|
193 | unsigned l = 0; |
---|
194 | while (a!=NULL) |
---|
195 | { |
---|
196 | pIter(a); |
---|
197 | l++; |
---|
198 | } |
---|
199 | return l; |
---|
200 | } |
---|
201 | |
---|
202 | // returns the length of a polynomial (numbers of monomials) and the last mon. |
---|
203 | // respect syzComp |
---|
204 | poly p_Last(const poly a, int &l, const ring r); |
---|
205 | |
---|
206 | /*----------------------------------------------------*/ |
---|
207 | |
---|
208 | void p_Norm(poly p1, const ring r); |
---|
209 | void p_Normalize(poly p,const ring r); |
---|
210 | void p_ProjectiveUnique(poly p,const ring r); |
---|
211 | |
---|
212 | void p_ContentForGB(poly p, const ring r); |
---|
213 | void p_Content(poly p, const ring r); |
---|
214 | void p_Content_n(poly p, number &c,const ring r); |
---|
215 | #if 1 |
---|
216 | // currently only used by Singular/janet |
---|
217 | void p_SimpleContent(poly p, int s, const ring r); |
---|
218 | number p_InitContent(poly ph, const ring r); |
---|
219 | #endif |
---|
220 | |
---|
221 | poly p_Cleardenom(poly p, const ring r); |
---|
222 | void p_Cleardenom_n(poly p, const ring r,number &c); |
---|
223 | //number p_GetAllDenom(poly ph, const ring r);// unused |
---|
224 | |
---|
225 | int p_Size( poly p, const ring r ); |
---|
226 | |
---|
227 | // homogenizes p by multiplying certain powers of the varnum-th variable |
---|
228 | poly p_Homogen (poly p, int varnum, const ring r); |
---|
229 | |
---|
230 | BOOLEAN p_IsHomogeneous (poly p, const ring r); |
---|
231 | BOOLEAN p_IsHomogeneousW (poly p, const intvec *w, const ring r); |
---|
232 | BOOLEAN p_IsHomogeneousW (poly p, const intvec *w, const intvec *module_w,const ring r); |
---|
233 | |
---|
234 | // Setm |
---|
235 | static inline void p_Setm(poly p, const ring r) |
---|
236 | { |
---|
237 | p_CheckRing2(r); |
---|
238 | r->p_Setm(p, r); |
---|
239 | } |
---|
240 | |
---|
241 | p_SetmProc p_GetSetmProc(const ring r); |
---|
242 | |
---|
243 | poly p_Subst(poly p, int n, poly e, const ring r); |
---|
244 | |
---|
245 | // TODO: |
---|
246 | #define p_SetmComp p_Setm |
---|
247 | |
---|
248 | // component |
---|
249 | static inline unsigned long p_SetComp(poly p, unsigned long c, ring r) |
---|
250 | { |
---|
251 | p_LmCheckPolyRing2(p, r); |
---|
252 | if (r->pCompIndex>=0) __p_GetComp(p,r) = c; |
---|
253 | return c; |
---|
254 | } |
---|
255 | // sets component of poly a to i |
---|
256 | static inline void p_SetCompP(poly p, int i, ring r) |
---|
257 | { |
---|
258 | if (p != NULL) |
---|
259 | { |
---|
260 | p_Test(p, r); |
---|
261 | if (rOrd_SetCompRequiresSetm(r)) |
---|
262 | { |
---|
263 | do |
---|
264 | { |
---|
265 | p_SetComp(p, i, r); |
---|
266 | p_SetmComp(p, r); |
---|
267 | pIter(p); |
---|
268 | } |
---|
269 | while (p != NULL); |
---|
270 | } |
---|
271 | else |
---|
272 | { |
---|
273 | do |
---|
274 | { |
---|
275 | p_SetComp(p, i, r); |
---|
276 | pIter(p); |
---|
277 | } |
---|
278 | while(p != NULL); |
---|
279 | } |
---|
280 | } |
---|
281 | } |
---|
282 | |
---|
283 | static inline void p_SetCompP(poly p, int i, ring lmRing, ring tailRing) |
---|
284 | { |
---|
285 | if (p != NULL) |
---|
286 | { |
---|
287 | p_SetComp(p, i, lmRing); |
---|
288 | p_SetmComp(p, lmRing); |
---|
289 | p_SetCompP(pNext(p), i, tailRing); |
---|
290 | } |
---|
291 | } |
---|
292 | |
---|
293 | // returns maximal column number in the modul element a (or 0) |
---|
294 | static inline long p_MaxComp(poly p, ring lmRing, ring tailRing) |
---|
295 | { |
---|
296 | long result,i; |
---|
297 | |
---|
298 | if(p==NULL) return 0; |
---|
299 | result = p_GetComp(p, lmRing); |
---|
300 | if (result != 0) |
---|
301 | { |
---|
302 | loop |
---|
303 | { |
---|
304 | pIter(p); |
---|
305 | if(p==NULL) break; |
---|
306 | i = p_GetComp(p, tailRing); |
---|
307 | if (i>result) result = i; |
---|
308 | } |
---|
309 | } |
---|
310 | return result; |
---|
311 | } |
---|
312 | |
---|
313 | static inline long p_MaxComp(poly p,ring lmRing) {return p_MaxComp(p,lmRing,lmRing);} |
---|
314 | |
---|
315 | static inline long p_MinComp(poly p, ring lmRing, ring tailRing) |
---|
316 | { |
---|
317 | long result,i; |
---|
318 | |
---|
319 | if(p==NULL) return 0; |
---|
320 | result = p_GetComp(p,lmRing); |
---|
321 | if (result != 0) |
---|
322 | { |
---|
323 | loop |
---|
324 | { |
---|
325 | pIter(p); |
---|
326 | if(p==NULL) break; |
---|
327 | i = p_GetComp(p,tailRing); |
---|
328 | if (i<result) result = i; |
---|
329 | } |
---|
330 | } |
---|
331 | return result; |
---|
332 | } |
---|
333 | |
---|
334 | static inline long p_MinComp(poly p,ring lmRing) {return p_MinComp(p,lmRing,lmRing);} |
---|
335 | |
---|
336 | |
---|
337 | static inline poly pReverse(poly p) |
---|
338 | { |
---|
339 | if (p == NULL || pNext(p) == NULL) return p; |
---|
340 | |
---|
341 | poly q = pNext(p), // == pNext(p) |
---|
342 | qn; |
---|
343 | pNext(p) = NULL; |
---|
344 | do |
---|
345 | { |
---|
346 | qn = pNext(q); |
---|
347 | pNext(q) = p; |
---|
348 | p = q; |
---|
349 | q = qn; |
---|
350 | } |
---|
351 | while (qn != NULL); |
---|
352 | return p; |
---|
353 | } |
---|
354 | void pEnlargeSet(poly**p, int length, int increment); |
---|
355 | |
---|
356 | |
---|
357 | /*************************************************************** |
---|
358 | * |
---|
359 | * I/O |
---|
360 | * |
---|
361 | ***************************************************************/ |
---|
362 | /// print p according to ShortOut in lmRing & tailRing |
---|
363 | void p_String0(poly p, ring lmRing, ring tailRing); |
---|
364 | char* p_String(poly p, ring lmRing, ring tailRing); |
---|
365 | void p_Write(poly p, ring lmRing, ring tailRing); |
---|
366 | void p_Write0(poly p, ring lmRing, ring tailRing); |
---|
367 | void p_wrp(poly p, ring lmRing, ring tailRing); |
---|
368 | |
---|
369 | /// print p in a short way, if possible |
---|
370 | void p_String0Short(const poly p, ring lmRing, ring tailRing); |
---|
371 | |
---|
372 | /// print p in a long way |
---|
373 | void p_String0Long(const poly p, ring lmRing, ring tailRing); |
---|
374 | |
---|
375 | |
---|
376 | /*************************************************************** |
---|
377 | * |
---|
378 | * Degree stuff -- see p_polys.cc for explainations |
---|
379 | * |
---|
380 | ***************************************************************/ |
---|
381 | |
---|
382 | static inline long p_FDeg(const poly p, const ring r) { return r->pFDeg(p,r); } |
---|
383 | static inline long p_LDeg(const poly p, int *l, const ring r) { return r->pLDeg(p,l,r); } |
---|
384 | |
---|
385 | long p_WFirstTotalDegree(poly p, ring r); |
---|
386 | long p_WTotaldegree(poly p, const ring r); |
---|
387 | long p_WDegree(poly p,const ring r); |
---|
388 | long pLDeg0(poly p,int *l, ring r); |
---|
389 | long pLDeg0c(poly p,int *l, ring r); |
---|
390 | long pLDegb(poly p,int *l, ring r); |
---|
391 | long pLDeg1(poly p,int *l, ring r); |
---|
392 | long pLDeg1c(poly p,int *l, ring r); |
---|
393 | long pLDeg1_Deg(poly p,int *l, ring r); |
---|
394 | long pLDeg1c_Deg(poly p,int *l, ring r); |
---|
395 | long pLDeg1_Totaldegree(poly p,int *l, ring r); |
---|
396 | long pLDeg1c_Totaldegree(poly p,int *l, ring r); |
---|
397 | long pLDeg1_WFirstTotalDegree(poly p,int *l, ring r); |
---|
398 | long pLDeg1c_WFirstTotalDegree(poly p,int *l, ring r); |
---|
399 | |
---|
400 | BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r); |
---|
401 | |
---|
402 | /// same as the usual p_EqualPolys for polys belonging to *equal* rings |
---|
403 | BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r1, const ring r2); |
---|
404 | |
---|
405 | long p_Deg(poly a, const ring r); |
---|
406 | |
---|
407 | |
---|
408 | /*************************************************************** |
---|
409 | * |
---|
410 | * Primitives for accessing and setting fields of a poly |
---|
411 | * |
---|
412 | ***************************************************************/ |
---|
413 | |
---|
414 | static inline number p_SetCoeff(poly p, number n, ring r) |
---|
415 | { |
---|
416 | p_LmCheckPolyRing2(p, r); |
---|
417 | n_Delete(&(p->coef), r->cf); |
---|
418 | (p)->coef=n; |
---|
419 | return n; |
---|
420 | } |
---|
421 | |
---|
422 | // order |
---|
423 | static inline long p_GetOrder(poly p, ring r) |
---|
424 | { |
---|
425 | p_LmCheckPolyRing2(p, r); |
---|
426 | if (r->typ==NULL) return ((p)->exp[r->pOrdIndex]); |
---|
427 | int i=0; |
---|
428 | loop |
---|
429 | { |
---|
430 | switch(r->typ[i].ord_typ) |
---|
431 | { |
---|
432 | case ro_am: |
---|
433 | case ro_wp_neg: |
---|
434 | return ((p->exp[r->pOrdIndex])-POLY_NEGWEIGHT_OFFSET); |
---|
435 | case ro_syzcomp: |
---|
436 | case ro_syz: |
---|
437 | case ro_cp: |
---|
438 | i++; |
---|
439 | break; |
---|
440 | //case ro_dp: |
---|
441 | //case ro_wp: |
---|
442 | default: |
---|
443 | return ((p)->exp[r->pOrdIndex]); |
---|
444 | } |
---|
445 | } |
---|
446 | } |
---|
447 | |
---|
448 | |
---|
449 | static inline unsigned long p_AddComp(poly p, unsigned long v, ring r) |
---|
450 | { |
---|
451 | p_LmCheckPolyRing2(p, r); |
---|
452 | pAssume2(rRing_has_Comp(r)); |
---|
453 | return __p_GetComp(p,r) += v; |
---|
454 | } |
---|
455 | static inline unsigned long p_SubComp(poly p, unsigned long v, ring r) |
---|
456 | { |
---|
457 | p_LmCheckPolyRing2(p, r); |
---|
458 | pAssume2(rRing_has_Comp(r)); |
---|
459 | _pPolyAssume2(__p_GetComp(p,r) >= v,p,r); |
---|
460 | return __p_GetComp(p,r) -= v; |
---|
461 | } |
---|
462 | |
---|
463 | #ifndef HAVE_EXPSIZES |
---|
464 | |
---|
465 | /// get a single variable exponent |
---|
466 | /// @Note: |
---|
467 | /// the integer VarOffset encodes: |
---|
468 | /// 1. the position of a variable in the exponent vector p->exp (lower 24 bits) |
---|
469 | /// 2. number of bits to shift to the right in the upper 8 bits (which takes at most 6 bits for 64 bit) |
---|
470 | /// Thus VarOffset always has 2 zero higher bits! |
---|
471 | static inline long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset) |
---|
472 | { |
---|
473 | pAssume2((VarOffset >> (24 + 6)) == 0); |
---|
474 | #if 0 |
---|
475 | int pos=(VarOffset & 0xffffff); |
---|
476 | int bitpos=(VarOffset >> 24); |
---|
477 | unsigned long exp=(p->exp[pos] >> bitmask) & iBitmask; |
---|
478 | return exp; |
---|
479 | #else |
---|
480 | return (long) |
---|
481 | ((p->exp[(VarOffset & 0xffffff)] >> (VarOffset >> 24)) |
---|
482 | & iBitmask); |
---|
483 | #endif |
---|
484 | } |
---|
485 | |
---|
486 | |
---|
487 | /// set a single variable exponent |
---|
488 | /// @Note: |
---|
489 | /// VarOffset encodes the position in p->exp @see p_GetExp |
---|
490 | static inline unsigned long p_SetExp(poly p, const unsigned long e, const unsigned long iBitmask, const int VarOffset) |
---|
491 | { |
---|
492 | pAssume2(e>=0); |
---|
493 | pAssume2(e<=iBitmask); |
---|
494 | pAssume2((VarOffset >> (24 + 6)) == 0); |
---|
495 | |
---|
496 | // shift e to the left: |
---|
497 | REGISTER int shift = VarOffset >> 24; |
---|
498 | unsigned long ee = e << shift /*(VarOffset >> 24)*/; |
---|
499 | // find the bits in the exponent vector |
---|
500 | REGISTER int offset = (VarOffset & 0xffffff); |
---|
501 | // clear the bits in the exponent vector: |
---|
502 | p->exp[offset] &= ~( iBitmask << shift ); |
---|
503 | // insert e with | |
---|
504 | p->exp[ offset ] |= ee; |
---|
505 | return e; |
---|
506 | } |
---|
507 | |
---|
508 | |
---|
509 | #else // #ifdef HAVE_EXPSIZES // EXPERIMENTAL!!! |
---|
510 | |
---|
511 | static inline unsigned long BitMask(unsigned long bitmask, int twobits) |
---|
512 | { |
---|
513 | // bitmask = 00000111111111111 |
---|
514 | // 0 must give bitmask! |
---|
515 | // 1, 2, 3 - anything like 00011..11 |
---|
516 | pAssume2((twobits >> 2) == 0); |
---|
517 | static const unsigned long _bitmasks[4] = {-1, 0x7fff, 0x7f, 0x3}; |
---|
518 | return bitmask & _bitmasks[twobits]; |
---|
519 | } |
---|
520 | |
---|
521 | |
---|
522 | /// @Note: we may add some more info (6 ) into VarOffset and thus encode |
---|
523 | static inline long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset) |
---|
524 | { |
---|
525 | int pos =(VarOffset & 0xffffff); |
---|
526 | int hbyte= (VarOffset >> 24); // the highest byte |
---|
527 | int bitpos = hbyte & 0x3f; // last 6 bits |
---|
528 | long bitmask = BitMask(iBitmask, hbyte >> 6); |
---|
529 | |
---|
530 | long exp=(p->exp[pos] >> bitpos) & bitmask; |
---|
531 | return exp; |
---|
532 | |
---|
533 | } |
---|
534 | |
---|
535 | static inline long p_SetExp(poly p, const long e, const unsigned long iBitmask, const int VarOffset) |
---|
536 | { |
---|
537 | pAssume2(e>=0); |
---|
538 | pAssume2(e <= BitMask(iBitmask, VarOffset >> 30)); |
---|
539 | |
---|
540 | // shift e to the left: |
---|
541 | REGISTER int hbyte = VarOffset >> 24; |
---|
542 | int bitmask = BitMask(iBitmask, hbyte >> 6); |
---|
543 | REGISTER int shift = hbyte & 0x3f; |
---|
544 | long ee = e << shift; |
---|
545 | // find the bits in the exponent vector |
---|
546 | REGISTER int offset = (VarOffset & 0xffffff); |
---|
547 | // clear the bits in the exponent vector: |
---|
548 | p->exp[offset] &= ~( bitmask << shift ); |
---|
549 | // insert e with | |
---|
550 | p->exp[ offset ] |= ee; |
---|
551 | return e; |
---|
552 | } |
---|
553 | |
---|
554 | #endif // #ifndef HAVE_EXPSIZES |
---|
555 | |
---|
556 | |
---|
557 | static inline long p_GetExp(const poly p, const ring r, const int VarOffset) |
---|
558 | { |
---|
559 | p_LmCheckPolyRing2(p, r); |
---|
560 | pAssume2(VarOffset != -1); |
---|
561 | return p_GetExp(p, r->bitmask, VarOffset); |
---|
562 | } |
---|
563 | |
---|
564 | static inline long p_SetExp(poly p, const long e, const ring r, const int VarOffset) |
---|
565 | { |
---|
566 | p_LmCheckPolyRing2(p, r); |
---|
567 | pAssume2(VarOffset != -1); |
---|
568 | return p_SetExp(p, e, r->bitmask, VarOffset); |
---|
569 | } |
---|
570 | |
---|
571 | |
---|
572 | |
---|
573 | /// get v^th exponent for a monomial |
---|
574 | static inline long p_GetExp(const poly p, const int v, const ring r) |
---|
575 | { |
---|
576 | p_LmCheckPolyRing2(p, r); |
---|
577 | pAssume2(v>0 && v <= r->N); |
---|
578 | pAssume2(r->VarOffset[v] != -1); |
---|
579 | return p_GetExp(p, r->bitmask, r->VarOffset[v]); |
---|
580 | } |
---|
581 | |
---|
582 | |
---|
583 | /// set v^th exponent for a monomial |
---|
584 | static inline long p_SetExp(poly p, const int v, const long e, const ring r) |
---|
585 | { |
---|
586 | p_LmCheckPolyRing2(p, r); |
---|
587 | pAssume2(v>0 && v <= r->N); |
---|
588 | pAssume2(r->VarOffset[v] != -1); |
---|
589 | return p_SetExp(p, e, r->bitmask, r->VarOffset[v]); |
---|
590 | } |
---|
591 | |
---|
592 | // the following should be implemented more efficiently |
---|
593 | static inline long p_IncrExp(poly p, int v, ring r) |
---|
594 | { |
---|
595 | p_LmCheckPolyRing2(p, r); |
---|
596 | int e = p_GetExp(p,v,r); |
---|
597 | e++; |
---|
598 | return p_SetExp(p,v,e,r); |
---|
599 | } |
---|
600 | static inline long p_DecrExp(poly p, int v, ring r) |
---|
601 | { |
---|
602 | p_LmCheckPolyRing2(p, r); |
---|
603 | int e = p_GetExp(p,v,r); |
---|
604 | pAssume2(e > 0); |
---|
605 | e--; |
---|
606 | return p_SetExp(p,v,e,r); |
---|
607 | } |
---|
608 | static inline long p_AddExp(poly p, int v, long ee, ring r) |
---|
609 | { |
---|
610 | p_LmCheckPolyRing2(p, r); |
---|
611 | int e = p_GetExp(p,v,r); |
---|
612 | e += ee; |
---|
613 | return p_SetExp(p,v,e,r); |
---|
614 | } |
---|
615 | static inline long p_SubExp(poly p, int v, long ee, ring r) |
---|
616 | { |
---|
617 | p_LmCheckPolyRing2(p, r); |
---|
618 | long e = p_GetExp(p,v,r); |
---|
619 | pAssume2(e >= ee); |
---|
620 | e -= ee; |
---|
621 | return p_SetExp(p,v,e,r); |
---|
622 | } |
---|
623 | static inline long p_MultExp(poly p, int v, long ee, ring r) |
---|
624 | { |
---|
625 | p_LmCheckPolyRing2(p, r); |
---|
626 | long e = p_GetExp(p,v,r); |
---|
627 | e *= ee; |
---|
628 | return p_SetExp(p,v,e,r); |
---|
629 | } |
---|
630 | |
---|
631 | static inline long p_GetExpSum(poly p1, poly p2, int i, ring r) |
---|
632 | { |
---|
633 | p_LmCheckPolyRing2(p1, r); |
---|
634 | p_LmCheckPolyRing2(p2, r); |
---|
635 | return p_GetExp(p1,i,r) + p_GetExp(p2,i,r); |
---|
636 | } |
---|
637 | static inline long p_GetExpDiff(poly p1, poly p2, int i, ring r) |
---|
638 | { |
---|
639 | return p_GetExp(p1,i,r) - p_GetExp(p2,i,r); |
---|
640 | } |
---|
641 | |
---|
642 | static inline int p_Comp_k_n(poly a, poly b, int k, ring r) |
---|
643 | { |
---|
644 | if ((a==NULL) || (b==NULL) ) return FALSE; |
---|
645 | p_LmCheckPolyRing2(a, r); |
---|
646 | p_LmCheckPolyRing2(b, r); |
---|
647 | pAssume2(k > 0 && k <= r->N); |
---|
648 | int i=k; |
---|
649 | for(;i<=r->N;i++) |
---|
650 | { |
---|
651 | if (p_GetExp(a,i,r) != p_GetExp(b,i,r)) return FALSE; |
---|
652 | // if (a->exp[(r->VarOffset[i] & 0xffffff)] != b->exp[(r->VarOffset[i] & 0xffffff)]) return FALSE; |
---|
653 | } |
---|
654 | return TRUE; |
---|
655 | } |
---|
656 | |
---|
657 | |
---|
658 | /*************************************************************** |
---|
659 | * |
---|
660 | * Allocation/Initalization/Deletion |
---|
661 | * |
---|
662 | ***************************************************************/ |
---|
663 | #if (OM_TRACK > 2) && defined(OM_TRACK_CUSTOM) |
---|
664 | static inline poly p_New(const ring r, omBin bin) |
---|
665 | #else |
---|
666 | static inline poly p_New(const ring /*r*/, omBin bin) |
---|
667 | #endif |
---|
668 | { |
---|
669 | p_CheckRing2(r); |
---|
670 | pAssume2(bin != NULL && omSizeWOfBin(r->PolyBin) == omSizeWOfBin(bin)); |
---|
671 | poly p; |
---|
672 | omTypeAllocBin(poly, p, bin); |
---|
673 | p_SetRingOfLm(p, r); |
---|
674 | return p; |
---|
675 | } |
---|
676 | |
---|
677 | static inline poly p_New(ring r) |
---|
678 | { |
---|
679 | return p_New(r, r->PolyBin); |
---|
680 | } |
---|
681 | |
---|
682 | #if (PDEBUG > 2) || defined(XALLOC_BIN) |
---|
683 | static inline void p_LmFree(poly p, ring r) |
---|
684 | #else |
---|
685 | static inline void p_LmFree(poly p, ring) |
---|
686 | #endif |
---|
687 | { |
---|
688 | p_LmCheckPolyRing2(p, r); |
---|
689 | #ifdef XALLOC_BIN |
---|
690 | omFreeBin(p,r->PolyBin); |
---|
691 | #else |
---|
692 | omFreeBinAddr(p); |
---|
693 | #endif |
---|
694 | } |
---|
695 | #if (PDEBUG > 2) || defined(XALLOC_BIN) |
---|
696 | static inline void p_LmFree(poly *p, ring r) |
---|
697 | #else |
---|
698 | static inline void p_LmFree(poly *p, ring) |
---|
699 | #endif |
---|
700 | { |
---|
701 | p_LmCheckPolyRing2(*p, r); |
---|
702 | poly h = *p; |
---|
703 | *p = pNext(h); |
---|
704 | #ifdef XALLOC_BIN |
---|
705 | omFreeBin(h,r->PolyBin); |
---|
706 | #else |
---|
707 | omFreeBinAddr(h); |
---|
708 | #endif |
---|
709 | } |
---|
710 | #if (PDEBUG > 2) || defined(XALLOC_BIN) |
---|
711 | static inline poly p_LmFreeAndNext(poly p, ring r) |
---|
712 | #else |
---|
713 | static inline poly p_LmFreeAndNext(poly p, ring) |
---|
714 | #endif |
---|
715 | { |
---|
716 | p_LmCheckPolyRing2(p, r); |
---|
717 | poly pnext = pNext(p); |
---|
718 | #ifdef XALLOC_BIN |
---|
719 | omFreeBin(p,r->PolyBin); |
---|
720 | #else |
---|
721 | omFreeBinAddr(p); |
---|
722 | #endif |
---|
723 | return pnext; |
---|
724 | } |
---|
725 | static inline void p_LmDelete(poly p, const ring r) |
---|
726 | { |
---|
727 | p_LmCheckPolyRing2(p, r); |
---|
728 | n_Delete(&pGetCoeff(p), r->cf); |
---|
729 | #ifdef XALLOC_BIN |
---|
730 | omFreeBin(p,r->PolyBin); |
---|
731 | #else |
---|
732 | omFreeBinAddr(p); |
---|
733 | #endif |
---|
734 | } |
---|
735 | static inline void p_LmDelete0(poly p, const ring r) |
---|
736 | { |
---|
737 | p_LmCheckPolyRing2(p, r); |
---|
738 | if (pGetCoeff(p)!=NULL) n_Delete(&pGetCoeff(p), r->cf); |
---|
739 | #ifdef XALLOC_BIN |
---|
740 | omFreeBin(p,r->PolyBin); |
---|
741 | #else |
---|
742 | omFreeBinAddr(p); |
---|
743 | #endif |
---|
744 | } |
---|
745 | static inline void p_LmDelete(poly *p, const ring r) |
---|
746 | { |
---|
747 | p_LmCheckPolyRing2(*p, r); |
---|
748 | poly h = *p; |
---|
749 | *p = pNext(h); |
---|
750 | n_Delete(&pGetCoeff(h), r->cf); |
---|
751 | #ifdef XALLOC_BIN |
---|
752 | omFreeBin(h,r->PolyBin); |
---|
753 | #else |
---|
754 | omFreeBinAddr(h); |
---|
755 | #endif |
---|
756 | } |
---|
757 | static inline poly p_LmDeleteAndNext(poly p, const ring r) |
---|
758 | { |
---|
759 | p_LmCheckPolyRing2(p, r); |
---|
760 | poly pnext = pNext(p); |
---|
761 | n_Delete(&pGetCoeff(p), r->cf); |
---|
762 | #ifdef XALLOC_BIN |
---|
763 | omFreeBin(p,r->PolyBin); |
---|
764 | #else |
---|
765 | omFreeBinAddr(p); |
---|
766 | #endif |
---|
767 | return pnext; |
---|
768 | } |
---|
769 | |
---|
770 | /*************************************************************** |
---|
771 | * |
---|
772 | * Misc routines |
---|
773 | * |
---|
774 | ***************************************************************/ |
---|
775 | |
---|
776 | /// return the maximal exponent of p in form of the maximal long var |
---|
777 | unsigned long p_GetMaxExpL(poly p, const ring r, unsigned long l_max = 0); |
---|
778 | |
---|
779 | /// return monomial r such that GetExp(r,i) is maximum of all |
---|
780 | /// monomials in p; coeff == 0, next == NULL, ord is not set |
---|
781 | poly p_GetMaxExpP(poly p, ring r); |
---|
782 | |
---|
783 | static inline unsigned long p_GetMaxExp(const unsigned long l, const ring r) |
---|
784 | { |
---|
785 | unsigned long bitmask = r->bitmask; |
---|
786 | unsigned long max = (l & bitmask); |
---|
787 | unsigned long j = r->ExpPerLong - 1; |
---|
788 | |
---|
789 | if (j > 0) |
---|
790 | { |
---|
791 | unsigned long i = r->BitsPerExp; |
---|
792 | long e; |
---|
793 | loop |
---|
794 | { |
---|
795 | e = ((l >> i) & bitmask); |
---|
796 | if ((unsigned long) e > max) |
---|
797 | max = e; |
---|
798 | j--; |
---|
799 | if (j==0) break; |
---|
800 | i += r->BitsPerExp; |
---|
801 | } |
---|
802 | } |
---|
803 | return max; |
---|
804 | } |
---|
805 | |
---|
806 | static inline unsigned long p_GetMaxExp(const poly p, const ring r) |
---|
807 | { |
---|
808 | return p_GetMaxExp(p_GetMaxExpL(p, r), r); |
---|
809 | } |
---|
810 | |
---|
811 | static inline unsigned long |
---|
812 | p_GetTotalDegree(const unsigned long l, const ring r, const int number_of_exps) |
---|
813 | { |
---|
814 | const unsigned long bitmask = r->bitmask; |
---|
815 | unsigned long sum = (l & bitmask); |
---|
816 | unsigned long j = number_of_exps - 1; |
---|
817 | |
---|
818 | if (j > 0) |
---|
819 | { |
---|
820 | unsigned long i = r->BitsPerExp; |
---|
821 | loop |
---|
822 | { |
---|
823 | sum += ((l >> i) & bitmask); |
---|
824 | j--; |
---|
825 | if (j==0) break; |
---|
826 | i += r->BitsPerExp; |
---|
827 | } |
---|
828 | } |
---|
829 | return sum; |
---|
830 | } |
---|
831 | |
---|
832 | /*************************************************************** |
---|
833 | * |
---|
834 | * Dispatcher to r->p_Procs, they do the tests/checks |
---|
835 | * |
---|
836 | ***************************************************************/ |
---|
837 | /// returns a copy of p (without any additional testing) |
---|
838 | static inline poly p_Copy_noCheck(poly p, const ring r) |
---|
839 | { |
---|
840 | /*assume(p!=NULL);*/ |
---|
841 | assume(r != NULL); |
---|
842 | assume(r->p_Procs != NULL); |
---|
843 | assume(r->p_Procs->p_Copy != NULL); |
---|
844 | return r->p_Procs->p_Copy(p, r); |
---|
845 | } |
---|
846 | |
---|
847 | /// returns a copy of p |
---|
848 | static inline poly p_Copy(poly p, const ring r) |
---|
849 | { |
---|
850 | if (p!=NULL) |
---|
851 | { |
---|
852 | p_Test(p,r); |
---|
853 | const poly pp = p_Copy_noCheck(p, r); |
---|
854 | p_Test(pp,r); |
---|
855 | return pp; |
---|
856 | } |
---|
857 | else |
---|
858 | return NULL; |
---|
859 | } |
---|
860 | |
---|
861 | /// copy the (leading) term of p |
---|
862 | static inline poly p_Head(const poly p, const ring r) |
---|
863 | { |
---|
864 | if (p == NULL) return NULL; |
---|
865 | p_LmCheckPolyRing1(p, r); |
---|
866 | poly np; |
---|
867 | omTypeAllocBin(poly, np, r->PolyBin); |
---|
868 | p_SetRingOfLm(np, r); |
---|
869 | memcpy(np->exp, p->exp, r->ExpL_Size*sizeof(long)); |
---|
870 | pNext(np) = NULL; |
---|
871 | pSetCoeff0(np, n_Copy(pGetCoeff(p), r->cf)); |
---|
872 | return np; |
---|
873 | } |
---|
874 | |
---|
875 | /// like p_Head, but allow NULL coeff |
---|
876 | poly p_Head0(const poly p, const ring r); |
---|
877 | |
---|
878 | /// like p_Head, but with coefficient 1 |
---|
879 | poly p_CopyPowerProduct(const poly p, const ring r); |
---|
880 | |
---|
881 | /// like p_Head, but with coefficient n |
---|
882 | poly p_CopyPowerProduct0(const poly p, const number n, const ring r); |
---|
883 | |
---|
884 | /// returns a copy of p with Lm(p) from lmRing and Tail(p) from tailRing |
---|
885 | static inline poly p_Copy(poly p, const ring lmRing, const ring tailRing) |
---|
886 | { |
---|
887 | if (p != NULL) |
---|
888 | { |
---|
889 | #ifndef PDEBUG |
---|
890 | if (tailRing == lmRing) |
---|
891 | return p_Copy_noCheck(p, tailRing); |
---|
892 | #endif |
---|
893 | poly pres = p_Head(p, lmRing); |
---|
894 | if (pNext(p)!=NULL) |
---|
895 | pNext(pres) = p_Copy_noCheck(pNext(p), tailRing); |
---|
896 | return pres; |
---|
897 | } |
---|
898 | else |
---|
899 | return NULL; |
---|
900 | } |
---|
901 | |
---|
902 | // deletes *p, and sets *p to NULL |
---|
903 | static inline void p_Delete(poly *p, const ring r) |
---|
904 | { |
---|
905 | assume( p!= NULL ); |
---|
906 | assume( r!= NULL ); |
---|
907 | if ((*p)!=NULL) r->p_Procs->p_Delete(p, r); |
---|
908 | } |
---|
909 | |
---|
910 | static inline void p_Delete(poly *p, const ring lmRing, const ring tailRing) |
---|
911 | { |
---|
912 | assume( p!= NULL ); |
---|
913 | if (*p != NULL) |
---|
914 | { |
---|
915 | #ifndef PDEBUG |
---|
916 | if (tailRing == lmRing) |
---|
917 | { |
---|
918 | p_Delete(p, tailRing); |
---|
919 | return; |
---|
920 | } |
---|
921 | #endif |
---|
922 | if (pNext(*p) != NULL) |
---|
923 | p_Delete(&pNext(*p), tailRing); |
---|
924 | p_LmDelete(p, lmRing); |
---|
925 | } |
---|
926 | } |
---|
927 | |
---|
928 | // copys monomials of p, allocates new monomials from bin, |
---|
929 | // deletes monomials of p |
---|
930 | static inline poly p_ShallowCopyDelete(poly p, const ring r, omBin bin) |
---|
931 | { |
---|
932 | p_LmCheckPolyRing2(p, r); |
---|
933 | pAssume2(omSizeWOfBin(r->PolyBin) == omSizeWOfBin(bin)); |
---|
934 | return r->p_Procs->p_ShallowCopyDelete(p, r, bin); |
---|
935 | } |
---|
936 | |
---|
937 | // returns p+q, destroys p and q |
---|
938 | static inline poly p_Add_q(poly p, poly q, const ring r) |
---|
939 | { |
---|
940 | assume( (p != q) || (p == NULL && q == NULL) ); |
---|
941 | if (q==NULL) return p; |
---|
942 | if (p==NULL) return q; |
---|
943 | int shorter; |
---|
944 | return r->p_Procs->p_Add_q(p, q, shorter, r); |
---|
945 | } |
---|
946 | |
---|
947 | /// like p_Add_q, except that if lp == pLength(lp) lq == pLength(lq) then lp == pLength(p+q) |
---|
948 | static inline poly p_Add_q(poly p, poly q, int &lp, int lq, const ring r) |
---|
949 | { |
---|
950 | assume( (p != q) || (p == NULL && q == NULL) ); |
---|
951 | if (q==NULL) return p; |
---|
952 | if (p==NULL) { lp=lq; return q; } |
---|
953 | int shorter; |
---|
954 | poly res = r->p_Procs->p_Add_q(p, q, shorter, r); |
---|
955 | lp += lq - shorter; |
---|
956 | return res; |
---|
957 | } |
---|
958 | |
---|
959 | // returns p*n, destroys p |
---|
960 | static inline poly p_Mult_nn(poly p, number n, const ring r) |
---|
961 | { |
---|
962 | if (p==NULL) return NULL; |
---|
963 | if (n_IsOne(n, r->cf)) |
---|
964 | return p; |
---|
965 | else if (n_IsZero(n, r->cf)) |
---|
966 | { |
---|
967 | p_Delete(&p, r); // NOTE: without p_Delete - memory leak! |
---|
968 | return NULL; |
---|
969 | } |
---|
970 | else |
---|
971 | return r->p_Procs->p_Mult_nn(p, n, r); |
---|
972 | } |
---|
973 | #define __p_Mult_nn(p,n,r) r->p_Procs->p_Mult_nn(p, n, r) |
---|
974 | |
---|
975 | static inline poly p_Mult_nn(poly p, number n, const ring lmRing, |
---|
976 | const ring tailRing) |
---|
977 | { |
---|
978 | assume(p!=NULL); |
---|
979 | #ifndef PDEBUG |
---|
980 | if (lmRing == tailRing) |
---|
981 | return p_Mult_nn(p, n, tailRing); |
---|
982 | #endif |
---|
983 | poly pnext = pNext(p); |
---|
984 | pNext(p) = NULL; |
---|
985 | p = lmRing->p_Procs->p_Mult_nn(p, n, lmRing); |
---|
986 | if (pnext!=NULL) |
---|
987 | { |
---|
988 | pNext(p) = tailRing->p_Procs->p_Mult_nn(pnext, n, tailRing); |
---|
989 | } |
---|
990 | return p; |
---|
991 | } |
---|
992 | |
---|
993 | // returns p*n, does not destroy p |
---|
994 | static inline poly pp_Mult_nn(poly p, number n, const ring r) |
---|
995 | { |
---|
996 | if (p==NULL) return NULL; |
---|
997 | if (n_IsOne(n, r->cf)) |
---|
998 | return p_Copy(p, r); |
---|
999 | else if (n_IsZero(n, r->cf)) |
---|
1000 | return NULL; |
---|
1001 | else |
---|
1002 | return r->p_Procs->pp_Mult_nn(p, n, r); |
---|
1003 | } |
---|
1004 | #define __pp_Mult_nn(p,n,r) r->p_Procs->pp_Mult_nn(p, n, r) |
---|
1005 | |
---|
1006 | // test if the monomial is a constant as a vector component |
---|
1007 | // i.e., test if all exponents are zero |
---|
1008 | static inline BOOLEAN p_LmIsConstantComp(const poly p, const ring r) |
---|
1009 | { |
---|
1010 | //p_LmCheckPolyRing(p, r); |
---|
1011 | int i = r->VarL_Size - 1; |
---|
1012 | |
---|
1013 | do |
---|
1014 | { |
---|
1015 | if (p->exp[r->VarL_Offset[i]] != 0) |
---|
1016 | return FALSE; |
---|
1017 | i--; |
---|
1018 | } |
---|
1019 | while (i >= 0); |
---|
1020 | return TRUE; |
---|
1021 | } |
---|
1022 | |
---|
1023 | // test if monomial is a constant, i.e. if all exponents and the component |
---|
1024 | // is zero |
---|
1025 | static inline BOOLEAN p_LmIsConstant(const poly p, const ring r) |
---|
1026 | { |
---|
1027 | if (p_LmIsConstantComp(p, r)) |
---|
1028 | return (p_GetComp(p, r) == 0); |
---|
1029 | return FALSE; |
---|
1030 | } |
---|
1031 | |
---|
1032 | // returns Copy(p)*m, does neither destroy p nor m |
---|
1033 | static inline poly pp_Mult_mm(poly p, poly m, const ring r) |
---|
1034 | { |
---|
1035 | if (p==NULL) return NULL; |
---|
1036 | if (p_LmIsConstant(m, r)) |
---|
1037 | return __pp_Mult_nn(p, pGetCoeff(m), r); |
---|
1038 | else |
---|
1039 | return r->p_Procs->pp_Mult_mm(p, m, r); |
---|
1040 | } |
---|
1041 | |
---|
1042 | // returns m*Copy(p), does neither destroy p nor m |
---|
1043 | static inline poly pp_mm_Mult(poly p, poly m, const ring r) |
---|
1044 | { |
---|
1045 | if (p==NULL) return NULL; |
---|
1046 | if (p_LmIsConstant(m, r)) |
---|
1047 | return __pp_Mult_nn(p, pGetCoeff(m), r); |
---|
1048 | else |
---|
1049 | return r->p_Procs->pp_mm_Mult(p, m, r); |
---|
1050 | } |
---|
1051 | |
---|
1052 | // returns p*m, destroys p, const: m |
---|
1053 | static inline poly p_Mult_mm(poly p, poly m, const ring r) |
---|
1054 | { |
---|
1055 | if (p==NULL) return NULL; |
---|
1056 | if (p_LmIsConstant(m, r)) |
---|
1057 | return __p_Mult_nn(p, pGetCoeff(m), r); |
---|
1058 | else |
---|
1059 | return r->p_Procs->p_Mult_mm(p, m, r); |
---|
1060 | } |
---|
1061 | |
---|
1062 | // returns m*p, destroys p, const: m |
---|
1063 | static inline poly p_mm_Mult(poly p, poly m, const ring r) |
---|
1064 | { |
---|
1065 | if (p==NULL) return NULL; |
---|
1066 | if (p_LmIsConstant(m, r)) |
---|
1067 | return __p_Mult_nn(p, pGetCoeff(m), r); |
---|
1068 | else |
---|
1069 | return r->p_Procs->p_mm_Mult(p, m, r); |
---|
1070 | } |
---|
1071 | |
---|
1072 | static inline poly p_Minus_mm_Mult_qq(poly p, const poly m, const poly q, int &lp, int lq, |
---|
1073 | const poly spNoether, const ring r) |
---|
1074 | { |
---|
1075 | int shorter; |
---|
1076 | const poly res = r->p_Procs->p_Minus_mm_Mult_qq(p, m, q, shorter, spNoether, r); |
---|
1077 | lp += lq - shorter; |
---|
1078 | // assume( lp == pLength(res) ); |
---|
1079 | return res; |
---|
1080 | } |
---|
1081 | |
---|
1082 | // return p - m*Copy(q), destroys p; const: p,m |
---|
1083 | static inline poly p_Minus_mm_Mult_qq(poly p, const poly m, const poly q, const ring r) |
---|
1084 | { |
---|
1085 | int shorter; |
---|
1086 | |
---|
1087 | return r->p_Procs->p_Minus_mm_Mult_qq(p, m, q, shorter, NULL, r); |
---|
1088 | } |
---|
1089 | |
---|
1090 | |
---|
1091 | // returns p*Coeff(m) for such monomials pm of p, for which m is divisble by pm |
---|
1092 | static inline poly pp_Mult_Coeff_mm_DivSelect(poly p, const poly m, const ring r) |
---|
1093 | { |
---|
1094 | int shorter; |
---|
1095 | return r->p_Procs->pp_Mult_Coeff_mm_DivSelect(p, m, shorter, r); |
---|
1096 | } |
---|
1097 | |
---|
1098 | // returns p*Coeff(m) for such monomials pm of p, for which m is divisble by pm |
---|
1099 | // if lp is length of p on input then lp is length of returned poly on output |
---|
1100 | static inline poly pp_Mult_Coeff_mm_DivSelect(poly p, int &lp, const poly m, const ring r) |
---|
1101 | { |
---|
1102 | int shorter; |
---|
1103 | poly pp = r->p_Procs->pp_Mult_Coeff_mm_DivSelect(p, m, shorter, r); |
---|
1104 | lp -= shorter; |
---|
1105 | return pp; |
---|
1106 | } |
---|
1107 | |
---|
1108 | // returns -p, destroys p |
---|
1109 | static inline poly p_Neg(poly p, const ring r) |
---|
1110 | { |
---|
1111 | return r->p_Procs->p_Neg(p, r); |
---|
1112 | } |
---|
1113 | |
---|
1114 | extern poly _p_Mult_q(poly p, poly q, const int copy, const ring r); |
---|
1115 | // returns p*q, destroys p and q |
---|
1116 | static inline poly p_Mult_q(poly p, poly q, const ring r) |
---|
1117 | { |
---|
1118 | assume( (p != q) || (p == NULL && q == NULL) ); |
---|
1119 | |
---|
1120 | if (p == NULL) |
---|
1121 | { |
---|
1122 | p_Delete(&q, r); |
---|
1123 | return NULL; |
---|
1124 | } |
---|
1125 | if (q == NULL) |
---|
1126 | { |
---|
1127 | p_Delete(&p, r); |
---|
1128 | return NULL; |
---|
1129 | } |
---|
1130 | |
---|
1131 | if (pNext(p) == NULL) |
---|
1132 | { |
---|
1133 | q = r->p_Procs->p_mm_Mult(q, p, r); |
---|
1134 | p_LmDelete(&p, r); |
---|
1135 | return q; |
---|
1136 | } |
---|
1137 | |
---|
1138 | if (pNext(q) == NULL) |
---|
1139 | { |
---|
1140 | p = r->p_Procs->p_Mult_mm(p, q, r); |
---|
1141 | p_LmDelete(&q, r); |
---|
1142 | return p; |
---|
1143 | } |
---|
1144 | #if defined(HAVE_PLURAL) || defined(HAVE_SHIFTBBA) |
---|
1145 | if (rIsNCRing(r)) |
---|
1146 | return _nc_p_Mult_q(p, q, r); |
---|
1147 | else |
---|
1148 | #endif |
---|
1149 | return _p_Mult_q(p, q, 0, r); |
---|
1150 | } |
---|
1151 | |
---|
1152 | // returns p*q, does neither destroy p nor q |
---|
1153 | static inline poly pp_Mult_qq(poly p, poly q, const ring r) |
---|
1154 | { |
---|
1155 | if (p == NULL || q == NULL) return NULL; |
---|
1156 | |
---|
1157 | if (pNext(p) == NULL) |
---|
1158 | { |
---|
1159 | return r->p_Procs->pp_mm_Mult(q, p, r); |
---|
1160 | } |
---|
1161 | |
---|
1162 | if (pNext(q) == NULL) |
---|
1163 | { |
---|
1164 | return r->p_Procs->pp_Mult_mm(p, q, r); |
---|
1165 | } |
---|
1166 | |
---|
1167 | poly qq = q; |
---|
1168 | if (p == q) |
---|
1169 | qq = p_Copy(q, r); |
---|
1170 | |
---|
1171 | poly res; |
---|
1172 | #if defined(HAVE_PLURAL) || defined(HAVE_SHIFTBBA) |
---|
1173 | if (rIsNCRing(r)) |
---|
1174 | res = _nc_pp_Mult_qq(p, qq, r); |
---|
1175 | else |
---|
1176 | #endif |
---|
1177 | res = _p_Mult_q(p, qq, 1, r); |
---|
1178 | |
---|
1179 | if (qq != q) |
---|
1180 | p_Delete(&qq, r); |
---|
1181 | return res; |
---|
1182 | } |
---|
1183 | |
---|
1184 | // returns p + m*q destroys p, const: q, m |
---|
1185 | static inline poly p_Plus_mm_Mult_qq(poly p, poly m, poly q, int &lp, int lq, |
---|
1186 | const ring r) |
---|
1187 | { |
---|
1188 | #ifdef HAVE_PLURAL |
---|
1189 | if (rIsPluralRing(r)) |
---|
1190 | return nc_p_Plus_mm_Mult_qq(p, m, q, lp, lq, r); |
---|
1191 | #endif |
---|
1192 | |
---|
1193 | // this should be implemented more efficiently |
---|
1194 | poly res; |
---|
1195 | int shorter; |
---|
1196 | number n_old = pGetCoeff(m); |
---|
1197 | number n_neg = n_Copy(n_old, r->cf); |
---|
1198 | n_neg = n_InpNeg(n_neg, r->cf); |
---|
1199 | pSetCoeff0(m, n_neg); |
---|
1200 | res = r->p_Procs->p_Minus_mm_Mult_qq(p, m, q, shorter, NULL, r); |
---|
1201 | lp = (lp + lq) - shorter; |
---|
1202 | pSetCoeff0(m, n_old); |
---|
1203 | n_Delete(&n_neg, r->cf); |
---|
1204 | return res; |
---|
1205 | } |
---|
1206 | |
---|
1207 | static inline poly p_Plus_mm_Mult_qq(poly p, poly m, poly q, const ring r) |
---|
1208 | { |
---|
1209 | int lp = 0, lq = 0; |
---|
1210 | return p_Plus_mm_Mult_qq(p, m, q, lp, lq, r); |
---|
1211 | } |
---|
1212 | |
---|
1213 | // returns merged p and q, assumes p and q have no monomials which are equal |
---|
1214 | static inline poly p_Merge_q(poly p, poly q, const ring r) |
---|
1215 | { |
---|
1216 | assume( (p != q) || (p == NULL && q == NULL) ); |
---|
1217 | return r->p_Procs->p_Merge_q(p, q, r); |
---|
1218 | } |
---|
1219 | |
---|
1220 | // like p_SortMerge, except that p may have equal monimals |
---|
1221 | static inline poly p_SortAdd(poly p, const ring r, BOOLEAN revert= FALSE) |
---|
1222 | { |
---|
1223 | if (revert) p = pReverse(p); |
---|
1224 | return sBucketSortAdd(p, r); |
---|
1225 | } |
---|
1226 | |
---|
1227 | // sorts p using bucket sort: returns sorted poly |
---|
1228 | // assumes that monomials of p are all different |
---|
1229 | // reverses it first, if revert == TRUE, use this if input p is "almost" sorted |
---|
1230 | // correctly |
---|
1231 | static inline poly p_SortMerge(poly p, const ring r, BOOLEAN revert= FALSE) |
---|
1232 | { |
---|
1233 | if (revert) p = pReverse(p); |
---|
1234 | return sBucketSortMerge(p, r); |
---|
1235 | } |
---|
1236 | |
---|
1237 | /*************************************************************** |
---|
1238 | * |
---|
1239 | * I/O |
---|
1240 | * |
---|
1241 | ***************************************************************/ |
---|
1242 | static inline char* p_String(poly p, ring p_ring) |
---|
1243 | { |
---|
1244 | return p_String(p, p_ring, p_ring); |
---|
1245 | } |
---|
1246 | static inline void p_String0(poly p, ring p_ring) |
---|
1247 | { |
---|
1248 | p_String0(p, p_ring, p_ring); |
---|
1249 | } |
---|
1250 | static inline void p_Write(poly p, ring p_ring) |
---|
1251 | { |
---|
1252 | p_Write(p, p_ring, p_ring); |
---|
1253 | } |
---|
1254 | static inline void p_Write0(poly p, ring p_ring) |
---|
1255 | { |
---|
1256 | p_Write0(p, p_ring, p_ring); |
---|
1257 | } |
---|
1258 | static inline void p_wrp(poly p, ring p_ring) |
---|
1259 | { |
---|
1260 | p_wrp(p, p_ring, p_ring); |
---|
1261 | } |
---|
1262 | |
---|
1263 | |
---|
1264 | #if PDEBUG > 0 |
---|
1265 | |
---|
1266 | #define _p_LmCmpAction(p, q, r, actionE, actionG, actionS) \ |
---|
1267 | do \ |
---|
1268 | { \ |
---|
1269 | int _cmp = p_LmCmp(p,q,r); \ |
---|
1270 | if (_cmp == 0) actionE; \ |
---|
1271 | if (_cmp == 1) actionG; \ |
---|
1272 | actionS; \ |
---|
1273 | } \ |
---|
1274 | while(0) |
---|
1275 | |
---|
1276 | #else |
---|
1277 | |
---|
1278 | #define _p_LmCmpAction(p, q, r, actionE, actionG, actionS) \ |
---|
1279 | p_MemCmp_LengthGeneral_OrdGeneral(p->exp, q->exp, r->CmpL_Size, r->ordsgn, \ |
---|
1280 | actionE, actionG, actionS) |
---|
1281 | |
---|
1282 | #endif |
---|
1283 | |
---|
1284 | #define pDivAssume(x) do {} while (0) |
---|
1285 | |
---|
1286 | |
---|
1287 | |
---|
1288 | /*************************************************************** |
---|
1289 | * |
---|
1290 | * Allocation/Initalization/Deletion |
---|
1291 | * |
---|
1292 | ***************************************************************/ |
---|
1293 | // adjustments for negative weights |
---|
1294 | static inline void p_MemAdd_NegWeightAdjust(poly p, const ring r) |
---|
1295 | { |
---|
1296 | if (r->NegWeightL_Offset != NULL) |
---|
1297 | { |
---|
1298 | for (int i=r->NegWeightL_Size-1; i>=0; i--) |
---|
1299 | { |
---|
1300 | p->exp[r->NegWeightL_Offset[i]] -= POLY_NEGWEIGHT_OFFSET; |
---|
1301 | } |
---|
1302 | } |
---|
1303 | } |
---|
1304 | static inline void p_MemSub_NegWeightAdjust(poly p, const ring r) |
---|
1305 | { |
---|
1306 | if (r->NegWeightL_Offset != NULL) |
---|
1307 | { |
---|
1308 | for (int i=r->NegWeightL_Size-1; i>=0; i--) |
---|
1309 | { |
---|
1310 | p->exp[r->NegWeightL_Offset[i]] += POLY_NEGWEIGHT_OFFSET; |
---|
1311 | } |
---|
1312 | } |
---|
1313 | } |
---|
1314 | // ExpVextor(d_p) = ExpVector(s_p) |
---|
1315 | static inline void p_ExpVectorCopy(poly d_p, poly s_p, const ring r) |
---|
1316 | { |
---|
1317 | p_LmCheckPolyRing1(d_p, r); |
---|
1318 | p_LmCheckPolyRing1(s_p, r); |
---|
1319 | memcpy(d_p->exp, s_p->exp, r->ExpL_Size*sizeof(long)); |
---|
1320 | } |
---|
1321 | |
---|
1322 | static inline poly p_Init(const ring r, omBin bin) |
---|
1323 | { |
---|
1324 | p_CheckRing1(r); |
---|
1325 | pAssume1(bin != NULL && omSizeWOfBin(r->PolyBin) == omSizeWOfBin(bin)); |
---|
1326 | poly p; |
---|
1327 | omTypeAlloc0Bin(poly, p, bin); |
---|
1328 | p_MemAdd_NegWeightAdjust(p, r); |
---|
1329 | p_SetRingOfLm(p, r); |
---|
1330 | return p; |
---|
1331 | } |
---|
1332 | static inline poly p_Init(const ring r) |
---|
1333 | { |
---|
1334 | return p_Init(r, r->PolyBin); |
---|
1335 | } |
---|
1336 | |
---|
1337 | static inline poly p_LmInit(poly p, const ring r) |
---|
1338 | { |
---|
1339 | p_LmCheckPolyRing1(p, r); |
---|
1340 | poly np; |
---|
1341 | omTypeAllocBin(poly, np, r->PolyBin); |
---|
1342 | p_SetRingOfLm(np, r); |
---|
1343 | memcpy(np->exp, p->exp, r->ExpL_Size*sizeof(long)); |
---|
1344 | pNext(np) = NULL; |
---|
1345 | pSetCoeff0(np, NULL); |
---|
1346 | return np; |
---|
1347 | } |
---|
1348 | static inline poly p_LmInit(poly s_p, const ring s_r, const ring d_r, omBin d_bin) |
---|
1349 | { |
---|
1350 | p_LmCheckPolyRing1(s_p, s_r); |
---|
1351 | p_CheckRing(d_r); |
---|
1352 | pAssume1(d_r->N <= s_r->N); |
---|
1353 | poly d_p = p_Init(d_r, d_bin); |
---|
1354 | for (unsigned i=d_r->N; i!=0; i--) |
---|
1355 | { |
---|
1356 | p_SetExp(d_p, i, p_GetExp(s_p, i,s_r), d_r); |
---|
1357 | } |
---|
1358 | if (rRing_has_Comp(d_r)) |
---|
1359 | { |
---|
1360 | p_SetComp(d_p, p_GetComp(s_p,s_r), d_r); |
---|
1361 | } |
---|
1362 | p_Setm(d_p, d_r); |
---|
1363 | return d_p; |
---|
1364 | } |
---|
1365 | static inline poly p_LmInit(poly s_p, const ring s_r, const ring d_r) |
---|
1366 | { |
---|
1367 | pAssume1(d_r != NULL); |
---|
1368 | return p_LmInit(s_p, s_r, d_r, d_r->PolyBin); |
---|
1369 | } |
---|
1370 | |
---|
1371 | // set all exponents l..k to 0, assume exp. k+1..n and 1..l-1 are in |
---|
1372 | // different blocks |
---|
1373 | // set coeff to 1 |
---|
1374 | static inline poly p_GetExp_k_n(poly p, int l, int k, const ring r) |
---|
1375 | { |
---|
1376 | if (p == NULL) return NULL; |
---|
1377 | p_LmCheckPolyRing1(p, r); |
---|
1378 | poly np; |
---|
1379 | omTypeAllocBin(poly, np, r->PolyBin); |
---|
1380 | p_SetRingOfLm(np, r); |
---|
1381 | memcpy(np->exp, p->exp, r->ExpL_Size*sizeof(long)); |
---|
1382 | pNext(np) = NULL; |
---|
1383 | pSetCoeff0(np, n_Init(1, r->cf)); |
---|
1384 | int i; |
---|
1385 | for(i=l;i<=k;i++) |
---|
1386 | { |
---|
1387 | //np->exp[(r->VarOffset[i] & 0xffffff)] =0; |
---|
1388 | p_SetExp(np,i,0,r); |
---|
1389 | } |
---|
1390 | p_Setm(np,r); |
---|
1391 | return np; |
---|
1392 | } |
---|
1393 | |
---|
1394 | // simialar to p_ShallowCopyDelete but does it only for leading monomial |
---|
1395 | static inline poly p_LmShallowCopyDelete(poly p, const ring r) |
---|
1396 | { |
---|
1397 | p_LmCheckPolyRing1(p, r); |
---|
1398 | pAssume1(omSizeWOfBin(bin) == omSizeWOfBin(r->PolyBin)); |
---|
1399 | poly new_p = p_New(r); |
---|
1400 | memcpy(new_p->exp, p->exp, r->ExpL_Size*sizeof(long)); |
---|
1401 | pSetCoeff0(new_p, pGetCoeff(p)); |
---|
1402 | pNext(new_p) = pNext(p); |
---|
1403 | omFreeBinAddr(p); |
---|
1404 | return new_p; |
---|
1405 | } |
---|
1406 | |
---|
1407 | /*************************************************************** |
---|
1408 | * |
---|
1409 | * Operation on ExpVectors |
---|
1410 | * |
---|
1411 | ***************************************************************/ |
---|
1412 | // ExpVector(p1) += ExpVector(p2) |
---|
1413 | static inline void p_ExpVectorAdd(poly p1, poly p2, const ring r) |
---|
1414 | { |
---|
1415 | p_LmCheckPolyRing1(p1, r); |
---|
1416 | p_LmCheckPolyRing1(p2, r); |
---|
1417 | #if PDEBUG >= 1 |
---|
1418 | for (int i=1; i<=r->N; i++) |
---|
1419 | pAssume1((unsigned long) (p_GetExp(p1, i, r) + p_GetExp(p2, i, r)) <= r->bitmask); |
---|
1420 | pAssume1(p_GetComp(p1, r) == 0 || p_GetComp(p2, r) == 0); |
---|
1421 | #endif |
---|
1422 | |
---|
1423 | p_MemAdd_LengthGeneral(p1->exp, p2->exp, r->ExpL_Size); |
---|
1424 | p_MemAdd_NegWeightAdjust(p1, r); |
---|
1425 | } |
---|
1426 | // ExpVector(pr) = ExpVector(p1) + ExpVector(p2) |
---|
1427 | static inline void p_ExpVectorSum(poly pr, poly p1, poly p2, const ring r) |
---|
1428 | { |
---|
1429 | p_LmCheckPolyRing1(p1, r); |
---|
1430 | p_LmCheckPolyRing1(p2, r); |
---|
1431 | p_LmCheckPolyRing1(pr, r); |
---|
1432 | #if PDEBUG >= 1 |
---|
1433 | for (int i=1; i<=r->N; i++) |
---|
1434 | pAssume1((unsigned long) (p_GetExp(p1, i, r) + p_GetExp(p2, i, r)) <= r->bitmask); |
---|
1435 | pAssume1(p_GetComp(p1, r) == 0 || p_GetComp(p2, r) == 0); |
---|
1436 | #endif |
---|
1437 | |
---|
1438 | p_MemSum_LengthGeneral(pr->exp, p1->exp, p2->exp, r->ExpL_Size); |
---|
1439 | p_MemAdd_NegWeightAdjust(pr, r); |
---|
1440 | } |
---|
1441 | // ExpVector(p1) -= ExpVector(p2) |
---|
1442 | static inline void p_ExpVectorSub(poly p1, poly p2, const ring r) |
---|
1443 | { |
---|
1444 | p_LmCheckPolyRing1(p1, r); |
---|
1445 | p_LmCheckPolyRing1(p2, r); |
---|
1446 | #if PDEBUG >= 1 |
---|
1447 | for (int i=1; i<=r->N; i++) |
---|
1448 | pAssume1(p_GetExp(p1, i, r) >= p_GetExp(p2, i, r)); |
---|
1449 | pAssume1(p_GetComp(p1, r) == 0 || p_GetComp(p2, r) == 0 || |
---|
1450 | p_GetComp(p1, r) == p_GetComp(p2, r)); |
---|
1451 | #endif |
---|
1452 | |
---|
1453 | p_MemSub_LengthGeneral(p1->exp, p2->exp, r->ExpL_Size); |
---|
1454 | p_MemSub_NegWeightAdjust(p1, r); |
---|
1455 | } |
---|
1456 | |
---|
1457 | // ExpVector(p1) += ExpVector(p2) - ExpVector(p3) |
---|
1458 | static inline void p_ExpVectorAddSub(poly p1, poly p2, poly p3, const ring r) |
---|
1459 | { |
---|
1460 | p_LmCheckPolyRing1(p1, r); |
---|
1461 | p_LmCheckPolyRing1(p2, r); |
---|
1462 | p_LmCheckPolyRing1(p3, r); |
---|
1463 | #if PDEBUG >= 1 |
---|
1464 | for (int i=1; i<=r->N; i++) |
---|
1465 | pAssume1(p_GetExp(p1, i, r) + p_GetExp(p2, i, r) >= p_GetExp(p3, i, r)); |
---|
1466 | pAssume1(p_GetComp(p1, r) == 0 || |
---|
1467 | (p_GetComp(p2, r) - p_GetComp(p3, r) == 0) || |
---|
1468 | (p_GetComp(p1, r) == p_GetComp(p2, r) - p_GetComp(p3, r))); |
---|
1469 | #endif |
---|
1470 | |
---|
1471 | p_MemAddSub_LengthGeneral(p1->exp, p2->exp, p3->exp, r->ExpL_Size); |
---|
1472 | // no need to adjust in case of NegWeights |
---|
1473 | } |
---|
1474 | |
---|
1475 | // ExpVector(pr) = ExpVector(p1) - ExpVector(p2) |
---|
1476 | static inline void p_ExpVectorDiff(poly pr, poly p1, poly p2, const ring r) |
---|
1477 | { |
---|
1478 | p_LmCheckPolyRing1(p1, r); |
---|
1479 | p_LmCheckPolyRing1(p2, r); |
---|
1480 | p_LmCheckPolyRing1(pr, r); |
---|
1481 | #if PDEBUG >= 2 |
---|
1482 | for (int i=1; i<=r->N; i++) |
---|
1483 | pAssume1(p_GetExp(p1, i, r) >= p_GetExp(p2, i, r)); |
---|
1484 | pAssume1(!rRing_has_Comp(r) || p_GetComp(p1, r) == p_GetComp(p2, r)); |
---|
1485 | #endif |
---|
1486 | |
---|
1487 | p_MemDiff_LengthGeneral(pr->exp, p1->exp, p2->exp, r->ExpL_Size); |
---|
1488 | p_MemSub_NegWeightAdjust(pr, r); |
---|
1489 | } |
---|
1490 | |
---|
1491 | static inline BOOLEAN p_ExpVectorEqual(poly p1, poly p2, const ring r) |
---|
1492 | { |
---|
1493 | p_LmCheckPolyRing1(p1, r); |
---|
1494 | p_LmCheckPolyRing1(p2, r); |
---|
1495 | |
---|
1496 | unsigned i = r->ExpL_Size; |
---|
1497 | unsigned long *ep = p1->exp; |
---|
1498 | unsigned long *eq = p2->exp; |
---|
1499 | |
---|
1500 | do |
---|
1501 | { |
---|
1502 | i--; |
---|
1503 | if (ep[i] != eq[i]) return FALSE; |
---|
1504 | } |
---|
1505 | while (i!=0); |
---|
1506 | return TRUE; |
---|
1507 | } |
---|
1508 | |
---|
1509 | static inline long p_Totaldegree(poly p, const ring r) |
---|
1510 | { |
---|
1511 | p_LmCheckPolyRing1(p, r); |
---|
1512 | unsigned long s = p_GetTotalDegree(p->exp[r->VarL_Offset[0]], |
---|
1513 | r, |
---|
1514 | r->ExpPerLong); |
---|
1515 | for (unsigned i=r->VarL_Size-1; i!=0; i--) |
---|
1516 | { |
---|
1517 | s += p_GetTotalDegree(p->exp[r->VarL_Offset[i]], r,r->ExpPerLong); |
---|
1518 | } |
---|
1519 | return (long)s; |
---|
1520 | } |
---|
1521 | |
---|
1522 | static inline void p_GetExpV(poly p, int *ev, const ring r) |
---|
1523 | { |
---|
1524 | p_LmCheckPolyRing1(p, r); |
---|
1525 | for (unsigned j = r->N; j!=0; j--) |
---|
1526 | ev[j] = p_GetExp(p, j, r); |
---|
1527 | |
---|
1528 | ev[0] = p_GetComp(p, r); |
---|
1529 | } |
---|
1530 | // p_GetExpVL is used in Singular,jl |
---|
1531 | static inline void p_GetExpVL(poly p, int64 *ev, const ring r) |
---|
1532 | { |
---|
1533 | p_LmCheckPolyRing1(p, r); |
---|
1534 | for (unsigned j = r->N; j!=0; j--) |
---|
1535 | ev[j-1] = p_GetExp(p, j, r); |
---|
1536 | } |
---|
1537 | // p_GetExpVLV is used in Singular,jl |
---|
1538 | static inline int64 p_GetExpVLV(poly p, int64 *ev, const ring r) |
---|
1539 | { |
---|
1540 | p_LmCheckPolyRing1(p, r); |
---|
1541 | for (unsigned j = r->N; j!=0; j--) |
---|
1542 | ev[j-1] = p_GetExp(p, j, r); |
---|
1543 | return (int64)p_GetComp(p,r); |
---|
1544 | } |
---|
1545 | // p_GetExpVL is used in Singular,jl |
---|
1546 | static inline void p_SetExpV(poly p, int *ev, const ring r) |
---|
1547 | { |
---|
1548 | p_LmCheckPolyRing1(p, r); |
---|
1549 | for (unsigned j = r->N; j!=0; j--) |
---|
1550 | p_SetExp(p, j, ev[j], r); |
---|
1551 | |
---|
1552 | if(ev[0]!=0) p_SetComp(p, ev[0],r); |
---|
1553 | p_Setm(p, r); |
---|
1554 | } |
---|
1555 | static inline void p_SetExpVL(poly p, int64 *ev, const ring r) |
---|
1556 | { |
---|
1557 | p_LmCheckPolyRing1(p, r); |
---|
1558 | for (unsigned j = r->N; j!=0; j--) |
---|
1559 | p_SetExp(p, j, ev[j-1], r); |
---|
1560 | p_SetComp(p, 0,r); |
---|
1561 | |
---|
1562 | p_Setm(p, r); |
---|
1563 | } |
---|
1564 | |
---|
1565 | // p_SetExpVLV is used in Singular,jl |
---|
1566 | static inline void p_SetExpVLV(poly p, int64 *ev, int64 comp, const ring r) |
---|
1567 | { |
---|
1568 | p_LmCheckPolyRing1(p, r); |
---|
1569 | for (unsigned j = r->N; j!=0; j--) |
---|
1570 | p_SetExp(p, j, ev[j-1], r); |
---|
1571 | p_SetComp(p, comp,r); |
---|
1572 | |
---|
1573 | p_Setm(p, r); |
---|
1574 | } |
---|
1575 | |
---|
1576 | /*************************************************************** |
---|
1577 | * |
---|
1578 | * Comparison w.r.t. monomial ordering |
---|
1579 | * |
---|
1580 | ***************************************************************/ |
---|
1581 | |
---|
1582 | static inline int p_LmCmp(poly p, poly q, const ring r) |
---|
1583 | { |
---|
1584 | p_LmCheckPolyRing1(p, r); |
---|
1585 | p_LmCheckPolyRing1(q, r); |
---|
1586 | |
---|
1587 | const unsigned long* _s1 = ((unsigned long*) p->exp); |
---|
1588 | const unsigned long* _s2 = ((unsigned long*) q->exp); |
---|
1589 | REGISTER unsigned long _v1; |
---|
1590 | REGISTER unsigned long _v2; |
---|
1591 | const unsigned long _l = r->CmpL_Size; |
---|
1592 | |
---|
1593 | REGISTER unsigned long _i=0; |
---|
1594 | |
---|
1595 | LengthGeneral_OrdGeneral_LoopTop: |
---|
1596 | _v1 = _s1[_i]; |
---|
1597 | _v2 = _s2[_i]; |
---|
1598 | if (_v1 == _v2) |
---|
1599 | { |
---|
1600 | _i++; |
---|
1601 | if (_i == _l) return 0; |
---|
1602 | goto LengthGeneral_OrdGeneral_LoopTop; |
---|
1603 | } |
---|
1604 | const long* _ordsgn = (long*) r->ordsgn; |
---|
1605 | #if 1 /* two variants*/ |
---|
1606 | if (_v1 > _v2) |
---|
1607 | { |
---|
1608 | return _ordsgn[_i]; |
---|
1609 | } |
---|
1610 | return -(_ordsgn[_i]); |
---|
1611 | #else |
---|
1612 | if (_v1 > _v2) |
---|
1613 | { |
---|
1614 | if (_ordsgn[_i] == 1) return 1; |
---|
1615 | return -1; |
---|
1616 | } |
---|
1617 | if (_ordsgn[_i] == 1) return -1; |
---|
1618 | return 1; |
---|
1619 | #endif |
---|
1620 | } |
---|
1621 | |
---|
1622 | // The coefficient will be compared in absolute value |
---|
1623 | static inline int p_LtCmp(poly p, poly q, const ring r) |
---|
1624 | { |
---|
1625 | int res = p_LmCmp(p,q,r); |
---|
1626 | if(res == 0) |
---|
1627 | { |
---|
1628 | if(p_GetCoeff(p,r) == NULL || p_GetCoeff(q,r) == NULL) |
---|
1629 | return res; |
---|
1630 | number pc = n_Copy(p_GetCoeff(p,r),r->cf); |
---|
1631 | number qc = n_Copy(p_GetCoeff(q,r),r->cf); |
---|
1632 | if(!n_GreaterZero(pc,r->cf)) |
---|
1633 | pc = n_InpNeg(pc,r->cf); |
---|
1634 | if(!n_GreaterZero(qc,r->cf)) |
---|
1635 | qc = n_InpNeg(qc,r->cf); |
---|
1636 | if(n_Greater(pc,qc,r->cf)) |
---|
1637 | res = 1; |
---|
1638 | else if(n_Greater(qc,pc,r->cf)) |
---|
1639 | res = -1; |
---|
1640 | else if(n_Equal(pc,qc,r->cf)) |
---|
1641 | res = 0; |
---|
1642 | n_Delete(&pc,r->cf); |
---|
1643 | n_Delete(&qc,r->cf); |
---|
1644 | } |
---|
1645 | return res; |
---|
1646 | } |
---|
1647 | |
---|
1648 | // The coefficient will be compared in absolute value |
---|
1649 | static inline int p_LtCmpNoAbs(poly p, poly q, const ring r) |
---|
1650 | { |
---|
1651 | int res = p_LmCmp(p,q,r); |
---|
1652 | if(res == 0) |
---|
1653 | { |
---|
1654 | if(p_GetCoeff(p,r) == NULL || p_GetCoeff(q,r) == NULL) |
---|
1655 | return res; |
---|
1656 | number pc = p_GetCoeff(p,r); |
---|
1657 | number qc = p_GetCoeff(q,r); |
---|
1658 | if(n_Greater(pc,qc,r->cf)) |
---|
1659 | res = 1; |
---|
1660 | if(n_Greater(qc,pc,r->cf)) |
---|
1661 | res = -1; |
---|
1662 | if(n_Equal(pc,qc,r->cf)) |
---|
1663 | res = 0; |
---|
1664 | } |
---|
1665 | return res; |
---|
1666 | } |
---|
1667 | |
---|
1668 | #ifdef HAVE_RINGS |
---|
1669 | // This is the equivalent of pLmCmp(p,q) != -currRing->OrdSgn for rings |
---|
1670 | // It is used in posInLRing and posInTRing |
---|
1671 | static inline int p_LtCmpOrdSgnDiffM(poly p, poly q, const ring r) |
---|
1672 | { |
---|
1673 | return(p_LtCmp(p,q,r) == r->OrdSgn); |
---|
1674 | } |
---|
1675 | #endif |
---|
1676 | |
---|
1677 | #ifdef HAVE_RINGS |
---|
1678 | // This is the equivalent of pLmCmp(p,q) != currRing->OrdSgn for rings |
---|
1679 | // It is used in posInLRing and posInTRing |
---|
1680 | static inline int p_LtCmpOrdSgnDiffP(poly p, poly q, const ring r) |
---|
1681 | { |
---|
1682 | if(r->OrdSgn == 1) |
---|
1683 | { |
---|
1684 | return(p_LmCmp(p,q,r) == -1); |
---|
1685 | } |
---|
1686 | else |
---|
1687 | { |
---|
1688 | return(p_LtCmp(p,q,r) != -1); |
---|
1689 | } |
---|
1690 | } |
---|
1691 | #endif |
---|
1692 | |
---|
1693 | #ifdef HAVE_RINGS |
---|
1694 | // This is the equivalent of pLmCmp(p,q) == -currRing->OrdSgn for rings |
---|
1695 | // It is used in posInLRing and posInTRing |
---|
1696 | static inline int p_LtCmpOrdSgnEqM(poly p, poly q, const ring r) |
---|
1697 | { |
---|
1698 | return(p_LtCmp(p,q,r) == -r->OrdSgn); |
---|
1699 | } |
---|
1700 | #endif |
---|
1701 | |
---|
1702 | #ifdef HAVE_RINGS |
---|
1703 | // This is the equivalent of pLmCmp(p,q) == currRing->OrdSgn for rings |
---|
1704 | // It is used in posInLRing and posInTRing |
---|
1705 | static inline int p_LtCmpOrdSgnEqP(poly p, poly q, const ring r) |
---|
1706 | { |
---|
1707 | return(p_LtCmp(p,q,r) == r->OrdSgn); |
---|
1708 | } |
---|
1709 | #endif |
---|
1710 | |
---|
1711 | /// returns TRUE if p1 is a skalar multiple of p2 |
---|
1712 | /// assume p1 != NULL and p2 != NULL |
---|
1713 | BOOLEAN p_ComparePolys(poly p1,poly p2, const ring r); |
---|
1714 | |
---|
1715 | |
---|
1716 | /*************************************************************** |
---|
1717 | * |
---|
1718 | * Comparisons: they are all done without regarding coeffs |
---|
1719 | * |
---|
1720 | ***************************************************************/ |
---|
1721 | #define p_LmCmpAction(p, q, r, actionE, actionG, actionS) \ |
---|
1722 | _p_LmCmpAction(p, q, r, actionE, actionG, actionS) |
---|
1723 | |
---|
1724 | // returns 1 if ExpVector(p)==ExpVector(q): does not compare numbers !! |
---|
1725 | #define p_LmEqual(p1, p2, r) p_ExpVectorEqual(p1, p2, r) |
---|
1726 | |
---|
1727 | // pCmp: args may be NULL |
---|
1728 | // returns: (p2==NULL ? 1 : (p1 == NULL ? -1 : p_LmCmp(p1, p2))) |
---|
1729 | static inline int p_Cmp(poly p1, poly p2, ring r) |
---|
1730 | { |
---|
1731 | if (p2==NULL) |
---|
1732 | { |
---|
1733 | if (p1==NULL) return 0; |
---|
1734 | return 1; |
---|
1735 | } |
---|
1736 | if (p1==NULL) |
---|
1737 | return -1; |
---|
1738 | return p_LmCmp(p1,p2,r); |
---|
1739 | } |
---|
1740 | |
---|
1741 | static inline int p_CmpPolys(poly p1, poly p2, ring r) |
---|
1742 | { |
---|
1743 | if (p2==NULL) |
---|
1744 | { |
---|
1745 | if (p1==NULL) return 0; |
---|
1746 | return 1; |
---|
1747 | } |
---|
1748 | if (p1==NULL) |
---|
1749 | return -1; |
---|
1750 | return p_ComparePolys(p1,p2,r); |
---|
1751 | } |
---|
1752 | |
---|
1753 | |
---|
1754 | /*************************************************************** |
---|
1755 | * |
---|
1756 | * divisibility |
---|
1757 | * |
---|
1758 | ***************************************************************/ |
---|
1759 | /// return: FALSE, if there exists i, such that a->exp[i] > b->exp[i] |
---|
1760 | /// TRUE, otherwise |
---|
1761 | /// (1) Consider long vars, instead of single exponents |
---|
1762 | /// (2) Clearly, if la > lb, then FALSE |
---|
1763 | /// (3) Suppose la <= lb, and consider first bits of single exponents in l: |
---|
1764 | /// if TRUE, then value of these bits is la ^ lb |
---|
1765 | /// if FALSE, then la-lb causes an "overflow" into one of those bits, i.e., |
---|
1766 | /// la ^ lb != la - lb |
---|
1767 | static inline BOOLEAN _p_LmDivisibleByNoComp(poly a, poly b, const ring r) |
---|
1768 | { |
---|
1769 | int i=r->VarL_Size - 1; |
---|
1770 | unsigned long divmask = r->divmask; |
---|
1771 | unsigned long la, lb; |
---|
1772 | |
---|
1773 | if (r->VarL_LowIndex >= 0) |
---|
1774 | { |
---|
1775 | i += r->VarL_LowIndex; |
---|
1776 | do |
---|
1777 | { |
---|
1778 | la = a->exp[i]; |
---|
1779 | lb = b->exp[i]; |
---|
1780 | if ((la > lb) || |
---|
1781 | (((la & divmask) ^ (lb & divmask)) != ((lb - la) & divmask))) |
---|
1782 | { |
---|
1783 | pDivAssume(p_DebugLmDivisibleByNoComp(a, b, r) == FALSE); |
---|
1784 | return FALSE; |
---|
1785 | } |
---|
1786 | i--; |
---|
1787 | } |
---|
1788 | while (i>=r->VarL_LowIndex); |
---|
1789 | } |
---|
1790 | else |
---|
1791 | { |
---|
1792 | do |
---|
1793 | { |
---|
1794 | la = a->exp[r->VarL_Offset[i]]; |
---|
1795 | lb = b->exp[r->VarL_Offset[i]]; |
---|
1796 | if ((la > lb) || |
---|
1797 | (((la & divmask) ^ (lb & divmask)) != ((lb - la) & divmask))) |
---|
1798 | { |
---|
1799 | pDivAssume(p_DebugLmDivisibleByNoComp(a, b, r) == FALSE); |
---|
1800 | return FALSE; |
---|
1801 | } |
---|
1802 | i--; |
---|
1803 | } |
---|
1804 | while (i>=0); |
---|
1805 | } |
---|
1806 | /*#ifdef HAVE_RINGS |
---|
1807 | pDivAssume(p_DebugLmDivisibleByNoComp(a, b, r) == n_DivBy(p_GetCoeff(b, r), p_GetCoeff(a, r), r->cf)); |
---|
1808 | return (!rField_is_Ring(r)) || n_DivBy(p_GetCoeff(b, r), p_GetCoeff(a, r), r->cf); |
---|
1809 | #else |
---|
1810 | */ |
---|
1811 | pDivAssume(p_DebugLmDivisibleByNoComp(a, b, r) == TRUE); |
---|
1812 | return TRUE; |
---|
1813 | //#endif |
---|
1814 | } |
---|
1815 | |
---|
1816 | static inline BOOLEAN _p_LmDivisibleByNoComp(poly a, const ring r_a, poly b, const ring r_b) |
---|
1817 | { |
---|
1818 | int i=r_a->N; |
---|
1819 | pAssume1(r_a->N == r_b->N); |
---|
1820 | |
---|
1821 | do |
---|
1822 | { |
---|
1823 | if (p_GetExp(a,i,r_a) > p_GetExp(b,i,r_b)) |
---|
1824 | { |
---|
1825 | return FALSE; |
---|
1826 | } |
---|
1827 | i--; |
---|
1828 | } |
---|
1829 | while (i); |
---|
1830 | /*#ifdef HAVE_RINGS |
---|
1831 | return n_DivBy(p_GetCoeff(b, r_b), p_GetCoeff(a, r_a), r_a->cf); |
---|
1832 | #else |
---|
1833 | */ |
---|
1834 | return TRUE; |
---|
1835 | //#endif |
---|
1836 | } |
---|
1837 | |
---|
1838 | #ifdef HAVE_RATGRING |
---|
1839 | static inline BOOLEAN _p_LmDivisibleByNoCompPart(poly a, const ring r_a, poly b, const ring r_b,const int start, const int end) |
---|
1840 | { |
---|
1841 | int i=end; |
---|
1842 | pAssume1(r_a->N == r_b->N); |
---|
1843 | |
---|
1844 | do |
---|
1845 | { |
---|
1846 | if (p_GetExp(a,i,r_a) > p_GetExp(b,i,r_b)) |
---|
1847 | return FALSE; |
---|
1848 | i--; |
---|
1849 | } |
---|
1850 | while (i>=start); |
---|
1851 | /*#ifdef HAVE_RINGS |
---|
1852 | return n_DivBy(p_GetCoeff(b, r_b), p_GetCoeff(a, r_a), r_a->cf); |
---|
1853 | #else |
---|
1854 | */ |
---|
1855 | return TRUE; |
---|
1856 | //#endif |
---|
1857 | } |
---|
1858 | static inline BOOLEAN _p_LmDivisibleByPart(poly a, const ring r_a, poly b, const ring r_b,const int start, const int end) |
---|
1859 | { |
---|
1860 | if (p_GetComp(a, r_a) == 0 || p_GetComp(a,r_a) == p_GetComp(b,r_b)) |
---|
1861 | return _p_LmDivisibleByNoCompPart(a, r_a, b, r_b,start,end); |
---|
1862 | return FALSE; |
---|
1863 | } |
---|
1864 | static inline BOOLEAN p_LmDivisibleByPart(poly a, poly b, const ring r,const int start, const int end) |
---|
1865 | { |
---|
1866 | p_LmCheckPolyRing1(b, r); |
---|
1867 | pIfThen1(a != NULL, p_LmCheckPolyRing1(b, r)); |
---|
1868 | if (p_GetComp(a, r) == 0 || p_GetComp(a,r) == p_GetComp(b,r)) |
---|
1869 | return _p_LmDivisibleByNoCompPart(a, r, b, r,start, end); |
---|
1870 | return FALSE; |
---|
1871 | } |
---|
1872 | #endif |
---|
1873 | static inline BOOLEAN _p_LmDivisibleBy(poly a, poly b, const ring r) |
---|
1874 | { |
---|
1875 | if (p_GetComp(a, r) == 0 || p_GetComp(a,r) == p_GetComp(b,r)) |
---|
1876 | return _p_LmDivisibleByNoComp(a, b, r); |
---|
1877 | return FALSE; |
---|
1878 | } |
---|
1879 | static inline BOOLEAN p_LmDivisibleByNoComp(poly a, poly b, const ring r) |
---|
1880 | { |
---|
1881 | p_LmCheckPolyRing1(a, r); |
---|
1882 | p_LmCheckPolyRing1(b, r); |
---|
1883 | return _p_LmDivisibleByNoComp(a, b, r); |
---|
1884 | } |
---|
1885 | |
---|
1886 | static inline BOOLEAN p_LmDivisibleByNoComp(poly a, const ring ra, poly b, const ring rb) |
---|
1887 | { |
---|
1888 | p_LmCheckPolyRing1(a, ra); |
---|
1889 | p_LmCheckPolyRing1(b, rb); |
---|
1890 | return _p_LmDivisibleByNoComp(a, ra, b, rb); |
---|
1891 | } |
---|
1892 | |
---|
1893 | static inline BOOLEAN p_LmDivisibleBy(poly a, poly b, const ring r) |
---|
1894 | { |
---|
1895 | p_LmCheckPolyRing1(b, r); |
---|
1896 | pIfThen1(a != NULL, p_LmCheckPolyRing1(b, r)); |
---|
1897 | if (p_GetComp(a, r) == 0 || p_GetComp(a,r) == p_GetComp(b,r)) |
---|
1898 | return _p_LmDivisibleByNoComp(a, b, r); |
---|
1899 | return FALSE; |
---|
1900 | } |
---|
1901 | |
---|
1902 | static inline BOOLEAN p_DivisibleBy(poly a, poly b, const ring r) |
---|
1903 | { |
---|
1904 | pIfThen1(b!=NULL, p_LmCheckPolyRing1(b, r)); |
---|
1905 | pIfThen1(a!=NULL, p_LmCheckPolyRing1(a, r)); |
---|
1906 | |
---|
1907 | if (a != NULL && (p_GetComp(a, r) == 0 || p_GetComp(a,r) == p_GetComp(b,r))) |
---|
1908 | return _p_LmDivisibleByNoComp(a,b,r); |
---|
1909 | return FALSE; |
---|
1910 | } |
---|
1911 | |
---|
1912 | static inline BOOLEAN p_LmShortDivisibleBy(poly a, unsigned long sev_a, |
---|
1913 | poly b, unsigned long not_sev_b, const ring r) |
---|
1914 | { |
---|
1915 | p_LmCheckPolyRing1(a, r); |
---|
1916 | p_LmCheckPolyRing1(b, r); |
---|
1917 | #ifndef PDIV_DEBUG |
---|
1918 | _pPolyAssume2(p_GetShortExpVector(a, r) == sev_a, a, r); |
---|
1919 | _pPolyAssume2(p_GetShortExpVector(b, r) == ~ not_sev_b, b, r); |
---|
1920 | |
---|
1921 | if (sev_a & not_sev_b) |
---|
1922 | { |
---|
1923 | pAssume1(p_LmDivisibleByNoComp(a, b, r) == FALSE); |
---|
1924 | return FALSE; |
---|
1925 | } |
---|
1926 | return p_LmDivisibleBy(a, b, r); |
---|
1927 | #else |
---|
1928 | return pDebugLmShortDivisibleBy(a, sev_a, r, b, not_sev_b, r); |
---|
1929 | #endif |
---|
1930 | } |
---|
1931 | |
---|
1932 | static inline BOOLEAN p_LmShortDivisibleByNoComp(poly a, unsigned long sev_a, |
---|
1933 | poly b, unsigned long not_sev_b, const ring r) |
---|
1934 | { |
---|
1935 | p_LmCheckPolyRing1(a, r); |
---|
1936 | p_LmCheckPolyRing1(b, r); |
---|
1937 | #ifndef PDIV_DEBUG |
---|
1938 | _pPolyAssume2(p_GetShortExpVector(a, r) == sev_a, a, r); |
---|
1939 | _pPolyAssume2(p_GetShortExpVector(b, r) == ~ not_sev_b, b, r); |
---|
1940 | |
---|
1941 | if (sev_a & not_sev_b) |
---|
1942 | { |
---|
1943 | pAssume1(p_LmDivisibleByNoComp(a, b, r) == FALSE); |
---|
1944 | return FALSE; |
---|
1945 | } |
---|
1946 | return p_LmDivisibleByNoComp(a, b, r); |
---|
1947 | #else |
---|
1948 | return pDebugLmShortDivisibleByNoComp(a, sev_a, r, b, not_sev_b, r); |
---|
1949 | #endif |
---|
1950 | } |
---|
1951 | |
---|
1952 | /*************************************************************** |
---|
1953 | * |
---|
1954 | * Misc things on Lm |
---|
1955 | * |
---|
1956 | ***************************************************************/ |
---|
1957 | |
---|
1958 | |
---|
1959 | /// like the respective p_LmIs* routines, except that p might be empty |
---|
1960 | static inline BOOLEAN p_IsConstantComp(const poly p, const ring r) |
---|
1961 | { |
---|
1962 | if (p == NULL) return TRUE; |
---|
1963 | return (pNext(p)==NULL) && p_LmIsConstantComp(p, r); |
---|
1964 | } |
---|
1965 | |
---|
1966 | static inline BOOLEAN p_IsConstant(const poly p, const ring r) |
---|
1967 | { |
---|
1968 | if (p == NULL) return TRUE; |
---|
1969 | return (pNext(p)==NULL) && p_LmIsConstant(p, r); |
---|
1970 | } |
---|
1971 | |
---|
1972 | /// either poly(1) or gen(k)?! |
---|
1973 | static inline BOOLEAN p_IsOne(const poly p, const ring R) |
---|
1974 | { |
---|
1975 | if (p == NULL) return FALSE; /* TODO check if 0 == 1 */ |
---|
1976 | p_Test(p, R); |
---|
1977 | return (p_IsConstant(p, R) && n_IsOne(p_GetCoeff(p, R), R->cf)); |
---|
1978 | } |
---|
1979 | |
---|
1980 | static inline BOOLEAN p_IsConstantPoly(const poly p, const ring r) |
---|
1981 | { |
---|
1982 | p_Test(p, r); |
---|
1983 | poly pp=p; |
---|
1984 | while(pp!=NULL) |
---|
1985 | { |
---|
1986 | if (! p_LmIsConstantComp(pp, r)) |
---|
1987 | return FALSE; |
---|
1988 | pIter(pp); |
---|
1989 | } |
---|
1990 | return TRUE; |
---|
1991 | } |
---|
1992 | |
---|
1993 | static inline BOOLEAN p_IsUnit(const poly p, const ring r) |
---|
1994 | { |
---|
1995 | if (p == NULL) return FALSE; |
---|
1996 | if (rField_is_Ring(r)) |
---|
1997 | return (p_LmIsConstant(p, r) && n_IsUnit(pGetCoeff(p),r->cf)); |
---|
1998 | return p_LmIsConstant(p, r); |
---|
1999 | } |
---|
2000 | |
---|
2001 | static inline BOOLEAN p_LmExpVectorAddIsOk(const poly p1, const poly p2, |
---|
2002 | const ring r) |
---|
2003 | { |
---|
2004 | p_LmCheckPolyRing(p1, r); |
---|
2005 | p_LmCheckPolyRing(p2, r); |
---|
2006 | unsigned long l1, l2, divmask = r->divmask; |
---|
2007 | int i; |
---|
2008 | |
---|
2009 | for (i=0; i<r->VarL_Size; i++) |
---|
2010 | { |
---|
2011 | l1 = p1->exp[r->VarL_Offset[i]]; |
---|
2012 | l2 = p2->exp[r->VarL_Offset[i]]; |
---|
2013 | // do the divisiblity trick |
---|
2014 | if ( (l1 > ULONG_MAX - l2) || |
---|
2015 | (((l1 & divmask) ^ (l2 & divmask)) != ((l1 + l2) & divmask))) |
---|
2016 | return FALSE; |
---|
2017 | } |
---|
2018 | return TRUE; |
---|
2019 | } |
---|
2020 | void p_Split(poly p, poly * r); /*p => IN(p), r => REST(p) */ |
---|
2021 | BOOLEAN p_HasNotCF(poly p1, poly p2, const ring r); |
---|
2022 | BOOLEAN p_HasNotCFRing(poly p1, poly p2, const ring r); |
---|
2023 | poly p_mInit(const char *s, BOOLEAN &ok, const ring r); /* monom s -> poly, interpreter */ |
---|
2024 | const char * p_Read(const char *s, poly &p,const ring r); /* monom -> poly */ |
---|
2025 | poly p_MDivide(poly a, poly b, const ring r); |
---|
2026 | poly p_DivideM(poly a, poly b, const ring r); |
---|
2027 | poly pp_DivideM(poly a, poly b, const ring r); |
---|
2028 | poly p_Div_nn(poly p, const number n, const ring r); |
---|
2029 | |
---|
2030 | // returns the LCM of the head terms of a and b in *m, does not p_Setm |
---|
2031 | void p_Lcm(const poly a, const poly b, poly m, const ring r); |
---|
2032 | // returns the LCM of the head terms of a and b, does p_Setm |
---|
2033 | poly p_Lcm(const poly a, const poly b, const ring r); |
---|
2034 | |
---|
2035 | #ifdef HAVE_RATGRING |
---|
2036 | poly p_LcmRat(const poly a, const poly b, const long lCompM, const ring r); |
---|
2037 | poly p_GetCoeffRat(poly p, int ishift, ring r); |
---|
2038 | void p_LmDeleteAndNextRat(poly *p, int ishift, ring r); |
---|
2039 | void p_ContentRat(poly &ph, const ring r); |
---|
2040 | #endif /* ifdef HAVE_RATGRING */ |
---|
2041 | |
---|
2042 | |
---|
2043 | poly p_Diff(poly a, int k, const ring r); |
---|
2044 | poly p_DiffOp(poly a, poly b,BOOLEAN multiply, const ring r); |
---|
2045 | int p_Weight(int c, const ring r); |
---|
2046 | |
---|
2047 | /// assumes that p and divisor are univariate polynomials in r, |
---|
2048 | /// mentioning the same variable; |
---|
2049 | /// assumes divisor != NULL; |
---|
2050 | /// p may be NULL; |
---|
2051 | /// assumes a global monomial ordering in r; |
---|
2052 | /// performs polynomial division of p by divisor: |
---|
2053 | /// - afterwards p contains the remainder of the division, i.e., |
---|
2054 | /// p_before = result * divisor + p_afterwards; |
---|
2055 | /// - if needResult == TRUE, then the method computes and returns 'result', |
---|
2056 | /// otherwise NULL is returned (This parametrization can be used when |
---|
2057 | /// one is only interested in the remainder of the division. In this |
---|
2058 | /// case, the method will be slightly faster.) |
---|
2059 | /// leaves divisor unmodified |
---|
2060 | poly p_PolyDiv(poly &p, const poly divisor, const BOOLEAN needResult, const ring r); |
---|
2061 | |
---|
2062 | /* syszygy stuff */ |
---|
2063 | BOOLEAN p_VectorHasUnitB(poly p, int * k, const ring r); |
---|
2064 | void p_VectorHasUnit(poly p, int * k, int * len, const ring r); |
---|
2065 | poly p_TakeOutComp1(poly * p, int k, const ring r); |
---|
2066 | // Splits *p into two polys: *q which consists of all monoms with |
---|
2067 | // component == comp and *p of all other monoms *lq == pLength(*q) |
---|
2068 | // On return all components pf *q == 0 |
---|
2069 | void p_TakeOutComp(poly *p, long comp, poly *q, int *lq, const ring r); |
---|
2070 | |
---|
2071 | // This is something weird -- Don't use it, unless you know what you are doing |
---|
2072 | poly p_TakeOutComp(poly * p, int k, const ring r); |
---|
2073 | |
---|
2074 | void p_DeleteComp(poly * p,int k, const ring r); |
---|
2075 | |
---|
2076 | /*-------------ring management:----------------------*/ |
---|
2077 | |
---|
2078 | // resets the pFDeg and pLDeg: if pLDeg is not given, it is |
---|
2079 | // set to currRing->pLDegOrig, i.e. to the respective LDegProc which |
---|
2080 | // only uses pFDeg (and not pDeg, or pTotalDegree, etc). |
---|
2081 | // If you use this, make sure your procs does not make any assumptions |
---|
2082 | // on ordering and/or OrdIndex -- otherwise they might return wrong results |
---|
2083 | // on strat->tailRing |
---|
2084 | void pSetDegProcs(ring r, pFDegProc new_FDeg, pLDegProc new_lDeg = NULL); |
---|
2085 | // restores pFDeg and pLDeg: |
---|
2086 | void pRestoreDegProcs(ring r, pFDegProc old_FDeg, pLDegProc old_lDeg); |
---|
2087 | |
---|
2088 | /*-------------pComp for syzygies:-------------------*/ |
---|
2089 | void p_SetModDeg(intvec *w, ring r); |
---|
2090 | |
---|
2091 | /*------------ Jet ----------------------------------*/ |
---|
2092 | poly pp_Jet(poly p, int m, const ring R); |
---|
2093 | poly p_Jet(poly p, int m,const ring R); |
---|
2094 | poly pp_JetW(poly p, int m, int *w, const ring R); |
---|
2095 | poly p_JetW(poly p, int m, int *w, const ring R); |
---|
2096 | |
---|
2097 | poly n_PermNumber(const number z, const int *par_perm, const int OldPar, const ring src, const ring dst); |
---|
2098 | |
---|
2099 | poly p_PermPoly (poly p, const int * perm,const ring OldRing, const ring dst, |
---|
2100 | nMapFunc nMap, const int *par_perm=NULL, int OldPar=0, |
---|
2101 | BOOLEAN use_mult=FALSE); |
---|
2102 | |
---|
2103 | /*----------------------------------------------------*/ |
---|
2104 | poly p_Series(int n,poly p,poly u, intvec *w, const ring R); |
---|
2105 | |
---|
2106 | /*----------------------------------------------------*/ |
---|
2107 | int p_Var(poly mi, const ring r); |
---|
2108 | /// the minimal index of used variables - 1 |
---|
2109 | int p_LowVar (poly p, const ring r); |
---|
2110 | |
---|
2111 | /*----------------------------------------------------*/ |
---|
2112 | /// shifts components of the vector p by i |
---|
2113 | void p_Shift (poly * p,int i, const ring r); |
---|
2114 | /*----------------------------------------------------*/ |
---|
2115 | |
---|
2116 | int p_Compare(const poly a, const poly b, const ring R); |
---|
2117 | |
---|
2118 | /// polynomial gcd for f=mon |
---|
2119 | poly p_GcdMon(poly f, poly g, const ring r); |
---|
2120 | |
---|
2121 | /// divide polynomial by monomial |
---|
2122 | poly p_Div_mm(poly p, const poly m, const ring r); |
---|
2123 | |
---|
2124 | |
---|
2125 | /// max exponent of variable x_i in p |
---|
2126 | int p_MaxExpPerVar(poly p, int i, const ring r); |
---|
2127 | #endif // P_POLYS_H |
---|
2128 | |
---|