1 | #ifndef POLYS_IMPL_H |
---|
2 | #define POLYS_IMPL_H |
---|
3 | /**************************************** |
---|
4 | * Computer Algebra System SINGULAR * |
---|
5 | ****************************************/ |
---|
6 | /* $Id: polys-impl.h,v 1.15 1998-01-24 17:22:06 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 | #define _pGetOrder(p) ((p)->Order) |
---|
279 | |
---|
280 | #if defined(PDEBUG) && PDEBUG != 0 |
---|
281 | extern Exponent_t pPDSetExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
282 | extern Exponent_t pPDGetExp(poly p, int v, char* f, int l); |
---|
283 | extern Exponent_t pPDIncrExp(poly p, int v, char* f, int l); |
---|
284 | extern Exponent_t pPDDecrExp(poly p, int v, char* f, int l); |
---|
285 | extern Exponent_t pPDAddExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
286 | extern Exponent_t pPDMultExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
287 | extern Exponent_t pPDSubExp(poly p, int v, Exponent_t e, char* f, int l); |
---|
288 | |
---|
289 | extern Exponent_t pPDRingSetExp(ring r,poly p,int v,Exponent_t e,char* f,int l); |
---|
290 | extern Exponent_t pPDRingGetExp(ring r,poly p, int v, char* f, int l); |
---|
291 | |
---|
292 | #define _pSetExp(p,v,e) pPDSetExp(p,v,e,__FILE__,__LINE__) |
---|
293 | #define _pGetExp(p,v) pPDGetExp(p,v,__FILE__,__LINE__) |
---|
294 | #define _pIncrExp(p,v) pPDIncrExp(p,v,__FILE__,__LINE__) |
---|
295 | #define _pDecrExp(p,v) pPDDecrExp(p,v,__FILE__,__LINE__) |
---|
296 | #define _pAddExp(p,i,v) pPDAddExp(p,i,v,__FILE__,__LINE__) |
---|
297 | #define _pSubExp(p,i,v) pPDSubExp(p,i,v,__FILE__,__LINE__) |
---|
298 | #define _pMultExp(p,i,v) pPDMultExp(p,i,v,__FILE__,__LINE__) |
---|
299 | |
---|
300 | #define _pRingSetExp(r,p,v,e) pPDRingSetExp(r,p,v,e,__FILE__,__LINE__) |
---|
301 | #define _pRingGetExp(r,p,v) pPDRingGetExp(r,p,v,__FILE__,__LINE__) |
---|
302 | |
---|
303 | #else // ! (defined(PDEBUG) && PDEBUG != 0) |
---|
304 | |
---|
305 | #define _pSetExp(p,v,e) (p)->exp[_pExpIndex(v)]=(e) |
---|
306 | #define _pGetExp(p,v) (p)->exp[_pExpIndex(v)] |
---|
307 | #define _pIncrExp(p,v) ((p)->exp[_pExpIndex(v)])++ |
---|
308 | #define _pDecrExp(p,v) ((p)->exp[_pExpIndex(v)])-- |
---|
309 | #define _pAddExp(p,i,v) ((p)->exp[_pExpIndex(i)]) += (v) |
---|
310 | #define _pSubExp(p,i,v) ((p)->exp[_pExpIndex(i)]) -= (v) |
---|
311 | #define _pMultExp(p,i,v) ((p)->exp[_pExpIndex(i)]) *= (v) |
---|
312 | |
---|
313 | #define _pRingSetExp(r,p,v,e) (p)->exp[_pRingExpIndex(r,v)]=(e) |
---|
314 | #define _pRingGetExp(r,p,v) (p)->exp[_pRingExpIndex(r,v)] |
---|
315 | |
---|
316 | #endif // defined(PDEBUG) && PDEBUG != 0 |
---|
317 | |
---|
318 | inline Exponent_t _pGetExpSum(poly p1, poly p2, int i) |
---|
319 | { |
---|
320 | int index = _pExpIndex(i); |
---|
321 | return p1->exp[index] + p2->exp[index]; |
---|
322 | } |
---|
323 | inline Exponent_t _pGetExpDiff(poly p1, poly p2, int i) |
---|
324 | { |
---|
325 | int index = _pExpIndex(i); |
---|
326 | return p1->exp[index] - p2->exp[index]; |
---|
327 | } |
---|
328 | |
---|
329 | #define _pSetComp(p,k) (p)->exp[pCompIndex] = (k) |
---|
330 | #define _pGetComp(p) (p)->exp[pCompIndex] |
---|
331 | #define _pIncrComp(p) (p)->exp[pCompIndex]++ |
---|
332 | #define _pDecrComp(p) (p)->exp[pCompIndex]-- |
---|
333 | #define _pAddComp(p,v) (p)->exp[pCompIndex] += (v) |
---|
334 | #define _pSubComp(p,v) (p)->exp[pCompIndex] -= (v) |
---|
335 | |
---|
336 | #ifdef COMP_FAST |
---|
337 | #define _pRingSetComp(r,p,k) (p)->exp[(r)->CompIndex] = (k) |
---|
338 | #define _pRingGetComp(r,p) (p)->exp[(r)->CompIndex] |
---|
339 | #else |
---|
340 | #define _pRingSetComp(r,p,k) _pSetComp(p,k) |
---|
341 | #define _pRingGetComp(r,p) _pGetComp(p) |
---|
342 | #endif |
---|
343 | |
---|
344 | |
---|
345 | /*************************************************************** |
---|
346 | * |
---|
347 | * Storage Managament Routines |
---|
348 | * |
---|
349 | ***************************************************************/ |
---|
350 | #ifdef MDEBUG |
---|
351 | |
---|
352 | poly pDBNew(char *f, int l); |
---|
353 | poly pDBInit(char * f,int l); |
---|
354 | |
---|
355 | void pDBDelete(poly * a, char * f, int l); |
---|
356 | void pDBDelete1(poly * a, char * f, int l); |
---|
357 | void pDBFree1(poly a, char * f, int l); |
---|
358 | |
---|
359 | poly pDBCopy(poly a, char *f, int l); |
---|
360 | poly pDBCopy1(poly a, char *f, int l); |
---|
361 | poly pDBHead(poly a, char *f, int l); |
---|
362 | poly pDBHead0(poly a, char *f, int l); |
---|
363 | poly pDBFetchCopy(ring r, poly a, char *f, int l); |
---|
364 | |
---|
365 | void pDBDelete(poly * a, char * f, int l); |
---|
366 | void pDBDelete1(poly * a, char * f, int l); |
---|
367 | void pDBFree1(poly a, char * f, int l); |
---|
368 | |
---|
369 | #define _pNew() pDBNew(__FILE__,__LINE__) |
---|
370 | #define _pInit() pDBInit(__FILE__,__LINE__) |
---|
371 | |
---|
372 | #define _pDelete(a) pDBDelete((a),__FILE__,__LINE__) |
---|
373 | #define _pDelete1(a) pDBDelete1((a),__FILE__,__LINE__) |
---|
374 | #define _pFree1(a) \ |
---|
375 | do \ |
---|
376 | { \ |
---|
377 | pDBFree1((a),__FILE__,__LINE__); \ |
---|
378 | (a)=NULL; \ |
---|
379 | } \ |
---|
380 | while(0) |
---|
381 | |
---|
382 | #define _pCopy(A) pDBCopy(A,__FILE__,__LINE__) |
---|
383 | #define _pCopy1(A) pDBCopy1(A, __FILE__,__LINE__) |
---|
384 | #define _pHead(A) pDBHead(A,__FILE__,__LINE__) |
---|
385 | #define _pHead0(A) pDBHead0(A, __FILE__,__LINE__) |
---|
386 | #ifdef COMP_FAST |
---|
387 | #define _pFetchCopy(r,A) pDBFetchCopy(r, A,__FILE__,__LINE__) |
---|
388 | #else |
---|
389 | #define _pFetchCopy(r,p) pOrdPolyInsertSetm(pCopy(p)) |
---|
390 | #endif |
---|
391 | |
---|
392 | #else // ! MDEBUG |
---|
393 | |
---|
394 | #define _pNew() (poly) mmAllocSpecialized() |
---|
395 | // #define _pNew() _pInit() |
---|
396 | |
---|
397 | #include <string.h> |
---|
398 | |
---|
399 | inline poly _pInit(void) |
---|
400 | { |
---|
401 | poly p=(poly)mmAllocSpecialized(); |
---|
402 | memset(p,0, pMonomSize); |
---|
403 | return p; |
---|
404 | } |
---|
405 | |
---|
406 | extern void _pDelete(poly * a); |
---|
407 | extern void _pDelete1(poly * a); |
---|
408 | #define _pFree1(a) mmFreeSpecialized((ADDRESS)a) |
---|
409 | |
---|
410 | extern poly _pCopy(poly a); |
---|
411 | extern poly _pCopy1(poly a); |
---|
412 | extern poly _pHead(poly a); |
---|
413 | extern poly _pHead0(poly a); |
---|
414 | #ifdef COMP_FAST |
---|
415 | extern poly _pFetchCopy(ring r,poly a); |
---|
416 | #else |
---|
417 | #define _pFetchCopy(r,p) pOrdPolyInsertSetm(pCopy(p)) |
---|
418 | #endif |
---|
419 | |
---|
420 | #endif // MDEBUG |
---|
421 | |
---|
422 | #define _pCopy2(p1, p2) memcpyW(p1, p2, pMonomSizeW) |
---|
423 | |
---|
424 | // Here is a handy Macro which disables inlining when run with |
---|
425 | // profiling and enables it otherwise |
---|
426 | |
---|
427 | #ifdef DO_PROFILE |
---|
428 | |
---|
429 | #ifndef POLYS_IMPL_CC |
---|
430 | #define DECLARE(type, arglist) type arglist; \ |
---|
431 | static type dummy_##arglist |
---|
432 | #else |
---|
433 | #define DECLARE(type, arglist) type arglist |
---|
434 | #endif // POLYS_IMPL_CC |
---|
435 | |
---|
436 | #else //! DO_PROFILE |
---|
437 | |
---|
438 | #define DECLARE(type, arglist ) inline type arglist |
---|
439 | |
---|
440 | #endif // DO_PROFILE |
---|
441 | |
---|
442 | |
---|
443 | /*************************************************************** |
---|
444 | * |
---|
445 | * Routines which work on vectors instead of single exponents |
---|
446 | * |
---|
447 | ***************************************************************/ |
---|
448 | |
---|
449 | #ifdef COMP_FAST |
---|
450 | |
---|
451 | // nice declaration isn't it ?? |
---|
452 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
453 | #define pMonAddFast(p1, p2) pDBMonAddFast(p1, p2, __FILE__, __LINE__) |
---|
454 | extern void pDBMonAddFast(poly p1, poly p2, char* f, int l); |
---|
455 | inline void _pMonAddFast(poly p1, poly p2) |
---|
456 | #else |
---|
457 | DECLARE(void, pMonAddFast(poly p1, poly p2)) |
---|
458 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
459 | { |
---|
460 | // OK -- this might be the only place where we are kind of quick and |
---|
461 | // dirty: the following only works correctly if all exponents are |
---|
462 | // positive and the sum of two exponents does not exceed |
---|
463 | // EXPONENT_MAX |
---|
464 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
465 | Exponent_t c2 = _pGetComp(p2); |
---|
466 | int i = pVariables1W; |
---|
467 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
468 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
469 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
470 | _pSetComp(p2, 0); |
---|
471 | #else |
---|
472 | int i = pVariables; |
---|
473 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
474 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
475 | #endif |
---|
476 | |
---|
477 | for (;;) |
---|
478 | { |
---|
479 | *s1 += *s2; |
---|
480 | i--; |
---|
481 | if (i==0) break; |
---|
482 | s1++; |
---|
483 | s2++; |
---|
484 | } |
---|
485 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
486 | // reset comp of p2 |
---|
487 | _pSetComp(p2, c2); |
---|
488 | #endif |
---|
489 | #ifdef TEST_MAC_ORDER |
---|
490 | if (bNoAdd) bSetm(p1);else |
---|
491 | #endif |
---|
492 | _pGetOrder(p1) += _pGetOrder(p2); |
---|
493 | } |
---|
494 | |
---|
495 | // Makes p1 a copy of p2 and adds on exponets of p3 |
---|
496 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
497 | #define _pCopyAddFast(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
498 | extern void pDBCopyAddFast(poly p1, poly p2, poly p3, char* f, int l); |
---|
499 | inline void __pCopyAddFast(poly p1, poly p2, poly p3) |
---|
500 | #else |
---|
501 | DECLARE(void, _pCopyAddFast(poly p1, poly p2, poly p3)) |
---|
502 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
503 | { |
---|
504 | p1->next = p2->next; |
---|
505 | p1->coef = p2->coef; |
---|
506 | // memset(p1, 0, pMonomSize); |
---|
507 | |
---|
508 | //memset(p1,0,pMonomSize); |
---|
509 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
510 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
511 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
512 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
513 | const unsigned long* const ub = s3 + pVariables1W; |
---|
514 | #else |
---|
515 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
516 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
517 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
518 | const Exponent_t* const ub = s3 + pVariables; |
---|
519 | // need to zero the "fill in" slots (i.e., empty exponents) |
---|
520 | #ifdef WORDS_BIGENDIAN |
---|
521 | *((unsigned long *) ((unsigned long*) p1) + pMonomSizeW -1) = 0; |
---|
522 | #else |
---|
523 | *((unsigned long *) p1->exp) = 0; |
---|
524 | #endif |
---|
525 | #endif |
---|
526 | |
---|
527 | for (;;) |
---|
528 | { |
---|
529 | *s1 = *s2 + *s3; |
---|
530 | s3++; |
---|
531 | if (s3 == ub) break; |
---|
532 | s1++; |
---|
533 | s2++; |
---|
534 | } |
---|
535 | // we first are supposed to do a copy from p2 to p1 -- therefore, |
---|
536 | // component of p1 is set to comp of p2 |
---|
537 | _pSetComp(p1, _pGetComp(p2)); |
---|
538 | #ifdef TEST_MAC_ORDER |
---|
539 | if (bNoAdd) bSetm(p1);else |
---|
540 | #endif |
---|
541 | _pGetOrder(p1) = _pGetOrder(p2) + _pGetOrder(p3); |
---|
542 | } |
---|
543 | |
---|
544 | // Similar to pCopyAddFast, except that we assume that the component |
---|
545 | // of p2 and p3 is zero component |
---|
546 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
547 | #define _pCopyAddFast1(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
548 | extern void pDBCopyAddFast(poly p1, poly p2, poly p3, char* f, int l); |
---|
549 | inline void __pCopyAddFast1(poly p1, poly p2, poly p3) |
---|
550 | #else |
---|
551 | DECLARE(void, _pCopyAddFast1(poly p1, poly p2, poly p3)) |
---|
552 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
553 | { |
---|
554 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
555 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
556 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
557 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
558 | const unsigned long* const ub = s3 + pVariables1W; |
---|
559 | #else |
---|
560 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
561 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
562 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
563 | const Exponent_t* const ub = s3 + pVariables; |
---|
564 | #ifdef WORDS_BIGENDIAN |
---|
565 | *((unsigned long *) ((unsigned long*) p1) + pMonomSizeW -1) = 0; |
---|
566 | #else |
---|
567 | *((unsigned long *) p1->exp) = 0; |
---|
568 | #endif |
---|
569 | #endif |
---|
570 | |
---|
571 | for (;;) |
---|
572 | { |
---|
573 | *s1 = *s2 + *s3; |
---|
574 | s3++; |
---|
575 | if (s3 == ub) break; |
---|
576 | s1++; |
---|
577 | s2++; |
---|
578 | } |
---|
579 | #ifdef TEST_MAC_ORDER |
---|
580 | if (bNoAdd) bSetm(p1);else |
---|
581 | #endif |
---|
582 | _pGetOrder(p1) = _pGetOrder(p2) + _pGetOrder(p3); |
---|
583 | } |
---|
584 | |
---|
585 | |
---|
586 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
587 | |
---|
588 | #if SIZEOF_LONG == 4 |
---|
589 | |
---|
590 | #if SIZEOF_EXPONENT == 1 |
---|
591 | #define P_DIV_MASK 0x80808080 |
---|
592 | #else // SIZEOF_EXPONENT == 2 |
---|
593 | #define P_DIV_MASK 0x80008000 |
---|
594 | #endif |
---|
595 | |
---|
596 | #else // SIZEOF_LONG == 8 |
---|
597 | |
---|
598 | #if SIZEOF_EXPONENT == 1 |
---|
599 | #define P_DIV_MASK 0x8080808080808080 |
---|
600 | #elif SIZEOF_EXPONENT == 2 |
---|
601 | #define P_DIV_MASK 0x8000800080008000 |
---|
602 | #else // SIZEOF_EXPONENT == 4 |
---|
603 | #define P_DIV_MASK 0x8000000080000000 |
---|
604 | #endif |
---|
605 | |
---|
606 | #endif |
---|
607 | |
---|
608 | DECLARE(BOOLEAN, __pDivisibleBy(poly a, poly b)) |
---|
609 | { |
---|
610 | const unsigned long* s1; |
---|
611 | const unsigned long* s2; |
---|
612 | const unsigned long* lb; |
---|
613 | |
---|
614 | #ifdef WORDS_BIGENDIAN |
---|
615 | lb = (unsigned long*) &(a->exp[0]); |
---|
616 | if (pVariables & ((SIZEOF_LONG / SIZEOF_EXPONENT) - 1)) |
---|
617 | { |
---|
618 | // now pVariables == pVariables1W, i.e. there are exponents in the |
---|
619 | // "first" word of exponentvector |
---|
620 | s1 = ((unsigned long*) a) + pMonomSizeW -1; |
---|
621 | s2 = ((unsigned long*) b) + pMonomSizeW -1; |
---|
622 | } |
---|
623 | else |
---|
624 | { |
---|
625 | // first exponent word has only component as significant field -- |
---|
626 | // Hence, do not bother |
---|
627 | s1 = ((unsigned long*) a) + pMonomSizeW -2; |
---|
628 | s2 = ((unsigned long*) b) + pMonomSizeW -2; |
---|
629 | } |
---|
630 | #else // !WORDS_BIGENDIAN |
---|
631 | lb = ((unsigned long*) a) + pMonomSizeW; |
---|
632 | if (pVariables & ((SIZEOF_LONG / SIZEOF_EXPONENT) - 1)) |
---|
633 | { |
---|
634 | s1 = (unsigned long*) &(a->exp[0]); |
---|
635 | s2 = (unsigned long*) &(b->exp[0]); |
---|
636 | } |
---|
637 | else |
---|
638 | { |
---|
639 | s1 = (unsigned long*) &(a->exp[0]) + 1; |
---|
640 | s2 = (unsigned long*) &(b->exp[0]) + 1; |
---|
641 | } |
---|
642 | #endif |
---|
643 | for (;;) |
---|
644 | { |
---|
645 | // O.K. -- and now comes a bit of magic. The following _really_ |
---|
646 | // works. Think about it! If you can prove it, please tell me, for I |
---|
647 | // did not bother to prove it formally (Hint: We can assume that our |
---|
648 | // exponents are always positive). |
---|
649 | if ((*s2 - *s1) & P_DIV_MASK) return FALSE; |
---|
650 | #ifdef WORDS_BIGENDIAN |
---|
651 | if (s1 == lb) return TRUE; |
---|
652 | s1--; |
---|
653 | s2--; |
---|
654 | #else |
---|
655 | s1++; |
---|
656 | if (s1 == lb) return TRUE; |
---|
657 | s2++; |
---|
658 | #endif |
---|
659 | } |
---|
660 | } |
---|
661 | |
---|
662 | #else // ! COMP_NO_EXP_VECTOR_OPS |
---|
663 | |
---|
664 | DECLARE(BOOLEAN, __pDivisibleBy(poly a, poly b)) |
---|
665 | { |
---|
666 | #ifdef WORDS_BIGENDIAN |
---|
667 | const Exponent_t* s1 = &(a->exp[pVarHighIndex]); |
---|
668 | const Exponent_t* s2 = &(b->exp[pVarHighIndex]); |
---|
669 | const Exponent_t* lb = s1 - pVariables; |
---|
670 | |
---|
671 | for (;;) |
---|
672 | { |
---|
673 | if (*s1 > *s2) return FALSE; |
---|
674 | s1--; |
---|
675 | if (s1 == lb) return TRUE; |
---|
676 | s2--; |
---|
677 | } |
---|
678 | |
---|
679 | #else // !WORDS_BIGENDIAN |
---|
680 | const Exponent_t* s1 = &(a->exp[pVarLowIndex]); |
---|
681 | const Exponent_t* s2 = &(b->exp[pVarLowIndex]); |
---|
682 | const Exponent_t* lb = s1 + pVariables; |
---|
683 | |
---|
684 | for (;;) |
---|
685 | { |
---|
686 | if (*s1 > *s2) return FALSE; |
---|
687 | s1++; |
---|
688 | if (s1 == lb) return TRUE; |
---|
689 | s2++; |
---|
690 | } |
---|
691 | #endif // WORDS_BIGENDIAN |
---|
692 | } |
---|
693 | |
---|
694 | #endif // COMP_NO_EXP_VECTOR_OPS |
---|
695 | |
---|
696 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
697 | #define _pDivisibleBy(a,b) pDBDivisibleBy(a, b, __FILE__, __LINE__) |
---|
698 | extern BOOLEAN pDBDivisibleBy(poly p1, poly p2, char* f, int l); |
---|
699 | inline BOOLEAN _pDivisibleBy_orig(poly a, poly b) |
---|
700 | #else |
---|
701 | inline BOOLEAN _pDivisibleBy(poly a, poly b) |
---|
702 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
703 | { |
---|
704 | if ((a!=NULL)&&((_pGetComp(a)==0) || (_pGetComp(a) == _pGetComp(b)))) |
---|
705 | { |
---|
706 | return __pDivisibleBy(a,b); |
---|
707 | } |
---|
708 | return FALSE; |
---|
709 | } |
---|
710 | |
---|
711 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
712 | #define _pDivisibleBy1(a,b) pDBDivisibleBy1(a, b, __FILE__, __LINE__) |
---|
713 | extern BOOLEAN pDBDivisibleBy1(poly p1, poly p2, char* f, int l); |
---|
714 | inline BOOLEAN _pDivisibleBy1_orig(poly a, poly b) |
---|
715 | #else |
---|
716 | inline BOOLEAN _pDivisibleBy1(poly a, poly b) |
---|
717 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
718 | { |
---|
719 | if (_pGetComp(a) == 0 || _pGetComp(a) == _pGetComp(b)) |
---|
720 | return __pDivisibleBy(a,b); |
---|
721 | return FALSE; |
---|
722 | } |
---|
723 | |
---|
724 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
725 | #define _pDivisibleBy2(a,b) pDBDivisibleBy2(a, b, __FILE__, __LINE__) |
---|
726 | extern BOOLEAN pDBDivisibleBy2(poly p1, poly p2, char* f, int l); |
---|
727 | #else |
---|
728 | #define _pDivisibleBy2(a,b) __pDivisibleBy(a,b) |
---|
729 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
730 | |
---|
731 | |
---|
732 | DECLARE(BOOLEAN, _pEqual(poly p1, poly p2)) |
---|
733 | { |
---|
734 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
735 | const long *s1 = (long*) &(p1->exp[0]); |
---|
736 | const long *s2 = (long*) &(p2->exp[0]); |
---|
737 | const long* const lb = s1 + pVariables1W; |
---|
738 | #else |
---|
739 | const Exponent_t *s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
740 | const Exponent_t *s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
741 | const Exponent_t* const lb = s1 + pVariables; |
---|
742 | if (_pGetComp(p1) != _pGetComp(p2)) return FALSE; |
---|
743 | #endif |
---|
744 | |
---|
745 | for(;;) |
---|
746 | { |
---|
747 | if (*s1 != *s2) return FALSE; |
---|
748 | s1++; |
---|
749 | if (s1 == lb) return TRUE; |
---|
750 | s2++; |
---|
751 | } |
---|
752 | } |
---|
753 | |
---|
754 | inline void _pGetExpV(poly p, Exponent_t *ev) |
---|
755 | { |
---|
756 | if (_pHasReverseExp) |
---|
757 | { |
---|
758 | for (int i = pVarLowIndex, j = pVariables; j; i++, j--) |
---|
759 | ev[j] = p->exp[i]; |
---|
760 | } |
---|
761 | else |
---|
762 | memcpy(&ev[1], &(p->exp[pVarLowIndex]), pVariables*sizeof(Exponent_t)); |
---|
763 | ev[0] = _pGetComp(p); |
---|
764 | } |
---|
765 | |
---|
766 | extern pSetmProc pSetm; |
---|
767 | inline void _pSetExpV(poly p, Exponent_t *ev) |
---|
768 | { |
---|
769 | if (_pHasReverseExp) |
---|
770 | { |
---|
771 | for (int i = pVarLowIndex, j = pVariables; j; i++, j--) |
---|
772 | p->exp[i] = ev[j]; |
---|
773 | } |
---|
774 | else |
---|
775 | memcpy(&(p->exp[pVarLowIndex]), &ev[1], pVariables*sizeof(Exponent_t)); |
---|
776 | _pSetComp(p, ev[0]); |
---|
777 | pSetm(p); |
---|
778 | } |
---|
779 | |
---|
780 | #else // ! COMP_FAST |
---|
781 | |
---|
782 | DECLARE(BOOLEAN, _pDivisibleBy(poly a, poly b)) |
---|
783 | { |
---|
784 | if ((a!=NULL)&&((a->exp[0]==0) || (a->exp[0] == b->exp[0]))) |
---|
785 | { |
---|
786 | int i=pVariables; |
---|
787 | short *e1=&(a->exp[1]); |
---|
788 | short *e2=&(b->exp[1]); |
---|
789 | if ((*e1) > (*e2)) return FALSE; |
---|
790 | do |
---|
791 | { |
---|
792 | i--; |
---|
793 | if (i == 0) return TRUE; |
---|
794 | e1++; |
---|
795 | e2++; |
---|
796 | } while ((*e1) <= (*e2)); |
---|
797 | } |
---|
798 | return FALSE; |
---|
799 | } |
---|
800 | |
---|
801 | #define _pDivisibleBy1(a,b) _pDivisibleBy(a,b) |
---|
802 | #define _pDivisibleBy2(a,b) _pDivisibleBy(a,b) |
---|
803 | |
---|
804 | #ifdef TEST_MAC_ORDER |
---|
805 | DECLARE(void, pMonAddFast(poly a, poly m)) |
---|
806 | { |
---|
807 | for(int ii =pVariables; ii; ii--) (a)->exp[ii] += (m)->exp[ii];\ |
---|
808 | if (bNoAdd) bSetm(a); else |
---|
809 | _pGetOrder(a) += _pGetOrder(m); |
---|
810 | } |
---|
811 | #else |
---|
812 | DECLARE(void, pMonAddFast(poly a, poly m)) |
---|
813 | { |
---|
814 | for(int ii =pVariables; ii; ii--) (a)->exp[ii] += (m)->exp[ii];\ |
---|
815 | _pGetOrder(a) += _pGetOrder(m); |
---|
816 | } |
---|
817 | #endif |
---|
818 | |
---|
819 | DECLARE(BOOLEAN, _pEqual(poly p1, poly p2)) |
---|
820 | { |
---|
821 | int i; |
---|
822 | short *e1=p1->exp; |
---|
823 | short *e2=p2->exp; |
---|
824 | |
---|
825 | if (p1->Order != p2->Order) return FALSE; |
---|
826 | for (i=pVariables; i>=0; i--,e1++,e2++) |
---|
827 | if (*e1 != *e2) return FALSE; |
---|
828 | return TRUE; |
---|
829 | } |
---|
830 | |
---|
831 | #define _pGetExpV(p,e) memcpy((e),(p)->exp,(pVariables+1)*sizeof(short)); |
---|
832 | //void pSetExpV(poly p, short * exp); |
---|
833 | #define _pSetExpV(p,e) {memcpy((p)->exp,(e),(pVariables+1)*sizeof(short));pSetm(p);} |
---|
834 | |
---|
835 | #endif // COMP_FAST |
---|
836 | |
---|
837 | /*************************************************************** |
---|
838 | * |
---|
839 | * Routines which implement low-level manipulations/operations on exponents |
---|
840 | * |
---|
841 | ***************************************************************/ |
---|
842 | |
---|
843 | DECLARE(int, __pExpQuerSum2(poly p, int from, int to)) |
---|
844 | { |
---|
845 | int j = p->exp[from]; |
---|
846 | int i = from + 1; |
---|
847 | |
---|
848 | for(;;) |
---|
849 | { |
---|
850 | if (i > to) return j; |
---|
851 | j += p->exp[i]; |
---|
852 | i++; |
---|
853 | } |
---|
854 | } |
---|
855 | |
---|
856 | #ifdef COMP_FAST |
---|
857 | #define _pExpQuerSum(p) __pExpQuerSum2(p, pVarLowIndex, pVarHighIndex) |
---|
858 | |
---|
859 | #define _pExpQuerSum1(p,to) \ |
---|
860 | (_pHasReverseExp ? \ |
---|
861 | __pExpQuerSum2(p, _pExpIndex(to), _pExpIndex(1)) : \ |
---|
862 | __pExpQuerSum2(p, _pExpIndex(1), _pExpIndex(to))) |
---|
863 | |
---|
864 | #define _pExpQuerSum2(p,from,to) \ |
---|
865 | (_pHasReverseExp ? \ |
---|
866 | __pExpQuerSum2(p, _pExpIndex(to), _pExpIndex(from)) : \ |
---|
867 | __pExpQuerSum2(p, _pExpIndex(from), _pExpIndex(to))) |
---|
868 | #else |
---|
869 | |
---|
870 | #define _pExpQuerSum(p) __pExpQuerSum2(p, 1, pVariables) |
---|
871 | #define _pExpQuerSum1(p, to) __pExpQuerSum2(p, 1, to) |
---|
872 | #define _pExpQuerSum2(p, from, to) __pExpQuerSum2(p, from, to) |
---|
873 | |
---|
874 | #endif |
---|
875 | |
---|
876 | /*************************************************************** |
---|
877 | * |
---|
878 | * Routines which implement macaulay ordering routines |
---|
879 | * |
---|
880 | ***************************************************************/ |
---|
881 | #ifdef TEST_MAC_ORDER |
---|
882 | |
---|
883 | DECLARE(void, _bSetm0(poly p)) |
---|
884 | { |
---|
885 | |
---|
886 | int i=1; |
---|
887 | int ord = -INT_MAX; |
---|
888 | Exponent_t *ep; |
---|
889 | |
---|
890 | if(_pHasReverseExp) |
---|
891 | { |
---|
892 | ep=&(p->exp[pVarHighIndex]); |
---|
893 | int *ip=bBinomials+(*ep); /*_pGetExp(p,1);*/ |
---|
894 | loop |
---|
895 | { |
---|
896 | ord += (*ip); |
---|
897 | if (i==pVariables) break; |
---|
898 | i++; |
---|
899 | //ip+=bHighdeg_1+_pGetExp(p,i); |
---|
900 | ep--; |
---|
901 | ip+=bHighdeg_1+(*ep); |
---|
902 | } |
---|
903 | } |
---|
904 | else |
---|
905 | { |
---|
906 | ep=&(p->exp[pVarLowIndex]); |
---|
907 | int *ip=bBinomials+(*ep); /*_pGetExp(p,1);*/ |
---|
908 | loop |
---|
909 | { |
---|
910 | ord += (*ip); |
---|
911 | if (i==pVariables) break; |
---|
912 | i++; |
---|
913 | //ip+=bHighdeg_1+_pGetExp(p,i); |
---|
914 | ep++; |
---|
915 | ip+=bHighdeg_1+(*ep); |
---|
916 | } |
---|
917 | } |
---|
918 | p->Order=ord; |
---|
919 | } |
---|
920 | |
---|
921 | DECLARE(void, _bSetm(poly p)) |
---|
922 | { |
---|
923 | int ord = _pExpQuerSum(p); |
---|
924 | |
---|
925 | if (ord<bHighdeg) |
---|
926 | _bSetm0(p); |
---|
927 | else |
---|
928 | p->Order=ord; |
---|
929 | } |
---|
930 | |
---|
931 | // ordering dp,c or c,dp, general case |
---|
932 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
933 | #define pbMonAddFast(p1, p2) pDBMonAddFast(p1, p2, __FILE__, __LINE__) |
---|
934 | extern void pbDBMonAddFast(poly p1, poly p2, char* f, int l); |
---|
935 | inline void _pbMonAddFast(poly p1, poly p2) |
---|
936 | #else |
---|
937 | #define pbMonAddFast(p1, p2) _pbMonAddFast(p1, p2) |
---|
938 | DECLARE(void, _pbMonAddFast(poly p1, poly p2)) |
---|
939 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
940 | { |
---|
941 | // OK -- this might be the only place where we are kind of quick and |
---|
942 | // dirty: the following only works correctly if all exponents are |
---|
943 | // positive and the sum of two exponents does not exceed |
---|
944 | // EXPONENT_MAX |
---|
945 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
946 | Exponent_t c2 = _pGetComp(p2); |
---|
947 | int i = pVariables1W; |
---|
948 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
949 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
950 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
951 | _pSetComp(p2, 0); |
---|
952 | #else |
---|
953 | int i = pVariables; |
---|
954 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
955 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
956 | #endif |
---|
957 | |
---|
958 | for (;;) |
---|
959 | { |
---|
960 | *s1 += *s2; |
---|
961 | i--; |
---|
962 | if (i==0) break; |
---|
963 | s1++; |
---|
964 | s2++; |
---|
965 | } |
---|
966 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
967 | // reset comp of p2 |
---|
968 | _pSetComp(p2, c2); |
---|
969 | #endif |
---|
970 | if ((_pGetOrder(p1)|_pGetOrder(p2))>0) // i.e. overflow of mac order |
---|
971 | _pGetOrder(p1) += _pGetOrder(p2); |
---|
972 | else |
---|
973 | _bSetm(p1); |
---|
974 | } |
---|
975 | |
---|
976 | // ordering dp,c or c,dp, below degree limit |
---|
977 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
978 | #define pbMonAddFast0(p1, p2) pbDBMonAddFast0(p1, p2, __FILE__, __LINE__) |
---|
979 | extern void pbDBMonAddFast0(poly p1, poly p2, char* f, int l); |
---|
980 | inline void _pbMonAddFast0(poly p1, poly p2) |
---|
981 | #else |
---|
982 | DECLARE(void, pbMonAddFast0(poly p1, poly p2)) |
---|
983 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
984 | { |
---|
985 | // OK -- this might be the only place where we are kind of quick and |
---|
986 | // dirty: the following only works correctly if all exponents are |
---|
987 | // positive and the sum of two exponents does not exceed |
---|
988 | // EXPONENT_MAX |
---|
989 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
990 | Exponent_t c2 = _pGetComp(p2); |
---|
991 | int i = pVariables1W; |
---|
992 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
993 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
994 | // set comp of p2 temporarily to 0, so that nothing is added to comp of p1 |
---|
995 | _pSetComp(p2, 0); |
---|
996 | #else |
---|
997 | int i = pVariables; |
---|
998 | Exponent_pt s1 = &(p1->exp[pVarLowIndex]); |
---|
999 | Exponent_pt s2 = &(p2->exp[pVarLowIndex]); |
---|
1000 | #endif |
---|
1001 | |
---|
1002 | for (;;) |
---|
1003 | { |
---|
1004 | *s1 += *s2; |
---|
1005 | i--; |
---|
1006 | if (i==0) break; |
---|
1007 | s1++; |
---|
1008 | s2++; |
---|
1009 | } |
---|
1010 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1011 | // reset comp of p2 |
---|
1012 | _pSetComp(p2, c2); |
---|
1013 | #endif |
---|
1014 | _bSetm0(p1); |
---|
1015 | } |
---|
1016 | |
---|
1017 | // ordering dp,c or c,dp, below degree limit |
---|
1018 | // Makes p1 a copy of p2 and adds on exponets of p3 |
---|
1019 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
1020 | #define _pbCopyAddFast0(p1, p2, p3) pDBCopyAddFast(p1, p2, p3, __FILE__, __LINE__) |
---|
1021 | inline void __pbCopyAddFast0(poly p1, poly p2, poly p3) |
---|
1022 | #else |
---|
1023 | DECLARE(void, _pbCopyAddFast0(poly p1, poly p2, poly p3)) |
---|
1024 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
1025 | { |
---|
1026 | p1->next = p2->next; |
---|
1027 | p1->coef = p2->coef; |
---|
1028 | |
---|
1029 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1030 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
1031 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
1032 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
1033 | const unsigned long* const ub = s3 + pVariables1W; |
---|
1034 | #else |
---|
1035 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
1036 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
1037 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
1038 | const Exponent_t* const ub = s3 + pVariables; |
---|
1039 | // need to zero the "fill in" slots (i.e., empty exponents) |
---|
1040 | #ifdef WORDS_BIGENDIAN |
---|
1041 | *((unsigned long*) p1 + pMonomSizeW -1) = 0; |
---|
1042 | #else |
---|
1043 | *((unsigned long *) p1->exp) = 0; |
---|
1044 | #endif |
---|
1045 | #endif |
---|
1046 | |
---|
1047 | for (;;) |
---|
1048 | { |
---|
1049 | *s1 = *s2 + *s3; |
---|
1050 | s3++; |
---|
1051 | if (s3 == ub) break; |
---|
1052 | s1++; |
---|
1053 | s2++; |
---|
1054 | } |
---|
1055 | // we first are supposed to do a copy from p2 to p1 -- therefore, |
---|
1056 | // component of p1 is set to comp of p2 |
---|
1057 | _pSetComp(p1, _pGetComp(p2)); |
---|
1058 | _bSetm0(p1); |
---|
1059 | } |
---|
1060 | |
---|
1061 | // Similar to pCopyAddFast, except that we assume that the component |
---|
1062 | // of p2 and p3 is zero component |
---|
1063 | #if defined(PDEBUG) && PDEBUG == 1 |
---|
1064 | #define _pbCopyAddFast10(p1, p2, p3) pbDBCopyAddFast0(p1, p2, p3, __FILE__, __LINE__) |
---|
1065 | extern void pbDBCopyAddFast0(poly p1, poly p2, poly p3, char* f, int l); |
---|
1066 | inline void __pbCopyAddFast10(poly p1, poly p2, poly p3) |
---|
1067 | #else |
---|
1068 | DECLARE(void, _pbCopyAddFast10(poly p1, poly p2, poly p3)) |
---|
1069 | #endif // defined(PDEBUG) && PDEBUG == 1 |
---|
1070 | { |
---|
1071 | p1->next = p2->next; |
---|
1072 | p1->coef = p2->coef; |
---|
1073 | |
---|
1074 | #ifndef COMP_NO_EXP_VECTOR_OPS |
---|
1075 | unsigned long* s1 = (unsigned long*) &(p1->exp[0]); |
---|
1076 | const unsigned long* s2 = (unsigned long*) &(p2->exp[0]); |
---|
1077 | const unsigned long* s3 = (unsigned long*) &(p3->exp[0]); |
---|
1078 | const unsigned long* const ub = s3 + pVariables1W; |
---|
1079 | #else |
---|
1080 | Exponent_t* s1 = (Exponent_t*) &(p1->exp[pVarLowIndex]); |
---|
1081 | const Exponent_t* s2 = (Exponent_t*) &(p2->exp[pVarLowIndex]); |
---|
1082 | const Exponent_t* s3 = (Exponent_t*) &(p3->exp[pVarLowIndex]); |
---|
1083 | const Exponent_t* const ub = s3 + pVariables; |
---|
1084 | #ifdef WORDS_BIGENDIAN |
---|
1085 | *((unsigned long*) p1 + pMonomSizeW -1) = 0; |
---|
1086 | #else |
---|
1087 | *((unsigned long *) p1->exp) = 0; |
---|
1088 | #endif |
---|
1089 | #endif |
---|
1090 | |
---|
1091 | for (;;) |
---|
1092 | { |
---|
1093 | *s1 = *s2 + *s3; |
---|
1094 | s3++; |
---|
1095 | if (s3 == ub) break; |
---|
1096 | s1++; |
---|
1097 | s2++; |
---|
1098 | } |
---|
1099 | _bSetm0(p1); |
---|
1100 | } |
---|
1101 | |
---|
1102 | #endif |
---|
1103 | |
---|
1104 | #endif // POLYS_IMPL_H |
---|
1105 | |
---|