1 | #ifndef POLYS_IMPL_H |
---|
2 | #define POLYS_IMPL_H |
---|
3 | /**************************************** |
---|
4 | * Computer Algebra System SINGULAR * |
---|
5 | ****************************************/ |
---|
6 | /* $Id: polys-impl.h,v 1.16 1998-01-27 18:51:21 Singular Exp $ */ |
---|
7 | |
---|
8 | /*************************************************************** |
---|
9 | * |
---|
10 | * File: polys-impl.h |
---|
11 | * Purpose: low-level definition and declarations for polys |
---|
12 | * |
---|
13 | * If you touch anything here, you better know what you are doing. |
---|
14 | * What is here should not be used directly from other routines -- the |
---|
15 | * encapsulations in polys.h should be used, instead. |
---|
16 | * |
---|
17 | ***************************************************************/ |
---|
18 | #include "tok.h" |
---|
19 | #include "structs.h" |
---|
20 | #include "mmemory.h" |
---|
21 | #include "binom.h" |
---|
22 | |
---|
23 | /*************************************************************** |
---|
24 | * |
---|
25 | * definition of the poly structure and its fields |
---|
26 | * |
---|
27 | ***************************************************************/ |
---|
28 | |
---|
29 | // EXPONENT_TYPE is determined by configure und defined in mod2.h |
---|
30 | typedef EXPONENT_TYPE Exponent_t; |
---|
31 | |
---|
32 | #define VARS (100) /*max. number of variables as constant*/ |
---|
33 | |
---|
34 | typedef Exponent_t monomial[VARS + 1]; |
---|
35 | typedef Exponent_t* Exponent_pt; |
---|
36 | |
---|
37 | typedef long Order_t; |
---|
38 | // make sure that exp is aligned |
---|
39 | struct spolyrec |
---|
40 | { |
---|
41 | poly next; |
---|
42 | number coef; |
---|
43 | Order_t Order; |
---|
44 | #ifdef TEST_MAC_DEBUG |
---|
45 | Order_t MOrder; |
---|
46 | #endif |
---|
47 | monomial exp; |
---|
48 | }; |
---|
49 | |
---|
50 | /*************************************************************** |
---|
51 | * MACROS CONTROLING MONOMIAL COMPARIONS: |
---|
52 | |
---|
53 | * COMP_TRADITIONAL |
---|
54 | Keeps the traditional comparison routines |
---|
55 | defined -- needed as long as their might be comparisons with |
---|
56 | negativ components. |
---|
57 | All the traditional routines are prefixed by t_ |
---|
58 | |
---|
59 | * COMP_FAST |
---|
60 | Implements monomial operations using the fast vector |
---|
61 | techniques and several other extensions which go along with that. |
---|
62 | Undefine in case there are problems. |
---|
63 | |
---|
64 | * COMP_STATISTIC |
---|
65 | Provides several routines for accumulating statistics on monomial |
---|
66 | comparisons and divisibility tests |
---|
67 | |
---|
68 | * COMP_DEBUG |
---|
69 | Turns on debugging of COMP_FAST by comparing the results of fast |
---|
70 | comparison with traditional comparison |
---|
71 | |
---|
72 | * COMP_NO_EXP_VECTOR_OPS |
---|
73 | Like COMP_FAST, except that it turns off "vector techniques" of |
---|
74 | monomial operations, i.e. does everything exponent-wise. |
---|
75 | ***************************************************************/ |
---|
76 | // #define COMP_FAST |
---|
77 | // #define COMP_DEBUG |
---|
78 | // #define COMP_NO_EXP_VECTOR_OPS |
---|
79 | #define COMP_TRADITIONAL |
---|
80 | |
---|
81 | #if defined(COMP_NO_EXP_VECTOR_OPS) && ! defined(COMP_FAST) |
---|
82 | #define COMP_FAST |
---|
83 | #endif |
---|
84 | |
---|
85 | #if defined(COMP_FAST) && ! defined(NDEBUG) |
---|
86 | #define COMP_DEBUG |
---|
87 | #endif |
---|
88 | |
---|
89 | // some relations between these flags |
---|
90 | #ifdef COMP_DEBUG |
---|
91 | #define COMP_TRADITIONAL |
---|
92 | #define COMP_FAST |
---|
93 | #undef COMP_PROFILE |
---|
94 | #undef COMP_STATISTICS |
---|
95 | #endif // COMP_DEBUG |
---|
96 | |
---|
97 | #ifdef COMP_STATISTICS |
---|
98 | #undef COMP_FAST |
---|
99 | #endif // COMP_STATISTICS |
---|
100 | |
---|
101 | // for the time being COMP_TRADITIONAL always has to be defined, since |
---|
102 | // traditional routines are needed in spolys.cc -- monomials with |
---|
103 | // negative exponents are compared there! |
---|
104 | #define COMP_TRADITIONAL |
---|
105 | |
---|
106 | /*************************************************************** |
---|
107 | * |
---|
108 | * variables used for storage management and monomial traversions |
---|
109 | * |
---|
110 | ***************************************************************/ |
---|
111 | |
---|
112 | // size of poly without exponents |
---|
113 | #ifdef TEST_MAC_DEBUG |
---|
114 | #define POLYSIZE (sizeof(poly) + sizeof(number) + 2*sizeof(Order_t)) |
---|
115 | #else |
---|
116 | #define POLYSIZE (sizeof(poly) + sizeof(number) + sizeof(Order_t)) |
---|
117 | #endif |
---|
118 | #define POLYSIZEW (POLYSIZE / sizeof(long)) |
---|
119 | // number of Variables |
---|
120 | extern int pVariables; |
---|
121 | // size of a monom in bytes - always a multiple of sizeof(void*) |
---|
122 | extern int pMonomSize; |
---|
123 | // size of a monom in units of sizeof(void*) -- i.e. in words |
---|
124 | extern int pMonomSizeW; |
---|
125 | #ifdef COMP_FAST |
---|
126 | // Ceiling((pVariables+1) / sizeof(void*)) == length of exp-vector in words |
---|
127 | extern int pVariables1W; |
---|
128 | // Ceiling((pVariables) / sizeof(void*)) |
---|
129 | extern int pVariablesW; |
---|
130 | extern int pCompIndex; |
---|
131 | extern int pVarOffset; |
---|
132 | extern int pVarLowIndex; |
---|
133 | extern int pVarHighIndex; |
---|
134 | #else |
---|
135 | #define pCompIndex 0 |
---|
136 | #endif |
---|
137 | |
---|
138 | /*************************************************************** |
---|
139 | * |
---|
140 | * Primitives for determening/setting the way exponents are arranged |
---|
141 | * |
---|
142 | ***************************************************************/ |
---|
143 | // And here is how we determine the way exponents are stored: |
---|
144 | // There are the following four possibilities: |
---|
145 | // |
---|
146 | // |
---|
147 | // BIGENDIAN -- lex order |
---|
148 | // e_1, e_2, ... , e_n,..,comp : pVarOffset = -1, |
---|
149 | // pCompIndex = pVariables + #(..) |
---|
150 | // pVarLowIndex = 0, |
---|
151 | // pVarHighIndex = pVariables-1 |
---|
152 | // BIGENDIAN -- rev lex order |
---|
153 | // e_n, ... , e_2, e_1,..,comp : pVarOffset = pVariables, |
---|
154 | // pCompIndex = pVariables + #(..) |
---|
155 | // pVarLowIndex = 0, |
---|
156 | // pVarHighIndex = pVariables-1 |
---|
157 | // LITTLEENDIAN -- rev lex order |
---|
158 | // comp,.., e_1, e_2, ... , e_n : pVarOffset = #(..), |
---|
159 | // pCompIndex = 0, |
---|
160 | // pVarLowIndex = 1 + #(..), |
---|
161 | // pVarHighIndex = #(..) + pVariables |
---|
162 | // LITTLEENDIAN -- lex order |
---|
163 | // comp,..,e_n, .... , e_2, e_1 : pVarOffset = pVariables + 1 + #(..) |
---|
164 | // pCompIndex = 0 |
---|
165 | // pVarLowIndex = 1 + #(..) |
---|
166 | // pVarHighIndex = #(..) + pVariables |
---|
167 | // |
---|
168 | // Furthermore, the size of the exponent vector is always a multiple |
---|
169 | // of the word size -- "empty exponents" (exactly #(..) ones) are |
---|
170 | // filled in between comp and first/last exponent -- i.e. comp and |
---|
171 | // first/last exponent might not be next to each other |
---|
172 | |
---|
173 | #ifdef COMP_FAST |
---|
174 | |
---|
175 | #ifdef WORDS_BIGENDIAN |
---|
176 | |
---|
177 | #define _pHasReverseExp (pVarOffset != -1) |
---|
178 | #define _pExpIndex(i) \ |
---|
179 | (pVarOffset == -1 ? (i) - 1 : pVarOffset - (i)) |
---|
180 | #define _pRingExpIndex(r, i) \ |
---|
181 | ((r)->VarOffset == -1 ? (i) - 1 : (r)->VarOffset - (i)) |
---|
182 | |
---|
183 | #else // ! WORDS_BIGENDIAN |
---|
184 | |
---|
185 | #define _pHasReverseExp (pVarOffset > (SIZEOF_LONG / SIZEOF_EXPONENT) - 1) |
---|
186 | #define _pExpIndex(i) \ |
---|
187 | (pVarOffset > (SIZEOF_LONG / SIZEOF_EXPONENT) - 1 ? \ |
---|
188 | pVarOffset - (i) : pVarOffset + (i)) |
---|
189 | #define _pRingExpIndex(r, i) \ |
---|
190 | ((r)->VarOffset > (SIZEOF_LONG / SIZEOF_EXPONENT) - 1 ? \ |
---|
191 | (r)->VarOffset - (i) : (r)->VarOffset + (i)) |
---|
192 | |
---|
193 | #endif // WORDS_BIGENDIAN |
---|
194 | |
---|
195 | inline void pGetVarIndicies_Lex(int nvars, |
---|
196 | int &VarOffset, int &VarCompIndex, |
---|
197 | int &VarLowIndex, int &VarHighIndex) |
---|
198 | { |
---|
199 | long temp = (nvars+1)*sizeof(Exponent_t); |
---|
200 | if ((temp % sizeof(long)) == 0) |
---|
201 | temp = temp / sizeof(long); |
---|
202 | else |
---|
203 | temp = (temp / sizeof(long)) + 1; // now temp == nvars1W |
---|
204 | #ifdef WORDS_BIGENDIAN |
---|
205 | VarCompIndex = temp * sizeof(long)/sizeof(Exponent_t) - 1; |
---|
206 | VarOffset = -1; |
---|
207 | VarLowIndex = 0; |
---|
208 | VarHighIndex = nvars - 1; |
---|
209 | #else // ! WORDS_BIGENDIAN |
---|
210 | VarHighIndex = temp * sizeof(long)/sizeof(Exponent_t) - 1; |
---|
211 | VarCompIndex = 0; |
---|
212 | VarOffset = VarHighIndex + 1; |
---|
213 | VarLowIndex = VarOffset - nvars; |
---|
214 | #endif // WORDS_BIGENDIAN |
---|
215 | } |
---|
216 | #define pSetVarIndicies_Lex(nvars) \ |
---|
217 | pGetVarIndicies_Lex(nvars,pVarOffset,pCompIndex,pVarLowIndex,pVarHighIndex) |
---|
218 | |
---|
219 | inline void pGetVarIndicies_RevLex(int nvars, |
---|
220 | int &VarOffset, int &VarCompIndex, |
---|
221 | int &VarLowIndex, int &VarHighIndex) |
---|
222 | { |
---|
223 | long temp = (nvars+1)*sizeof(Exponent_t); |
---|
224 | if ((temp % sizeof(long)) == 0) |
---|
225 | temp = temp / sizeof(long); |
---|
226 | else |
---|
227 | temp = (temp / sizeof(long)) + 1; |
---|
228 | #ifdef WORDS_BIGENDIAN |
---|
229 | VarCompIndex = temp * sizeof(long)/sizeof(Exponent_t) - 1; |
---|
230 | VarOffset = nvars; |
---|
231 | VarLowIndex = 0; |
---|
232 | VarHighIndex = nvars-1; |
---|
233 | #else // ! WORDS_BIGENDIAN |
---|
234 | // comp, ..., e_1, e_2, ... , e_n |
---|
235 | VarHighIndex = temp * sizeof(long)/sizeof(Exponent_t) - 1; |
---|
236 | VarCompIndex = 0; |
---|
237 | VarLowIndex = VarHighIndex - nvars + 1; |
---|
238 | VarOffset = VarLowIndex - 1; |
---|
239 | #endif // WORDS_BIGENDIAN |
---|
240 | } |
---|
241 | #define pSetVarIndicies_RevLex(nvars) \ |
---|
242 | pGetVarIndicies_RevLex(nvars,pVarOffset,pCompIndex,pVarLowIndex,pVarHighIndex) |
---|
243 | |
---|
244 | // The default settings: |
---|
245 | inline void pGetVarIndicies(int nvars, |
---|
246 | int &VarOffset, int &VarCompIndex, |
---|
247 | int &VarLowIndex, int &VarHighIndex) |
---|
248 | { |
---|
249 | pGetVarIndicies_Lex(nvars,VarOffset,VarCompIndex,VarLowIndex,VarHighIndex); |
---|
250 | } |
---|
251 | |
---|
252 | // gets var indicies w.r.t. the ring r |
---|
253 | extern void pGetVarIndicies(ring r, |
---|
254 | int &VarOffset, int &VarCompIndex, |
---|
255 | int &VarLowIndex, int &VarHighIndex); |
---|
256 | |
---|
257 | #define pSetVarIndicies(nvars) \ |
---|
258 | pGetVarIndicies(nvars, pVarOffset, pCompIndex, pVarLowIndex, pVarHighIndex) |
---|
259 | |
---|
260 | |
---|
261 | #else // ! COMP_FAST |
---|
262 | #define _pExpIndex(i) (i) |
---|
263 | #define _pRingExpIndex(r,i) (i) |
---|
264 | #endif // COMP_FAST |
---|
265 | |
---|
266 | /*************************************************************** |
---|
267 | * |
---|
268 | * Primitives for accessing and seeting fields of a poly |
---|
269 | * |
---|
270 | ***************************************************************/ |
---|
271 | #define _pNext(p) ((p)->next) |
---|
272 | #define _pIter(p) ((p) = (p)->next) |
---|
273 | |
---|
274 | #define _pGetCoeff(p) ((p)->coef) |
---|
275 | #define _pSetCoeff(p,n) {nDelete(&((p)->coef));(p)->coef=n;} |
---|
276 | #define _pSetCoeff0(p,n) (p)->coef=n |
---|
277 | |
---|
278 | #ifndef TEST_MAC_DEBUG |
---|
279 | #define _pGetOrder(p) ((p)->Order) |
---|
280 | #else |
---|
281 | #define _pGetOrder(p) ((p)->MOrder) |
---|
282 | #endif |
---|
283 | |
---|
284 | #if defined(PDEBUG) && PDEBUG != 0 |
---|
285 | extern Exponent_t pPDSetExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
286 | extern Exponent_t pPDGetExp(poly p, int v, char* f, int l); |
---|
287 | extern Exponent_t pPDIncrExp(poly p, int v, char* f, int l); |
---|
288 | extern Exponent_t pPDDecrExp(poly p, int v, char* f, int l); |
---|
289 | extern Exponent_t pPDAddExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
290 | extern Exponent_t pPDMultExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
291 | extern Exponent_t pPDSubExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
292 | |
---|
293 | extern Exponent_t pPDRingSetExp(ring r,poly p,int v,Exponent_t e,char* f,int l); |
---|
294 | extern Exponent_t pPDRingGetExp(ring r,poly p, int v, char* f, int l); |
---|
295 | |
---|
296 | #define _pSetExp(p,v,e) pPDSetExp(p,v,e,__FILE__,__LINE__) |
---|
297 | #define _pGetExp(p,v) pPDGetExp(p,v,__FILE__,__LINE__) |
---|
298 | #define _pIncrExp(p,v) pPDIncrExp(p,v,__FILE__,__LINE__) |
---|
299 | #define _pDecrExp(p,v) pPDDecrExp(p,v,__FILE__,__LINE__) |
---|
300 | #define _pAddExp(p,i,v) pPDAddExp(p,i,v,__FILE__,__LINE__) |
---|
301 | #define _pSubExp(p,i,v) pPDSubExp(p,i,v,__FILE__,__LINE__) |
---|
302 | #define _pMultExp(p,i,v) pPDMultExp(p,i,v,__FILE__,__LINE__) |
---|
303 | |
---|
304 | #define _pRingSetExp(r,p,v,e) pPDRingSetExp(r,p,v,e,__FILE__,__LINE__) |
---|
305 | #define _pRingGetExp(r,p,v) pPDRingGetExp(r,p,v,__FILE__,__LINE__) |
---|
306 | |
---|
307 | #else // ! (defined(PDEBUG) && PDEBUG != 0) |
---|
308 | |
---|
309 | #define _pSetExp(p,v,e) (p)->exp[_pExpIndex(v)]=(e) |
---|
310 | #define _pGetExp(p,v) (p)->exp[_pExpIndex(v)] |
---|
311 | #define _pIncrExp(p,v) ((p)->exp[_pExpIndex(v)])++ |
---|
312 | #define _pDecrExp(p,v) ((p)->exp[_pExpIndex(v)])-- |
---|
313 | #define _pAddExp(p,i,v) ((p)->exp[_pExpIndex(i)]) += (v) |
---|
314 | #define _pSubExp(p,i,v) ((p)->exp[_pExpIndex(i)]) -= (v) |
---|
315 | #define _pMultExp(p,i,v) ((p)->exp[_pExpIndex(i)]) *= (v) |
---|
316 | |
---|
317 | #define _pRingSetExp(r,p,v,e) (p)->exp[_pRingExpIndex(r,v)]=(e) |
---|
318 | #define _pRingGetExp(r,p,v) (p)->exp[_pRingExpIndex(r,v)] |
---|
319 | |
---|
320 | #endif // defined(PDEBUG) && PDEBUG != 0 |
---|
321 | |
---|
322 | inline Exponent_t _pGetExpSum(poly p1, poly p2, int i) |
---|
323 | { |
---|
324 | int index = _pExpIndex(i); |
---|
325 | return p1->exp[index] + p2->exp[index]; |
---|
326 | } |
---|
327 | inline Exponent_t _pGetExpDiff(poly p1, poly p2, int i) |
---|
328 | { |
---|
329 | int index = _pExpIndex(i); |
---|
330 | return p1->exp[index] - p2->exp[index]; |
---|
331 | } |
---|
332 | |
---|
333 | #define _pSetComp(p,k) (p)->exp[pCompIndex] = (k) |
---|
334 | #define _pGetComp(p) (p)->exp[pCompIndex] |
---|
335 | #define _pIncrComp(p) (p)->exp[pCompIndex]++ |
---|
336 | #define _pDecrComp(p) (p)->exp[pCompIndex]-- |
---|
337 | #define _pAddComp(p,v) (p)->exp[pCompIndex] += (v) |
---|
338 | #define _pSubComp(p,v) (p)->exp[pCompIndex] -= (v) |
---|
339 | |
---|
340 | #ifdef COMP_FAST |
---|
341 | #define _pRingSetComp(r,p,k) (p)->exp[(r)->CompIndex] = (k) |
---|
342 | #define _pRingGetComp(r,p) (p)->exp[(r)->CompIndex] |
---|
343 | #else |
---|
344 | #define _pRingSetComp(r,p,k) _pSetComp(p,k) |
---|
345 | #define _pRingGetComp(r,p) _pGetComp(p) |
---|
346 | #endif |
---|
347 | |
---|
348 | |
---|
349 | /*************************************************************** |
---|
350 | * |
---|
351 | * Storage Managament Routines |
---|
352 | * |
---|
353 | ***************************************************************/ |
---|
354 | #ifdef PDEBUG |
---|
355 | |
---|
356 | poly pDBNew(char *f, int l); |
---|
357 | poly pDBInit(char * f,int l); |
---|
358 | |
---|
359 | void pDBDelete(poly * a, char * f, int l); |
---|
360 | void pDBDelete1(poly * a, char * f, int l); |
---|
361 | void pDBFree1(poly a, char * f, int l); |
---|
362 | |
---|
363 | poly pDBCopy(poly a, char *f, int l); |
---|
364 | poly pDBCopy1(poly a, char *f, int l); |
---|
365 | poly pDBHead(poly a, char *f, int l); |
---|
366 | poly pDBHead0(poly a, char *f, int l); |
---|
367 | poly pDBFetchCopy(ring r, poly a, char *f, int l); |
---|
368 | |
---|
369 | void pDBDelete(poly * a, char * f, int l); |
---|
370 | void pDBDelete1(poly * a, char * f, int l); |
---|
371 | void pDBFree1(poly a, char * f, int l); |
---|
372 | |
---|
373 | #define _pNew() pDBNew(__FILE__,__LINE__) |
---|
374 | #define _pInit() pDBInit(__FILE__,__LINE__) |
---|
375 | |
---|
376 | #define _pDelete(a) pDBDelete((a),__FILE__,__LINE__) |
---|
377 | #define _pDelete1(a) pDBDelete1((a),__FILE__,__LINE__) |
---|
378 | #define _pFree1(a) \ |
---|
379 | do \ |
---|
380 | { \ |
---|
381 | pDBFree1((a),__FILE__,__LINE__); \ |
---|
382 | (a)=NULL; \ |
---|
383 | } \ |
---|
384 | while(0) |
---|
385 | |
---|
386 | #define _pCopy(A) pDBCopy(A,__FILE__,__LINE__) |
---|
387 | #define _pCopy1(A) pDBCopy1(A, __FILE__,__LINE__) |
---|
388 | #define _pHead(A) pDBHead(A,__FILE__,__LINE__) |
---|
389 | #define _pHead0(A) pDBHead0(A, __FILE__,__LINE__) |
---|
390 | #ifdef COMP_FAST |
---|
391 | #define _pFetchCopy(r,A) pDBFetchCopy(r, A,__FILE__,__LINE__) |
---|
392 | #else |
---|
393 | #define _pFetchCopy(r,p) pOrdPolyInsertSetm(pCopy(p)) |
---|
394 | #endif |
---|
395 | |
---|
396 | #else // ! PDEBUG |
---|
397 | |
---|
398 | #ifdef MDEBUG |
---|
399 | #define _pNew() (poly) mmDBAllocSpecialized(__FILE__,__LINE__) |
---|
400 | #else |
---|
401 | #define _pNew() (poly) mmAllocSpecialized() |
---|
402 | #endif |
---|
403 | // #define _pNew() _pInit() |
---|
404 | |
---|
405 | #include <string.h> |
---|
406 | |
---|
407 | inline poly _pInit(void) |
---|
408 | { |
---|
409 | #ifdef MDEBUG |
---|
410 | poly p=(poly)mmDBAllocSpecialized(__FILE__,__LINE__); |
---|
411 | #else |
---|
412 | poly p=(poly)mmAllocSpecialized(); |
---|
413 | #endif |
---|
414 | memset(p,0, pMonomSize); |
---|
415 | return p; |
---|
416 | } |
---|
417 | |
---|
418 | extern void _pDelete(poly * a); |
---|
419 | extern void _pDelete1(poly * a); |
---|
420 | #ifdef MDEBUG |
---|
421 | #define _pFree1(a) mmDBFreeSpecialized((ADDRESS)a,__FILE__,__LINE__) |
---|
422 | #else |
---|
423 | #define _pFree1(a) mmFreeSpecialized((ADDRESS)a) |
---|
424 | #endif |
---|
425 | |
---|
426 | extern poly _pCopy(poly a); |
---|
427 | extern poly _pCopy1(poly a); |
---|
428 | extern poly _pHead(poly a); |
---|
429 | extern poly _pHead0(poly a); |
---|
430 | #ifdef COMP_FAST |
---|
431 | extern poly _pFetchCopy(ring r,poly a); |
---|
432 | #else |
---|
433 | #define _pFetchCopy(r,p) pOrdPolyInsertSetm(pCopy(p)) |
---|
434 | #endif |
---|
435 | |
---|
436 | #endif // PDEBUG |
---|
437 | |
---|
438 | #define _pCopy2(p1, p2) memcpyW(p1, p2, pMonomSizeW) |
---|
439 | |
---|
440 | // Here is a handy Macro which disables inlining when run with |
---|
441 | // profiling and enables it otherwise |
---|
442 | |
---|
443 | #ifdef DO_PROFILE |
---|
444 | |
---|
445 | #ifndef POLYS_IMPL_CC |
---|
446 | #define DECLARE(type, arglist) type arglist; \ |
---|
447 | static type dummy_##arglist |
---|
448 | #else |
---|
449 | #define DECLARE(type, arglist) type arglist |
---|
450 | #endif // POLYS_IMPL_CC |
---|
451 | |
---|
452 | #else //! DO_PROFILE |
---|
453 | |
---|
454 | #define DECLARE(type, arglist ) inline type arglist |
---|
455 | |
---|
456 | #endif // DO_PROFILE |
---|
457 | |
---|
458 | |
---|
459 | /*************************************************************** |
---|
460 | * |
---|
461 | * Routines which work on vectors instead of single exponents |
---|
462 | * |
---|
463 | ***************************************************************/ |
---|
464 | |
---|
465 | #ifdef COMP_FAST |
---|
466 | |
---|
467 | // nice declaration isn't it ?? |
---|
468 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
469 | #define pMonAddFast(p1, p2) pDBMonAddFast(p1, p2, __FILE__, __LINE__) |
---|
470 | extern void pDBMonAddFast(poly p1, poly p2, char* f, int l); |
---|
471 | inline void _pMonAddFast(poly p1, poly p2) |
---|
472 | #else |
---|
473 | DECLARE(void, pMonAddFast(poly p1, poly p2)) |
---|
474 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
475 | { |
---|
476 | // OK -- this might be the only place where we are kind of quick and |
---|
477 | // dirty: the following only works correctly if all exponents are |
---|
478 | // positive and the sum of two exponents does not exceed |
---|
479 | // EXPONENT_MAX |
---|
480 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
481 | Exponent_t c2 = _pGetComp(p2); |
---|
482 | int i = pVariables1W; |
---|
483 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
484 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
485 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
486 | _pSetComp(p2, 0); |
---|
487 | #else |
---|
488 | int i = pVariables; |
---|
489 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
490 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
491 | #endif |
---|
492 | |
---|
493 | for (;;) |
---|
494 | { |
---|
495 | *s1 += *s2; |
---|
496 | i--; |
---|
497 | if (i==0) break; |
---|
498 | s1++; |
---|
499 | s2++; |
---|
500 | } |
---|
501 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
502 | // reset comp of p2 |
---|
503 | _pSetComp(p2, c2); |
---|
504 | #endif |
---|
505 | #ifdef TEST_MAC_ORDER |
---|
506 | if (bNoAdd) bSetm(p1);else |
---|
507 | #endif |
---|
508 | _pGetOrder(p1) += _pGetOrder(p2); |
---|
509 | } |
---|
510 | |
---|
511 | // Makes p1 a copy of p2 and adds on exponets of p3 |
---|
512 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
513 | #define _pCopyAddFast(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
514 | extern void pDBCopyAddFast(poly p1, poly p2, poly p3, char* f, int l); |
---|
515 | inline void __pCopyAddFast(poly p1, poly p2, poly p3) |
---|
516 | #else |
---|
517 | DECLARE(void, _pCopyAddFast(poly p1, poly p2, poly p3)) |
---|
518 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
519 | { |
---|
520 | p1->next = p2->next; |
---|
521 | p1->coef = p2->coef; |
---|
522 | // memset(p1, 0, pMonomSize); |
---|
523 | |
---|
524 | //memset(p1,0,pMonomSize); |
---|
525 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
526 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
527 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
528 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
529 | const unsigned long* const ub = s3 + pVariables1W; |
---|
530 | #else |
---|
531 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
532 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
533 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
534 | const Exponent_t* const ub = s3 + pVariables; |
---|
535 | // need to zero the "fill in" slots (i.e., empty exponents) |
---|
536 | #ifdef WORDS_BIGENDIAN |
---|
537 | *((unsigned long *) ((unsigned long*) p1) + pMonomSizeW -1) = 0; |
---|
538 | #else |
---|
539 | *((unsigned long *) p1->exp) = 0; |
---|
540 | #endif |
---|
541 | #endif |
---|
542 | |
---|
543 | for (;;) |
---|
544 | { |
---|
545 | *s1 = *s2 + *s3; |
---|
546 | s3++; |
---|
547 | if (s3 == ub) break; |
---|
548 | s1++; |
---|
549 | s2++; |
---|
550 | } |
---|
551 | // we first are supposed to do a copy from p2 to p1 -- therefore, |
---|
552 | // component of p1 is set to comp of p2 |
---|
553 | _pSetComp(p1, _pGetComp(p2)); |
---|
554 | #ifdef TEST_MAC_ORDER |
---|
555 | if (bNoAdd) bSetm(p1);else |
---|
556 | #endif |
---|
557 | _pGetOrder(p1) = _pGetOrder(p2) + _pGetOrder(p3); |
---|
558 | } |
---|
559 | |
---|
560 | // Similar to pCopyAddFast, except that we assume that the component |
---|
561 | // of p2 and p3 is zero component |
---|
562 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
563 | #define _pCopyAddFast1(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
564 | extern void pDBCopyAddFast(poly p1, poly p2, poly p3, char* f, int l); |
---|
565 | inline void __pCopyAddFast1(poly p1, poly p2, poly p3) |
---|
566 | #else |
---|
567 | DECLARE(void, _pCopyAddFast1(poly p1, poly p2, poly p3)) |
---|
568 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
569 | { |
---|
570 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
571 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
572 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
573 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
574 | const unsigned long* const ub = s3 + pVariables1W; |
---|
575 | #else |
---|
576 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
577 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
578 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
579 | const Exponent_t* const ub = s3 + pVariables; |
---|
580 | #ifdef WORDS_BIGENDIAN |
---|
581 | *((unsigned long *) ((unsigned long*) p1) + pMonomSizeW -1) = 0; |
---|
582 | #else |
---|
583 | *((unsigned long *) p1->exp) = 0; |
---|
584 | #endif |
---|
585 | #endif |
---|
586 | |
---|
587 | for (;;) |
---|
588 | { |
---|
589 | *s1 = *s2 + *s3; |
---|
590 | s3++; |
---|
591 | if (s3 == ub) break; |
---|
592 | s1++; |
---|
593 | s2++; |
---|
594 | } |
---|
595 | #ifdef TEST_MAC_ORDER |
---|
596 | if (bNoAdd) bSetm(p1);else |
---|
597 | #endif |
---|
598 | _pGetOrder(p1) = _pGetOrder(p2) + _pGetOrder(p3); |
---|
599 | } |
---|
600 | |
---|
601 | |
---|
602 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
603 | |
---|
604 | #if SIZEOF_LONG == 4 |
---|
605 | |
---|
606 | #if SIZEOF_EXPONENT == 1 |
---|
607 | #define P_DIV_MASK 0x80808080 |
---|
608 | #else // SIZEOF_EXPONENT == 2 |
---|
609 | #define P_DIV_MASK 0x80008000 |
---|
610 | #endif |
---|
611 | |
---|
612 | #else // SIZEOF_LONG == 8 |
---|
613 | |
---|
614 | #if SIZEOF_EXPONENT == 1 |
---|
615 | #define P_DIV_MASK 0x8080808080808080 |
---|
616 | #elif SIZEOF_EXPONENT == 2 |
---|
617 | #define P_DIV_MASK 0x8000800080008000 |
---|
618 | #else // SIZEOF_EXPONENT == 4 |
---|
619 | #define P_DIV_MASK 0x8000000080000000 |
---|
620 | #endif |
---|
621 | |
---|
622 | #endif |
---|
623 | |
---|
624 | DECLARE(BOOLEAN, __pDivisibleBy(poly a, poly b)) |
---|
625 | { |
---|
626 | const unsigned long* s1; |
---|
627 | const unsigned long* s2; |
---|
628 | const unsigned long* lb; |
---|
629 | |
---|
630 | #ifdef WORDS_BIGENDIAN |
---|
631 | lb = (unsigned long*) &(a->exp[0]); |
---|
632 | if (pVariables & ((SIZEOF_LONG / SIZEOF_EXPONENT) - 1)) |
---|
633 | { |
---|
634 | // now pVariables == pVariables1W, i.e. there are exponents in the |
---|
635 | // "first" word of exponentvector |
---|
636 | s1 = ((unsigned long*) a) + pMonomSizeW -1; |
---|
637 | s2 = ((unsigned long*) b) + pMonomSizeW -1; |
---|
638 | } |
---|
639 | else |
---|
640 | { |
---|
641 | // first exponent word has only component as significant field -- |
---|
642 | // Hence, do not bother |
---|
643 | s1 = ((unsigned long*) a) + pMonomSizeW -2; |
---|
644 | s2 = ((unsigned long*) b) + pMonomSizeW -2; |
---|
645 | } |
---|
646 | #else // !WORDS_BIGENDIAN |
---|
647 | lb = ((unsigned long*) a) + pMonomSizeW; |
---|
648 | if (pVariables & ((SIZEOF_LONG / SIZEOF_EXPONENT) - 1)) |
---|
649 | { |
---|
650 | s1 = (unsigned long*) &(a->exp[0]); |
---|
651 | s2 = (unsigned long*) &(b->exp[0]); |
---|
652 | } |
---|
653 | else |
---|
654 | { |
---|
655 | s1 = (unsigned long*) &(a->exp[0]) + 1; |
---|
656 | s2 = (unsigned long*) &(b->exp[0]) + 1; |
---|
657 | } |
---|
658 | #endif |
---|
659 | for (;;) |
---|
660 | { |
---|
661 | // O.K. -- and now comes a bit of magic. The following _really_ |
---|
662 | // works. Think about it! If you can prove it, please tell me, for I |
---|
663 | // did not bother to prove it formally (Hint: We can assume that our |
---|
664 | // exponents are always positive). |
---|
665 | if ((*s2 - *s1) & P_DIV_MASK) return FALSE; |
---|
666 | #ifdef WORDS_BIGENDIAN |
---|
667 | if (s1 == lb) return TRUE; |
---|
668 | s1--; |
---|
669 | s2--; |
---|
670 | #else |
---|
671 | s1++; |
---|
672 | if (s1 == lb) return TRUE; |
---|
673 | s2++; |
---|
674 | #endif |
---|
675 | } |
---|
676 | } |
---|
677 | |
---|
678 | #else // ! COMP_NO_EXP_VECTOR_OPS |
---|
679 | |
---|
680 | DECLARE(BOOLEAN, __pDivisibleBy(poly a, poly b)) |
---|
681 | { |
---|
682 | #ifdef WORDS_BIGENDIAN |
---|
683 | const Exponent_t* s1 = &(a->exp[pVarHighIndex]); |
---|
684 | const Exponent_t* s2 = &(b->exp[pVarHighIndex]); |
---|
685 | const Exponent_t* lb = s1 - pVariables; |
---|
686 | |
---|
687 | for (;;) |
---|
688 | { |
---|
689 | if (*s1 > *s2) return FALSE; |
---|
690 | s1--; |
---|
691 | if (s1 == lb) return TRUE; |
---|
692 | s2--; |
---|
693 | } |
---|
694 | |
---|
695 | #else // !WORDS_BIGENDIAN |
---|
696 | const Exponent_t* s1 = &(a->exp[pVarLowIndex]); |
---|
697 | const Exponent_t* s2 = &(b->exp[pVarLowIndex]); |
---|
698 | const Exponent_t* lb = s1 + pVariables; |
---|
699 | |
---|
700 | for (;;) |
---|
701 | { |
---|
702 | if (*s1 > *s2) return FALSE; |
---|
703 | s1++; |
---|
704 | if (s1 == lb) return TRUE; |
---|
705 | s2++; |
---|
706 | } |
---|
707 | #endif // WORDS_BIGENDIAN |
---|
708 | } |
---|
709 | |
---|
710 | #endif // COMP_NO_EXP_VECTOR_OPS |
---|
711 | |
---|
712 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
713 | #define _pDivisibleBy(a,b) pDBDivisibleBy(a, b, __FILE__, __LINE__) |
---|
714 | extern BOOLEAN pDBDivisibleBy(poly p1, poly p2, char* f, int l); |
---|
715 | inline BOOLEAN _pDivisibleBy_orig(poly a, poly b) |
---|
716 | #else |
---|
717 | inline BOOLEAN _pDivisibleBy(poly a, poly b) |
---|
718 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
719 | { |
---|
720 | if ((a!=NULL)&&((_pGetComp(a)==0) || (_pGetComp(a) == _pGetComp(b)))) |
---|
721 | { |
---|
722 | return __pDivisibleBy(a,b); |
---|
723 | } |
---|
724 | return FALSE; |
---|
725 | } |
---|
726 | |
---|
727 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
728 | #define _pDivisibleBy1(a,b) pDBDivisibleBy1(a, b, __FILE__, __LINE__) |
---|
729 | extern BOOLEAN pDBDivisibleBy1(poly p1, poly p2, char* f, int l); |
---|
730 | inline BOOLEAN _pDivisibleBy1_orig(poly a, poly b) |
---|
731 | #else |
---|
732 | inline BOOLEAN _pDivisibleBy1(poly a, poly b) |
---|
733 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
734 | { |
---|
735 | if (_pGetComp(a) == 0 || _pGetComp(a) == _pGetComp(b)) |
---|
736 | return __pDivisibleBy(a,b); |
---|
737 | return FALSE; |
---|
738 | } |
---|
739 | |
---|
740 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
741 | #define _pDivisibleBy2(a,b) pDBDivisibleBy2(a, b, __FILE__, __LINE__) |
---|
742 | extern BOOLEAN pDBDivisibleBy2(poly p1, poly p2, char* f, int l); |
---|
743 | #else |
---|
744 | #define _pDivisibleBy2(a,b) __pDivisibleBy(a,b) |
---|
745 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
746 | |
---|
747 | |
---|
748 | DECLARE(BOOLEAN, _pEqual(poly p1, poly p2)) |
---|
749 | { |
---|
750 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
751 | const long *s1 = (long*) &(p1->exp[0]); |
---|
752 | const long *s2 = (long*) &(p2->exp[0]); |
---|
753 | const long* const lb = s1 + pVariables1W; |
---|
754 | #else |
---|
755 | const Exponent_t *s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
756 | const Exponent_t *s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
757 | const Exponent_t* const lb = s1 + pVariables; |
---|
758 | if (_pGetComp(p1) != _pGetComp(p2)) return FALSE; |
---|
759 | #endif |
---|
760 | |
---|
761 | for(;;) |
---|
762 | { |
---|
763 | if (*s1 != *s2) return FALSE; |
---|
764 | s1++; |
---|
765 | if (s1 == lb) return TRUE; |
---|
766 | s2++; |
---|
767 | } |
---|
768 | } |
---|
769 | |
---|
770 | inline void _pGetExpV(poly p, Exponent_t *ev) |
---|
771 | { |
---|
772 | if (_pHasReverseExp) |
---|
773 | { |
---|
774 | for (int i = pVarLowIndex, j = pVariables; j; i++, j--) |
---|
775 | ev[j] = p->exp[i]; |
---|
776 | } |
---|
777 | else |
---|
778 | memcpy(&ev[1], &(p->exp[pVarLowIndex]), pVariables*sizeof(Exponent_t)); |
---|
779 | ev[0] = _pGetComp(p); |
---|
780 | } |
---|
781 | |
---|
782 | extern pSetmProc pSetm; |
---|
783 | inline void _pSetExpV(poly p, Exponent_t *ev) |
---|
784 | { |
---|
785 | if (_pHasReverseExp) |
---|
786 | { |
---|
787 | for (int i = pVarLowIndex, j = pVariables; j; i++, j--) |
---|
788 | p->exp[i] = ev[j]; |
---|
789 | } |
---|
790 | else |
---|
791 | memcpy(&(p->exp[pVarLowIndex]), &ev[1], pVariables*sizeof(Exponent_t)); |
---|
792 | _pSetComp(p, ev[0]); |
---|
793 | pSetm(p); |
---|
794 | } |
---|
795 | |
---|
796 | #else // ! COMP_FAST |
---|
797 | |
---|
798 | DECLARE(BOOLEAN, _pDivisibleBy(poly a, poly b)) |
---|
799 | { |
---|
800 | if ((a!=NULL)&&((a->exp[0]==0) || (a->exp[0] == b->exp[0]))) |
---|
801 | { |
---|
802 | int i=pVariables; |
---|
803 | short *e1=&(a->exp[1]); |
---|
804 | short *e2=&(b->exp[1]); |
---|
805 | if ((*e1) > (*e2)) return FALSE; |
---|
806 | do |
---|
807 | { |
---|
808 | i--; |
---|
809 | if (i == 0) return TRUE; |
---|
810 | e1++; |
---|
811 | e2++; |
---|
812 | } while ((*e1) <= (*e2)); |
---|
813 | } |
---|
814 | return FALSE; |
---|
815 | } |
---|
816 | |
---|
817 | #define _pDivisibleBy1(a,b) _pDivisibleBy(a,b) |
---|
818 | #define _pDivisibleBy2(a,b) _pDivisibleBy(a,b) |
---|
819 | |
---|
820 | #ifdef TEST_MAC_ORDER |
---|
821 | DECLARE(void, pMonAddFast(poly a, poly m)) |
---|
822 | { |
---|
823 | for(int ii =pVariables; ii; ii--) (a)->exp[ii] += (m)->exp[ii];\ |
---|
824 | if (bNoAdd) bSetm(a); else |
---|
825 | _pGetOrder(a) += _pGetOrder(m); |
---|
826 | } |
---|
827 | #else |
---|
828 | DECLARE(void, pMonAddFast(poly a, poly m)) |
---|
829 | { |
---|
830 | for(int ii =pVariables; ii; ii--) (a)->exp[ii] += (m)->exp[ii];\ |
---|
831 | _pGetOrder(a) += _pGetOrder(m); |
---|
832 | } |
---|
833 | #endif |
---|
834 | |
---|
835 | DECLARE(BOOLEAN, _pEqual(poly p1, poly p2)) |
---|
836 | { |
---|
837 | int i; |
---|
838 | short *e1=p1->exp; |
---|
839 | short *e2=p2->exp; |
---|
840 | |
---|
841 | if (p1->Order != p2->Order) return FALSE; |
---|
842 | for (i=pVariables; i>=0; i--,e1++,e2++) |
---|
843 | if (*e1 != *e2) return FALSE; |
---|
844 | return TRUE; |
---|
845 | } |
---|
846 | |
---|
847 | #define _pGetExpV(p,e) memcpy((e),(p)->exp,(pVariables+1)*sizeof(short)); |
---|
848 | //void pSetExpV(poly p, short * exp); |
---|
849 | #define _pSetExpV(p,e) {memcpy((p)->exp,(e),(pVariables+1)*sizeof(short));pSetm(p);} |
---|
850 | |
---|
851 | #endif // COMP_FAST |
---|
852 | |
---|
853 | /*************************************************************** |
---|
854 | * |
---|
855 | * Routines which implement low-level manipulations/operations on exponents |
---|
856 | * |
---|
857 | ***************************************************************/ |
---|
858 | |
---|
859 | DECLARE(int, __pExpQuerSum2(poly p, int from, int to)) |
---|
860 | { |
---|
861 | int j = p->exp[from]; |
---|
862 | int i = from + 1; |
---|
863 | |
---|
864 | for(;;) |
---|
865 | { |
---|
866 | if (i > to) return j; |
---|
867 | j += p->exp[i]; |
---|
868 | i++; |
---|
869 | } |
---|
870 | } |
---|
871 | |
---|
872 | #ifdef COMP_FAST |
---|
873 | #define _pExpQuerSum(p) __pExpQuerSum2(p, pVarLowIndex, pVarHighIndex) |
---|
874 | |
---|
875 | #define _pExpQuerSum1(p,to) \ |
---|
876 | (_pHasReverseExp ? \ |
---|
877 | __pExpQuerSum2(p, _pExpIndex(to), _pExpIndex(1)) : \ |
---|
878 | __pExpQuerSum2(p, _pExpIndex(1), _pExpIndex(to))) |
---|
879 | |
---|
880 | #define _pExpQuerSum2(p,from,to) \ |
---|
881 | (_pHasReverseExp ? \ |
---|
882 | __pExpQuerSum2(p, _pExpIndex(to), _pExpIndex(from)) : \ |
---|
883 | __pExpQuerSum2(p, _pExpIndex(from), _pExpIndex(to))) |
---|
884 | #else |
---|
885 | |
---|
886 | #define _pExpQuerSum(p) __pExpQuerSum2(p, 1, pVariables) |
---|
887 | #define _pExpQuerSum1(p, to) __pExpQuerSum2(p, 1, to) |
---|
888 | #define _pExpQuerSum2(p, from, to) __pExpQuerSum2(p, from, to) |
---|
889 | |
---|
890 | #endif |
---|
891 | |
---|
892 | /*************************************************************** |
---|
893 | * |
---|
894 | * Routines which implement macaulay ordering routines |
---|
895 | * |
---|
896 | ***************************************************************/ |
---|
897 | #ifdef TEST_MAC_ORDER |
---|
898 | |
---|
899 | DECLARE(void, _bSetm0(poly p)) |
---|
900 | { |
---|
901 | |
---|
902 | int i=1; |
---|
903 | int ord = -INT_MAX; |
---|
904 | Exponent_t *ep; |
---|
905 | |
---|
906 | if(_pHasReverseExp) |
---|
907 | { |
---|
908 | ep=&(p->exp[pVarHighIndex]); |
---|
909 | int *ip=bBinomials+(*ep); /*_pGetExp(p,1);*/ |
---|
910 | loop |
---|
911 | { |
---|
912 | ord += (*ip); |
---|
913 | if (i==pVariables) break; |
---|
914 | i++; |
---|
915 | //ip+=bHighdeg_1+_pGetExp(p,i); |
---|
916 | ep--; |
---|
917 | ip+=bHighdeg_1+(*ep); |
---|
918 | } |
---|
919 | #if 0 |
---|
920 | int *ip=bBinomials+_pGetExp(p,1); |
---|
921 | loop |
---|
922 | { |
---|
923 | ord += (*ip); |
---|
924 | if (i==pVariables) break; |
---|
925 | i++; |
---|
926 | ip+=bHighdeg_1+_pGetExp(p,i); |
---|
927 | } |
---|
928 | #endif |
---|
929 | } |
---|
930 | else |
---|
931 | { |
---|
932 | ep=&(p->exp[pVarLowIndex]); |
---|
933 | int *ip=bBinomials+(*ep); /*_pGetExp(p,1);*/ |
---|
934 | loop |
---|
935 | { |
---|
936 | ord += (*ip); |
---|
937 | if (i==pVariables) break; |
---|
938 | i++; |
---|
939 | //ip+=bHighdeg_1+_pGetExp(p,i); |
---|
940 | ep++; |
---|
941 | ip+=bHighdeg_1+(*ep); |
---|
942 | } |
---|
943 | #if 0 |
---|
944 | int *ip=bBinomials+_pGetExp(p,1); |
---|
945 | loop |
---|
946 | { |
---|
947 | ord += (*ip); |
---|
948 | if (i==pVariables) break; |
---|
949 | i++; |
---|
950 | ip+=bHighdeg_1+_pGetExp(p,i); |
---|
951 | } |
---|
952 | #endif |
---|
953 | } |
---|
954 | p->Order=ord; |
---|
955 | } |
---|
956 | |
---|
957 | DECLARE(void, _bSetm(poly p)) |
---|
958 | { |
---|
959 | int ord = _pExpQuerSum(p); |
---|
960 | |
---|
961 | if (ord<bHighdeg) |
---|
962 | _bSetm0(p); |
---|
963 | else |
---|
964 | p->Order=ord; |
---|
965 | } |
---|
966 | |
---|
967 | // ordering dp,c or c,dp, general case |
---|
968 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
969 | #define pbMonAddFast(p1, p2) pDBMonAddFast(p1, p2, __FILE__, __LINE__) |
---|
970 | extern void pbDBMonAddFast(poly p1, poly p2, char* f, int l); |
---|
971 | inline void _pbMonAddFast(poly p1, poly p2) |
---|
972 | #else |
---|
973 | #define pbMonAddFast(p1, p2) _pbMonAddFast(p1, p2) |
---|
974 | DECLARE(void, _pbMonAddFast(poly p1, poly p2)) |
---|
975 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
976 | { |
---|
977 | // OK -- this might be the only place where we are kind of quick and |
---|
978 | // dirty: the following only works correctly if all exponents are |
---|
979 | // positive and the sum of two exponents does not exceed |
---|
980 | // EXPONENT_MAX |
---|
981 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
982 | Exponent_t c2 = _pGetComp(p2); |
---|
983 | int i = pVariables1W; |
---|
984 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
985 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
986 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
987 | _pSetComp(p2, 0); |
---|
988 | #else |
---|
989 | int i = pVariables; |
---|
990 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
991 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
992 | #endif |
---|
993 | |
---|
994 | for (;;) |
---|
995 | { |
---|
996 | *s1 += *s2; |
---|
997 | i--; |
---|
998 | if (i==0) break; |
---|
999 | s1++; |
---|
1000 | s2++; |
---|
1001 | } |
---|
1002 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1003 | // reset comp of p2 |
---|
1004 | _pSetComp(p2, c2); |
---|
1005 | #endif |
---|
1006 | if ((_pGetOrder(p1)|_pGetOrder(p2))>0) // i.e. overflow of mac order |
---|
1007 | _pGetOrder(p1) += _pGetOrder(p2); |
---|
1008 | else |
---|
1009 | _bSetm(p1); |
---|
1010 | } |
---|
1011 | |
---|
1012 | // ordering dp,c or c,dp, below degree limit |
---|
1013 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
1014 | #define pbMonAddFast0(p1, p2) pbDBMonAddFast0(p1, p2, __FILE__, __LINE__) |
---|
1015 | extern void pbDBMonAddFast0(poly p1, poly p2, char* f, int l); |
---|
1016 | inline void _pbMonAddFast0(poly p1, poly p2) |
---|
1017 | #else |
---|
1018 | DECLARE(void, pbMonAddFast0(poly p1, poly p2)) |
---|
1019 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
1020 | { |
---|
1021 | // OK -- this might be the only place where we are kind of quick and |
---|
1022 | // dirty: the following only works correctly if all exponents are |
---|
1023 | // positive and the sum of two exponents does not exceed |
---|
1024 | // EXPONENT_MAX |
---|
1025 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1026 | Exponent_t c2 = _pGetComp(p2); |
---|
1027 | int i = pVariables1W; |
---|
1028 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
1029 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
1030 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
1031 | _pSetComp(p2, 0); |
---|
1032 | #else |
---|
1033 | int i = pVariables; |
---|
1034 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
1035 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
1036 | #endif |
---|
1037 | |
---|
1038 | for (;;) |
---|
1039 | { |
---|
1040 | *s1 += *s2; |
---|
1041 | i--; |
---|
1042 | if (i==0) break; |
---|
1043 | s1++; |
---|
1044 | s2++; |
---|
1045 | } |
---|
1046 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1047 | // reset comp of p2 |
---|
1048 | _pSetComp(p2, c2); |
---|
1049 | #endif |
---|
1050 | _bSetm0(p1); |
---|
1051 | } |
---|
1052 | |
---|
1053 | // ordering dp,c or c,dp, below degree limit |
---|
1054 | // Makes p1 a copy of p2 and adds on exponets of p3 |
---|
1055 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
1056 | #define _pbCopyAddFast0(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
1057 | inline void __pbCopyAddFast0(poly p1, poly p2, poly p3) |
---|
1058 | #else |
---|
1059 | DECLARE(void, _pbCopyAddFast0(poly p1, poly p2, poly p3)) |
---|
1060 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
1061 | { |
---|
1062 | p1->next = p2->next; |
---|
1063 | p1->coef = p2->coef; |
---|
1064 | |
---|
1065 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1066 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
1067 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
1068 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
1069 | const unsigned long* const ub = s3 + pVariables1W; |
---|
1070 | #else |
---|
1071 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
1072 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
1073 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
1074 | const Exponent_t* const ub = s3 + pVariables; |
---|
1075 | // need to zero the "fill in" slots (i.e., empty exponents) |
---|
1076 | #ifdef WORDS_BIGENDIAN |
---|
1077 | *((unsigned long*) p1 + pMonomSizeW -1) = 0; |
---|
1078 | #else |
---|
1079 | *((unsigned long *) p1->exp) = 0; |
---|
1080 | #endif |
---|
1081 | #endif |
---|
1082 | |
---|
1083 | for (;;) |
---|
1084 | { |
---|
1085 | *s1 = *s2 + *s3; |
---|
1086 | s3++; |
---|
1087 | if (s3 == ub) break; |
---|
1088 | s1++; |
---|
1089 | s2++; |
---|
1090 | } |
---|
1091 | // we first are supposed to do a copy from p2 to p1 -- therefore, |
---|
1092 | // component of p1 is set to comp of p2 |
---|
1093 | _pSetComp(p1, _pGetComp(p2)); |
---|
1094 | _bSetm0(p1); |
---|
1095 | } |
---|
1096 | |
---|
1097 | // Similar to pCopyAddFast, except that we assume that the component |
---|
1098 | // of p2 and p3 is zero component |
---|
1099 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
1100 | #define _pbCopyAddFast10(p1, p2, p3) pbDBCopyAddFast0(p1, p2, p3, __FILE__, __LINE__) |
---|
1101 | extern void pbDBCopyAddFast0(poly p1, poly p2, poly p3, char* f, int l); |
---|
1102 | inline void __pbCopyAddFast10(poly p1, poly p2, poly p3) |
---|
1103 | #else |
---|
1104 | DECLARE(void, _pbCopyAddFast10(poly p1, poly p2, poly p3)) |
---|
1105 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
1106 | { |
---|
1107 | p1->next = p2->next; |
---|
1108 | p1->coef = p2->coef; |
---|
1109 | |
---|
1110 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1111 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
1112 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
1113 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
1114 | const unsigned long* const ub = s3 + pVariables1W; |
---|
1115 | #else |
---|
1116 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
1117 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
1118 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
1119 | const Exponent_t* const ub = s3 + pVariables; |
---|
1120 | #ifdef WORDS_BIGENDIAN |
---|
1121 | *((unsigned long*) p1 + pMonomSizeW -1) = 0; |
---|
1122 | #else |
---|
1123 | *((unsigned long *) p1->exp) = 0; |
---|
1124 | #endif |
---|
1125 | #endif |
---|
1126 | |
---|
1127 | for (;;) |
---|
1128 | { |
---|
1129 | *s1 = *s2 + *s3; |
---|
1130 | s3++; |
---|
1131 | if (s3 == ub) break; |
---|
1132 | s1++; |
---|
1133 | s2++; |
---|
1134 | } |
---|
1135 | _bSetm0(p1); |
---|
1136 | } |
---|
1137 | |
---|
1138 | #endif |
---|
1139 | |
---|
1140 | #endif // POLYS_IMPL_H |
---|
1141 | |
---|