![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to nd.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / engine|
Return to nd.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / engine |
1.22 ! noro 1: /* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.21 2003/08/01 08:39:54 noro Exp $ */
1.2 noro 2:
1.1 noro 3: #include "ca.h"
4: #include "inline.h"
5:
6: #if defined(__GNUC__)
7: #define INLINE inline
8: #elif defined(VISUAL)
9: #define INLINE __inline
10: #else
11: #define INLINE
12: #endif
13:
14: #define REDTAB_LEN 32003
15:
16: typedef struct oPGeoBucket {
17: int m;
18: struct oND *body[32];
19: } *PGeoBucket;
20:
21: typedef struct oND {
22: struct oNM *body;
23: int nv;
24: int sugar;
25: } *ND;
26:
1.3 noro 27: typedef struct oNDV {
28: struct oNMV *body;
29: int nv;
30: int sugar;
31: int len;
32: } *NDV;
33:
1.1 noro 34: typedef struct oNM {
35: struct oNM *next;
1.14 noro 36: union {
37: int m;
38: Q z;
39: } c;
1.1 noro 40: int td;
41: unsigned int dl[1];
42: } *NM;
43:
1.3 noro 44: typedef struct oNMV {
1.14 noro 45: union {
46: int m;
47: Q z;
48: } c;
1.3 noro 49: int td;
50: unsigned int dl[1];
51: } *NMV;
52:
1.13 noro 53: typedef struct oRHist {
54: struct oRHist *next;
55: int index;
1.20 noro 56: int td,sugar;
1.13 noro 57: unsigned int dl[1];
58: } *RHist;
59:
1.1 noro 60: typedef struct oND_pairs {
61: struct oND_pairs *next;
62: int i1,i2;
63: int td,sugar;
64: unsigned int lcm[1];
65: } *ND_pairs;
66:
1.21 noro 67: double nd_scale=2;
1.1 noro 68: static unsigned int **nd_bound;
1.19 noro 69: int nd_nvar;
1.1 noro 70: int is_rlex;
71: int nd_epw,nd_bpe,nd_wpd;
72: unsigned int nd_mask[32];
73: unsigned int nd_mask0,nd_mask1;
1.20 noro 74:
1.1 noro 75: NM _nm_free_list;
76: ND _nd_free_list;
77: ND_pairs _ndp_free_list;
1.20 noro 78:
79: static NDV *nd_ps;
80: static NDV *nd_psq;
81: int *nd_psl;
82: RHist *nd_psh;
83: int nd_psn,nd_pslen;
84:
1.13 noro 85: RHist *nd_red;
1.1 noro 86: int nd_red_len;
87:
88: int nd_found,nd_create,nd_notfirst;
1.13 noro 89: int nm_adv;
1.20 noro 90: int nmv_adv;
91: int nmv_len;
92: NDV ndv_red;
1.1 noro 93:
1.20 noro 94: extern int Top,Reverse;
1.1 noro 95:
96: #define HTD(d) ((d)->body->td)
97: #define HDL(d) ((d)->body->dl)
1.14 noro 98: #define HCM(d) ((d)->body->c.m)
1.16 noro 99: #define HCQ(d) ((d)->body->c.z)
1.14 noro 100: #define CM(a) ((a)->c.m)
1.16 noro 101: #define CQ(a) ((a)->c.z)
1.14 noro 102: #define DL(a) ((a)->dl)
103: #define TD(a) ((a)->td)
104: #define SG(a) ((a)->sugar)
105: #define LEN(a) ((a)->len)
1.1 noro 106:
1.20 noro 107: #define NM_ADV(m) (m = (NM)(((char *)m)+nm_adv))
1.15 noro 108: #define NEWRHist(r) \
109: ((r)=(RHist)MALLOC(sizeof(struct oRHist)+(nd_wpd-1)*sizeof(unsigned int)))
1.1 noro 110: #define NEWND_pairs(m) if(!_ndp_free_list)_NDP_alloc(); (m)=_ndp_free_list; _ndp_free_list = NEXT(_ndp_free_list)
111: #define NEWNM(m) if(!_nm_free_list)_NM_alloc(); (m)=_nm_free_list; _nm_free_list = NEXT(_nm_free_list)
112: #define MKND(n,m,d) if(!_nd_free_list)_ND_alloc(); (d)=_nd_free_list; _nd_free_list = (ND)BDY(_nd_free_list); (d)->nv=(n); BDY(d)=(m)
113:
1.13 noro 114: #define NEXTRHist(r,c) \
115: if(!(r)){NEWRHist(r);(c)=(r);}else{NEWRHist(NEXT(c));(c)=NEXT(c);}
1.1 noro 116: #define NEXTNM(r,c) \
117: if(!(r)){NEWNM(r);(c)=(r);}else{NEWNM(NEXT(c));(c)=NEXT(c);}
118: #define NEXTNM2(r,c,s) \
119: if(!(r)){(c)=(r)=(s);}else{NEXT(c)=(s);(c)=(s);}
120: #define FREENM(m) NEXT(m)=_nm_free_list; _nm_free_list=(m)
121: #define FREENDP(m) NEXT(m)=_ndp_free_list; _ndp_free_list=(m)
122: #define FREEND(m) BDY(m)=(NM)_nd_free_list; _nd_free_list=(m)
123:
124: #define NEXTND_pairs(r,c) \
125: if(!(r)){NEWND_pairs(r);(c)=(r);}else{NEWND_pairs(NEXT(c));(c)=NEXT(c);}
126:
1.20 noro 127: void nd_removecont(int mod,ND p);
1.21 noro 128: void nd_removecont2(ND p1,ND p2);
1.20 noro 129: void ndv_removecont(int mod,NDV p);
1.16 noro 130: void ndv_mul_c_q(NDV p,Q mul);
131: void nd_mul_c_q(ND p,Q mul);
132:
1.20 noro 133: void GC_gcollect();
134: NODE append_one(NODE,int);
135:
1.21 noro 136: void removecont_array(Q *c,int n);
1.1 noro 137: ND_pairs crit_B( ND_pairs d, int s );
138: void nd_gr(LIST f,LIST v,int m,struct order_spec *ord,LIST *rp);
1.20 noro 139: void nd_gr_trace(LIST f,LIST v,int m,struct order_spec *ord,LIST *rp);
1.16 noro 140: NODE nd_setup(int mod,NODE f);
141: int nd_newps(int mod,ND a);
1.20 noro 142: int nd_newps_trace(int mod,ND nf,ND nfq);
1.1 noro 143: ND_pairs nd_minp( ND_pairs d, ND_pairs *prest );
144: NODE update_base(NODE nd,int ndp);
145: static ND_pairs equivalent_pairs( ND_pairs d1, ND_pairs *prest );
146: int crit_2( int dp1, int dp2 );
147: ND_pairs crit_F( ND_pairs d1 );
148: ND_pairs crit_M( ND_pairs d1 );
149: ND_pairs nd_newpairs( NODE g, int t );
150: ND_pairs update_pairs( ND_pairs d, NODE /* of index */ g, int t);
1.16 noro 151: NODE nd_gb(int m,NODE f);
1.20 noro 152: NODE nd_gb_trace(int m,NODE f);
1.1 noro 153: void nd_free_private_storage();
154: void _NM_alloc();
155: void _ND_alloc();
156: int ndl_td(unsigned int *d);
1.19 noro 157: ND nd_add(int mod,ND p1,ND p2);
1.17 noro 158: ND nd_add_q(ND p1,ND p2);
159: ND nd_mul_nm(int mod,ND p,NM m0);
160: ND nd_mul_ind_nm(int mod,int index,NM m0);
1.16 noro 161: int nd_sp(int mod,ND_pairs p,ND *nf);
1.6 noro 162: int nd_find_reducer(ND g);
1.16 noro 163: int nd_nf(int mod,ND g,int full,ND *nf);
1.1 noro 164: ND nd_reduce(ND p1,ND p2);
165: ND nd_reduce_special(ND p1,ND p2);
166: void nd_free(ND p);
167: void ndl_print(unsigned int *dl);
168: void nd_print(ND p);
1.16 noro 169: void nd_print_q(ND p);
170: void ndv_print(NDV p);
171: void ndv_print_q(NDV p);
1.1 noro 172: void ndp_print(ND_pairs d);
173: int nd_length(ND p);
1.19 noro 174: void nd_mul_c(int mod,ND p,int mul);
1.1 noro 175: void nd_free_redlist();
176: void nd_append_red(unsigned int *d,int td,int i);
177: unsigned int *nd_compute_bound(ND p);
1.5 noro 178: unsigned int *dp_compute_bound(DP p);
1.20 noro 179: ND_pairs nd_reconstruct(int mod,int trace,ND_pairs ndp);
1.1 noro 180: void nd_setup_parameters();
1.11 noro 181: void nd_realloc(ND p,int obpe);
1.6 noro 182: ND nd_copy(ND p);
1.1 noro 183: void ndl_dup(int obpe,unsigned int *d,unsigned int *r);
1.4 noro 184:
185: #define NMV_ADV(m) (m = (NMV)(((char *)m)+nmv_adv))
186: #define NEWNDV(d) ((d)=(NDV)MALLOC(sizeof(struct oNDV)))
1.14 noro 187: #define MKNDV(n,m,l,d) NEWNDV(d); NV(d)=(n); BDY(d)=(m); LEN(d) = l;
1.19 noro 188: void ndv_mul_c(int mod,NDV p,int mul);
189: ND ndv_add(int mod,ND p1,NDV p2);
1.16 noro 190: ND ndv_add_q(ND p1,NDV p2);
191: NDV ndtondv(int mod,ND p);
192: void ndv_mul_nm(int mod,NDV pv,NM m,NDV r);
193: ND ndv_mul_nm_create(int mod,NDV p,NM m0);
1.11 noro 194: void ndv_realloc(NDV p,int obpe,int oadv);
1.16 noro 195: NDV dptondv(int,DP);
196: DP ndvtodp(int,NDV);
1.17 noro 197: ND dptond(int,DP);
198: DP ndtodp(int,ND);
1.1 noro 199:
200: void nd_free_private_storage()
201: {
202: _nd_free_list = 0;
203: _nm_free_list = 0;
1.5 noro 204: _ndp_free_list = 0;
1.13 noro 205: bzero(nd_red,sizeof(REDTAB_LEN*sizeof(RHist)));
1.1 noro 206: GC_gcollect();
207: }
208:
209: void _NM_alloc()
210: {
211: NM p;
212: int i;
213:
1.11 noro 214: for ( i = 0; i < 1024; i++ ) {
1.1 noro 215: p = (NM)GC_malloc(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int));
216: p->next = _nm_free_list; _nm_free_list = p;
217: }
218: }
219:
220: void _ND_alloc()
221: {
222: ND p;
223: int i;
224:
225: for ( i = 0; i < 1024; i++ ) {
226: p = (ND)GC_malloc(sizeof(struct oND));
227: p->body = (NM)_nd_free_list; _nd_free_list = p;
228: }
229: }
230:
231: void _NDP_alloc()
232: {
233: ND_pairs p;
234: int i;
235:
1.11 noro 236: for ( i = 0; i < 1024; i++ ) {
1.1 noro 237: p = (ND_pairs)GC_malloc(sizeof(struct oND_pairs)
238: +(nd_wpd-1)*sizeof(unsigned int));
239: p->next = _ndp_free_list; _ndp_free_list = p;
240: }
241: }
242:
243: INLINE nd_length(ND p)
244: {
245: NM m;
246: int i;
247:
248: if ( !p )
249: return 0;
250: else {
251: for ( i = 0, m = BDY(p); m; m = NEXT(m), i++ );
252: return i;
253: }
254: }
255:
256: INLINE int ndl_reducible(unsigned int *d1,unsigned int *d2)
257: {
258: unsigned int u1,u2;
259: int i,j;
260:
261: switch ( nd_bpe ) {
262: case 4:
263: for ( i = 0; i < nd_wpd; i++ ) {
264: u1 = d1[i]; u2 = d2[i];
265: if ( (u1&0xf0000000) < (u2&0xf0000000) ) return 0;
266: if ( (u1&0xf000000) < (u2&0xf000000) ) return 0;
267: if ( (u1&0xf00000) < (u2&0xf00000) ) return 0;
268: if ( (u1&0xf0000) < (u2&0xf0000) ) return 0;
269: if ( (u1&0xf000) < (u2&0xf000) ) return 0;
270: if ( (u1&0xf00) < (u2&0xf00) ) return 0;
271: if ( (u1&0xf0) < (u2&0xf0) ) return 0;
272: if ( (u1&0xf) < (u2&0xf) ) return 0;
273: }
274: return 1;
275: break;
276: case 6:
277: for ( i = 0; i < nd_wpd; i++ ) {
278: u1 = d1[i]; u2 = d2[i];
279: if ( (u1&0x3f000000) < (u2&0x3f000000) ) return 0;
280: if ( (u1&0xfc0000) < (u2&0xfc0000) ) return 0;
281: if ( (u1&0x3f000) < (u2&0x3f000) ) return 0;
282: if ( (u1&0xfc0) < (u2&0xfc0) ) return 0;
283: if ( (u1&0x3f) < (u2&0x3f) ) return 0;
284: }
285: return 1;
286: break;
287: case 8:
288: for ( i = 0; i < nd_wpd; i++ ) {
289: u1 = d1[i]; u2 = d2[i];
290: if ( (u1&0xff000000) < (u2&0xff000000) ) return 0;
291: if ( (u1&0xff0000) < (u2&0xff0000) ) return 0;
292: if ( (u1&0xff00) < (u2&0xff00) ) return 0;
293: if ( (u1&0xff) < (u2&0xff) ) return 0;
294: }
295: return 1;
296: break;
297: case 16:
298: for ( i = 0; i < nd_wpd; i++ ) {
299: u1 = d1[i]; u2 = d2[i];
300: if ( (u1&0xffff0000) < (u2&0xffff0000) ) return 0;
301: if ( (u1&0xffff) < (u2&0xffff) ) return 0;
302: }
303: return 1;
304: break;
305: case 32:
306: for ( i = 0; i < nd_wpd; i++ )
307: if ( d1[i] < d2[i] ) return 0;
308: return 1;
309: break;
310: default:
311: for ( i = 0; i < nd_wpd; i++ ) {
312: u1 = d1[i]; u2 = d2[i];
313: for ( j = 0; j < nd_epw; j++ )
314: if ( (u1&nd_mask[j]) < (u2&nd_mask[j]) ) return 0;
315: }
316: return 1;
317: }
318: }
319:
320: void ndl_lcm(unsigned int *d1,unsigned *d2,unsigned int *d)
321: {
322: unsigned int t1,t2,u,u1,u2;
323: int i,j;
324:
325: switch ( nd_bpe ) {
326: case 4:
327: for ( i = 0; i < nd_wpd; i++ ) {
328: u1 = d1[i]; u2 = d2[i];
329: t1 = (u1&0xf0000000); t2 = (u2&0xf0000000); u = t1>t2?t1:t2;
330: t1 = (u1&0xf000000); t2 = (u2&0xf000000); u |= t1>t2?t1:t2;
331: t1 = (u1&0xf00000); t2 = (u2&0xf00000); u |= t1>t2?t1:t2;
332: t1 = (u1&0xf0000); t2 = (u2&0xf0000); u |= t1>t2?t1:t2;
333: t1 = (u1&0xf000); t2 = (u2&0xf000); u |= t1>t2?t1:t2;
334: t1 = (u1&0xf00); t2 = (u2&0xf00); u |= t1>t2?t1:t2;
335: t1 = (u1&0xf0); t2 = (u2&0xf0); u |= t1>t2?t1:t2;
336: t1 = (u1&0xf); t2 = (u2&0xf); u |= t1>t2?t1:t2;
337: d[i] = u;
338: }
339: break;
340: case 6:
341: for ( i = 0; i < nd_wpd; i++ ) {
342: u1 = d1[i]; u2 = d2[i];
343: t1 = (u1&0x3f000000); t2 = (u2&0x3f000000); u = t1>t2?t1:t2;
344: t1 = (u1&0xfc0000); t2 = (u2&0xfc0000); u |= t1>t2?t1:t2;
345: t1 = (u1&0x3f000); t2 = (u2&0x3f000); u |= t1>t2?t1:t2;
346: t1 = (u1&0xfc0); t2 = (u2&0xfc0); u |= t1>t2?t1:t2;
347: t1 = (u1&0x3f); t2 = (u2&0x3f); u |= t1>t2?t1:t2;
348: d[i] = u;
349: }
350: break;
351: case 8:
352: for ( i = 0; i < nd_wpd; i++ ) {
353: u1 = d1[i]; u2 = d2[i];
354: t1 = (u1&0xff000000); t2 = (u2&0xff000000); u = t1>t2?t1:t2;
355: t1 = (u1&0xff0000); t2 = (u2&0xff0000); u |= t1>t2?t1:t2;
356: t1 = (u1&0xff00); t2 = (u2&0xff00); u |= t1>t2?t1:t2;
357: t1 = (u1&0xff); t2 = (u2&0xff); u |= t1>t2?t1:t2;
358: d[i] = u;
359: }
360: break;
361: case 16:
362: for ( i = 0; i < nd_wpd; i++ ) {
363: u1 = d1[i]; u2 = d2[i];
364: t1 = (u1&0xffff0000); t2 = (u2&0xffff0000); u = t1>t2?t1:t2;
365: t1 = (u1&0xffff); t2 = (u2&0xffff); u |= t1>t2?t1:t2;
366: d[i] = u;
367: }
368: break;
369: case 32:
370: for ( i = 0; i < nd_wpd; i++ ) {
371: u1 = d1[i]; u2 = d2[i];
372: d[i] = u1>u2?u1:u2;
373: }
374: break;
375: default:
376: for ( i = 0; i < nd_wpd; i++ ) {
377: u1 = d1[i]; u2 = d2[i];
378: for ( j = 0, u = 0; j < nd_epw; j++ ) {
379: t1 = (u1&nd_mask[j]); t2 = (u2&nd_mask[j]); u |= t1>t2?t1:t2;
380: }
381: d[i] = u;
382: }
383: break;
384: }
385: }
386:
387: int ndl_td(unsigned int *d)
388: {
389: unsigned int t,u;
390: int i,j;
391:
392: for ( t = 0, i = 0; i < nd_wpd; i++ ) {
393: u = d[i];
394: for ( j = 0; j < nd_epw; j++, u>>=nd_bpe )
395: t += (u&nd_mask0);
396: }
397: return t;
398: }
399:
400: INLINE int ndl_compare(unsigned int *d1,unsigned int *d2)
401: {
402: int i;
403:
404: for ( i = 0; i < nd_wpd; i++, d1++, d2++ )
405: if ( *d1 > *d2 )
406: return is_rlex ? -1 : 1;
407: else if ( *d1 < *d2 )
408: return is_rlex ? 1 : -1;
409: return 0;
410: }
411:
412: INLINE int ndl_equal(unsigned int *d1,unsigned int *d2)
413: {
414: int i;
415:
416: for ( i = 0; i < nd_wpd; i++ )
417: if ( d1[i] != d2[i] )
418: return 0;
419: return 1;
420: }
421:
1.6 noro 422: INLINE void ndl_copy(unsigned int *d1,unsigned int *d2)
423: {
424: int i;
425:
426: switch ( nd_wpd ) {
427: case 1:
428: d2[0] = d1[0];
429: break;
430: case 2:
431: d2[0] = d1[0];
432: d2[1] = d1[1];
433: break;
434: default:
435: for ( i = 0; i < nd_wpd; i++ )
436: d2[i] = d1[i];
437: break;
438: }
439: }
440:
1.1 noro 441: INLINE void ndl_add(unsigned int *d1,unsigned int *d2,unsigned int *d)
442: {
443: int i;
444:
1.6 noro 445: switch ( nd_wpd ) {
446: case 1:
447: d[0] = d1[0]+d2[0];
448: break;
449: case 2:
450: d[0] = d1[0]+d2[0];
451: d[1] = d1[1]+d2[1];
452: break;
453: default:
454: for ( i = 0; i < nd_wpd; i++ )
455: d[i] = d1[i]+d2[i];
456: break;
457: }
458: }
459:
460: INLINE void ndl_add2(unsigned int *d1,unsigned int *d2)
461: {
462: int i;
463:
464: switch ( nd_wpd ) {
465: case 1:
466: d2[0] += d1[0];
467: break;
468: case 2:
469: d2[0] += d1[0];
470: d2[1] += d1[1];
471: break;
472: default:
473: for ( i = 0; i < nd_wpd; i++ )
474: d2[i] += d1[i];
475: break;
1.1 noro 476: }
477: }
478:
479: void ndl_sub(unsigned int *d1,unsigned int *d2,unsigned int *d)
480: {
481: int i;
482:
483: for ( i = 0; i < nd_wpd; i++ )
484: d[i] = d1[i]-d2[i];
485: }
486:
487: int ndl_disjoint(unsigned int *d1,unsigned int *d2)
488: {
489: unsigned int t1,t2,u,u1,u2;
490: int i,j;
491:
492: switch ( nd_bpe ) {
493: case 4:
494: for ( i = 0; i < nd_wpd; i++ ) {
495: u1 = d1[i]; u2 = d2[i];
496: t1 = u1&0xf0000000; t2 = u2&0xf0000000; if ( t1&&t2 ) return 0;
497: t1 = u1&0xf000000; t2 = u2&0xf000000; if ( t1&&t2 ) return 0;
498: t1 = u1&0xf00000; t2 = u2&0xf00000; if ( t1&&t2 ) return 0;
499: t1 = u1&0xf0000; t2 = u2&0xf0000; if ( t1&&t2 ) return 0;
500: t1 = u1&0xf000; t2 = u2&0xf000; if ( t1&&t2 ) return 0;
501: t1 = u1&0xf00; t2 = u2&0xf00; if ( t1&&t2 ) return 0;
502: t1 = u1&0xf0; t2 = u2&0xf0; if ( t1&&t2 ) return 0;
503: t1 = u1&0xf; t2 = u2&0xf; if ( t1&&t2 ) return 0;
504: }
505: return 1;
506: break;
507: case 6:
508: for ( i = 0; i < nd_wpd; i++ ) {
509: u1 = d1[i]; u2 = d2[i];
510: t1 = u1&0x3f000000; t2 = u2&0x3f000000; if ( t1&&t2 ) return 0;
511: t1 = u1&0xfc0000; t2 = u2&0xfc0000; if ( t1&&t2 ) return 0;
512: t1 = u1&0x3f000; t2 = u2&0x3f000; if ( t1&&t2 ) return 0;
513: t1 = u1&0xfc0; t2 = u2&0xfc0; if ( t1&&t2 ) return 0;
514: t1 = u1&0x3f; t2 = u2&0x3f; if ( t1&&t2 ) return 0;
515: }
516: return 1;
517: break;
518: case 8:
519: for ( i = 0; i < nd_wpd; i++ ) {
520: u1 = d1[i]; u2 = d2[i];
521: t1 = u1&0xff000000; t2 = u2&0xff000000; if ( t1&&t2 ) return 0;
522: t1 = u1&0xff0000; t2 = u2&0xff0000; if ( t1&&t2 ) return 0;
523: t1 = u1&0xff00; t2 = u2&0xff00; if ( t1&&t2 ) return 0;
524: t1 = u1&0xff; t2 = u2&0xff; if ( t1&&t2 ) return 0;
525: }
526: return 1;
527: break;
528: case 16:
529: for ( i = 0; i < nd_wpd; i++ ) {
530: u1 = d1[i]; u2 = d2[i];
531: t1 = u1&0xffff0000; t2 = u2&0xffff0000; if ( t1&&t2 ) return 0;
532: t1 = u1&0xffff; t2 = u2&0xffff; if ( t1&&t2 ) return 0;
533: }
534: return 1;
535: break;
536: case 32:
537: for ( i = 0; i < nd_wpd; i++ )
538: if ( d1[i] && d2[i] ) return 0;
539: return 1;
540: break;
541: default:
542: for ( i = 0; i < nd_wpd; i++ ) {
543: u1 = d1[i]; u2 = d2[i];
544: for ( j = 0; j < nd_epw; j++ ) {
545: if ( (u1&nd_mask0) && (u2&nd_mask0) ) return 0;
546: u1 >>= nd_bpe; u2 >>= nd_bpe;
547: }
548: }
549: return 1;
550: break;
551: }
552: }
553:
1.5 noro 554: int ndl_check_bound2(int index,unsigned int *d2)
1.1 noro 555: {
1.5 noro 556: unsigned int u2;
557: unsigned int *d1;
558: int i,j,ind,k;
1.1 noro 559:
1.5 noro 560: d1 = nd_bound[index];
561: ind = 0;
562: switch ( nd_bpe ) {
563: case 4:
564: for ( i = 0; i < nd_wpd; i++ ) {
565: u2 = d2[i];
566: if ( d1[ind++]+((u2>>28)&0xf) >= 0x10 ) return 1;
567: if ( d1[ind++]+((u2>>24)&0xf) >= 0x10 ) return 1;
568: if ( d1[ind++]+((u2>>20)&0xf) >= 0x10 ) return 1;
569: if ( d1[ind++]+((u2>>16)&0xf) >= 0x10 ) return 1;
570: if ( d1[ind++]+((u2>>12)&0xf) >= 0x10 ) return 1;
571: if ( d1[ind++]+((u2>>8)&0xf) >= 0x10 ) return 1;
572: if ( d1[ind++]+((u2>>4)&0xf) >= 0x10 ) return 1;
573: if ( d1[ind++]+(u2&0xf) >= 0x10 ) return 1;
574: }
575: return 0;
576: break;
577: case 6:
578: for ( i = 0; i < nd_wpd; i++ ) {
579: u2 = d2[i];
580: if ( d1[ind++]+((u2>>24)&0x3f) >= 0x40 ) return 1;
581: if ( d1[ind++]+((u2>>18)&0x3f) >= 0x40 ) return 1;
582: if ( d1[ind++]+((u2>>12)&0x3f) >= 0x40 ) return 1;
583: if ( d1[ind++]+((u2>>6)&0x3f) >= 0x40 ) return 1;
584: if ( d1[ind++]+(u2&0x3f) >= 0x40 ) return 1;
585: }
586: return 0;
587: break;
588: case 8:
589: for ( i = 0; i < nd_wpd; i++ ) {
590: u2 = d2[i];
591: if ( d1[ind++]+((u2>>24)&0xff) >= 0x100 ) return 1;
592: if ( d1[ind++]+((u2>>16)&0xff) >= 0x100 ) return 1;
593: if ( d1[ind++]+((u2>>8)&0xff) >= 0x100 ) return 1;
594: if ( d1[ind++]+(u2&0xff) >= 0x100 ) return 1;
595: }
596: return 0;
597: break;
598: case 16:
599: for ( i = 0; i < nd_wpd; i++ ) {
600: u2 = d2[i];
601: if ( d1[ind++]+((u2>>16)&0xffff) > 0x10000 ) return 1;
602: if ( d1[ind++]+(u2&0xffff) > 0x10000 ) return 1;
603: }
604: return 0;
605: break;
606: case 32:
607: for ( i = 0; i < nd_wpd; i++ )
608: if ( d1[i]+d2[i]<d1[i] ) return 1;
609: return 0;
610: break;
611: default:
612: for ( i = 0; i < nd_wpd; i++ ) {
613: u2 = d2[i];
614: k = (nd_epw-1)*nd_bpe;
615: for ( j = 0; j < nd_epw; j++, k -= nd_bpe )
616: if ( d1[ind++]+((u2>>k)&nd_mask0) > nd_mask0 ) return 1;
617: }
618: return 0;
619: break;
620: }
1.1 noro 621: }
622:
1.6 noro 623: INLINE int ndl_hash_value(int td,unsigned int *d)
1.1 noro 624: {
625: int i;
626: int r;
627:
628: r = td;
629: for ( i = 0; i < nd_wpd; i++ )
630: r = ((r<<16)+d[i])%REDTAB_LEN;
631: return r;
632: }
633:
1.9 noro 634: INLINE int nd_find_reducer(ND g)
1.1 noro 635: {
1.13 noro 636: RHist r;
1.6 noro 637: int d,k,i;
1.1 noro 638:
639: d = ndl_hash_value(HTD(g),HDL(g));
1.13 noro 640: for ( r = nd_red[d], k = 0; r; r = NEXT(r), k++ ) {
1.14 noro 641: if ( HTD(g) == TD(r) && ndl_equal(HDL(g),DL(r)) ) {
1.1 noro 642: if ( k > 0 ) nd_notfirst++;
643: nd_found++;
1.13 noro 644: return r->index;
1.1 noro 645: }
646: }
647:
1.13 noro 648: if ( Reverse )
649: for ( i = nd_psn-1; i >= 0; i-- ) {
650: r = nd_psh[i];
1.14 noro 651: if ( HTD(g) >= TD(r) && ndl_reducible(HDL(g),DL(r)) ) {
1.13 noro 652: nd_create++;
653: nd_append_red(HDL(g),HTD(g),i);
654: return i;
655: }
656: }
657: else
658: for ( i = 0; i < nd_psn; i++ ) {
659: r = nd_psh[i];
1.14 noro 660: if ( HTD(g) >= TD(r) && ndl_reducible(HDL(g),DL(r)) ) {
1.13 noro 661: nd_create++;
662: nd_append_red(HDL(g),HTD(g),i);
663: return i;
664: }
1.1 noro 665: }
1.6 noro 666: return -1;
1.1 noro 667: }
668:
1.19 noro 669: ND nd_add(int mod,ND p1,ND p2)
1.1 noro 670: {
671: int n,c;
672: int t;
673: ND r;
674: NM m1,m2,mr0,mr,s;
675:
676: if ( !p1 )
677: return p2;
678: else if ( !p2 )
679: return p1;
680: else {
681: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) {
1.14 noro 682: if ( TD(m1) > TD(m2) )
1.1 noro 683: c = 1;
1.14 noro 684: else if ( TD(m1) < TD(m2) )
1.1 noro 685: c = -1;
686: else
1.14 noro 687: c = ndl_compare(DL(m1),DL(m2));
1.1 noro 688: switch ( c ) {
689: case 0:
1.19 noro 690: t = ((CM(m1))+(CM(m2))) - mod;
1.1 noro 691: if ( t < 0 )
1.19 noro 692: t += mod;
1.1 noro 693: s = m1; m1 = NEXT(m1);
694: if ( t ) {
1.14 noro 695: NEXTNM2(mr0,mr,s); CM(mr) = (t);
1.1 noro 696: } else {
697: FREENM(s);
698: }
699: s = m2; m2 = NEXT(m2); FREENM(s);
700: break;
701: case 1:
702: s = m1; m1 = NEXT(m1); NEXTNM2(mr0,mr,s);
703: break;
704: case -1:
705: s = m2; m2 = NEXT(m2); NEXTNM2(mr0,mr,s);
706: break;
707: }
708: }
709: if ( !mr0 )
710: if ( m1 )
711: mr0 = m1;
712: else if ( m2 )
713: mr0 = m2;
714: else
715: return 0;
716: else if ( m1 )
717: NEXT(mr) = m1;
718: else if ( m2 )
719: NEXT(mr) = m2;
720: else
721: NEXT(mr) = 0;
722: BDY(p1) = mr0;
1.14 noro 723: SG(p1) = MAX(SG(p1),SG(p2));
1.1 noro 724: FREEND(p2);
725: return p1;
726: }
727: }
728:
1.17 noro 729: ND nd_add_q(ND p1,ND p2)
730: {
731: int n,c;
732: ND r;
733: NM m1,m2,mr0,mr,s;
734: Q t;
735:
736: if ( !p1 )
737: return p2;
738: else if ( !p2 )
739: return p1;
740: else {
741: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) {
742: if ( TD(m1) > TD(m2) )
743: c = 1;
744: else if ( TD(m1) < TD(m2) )
745: c = -1;
746: else
747: c = ndl_compare(DL(m1),DL(m2));
748: switch ( c ) {
749: case 0:
750: addq(CQ(m1),CQ(m2),&t);
751: s = m1; m1 = NEXT(m1);
752: if ( t ) {
753: NEXTNM2(mr0,mr,s); CQ(mr) = (t);
754: } else {
755: FREENM(s);
756: }
757: s = m2; m2 = NEXT(m2); FREENM(s);
758: break;
759: case 1:
760: s = m1; m1 = NEXT(m1); NEXTNM2(mr0,mr,s);
761: break;
762: case -1:
763: s = m2; m2 = NEXT(m2); NEXTNM2(mr0,mr,s);
764: break;
765: }
766: }
767: if ( !mr0 )
768: if ( m1 )
769: mr0 = m1;
770: else if ( m2 )
771: mr0 = m2;
772: else
773: return 0;
774: else if ( m1 )
775: NEXT(mr) = m1;
776: else if ( m2 )
777: NEXT(mr) = m2;
778: else
779: NEXT(mr) = 0;
780: BDY(p1) = mr0;
781: SG(p1) = MAX(SG(p1),SG(p2));
782: FREEND(p2);
783: return p1;
784: }
785: }
786:
1.1 noro 787: #if 1
788: /* ret=1 : success, ret=0 : overflow */
1.16 noro 789: int nd_nf(int mod,ND g,int full,ND *rp)
1.1 noro 790: {
1.11 noro 791: ND d;
1.1 noro 792: NM m,mrd,tail;
1.7 noro 793: NM mul;
1.10 noro 794: int n,sugar,psugar,sugar0,stat,index;
1.6 noro 795: int c,c1,c2;
1.17 noro 796: RHist h;
1.11 noro 797: NDV p,red;
1.16 noro 798: Q cg,cred,gcd;
1.21 noro 799: double hmag;
1.1 noro 800:
801: if ( !g ) {
802: *rp = 0;
803: return 1;
804: }
1.21 noro 805: if ( !mod )
806: hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;
807:
1.14 noro 808: sugar0 = sugar = SG(g);
1.1 noro 809: n = NV(g);
1.7 noro 810: mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int));
1.1 noro 811: for ( d = 0; g; ) {
1.6 noro 812: index = nd_find_reducer(g);
813: if ( index >= 0 ) {
1.17 noro 814: h = nd_psh[index];
815: ndl_sub(HDL(g),DL(h),DL(mul));
816: TD(mul) = HTD(g)-TD(h);
1.10 noro 817: #if 0
1.14 noro 818: if ( d && (SG(p)+TD(mul)) > sugar ) {
1.10 noro 819: goto afo;
820: }
821: #endif
1.14 noro 822: if ( ndl_check_bound2(index,DL(mul)) ) {
1.6 noro 823: nd_free(g); nd_free(d);
824: return 0;
825: }
1.16 noro 826: if ( mod ) {
1.17 noro 827: p = nd_ps[index];
1.19 noro 828: c1 = invm(HCM(p),mod); c2 = mod-HCM(g);
829: DMAR(c1,c2,0,mod,c); CM(mul) = c;
1.16 noro 830: } else {
1.20 noro 831: p = nd_psq[index];
1.17 noro 832: igcd_cofactor(HCQ(g),HCQ(p),&gcd,&cg,&cred);
1.16 noro 833: chsgnq(cg,&CQ(mul));
1.20 noro 834: nd_mul_c_q(d,cred); nd_mul_c_q(g,cred);
1.16 noro 835: }
1.20 noro 836: ndv_mul_nm(mod,p,mul,ndv_red);
837: g = ndv_add(mod,g,ndv_red);
1.14 noro 838: sugar = MAX(sugar,SG(ndv_red));
1.22 ! noro 839: if ( !mod && hmag && g && ((double)(p_mag((P)HCQ(g))) > hmag) ) {
1.21 noro 840: nd_removecont2(d,g);
841: hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;
842: }
1.1 noro 843: } else if ( !full ) {
844: *rp = g;
845: return 1;
846: } else {
1.10 noro 847: afo:
1.1 noro 848: m = BDY(g);
849: if ( NEXT(m) ) {
850: BDY(g) = NEXT(m); NEXT(m) = 0;
851: } else {
852: FREEND(g); g = 0;
853: }
854: if ( d ) {
855: NEXT(tail)=m;
856: tail=m;
857: } else {
858: MKND(n,m,d);
859: tail = BDY(d);
860: }
861: }
862: }
863: if ( d )
1.14 noro 864: SG(d) = sugar;
1.1 noro 865: *rp = d;
866: return 1;
867: }
868: #else
869:
870: ND nd_remove_head(ND p)
871: {
872: NM m;
873:
874: m = BDY(p);
875: if ( !NEXT(m) ) {
876: FREEND(p);
877: p = 0;
878: } else
879: BDY(p) = NEXT(m);
880: FREENM(m);
881: return p;
882: }
883:
884: PGeoBucket create_pbucket()
885: {
886: PGeoBucket g;
887:
888: g = CALLOC(1,sizeof(struct oPGeoBucket));
889: g->m = -1;
890: return g;
891: }
892:
1.19 noro 893: void add_pbucket(int mod,PGeoBucket g,ND d)
1.1 noro 894: {
895: int l,k,m;
896:
897: l = nd_length(d);
898: for ( k = 0, m = 1; l > m; k++, m <<= 2 );
899: /* 4^(k-1) < l <= 4^k */
1.19 noro 900: d = nd_add(mod,g->body[k],d);
1.1 noro 901: for ( ; d && nd_length(d) > 1<<(2*k); k++ ) {
902: g->body[k] = 0;
1.19 noro 903: d = nd_add(int mod,g->body[k+1],d);
1.1 noro 904: }
905: g->body[k] = d;
906: g->m = MAX(g->m,k);
907: }
908:
1.19 noro 909: int head_pbucket(int mod,PGeoBucket g)
1.1 noro 910: {
911: int j,i,c,k,nv,sum;
912: unsigned int *di,*dj;
913: ND gi,gj;
914:
915: k = g->m;
916: while ( 1 ) {
917: j = -1;
918: for ( i = 0; i <= k; i++ ) {
919: if ( !(gi = g->body[i]) )
920: continue;
921: if ( j < 0 ) {
922: j = i;
923: gj = g->body[j];
924: dj = HDL(gj);
1.14 noro 925: sum = HCM(gj);
1.1 noro 926: } else {
927: di = HDL(gi);
928: nv = NV(gi);
929: if ( HTD(gi) > HTD(gj) )
930: c = 1;
931: else if ( HTD(gi) < HTD(gj) )
932: c = -1;
933: else
934: c = ndl_compare(di,dj);
935: if ( c > 0 ) {
936: if ( sum )
1.14 noro 937: HCM(gj) = sum;
1.1 noro 938: else
939: g->body[j] = nd_remove_head(gj);
940: j = i;
941: gj = g->body[j];
942: dj = HDL(gj);
1.14 noro 943: sum = HCM(gj);
1.1 noro 944: } else if ( c == 0 ) {
1.19 noro 945: sum = sum+HCM(gi)-mod;
1.1 noro 946: if ( sum < 0 )
1.19 noro 947: sum += mod;
1.1 noro 948: g->body[i] = nd_remove_head(gi);
949: }
950: }
951: }
952: if ( j < 0 )
953: return -1;
954: else if ( sum ) {
1.14 noro 955: HCM(gj) = sum;
1.1 noro 956: return j;
957: } else
958: g->body[j] = nd_remove_head(gj);
959: }
960: }
961:
962: ND normalize_pbucket(PGeoBucket g)
963: {
964: int i;
965: ND r,t;
966:
967: r = 0;
968: for ( i = 0; i <= g->m; i++ )
1.19 noro 969: r = nd_add(mod,r,g->body[i]);
1.1 noro 970: return r;
971: }
972:
973: ND nd_nf(ND g,int full)
974: {
975: ND u,p,d,red;
976: NODE l;
977: NM m,mrd;
978: int sugar,psugar,n,h_reducible,h;
979: PGeoBucket bucket;
980:
981: if ( !g ) {
982: return 0;
983: }
1.14 noro 984: sugar = SG(g);
985: n = NV(g);
1.1 noro 986: bucket = create_pbucket();
987: add_pbucket(bucket,g);
988: d = 0;
989: while ( 1 ) {
990: h = head_pbucket(bucket);
991: if ( h < 0 ) {
992: if ( d )
1.14 noro 993: SG(d) = sugar;
1.1 noro 994: return d;
995: }
996: g = bucket->body[h];
997: red = nd_find_reducer(g);
998: if ( red ) {
999: bucket->body[h] = nd_remove_head(g);
1000: red = nd_remove_head(red);
1001: add_pbucket(bucket,red);
1.14 noro 1002: sugar = MAX(sugar,SG(red));
1.1 noro 1003: } else if ( !full ) {
1004: g = normalize_pbucket(bucket);
1005: if ( g )
1.14 noro 1006: SG(g) = sugar;
1.1 noro 1007: return g;
1008: } else {
1009: m = BDY(g);
1010: if ( NEXT(m) ) {
1011: BDY(g) = NEXT(m); NEXT(m) = 0;
1012: } else {
1013: FREEND(g); g = 0;
1014: }
1015: bucket->body[h] = g;
1016: NEXT(m) = 0;
1017: if ( d ) {
1018: for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) );
1019: NEXT(mrd) = m;
1020: } else {
1021: MKND(n,m,d);
1022: }
1023: }
1024: }
1025: }
1026: #endif
1027:
1.16 noro 1028: NODE nd_gb(int m,NODE f)
1.1 noro 1029: {
1030: int i,nh,sugar,stat;
1031: NODE r,g,gall;
1032: ND_pairs d;
1033: ND_pairs l;
1034: ND h,nf;
1035:
1036: for ( gall = g = 0, d = 0, r = f; r; r = NEXT(r) ) {
1037: i = (int)BDY(r);
1038: d = update_pairs(d,g,i);
1039: g = update_base(g,i);
1040: gall = append_one(gall,i);
1041: }
1042: sugar = 0;
1043: while ( d ) {
1044: again:
1045: l = nd_minp(d,&d);
1.14 noro 1046: if ( SG(l) != sugar ) {
1047: sugar = SG(l);
1.1 noro 1048: fprintf(asir_out,"%d",sugar);
1049: }
1.16 noro 1050: stat = nd_sp(m,l,&h);
1.1 noro 1051: if ( !stat ) {
1052: NEXT(l) = d; d = l;
1.20 noro 1053: d = nd_reconstruct(m,0,d);
1.1 noro 1054: goto again;
1055: }
1.16 noro 1056: stat = nd_nf(m,h,!Top,&nf);
1.1 noro 1057: if ( !stat ) {
1058: NEXT(l) = d; d = l;
1.20 noro 1059: d = nd_reconstruct(m,0,d);
1.1 noro 1060: goto again;
1061: } else if ( nf ) {
1062: printf("+"); fflush(stdout);
1.16 noro 1063: nh = nd_newps(m,nf);
1.1 noro 1064: d = update_pairs(d,g,nh);
1065: g = update_base(g,nh);
1066: gall = append_one(gall,nh);
1067: FREENDP(l);
1068: } else {
1069: printf("."); fflush(stdout);
1070: FREENDP(l);
1071: }
1072: }
1073: return g;
1074: }
1075:
1.20 noro 1076: NODE nd_gb_trace(int m,NODE f)
1077: {
1078: int i,nh,sugar,stat;
1079: NODE r,g,gall;
1080: ND_pairs d;
1081: ND_pairs l;
1082: ND h,nf,nfq;
1083:
1084: for ( gall = g = 0, d = 0, r = f; r; r = NEXT(r) ) {
1085: i = (int)BDY(r);
1086: d = update_pairs(d,g,i);
1087: g = update_base(g,i);
1088: gall = append_one(gall,i);
1089: }
1090: sugar = 0;
1091: while ( d ) {
1092: again:
1093: l = nd_minp(d,&d);
1094: if ( SG(l) != sugar ) {
1095: sugar = SG(l);
1096: fprintf(asir_out,"%d",sugar);
1097: }
1098: stat = nd_sp(m,l,&h);
1099: if ( !stat ) {
1100: NEXT(l) = d; d = l;
1101: d = nd_reconstruct(m,1,d);
1102: goto again;
1103: }
1104: stat = nd_nf(m,h,!Top,&nf);
1105: if ( !stat ) {
1106: NEXT(l) = d; d = l;
1107: d = nd_reconstruct(m,1,d);
1108: goto again;
1109: } else if ( nf ) {
1110: /* overflow does not occur */
1111: nd_sp(0,l,&h);
1112: nd_nf(0,h,!Top,&nfq);
1113: if ( nfq ) {
1114: printf("+"); fflush(stdout);
1115: nh = nd_newps_trace(m,nf,nfq);
1116: d = update_pairs(d,g,nh);
1117: g = update_base(g,nh);
1118: gall = append_one(gall,nh);
1119: } else {
1120: printf("*"); fflush(stdout);
1121: }
1122: } else {
1123: printf("."); fflush(stdout);
1124: }
1125: FREENDP(l);
1126: }
1127: return g;
1128: }
1129:
1.1 noro 1130: ND_pairs update_pairs( ND_pairs d, NODE /* of index */ g, int t)
1131: {
1132: ND_pairs d1,nd,cur,head,prev,remove;
1133:
1134: if ( !g ) return d;
1135: d = crit_B(d,t);
1136: d1 = nd_newpairs(g,t);
1137: d1 = crit_M(d1);
1138: d1 = crit_F(d1);
1139: prev = 0; cur = head = d1;
1140: while ( cur ) {
1141: if ( crit_2( cur->i1,cur->i2 ) ) {
1142: remove = cur;
1143: if ( !prev ) {
1144: head = cur = NEXT(cur);
1145: } else {
1146: cur = NEXT(prev) = NEXT(cur);
1147: }
1148: FREENDP(remove);
1149: } else {
1150: prev = cur;
1151: cur = NEXT(cur);
1152: }
1153: }
1154: if ( !d )
1155: return head;
1156: else {
1157: nd = d;
1158: while ( NEXT(nd) )
1159: nd = NEXT(nd);
1160: NEXT(nd) = head;
1161: return d;
1162: }
1163: }
1164:
1165: ND_pairs nd_newpairs( NODE g, int t )
1166: {
1167: NODE h;
1168: unsigned int *dl;
1169: int td,ts,s;
1170: ND_pairs r,r0;
1171:
1.20 noro 1172: dl = DL(nd_psh[t]);
1173: td = TD(nd_psh[t]);
1174: ts = SG(nd_psh[t]) - td;
1.1 noro 1175: for ( r0 = 0, h = g; h; h = NEXT(h) ) {
1176: NEXTND_pairs(r0,r);
1177: r->i1 = (int)BDY(h);
1178: r->i2 = t;
1.20 noro 1179: ndl_lcm(DL(nd_psh[r->i1]),dl,r->lcm);
1.14 noro 1180: TD(r) = ndl_td(r->lcm);
1.20 noro 1181: s = SG(nd_psh[r->i1])-TD(nd_psh[r->i1]);
1.14 noro 1182: SG(r) = MAX(s,ts) + TD(r);
1.1 noro 1183: }
1184: NEXT(r) = 0;
1185: return r0;
1186: }
1187:
1188: ND_pairs crit_B( ND_pairs d, int s )
1189: {
1190: ND_pairs cur,head,prev,remove;
1191: unsigned int *t,*tl,*lcm;
1192: int td,tdl;
1193:
1194: if ( !d ) return 0;
1.20 noro 1195: t = DL(nd_psh[s]);
1.1 noro 1196: prev = 0;
1197: head = cur = d;
1198: lcm = (unsigned int *)ALLOCA(nd_wpd*sizeof(unsigned int));
1199: while ( cur ) {
1200: tl = cur->lcm;
1201: if ( ndl_reducible(tl,t)
1.20 noro 1202: && (ndl_lcm(DL(nd_psh[cur->i1]),t,lcm),!ndl_equal(lcm,tl))
1203: && (ndl_lcm(DL(nd_psh[cur->i2]),t,lcm),!ndl_equal(lcm,tl)) ) {
1.1 noro 1204: remove = cur;
1205: if ( !prev ) {
1206: head = cur = NEXT(cur);
1207: } else {
1208: cur = NEXT(prev) = NEXT(cur);
1209: }
1210: FREENDP(remove);
1211: } else {
1212: prev = cur;
1213: cur = NEXT(cur);
1214: }
1215: }
1216: return head;
1217: }
1218:
1219: ND_pairs crit_M( ND_pairs d1 )
1220: {
1221: ND_pairs e,d2,d3,dd,p;
1222: unsigned int *id,*jd;
1223: int itd,jtd;
1224:
1225: for ( dd = 0, e = d1; e; e = d3 ) {
1226: if ( !(d2 = NEXT(e)) ) {
1227: NEXT(e) = dd;
1228: return e;
1229: }
1230: id = e->lcm;
1.14 noro 1231: itd = TD(e);
1.1 noro 1232: for ( d3 = 0; d2; d2 = p ) {
1233: p = NEXT(d2),
1234: jd = d2->lcm;
1.14 noro 1235: jtd = TD(d2);
1.1 noro 1236: if ( jtd == itd )
1237: if ( id == jd );
1238: else if ( ndl_reducible(jd,id) ) continue;
1239: else if ( ndl_reducible(id,jd) ) goto delit;
1240: else ;
1241: else if ( jtd > itd )
1242: if ( ndl_reducible(jd,id) ) continue;
1243: else ;
1244: else if ( ndl_reducible(id,jd ) ) goto delit;
1245: NEXT(d2) = d3;
1246: d3 = d2;
1247: }
1248: NEXT(e) = dd;
1249: dd = e;
1250: continue;
1251: /**/
1252: delit: NEXT(d2) = d3;
1253: d3 = d2;
1254: for ( ; p; p = d2 ) {
1255: d2 = NEXT(p);
1256: NEXT(p) = d3;
1257: d3 = p;
1258: }
1259: FREENDP(e);
1260: }
1261: return dd;
1262: }
1263:
1264: ND_pairs crit_F( ND_pairs d1 )
1265: {
1266: ND_pairs rest, head,remove;
1267: ND_pairs last, p, r, w;
1268: int s;
1269:
1270: for ( head = last = 0, p = d1; NEXT(p); ) {
1271: r = w = equivalent_pairs(p,&rest);
1.14 noro 1272: s = SG(r);
1.1 noro 1273: w = NEXT(w);
1274: while ( w ) {
1275: if ( crit_2(w->i1,w->i2) ) {
1276: r = w;
1277: w = NEXT(w);
1278: while ( w ) {
1279: remove = w;
1280: w = NEXT(w);
1281: FREENDP(remove);
1282: }
1283: break;
1.14 noro 1284: } else if ( SG(w) < s ) {
1.1 noro 1285: FREENDP(r);
1286: r = w;
1.14 noro 1287: s = SG(r);
1.1 noro 1288: w = NEXT(w);
1289: } else {
1290: remove = w;
1291: w = NEXT(w);
1292: FREENDP(remove);
1293: }
1294: }
1295: if ( last ) NEXT(last) = r;
1296: else head = r;
1297: NEXT(last = r) = 0;
1298: p = rest;
1299: if ( !p ) return head;
1300: }
1301: if ( !last ) return p;
1302: NEXT(last) = p;
1303: return head;
1304: }
1305:
1306: int crit_2( int dp1, int dp2 )
1307: {
1.20 noro 1308: return ndl_disjoint(DL(nd_psh[dp1]),DL(nd_psh[dp2]));
1.1 noro 1309: }
1310:
1311: static ND_pairs equivalent_pairs( ND_pairs d1, ND_pairs *prest )
1312: {
1313: ND_pairs w,p,r,s;
1314: unsigned int *d;
1315: int td;
1316:
1317: w = d1;
1318: d = w->lcm;
1.14 noro 1319: td = TD(w);
1.1 noro 1320: s = NEXT(w);
1321: NEXT(w) = 0;
1322: for ( r = 0; s; s = p ) {
1323: p = NEXT(s);
1.14 noro 1324: if ( td == TD(s) && ndl_equal(d,s->lcm) ) {
1.1 noro 1325: NEXT(s) = w;
1326: w = s;
1327: } else {
1328: NEXT(s) = r;
1329: r = s;
1330: }
1331: }
1332: *prest = r;
1333: return w;
1334: }
1335:
1336: NODE update_base(NODE nd,int ndp)
1337: {
1338: unsigned int *dl, *dln;
1339: NODE last, p, head;
1340: int td,tdn;
1341:
1.20 noro 1342: dl = DL(nd_psh[ndp]);
1343: td = TD(nd_psh[ndp]);
1.1 noro 1344: for ( head = last = 0, p = nd; p; ) {
1.20 noro 1345: dln = DL(nd_psh[(int)BDY(p)]);
1346: tdn = TD(nd_psh[(int)BDY(p)]);
1.1 noro 1347: if ( tdn >= td && ndl_reducible( dln, dl ) ) {
1348: p = NEXT(p);
1349: if ( last ) NEXT(last) = p;
1350: } else {
1351: if ( !last ) head = p;
1352: p = NEXT(last = p);
1353: }
1354: }
1355: head = append_one(head,ndp);
1356: return head;
1357: }
1358:
1359: ND_pairs nd_minp( ND_pairs d, ND_pairs *prest )
1360: {
1361: ND_pairs m,ml,p,l;
1362: unsigned int *lcm;
1363: int s,td,len,tlen,c;
1364:
1365: if ( !(p = NEXT(m = d)) ) {
1366: *prest = p;
1367: NEXT(m) = 0;
1368: return m;
1369: }
1370: lcm = m->lcm;
1.14 noro 1371: s = SG(m);
1372: td = TD(m);
1.6 noro 1373: len = nd_psl[m->i1]+nd_psl[m->i2];
1.1 noro 1374: for ( ml = 0, l = m; p; p = NEXT(l = p) ) {
1.14 noro 1375: if (SG(p) < s)
1.1 noro 1376: goto find;
1.14 noro 1377: else if ( SG(p) == s ) {
1378: if ( TD(p) < td )
1.1 noro 1379: goto find;
1.14 noro 1380: else if ( TD(p) == td ) {
1.1 noro 1381: c = ndl_compare(p->lcm,lcm);
1382: if ( c < 0 )
1383: goto find;
1.10 noro 1384: #if 0
1.1 noro 1385: else if ( c == 0 ) {
1.6 noro 1386: tlen = nd_psl[p->i1]+nd_psl[p->i2];
1.1 noro 1387: if ( tlen < len )
1388: goto find;
1389: }
1.10 noro 1390: #endif
1.1 noro 1391: }
1392: }
1393: continue;
1394: find:
1395: ml = l;
1396: m = p;
1397: lcm = m->lcm;
1.14 noro 1398: s = SG(m);
1399: td = TD(m);
1.1 noro 1400: len = tlen;
1401: }
1402: if ( !ml ) *prest = NEXT(m);
1403: else {
1404: NEXT(ml) = NEXT(m);
1405: *prest = d;
1406: }
1407: NEXT(m) = 0;
1408: return m;
1409: }
1410:
1.16 noro 1411: int nd_newps(int mod,ND a)
1.1 noro 1412: {
1.3 noro 1413: int len;
1.13 noro 1414: RHist r;
1.20 noro 1415: NDV b;
1.3 noro 1416:
1.1 noro 1417: if ( nd_psn == nd_pslen ) {
1418: nd_pslen *= 2;
1.6 noro 1419: nd_psl = (int *)REALLOC((char *)nd_psl,nd_pslen*sizeof(int));
1.11 noro 1420: nd_ps = (NDV *)REALLOC((char *)nd_ps,nd_pslen*sizeof(NDV));
1.20 noro 1421: nd_psq = (NDV *)REALLOC((char *)nd_psq,nd_pslen*sizeof(NDV));
1.13 noro 1422: nd_psh = (RHist *)REALLOC((char *)nd_psh,nd_pslen*sizeof(RHist));
1.1 noro 1423: nd_bound = (unsigned int **)
1424: REALLOC((char *)nd_bound,nd_pslen*sizeof(unsigned int *));
1425: }
1.20 noro 1426: nd_removecont(mod,a);
1427: nd_bound[nd_psn] = nd_compute_bound(a);
1428: NEWRHist(r); SG(r) = SG(a); TD(r) = HTD(a); ndl_copy(HDL(a),DL(r));
1429: nd_psh[nd_psn] = r;
1430: b = ndtondv(mod,a);
1431: len = LEN(b);
1.16 noro 1432: if ( mod )
1.20 noro 1433: nd_ps[nd_psn] = b;
1.16 noro 1434: else
1.20 noro 1435: nd_psq[nd_psn] = b;
1436: nd_psl[nd_psn] = len;
1.13 noro 1437: nd_free(a);
1.20 noro 1438: if ( len > nmv_len ) {
1439: nmv_len = 2*len;
1440: BDY(ndv_red) = (NMV)REALLOC(BDY(ndv_red),nmv_len*nmv_adv);
1441: }
1442: return nd_psn++;
1443: }
1444:
1445: int nd_newps_trace(int mod,ND nf,ND nfq)
1446: {
1447: int len;
1448: RHist r;
1449: NDV b;
1450:
1451: if ( nd_psn == nd_pslen ) {
1452: nd_pslen *= 2;
1453: nd_psl = (int *)REALLOC((char *)nd_psl,nd_pslen*sizeof(int));
1454: nd_ps = (NDV *)REALLOC((char *)nd_ps,nd_pslen*sizeof(NDV));
1455: nd_psq = (NDV *)REALLOC((char *)nd_psq,nd_pslen*sizeof(NDV));
1456: nd_psh = (RHist *)REALLOC((char *)nd_psh,nd_pslen*sizeof(RHist));
1457: nd_bound = (unsigned int **)
1458: REALLOC((char *)nd_bound,nd_pslen*sizeof(unsigned int *));
1459: }
1460: nd_removecont(mod,nf);
1461: nd_ps[nd_psn] = ndtondv(mod,nf);
1462:
1463: nd_removecont(0,nfq);
1464: nd_psq[nd_psn] = ndtondv(0,nfq);
1465:
1466: nd_bound[nd_psn] = nd_compute_bound(nfq);
1467: NEWRHist(r); SG(r) = SG(nf); TD(r) = HTD(nf); ndl_copy(HDL(nf),DL(r));
1468: nd_psh[nd_psn] = r;
1469:
1470: len = LEN(nd_psq[nd_psn]);
1471: nd_psl[nd_psn] = len;
1472:
1473: nd_free(nf); nd_free(nfq);
1.3 noro 1474: if ( len > nmv_len ) {
1475: nmv_len = 2*len;
1476: BDY(ndv_red) = (NMV)REALLOC(BDY(ndv_red),nmv_len*nmv_adv);
1477: }
1.1 noro 1478: return nd_psn++;
1479: }
1480:
1.16 noro 1481: NODE nd_setup(int mod,NODE f)
1.1 noro 1482: {
1.5 noro 1483: int i,j,td,len,max;
1.1 noro 1484: NODE s,s0,f0;
1.5 noro 1485: unsigned int *d;
1.13 noro 1486: RHist r;
1.20 noro 1487: NDV a;
1.11 noro 1488:
1489: nd_found = 0; nd_notfirst = 0; nd_create = 0;
1.1 noro 1490:
1491: nd_psn = length(f); nd_pslen = 2*nd_psn;
1.6 noro 1492: nd_psl = (int *)MALLOC(nd_pslen*sizeof(int));
1.11 noro 1493: nd_ps = (NDV *)MALLOC(nd_pslen*sizeof(NDV));
1.20 noro 1494: nd_psq = (NDV *)MALLOC(nd_pslen*sizeof(NDV));
1.13 noro 1495: nd_psh = (RHist *)MALLOC(nd_pslen*sizeof(RHist));
1.1 noro 1496: nd_bound = (unsigned int **)MALLOC(nd_pslen*sizeof(unsigned int *));
1.5 noro 1497: for ( max = 0, i = 0, s = f; i < nd_psn; i++, s = NEXT(s) ) {
1498: nd_bound[i] = d = dp_compute_bound((DP)BDY(s));
1499: for ( j = 0; j < nd_nvar; j++ )
1500: max = MAX(d[j],max);
1501: }
1.11 noro 1502: if ( !nd_red )
1.13 noro 1503: nd_red = (RHist *)MALLOC(REDTAB_LEN*sizeof(RHist));
1504: bzero(nd_red,REDTAB_LEN*sizeof(RHist));
1.5 noro 1505:
1506: if ( max < 2 )
1507: nd_bpe = 2;
1508: else if ( max < 4 )
1509: nd_bpe = 4;
1510: else if ( max < 64 )
1511: nd_bpe = 6;
1512: else if ( max < 256 )
1513: nd_bpe = 8;
1514: else if ( max < 65536 )
1515: nd_bpe = 16;
1516: else
1517: nd_bpe = 32;
1.13 noro 1518:
1.1 noro 1519: nd_setup_parameters();
1520: nd_free_private_storage();
1.3 noro 1521: len = 0;
1.1 noro 1522: for ( i = 0; i < nd_psn; i++, f = NEXT(f) ) {
1.20 noro 1523: NEWRHist(r);
1524: a = dptondv(mod,(DP)BDY(f));
1525: ndv_removecont(mod,a);
1526: len = MAX(len,LEN(a));
1527: SG(r) = HTD(a); TD(r) = HTD(a); ndl_copy(HDL(a),DL(r));
1.16 noro 1528: if ( mod )
1.20 noro 1529: nd_ps[i] = a;
1.16 noro 1530: else
1.20 noro 1531: nd_psq[i] = a;
1.13 noro 1532: nd_psh[i] = r;
1.1 noro 1533: }
1.3 noro 1534: nmv_len = 16*len;
1535: NEWNDV(ndv_red);
1.16 noro 1536: if ( mod )
1537: BDY(ndv_red) = (NMV)MALLOC_ATOMIC(nmv_len*nmv_adv);
1538: else
1539: BDY(ndv_red) = (NMV)MALLOC(nmv_len*nmv_adv);
1.1 noro 1540: for ( s0 = 0, i = 0; i < nd_psn; i++ ) {
1541: NEXTNODE(s0,s); BDY(s) = (pointer)i;
1542: }
1543: if ( s0 ) NEXT(s) = 0;
1544: return s0;
1545: }
1546:
1.20 noro 1547: NODE nd_setup_trace(int mod,NODE f)
1548: {
1549: int i,j,td,len,max;
1550: NODE s,s0,f0;
1551: unsigned int *d;
1552: RHist r;
1553: NDV a;
1554:
1555: nd_found = 0; nd_notfirst = 0; nd_create = 0;
1556:
1557: nd_psn = length(f); nd_pslen = 2*nd_psn;
1558: nd_psl = (int *)MALLOC(nd_pslen*sizeof(int));
1559: nd_ps = (NDV *)MALLOC(nd_pslen*sizeof(NDV));
1560: nd_psq = (NDV *)MALLOC(nd_pslen*sizeof(NDV));
1561: nd_psh = (RHist *)MALLOC(nd_pslen*sizeof(RHist));
1562: nd_bound = (unsigned int **)MALLOC(nd_pslen*sizeof(unsigned int *));
1563: for ( max = 0, i = 0, s = f; i < nd_psn; i++, s = NEXT(s) ) {
1564: nd_bound[i] = d = dp_compute_bound((DP)BDY(s));
1565: for ( j = 0; j < nd_nvar; j++ )
1566: max = MAX(d[j],max);
1567: }
1568: if ( !nd_red )
1569: nd_red = (RHist *)MALLOC(REDTAB_LEN*sizeof(RHist));
1570: bzero(nd_red,REDTAB_LEN*sizeof(RHist));
1571:
1572: if ( max < 2 )
1573: nd_bpe = 2;
1574: else if ( max < 4 )
1575: nd_bpe = 4;
1576: else if ( max < 64 )
1577: nd_bpe = 6;
1578: else if ( max < 256 )
1579: nd_bpe = 8;
1580: else if ( max < 65536 )
1581: nd_bpe = 16;
1582: else
1583: nd_bpe = 32;
1584:
1585: nd_setup_parameters();
1586: nd_free_private_storage();
1587: len = 0;
1588: for ( i = 0; i < nd_psn; i++, f = NEXT(f) ) {
1589: a = dptondv(mod,(DP)BDY(f)); ndv_removecont(mod,a); nd_ps[i] = a;
1590: a = dptondv(0,(DP)BDY(f)); ndv_removecont(0,a); nd_psq[i] = a;
1591: NEWRHist(r);
1592: len = MAX(len,LEN(a));
1593: SG(r) = HTD(a); TD(r) = HTD(a); ndl_copy(HDL(a),DL(r));
1594: nd_psh[i] = r;
1595: }
1596: nmv_len = 16*len;
1597: NEWNDV(ndv_red);
1598: BDY(ndv_red) = (NMV)MALLOC(nmv_len*nmv_adv);
1599: for ( s0 = 0, i = 0; i < nd_psn; i++ ) {
1600: NEXTNODE(s0,s); BDY(s) = (pointer)i;
1601: }
1602: if ( s0 ) NEXT(s) = 0;
1603: return s0;
1604: }
1605:
1.1 noro 1606: void nd_gr(LIST f,LIST v,int m,struct order_spec *ord,LIST *rp)
1607: {
1608: struct order_spec ord1;
1609: VL fv,vv,vc;
1610: NODE fd,fd0,r,r0,t,x,s,xx;
1611: DP a,b,c;
1612:
1613: get_vars((Obj)f,&fv); pltovl(v,&vv);
1614: nd_nvar = length(vv);
1615: if ( ord->id )
1616: error("nd_gr : unsupported order");
1617: switch ( ord->ord.simple ) {
1618: case 0:
1619: is_rlex = 1;
1620: break;
1621: case 1:
1622: is_rlex = 0;
1623: break;
1624: default:
1625: error("nd_gr : unsupported order");
1626: }
1627: initd(ord);
1628: for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
1629: ptod(CO,vv,(P)BDY(t),&b);
1.20 noro 1630: NEXTNODE(fd0,fd); BDY(fd) = (pointer)b;
1.1 noro 1631: }
1632: if ( fd0 ) NEXT(fd) = 0;
1.16 noro 1633: s = nd_setup(m,fd0);
1634: x = nd_gb(m,s);
1.1 noro 1635: #if 0
1636: x = nd_reduceall(x,m);
1637: #endif
1638: for ( r0 = 0; x; x = NEXT(x) ) {
1639: NEXTNODE(r0,r);
1.20 noro 1640: if ( m ) {
1641: a = ndvtodp(m,nd_ps[(int)BDY(x)]);
1.16 noro 1642: _dtop_mod(CO,vv,a,(P *)&BDY(r));
1.20 noro 1643: } else {
1644: a = ndvtodp(m,nd_psq[(int)BDY(x)]);
1.16 noro 1645: dtop(CO,vv,a,(P *)&BDY(r));
1.20 noro 1646: }
1647: }
1648: if ( r0 ) NEXT(r) = 0;
1649: MKLIST(*rp,r0);
1650: fprintf(asir_out,"found=%d,notfirst=%d,create=%d\n",
1651: nd_found,nd_notfirst,nd_create);
1652: }
1653:
1654: void nd_gr_trace(LIST f,LIST v,int m,struct order_spec *ord,LIST *rp)
1655: {
1656: struct order_spec ord1;
1657: VL fv,vv,vc;
1658: NODE fd,fd0,r,r0,t,x,s,xx;
1659: DP a,b,c;
1660:
1661: get_vars((Obj)f,&fv); pltovl(v,&vv);
1662: nd_nvar = length(vv);
1663: if ( ord->id )
1664: error("nd_gr : unsupported order");
1665: switch ( ord->ord.simple ) {
1666: case 0:
1667: is_rlex = 1;
1668: break;
1669: case 1:
1670: is_rlex = 0;
1671: break;
1672: default:
1673: error("nd_gr : unsupported order");
1674: }
1675: initd(ord);
1676: for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
1677: ptod(CO,vv,(P)BDY(t),&c);
1678: if ( c ) {
1679: NEXTNODE(fd0,fd); BDY(fd) = (pointer)c;
1680: }
1681: }
1682: if ( fd0 ) NEXT(fd) = 0;
1683: /* setup over GF(m) */
1684: s = nd_setup_trace(m,fd0);
1685: x = nd_gb_trace(m,s);
1686: #if 0
1687: x = nd_reduceall(x,m);
1688: #endif
1689: for ( r0 = 0; x; x = NEXT(x) ) {
1690: NEXTNODE(r0,r);
1691: a = ndvtodp(0,nd_psq[(int)BDY(x)]);
1692: dtop(CO,vv,a,(P *)&BDY(r));
1.1 noro 1693: }
1694: if ( r0 ) NEXT(r) = 0;
1695: MKLIST(*rp,r0);
1696: fprintf(asir_out,"found=%d,notfirst=%d,create=%d\n",
1697: nd_found,nd_notfirst,nd_create);
1698: }
1699:
1700: void dltondl(int n,DL dl,unsigned int *r)
1701: {
1702: unsigned int *d;
1703: int i;
1704:
1705: d = dl->d;
1706: bzero(r,nd_wpd*sizeof(unsigned int));
1707: if ( is_rlex )
1708: for ( i = 0; i < n; i++ )
1709: r[(n-1-i)/nd_epw] |= (d[i]<<((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe));
1710: else
1711: for ( i = 0; i < n; i++ )
1712: r[i/nd_epw] |= d[i]<<((nd_epw-(i%nd_epw)-1)*nd_bpe);
1713: }
1714:
1715: DL ndltodl(int n,int td,unsigned int *ndl)
1716: {
1717: DL dl;
1718: int *d;
1719: int i;
1720:
1721: NEWDL(dl,n);
1.14 noro 1722: TD(dl) = td;
1.1 noro 1723: d = dl->d;
1724: if ( is_rlex )
1725: for ( i = 0; i < n; i++ )
1726: d[i] = (ndl[(n-1-i)/nd_epw]>>((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe))
1727: &((1<<nd_bpe)-1);
1728: else
1729: for ( i = 0; i < n; i++ )
1730: d[i] = (ndl[i/nd_epw]>>((nd_epw-(i%nd_epw)-1)*nd_bpe))
1731: &((1<<nd_bpe)-1);
1732: return dl;
1733: }
1734:
1.17 noro 1735: ND dptond(int mod,DP p)
1.1 noro 1736: {
1737: ND d;
1738: NM m0,m;
1739: MP t;
1740: int n;
1741:
1742: if ( !p )
1743: return 0;
1744: n = NV(p);
1745: m0 = 0;
1746: for ( t = BDY(p); t; t = NEXT(t) ) {
1747: NEXTNM(m0,m);
1.17 noro 1748: if ( mod )
1749: CM(m) = ITOS(C(t));
1750: else
1751: CQ(m) = (Q)C(t);
1.14 noro 1752: TD(m) = TD(DL(t));
1753: dltondl(n,DL(t),DL(m));
1.1 noro 1754: }
1755: NEXT(m) = 0;
1756: MKND(n,m0,d);
1.14 noro 1757: NV(d) = n;
1758: SG(d) = SG(p);
1.1 noro 1759: return d;
1760: }
1761:
1.17 noro 1762: DP ndtodp(int mod,ND p)
1.1 noro 1763: {
1764: DP d;
1765: MP m0,m;
1766: NM t;
1767: int n;
1768:
1769: if ( !p )
1770: return 0;
1771: n = NV(p);
1772: m0 = 0;
1773: for ( t = BDY(p); t; t = NEXT(t) ) {
1774: NEXTMP(m0,m);
1.17 noro 1775: if ( mod )
1776: C(m) = STOI(CM(t));
1777: else
1778: C(m) = (P)CQ(t);
1.14 noro 1779: DL(m) = ndltodl(n,TD(t),DL(t));
1.1 noro 1780: }
1781: NEXT(m) = 0;
1782: MKDP(n,m0,d);
1.14 noro 1783: SG(d) = SG(p);
1.1 noro 1784: return d;
1785: }
1786:
1787: void ndl_print(unsigned int *dl)
1788: {
1789: int n;
1790: int i;
1791:
1792: n = nd_nvar;
1793: printf("<<");
1794: if ( is_rlex )
1795: for ( i = 0; i < n; i++ )
1796: printf(i==n-1?"%d":"%d,",
1797: (dl[(n-1-i)/nd_epw]>>((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe))
1798: &((1<<nd_bpe)-1));
1799: else
1800: for ( i = 0; i < n; i++ )
1801: printf(i==n-1?"%d":"%d,",
1802: (dl[i/nd_epw]>>((nd_epw-(i%nd_epw)-1)*nd_bpe))
1803: &((1<<nd_bpe)-1));
1804: printf(">>");
1805: }
1806:
1807: void nd_print(ND p)
1808: {
1809: NM m;
1810:
1811: if ( !p )
1812: printf("0\n");
1813: else {
1814: for ( m = BDY(p); m; m = NEXT(m) ) {
1.14 noro 1815: printf("+%d*",CM(m));
1816: ndl_print(DL(m));
1.1 noro 1817: }
1818: printf("\n");
1819: }
1820: }
1821:
1.16 noro 1822: void nd_print_q(ND p)
1823: {
1824: NM m;
1825:
1826: if ( !p )
1827: printf("0\n");
1828: else {
1829: for ( m = BDY(p); m; m = NEXT(m) ) {
1830: printf("+");
1831: printexpr(CO,CQ(m));
1832: printf("*");
1833: ndl_print(DL(m));
1834: }
1835: printf("\n");
1836: }
1837: }
1838:
1.1 noro 1839: void ndp_print(ND_pairs d)
1840: {
1841: ND_pairs t;
1842:
1843: for ( t = d; t; t = NEXT(t) ) {
1844: printf("%d,%d ",t->i1,t->i2);
1845: }
1846: printf("\n");
1847: }
1848:
1.20 noro 1849: void nd_removecont(int mod,ND p)
1.16 noro 1850: {
1851: int i,n;
1852: Q *w;
1853: Q dvr,t;
1854: NM m;
1.21 noro 1855: struct oVECT v;
1856: N q,r;
1.16 noro 1857:
1.20 noro 1858: if ( mod )
1859: nd_mul_c(mod,p,invm(HCM(p),mod));
1860: else {
1861: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
1862: w = (Q *)ALLOCA(n*sizeof(Q));
1.21 noro 1863: v.len = n;
1864: v.body = (pointer *)w;
1.20 noro 1865: for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ )
1866: w[i] = CQ(m);
1.21 noro 1867: removecont_array(w,n);
1868: for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) CQ(m) = w[i];
1.16 noro 1869: }
1870: }
1871:
1.21 noro 1872: void nd_removecont2(ND p1,ND p2)
1873: {
1874: int i,n1,n2,n;
1875: Q *w;
1876: Q dvr,t;
1877: NM m;
1878: struct oVECT v;
1879: N q,r;
1880:
1881: if ( !p1 ) {
1882: nd_removecont(0,p2); return;
1883: } else if ( !p2 ) {
1884: nd_removecont(0,p1); return;
1885: }
1886: n1 = nd_length(p1);
1887: n2 = nd_length(p2);
1888: n = n1+n2;
1889: w = (Q *)ALLOCA(n*sizeof(Q));
1890: v.len = n;
1891: v.body = (pointer *)w;
1892: for ( m = BDY(p1), i = 0; i < n1; m = NEXT(m), i++ )
1893: w[i] = CQ(m);
1894: for ( m = BDY(p2); i < n; m = NEXT(m), i++ )
1895: w[i] = CQ(m);
1896: removecont_array(w,n);
1897: for ( m = BDY(p1), i = 0; i < n1; m = NEXT(m), i++ ) CQ(m) = w[i];
1898: for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) CQ(m) = w[i];
1899: }
1900:
1.20 noro 1901: void ndv_removecont(int mod,NDV p)
1.16 noro 1902: {
1903: int i,len;
1904: Q *w;
1905: Q dvr,t;
1906: NMV m;
1907:
1.20 noro 1908: if ( mod )
1909: ndv_mul_c(mod,p,invm(HCM(p),mod));
1910: else {
1911: len = p->len;
1912: w = (Q *)ALLOCA(len*sizeof(Q));
1913: for ( m = BDY(p), i = 0; i < len; NMV_ADV(m), i++ )
1914: w[i] = CQ(m);
1915: sortbynm(w,len);
1916: qltozl(w,len,&dvr);
1917: for ( m = BDY(p), i = 0; i < len; NMV_ADV(m), i++ ) {
1918: divq(CQ(m),dvr,&t); CQ(m) = t;
1919: }
1.16 noro 1920: }
1.21 noro 1921: }
1922:
1923: void removecont_array(Q *c,int n)
1924: {
1925: struct oVECT v;
1926: Q d0,d1,a,u,u1,gcd;
1927: int i;
1928: N qn,rn,gn;
1929: Q *q,*r;
1930:
1931: q = (Q *)ALLOCA(n*sizeof(Q));
1932: r = (Q *)ALLOCA(n*sizeof(Q));
1933: v.id = O_VECT; v.len = n; v.body = (pointer *)c;
1934: igcdv_estimate(&v,&d0);
1935: for ( i = 0; i < n; i++ ) {
1936: divn(NM(c[i]),NM(d0),&qn,&rn);
1937: NTOQ(qn,SGN(c[i])*SGN(d0),q[i]);
1938: NTOQ(rn,SGN(c[i]),r[i]);
1939: }
1940: for ( i = 0; i < n; i++ )
1941: if ( r[i] )
1942: break;
1943: if ( i < n ) {
1944: v.id = O_VECT; v.len = n; v.body = (pointer *)r;
1945: igcdv(&v,&d1);
1946: gcdn(NM(d0),NM(d1),&gn); NTOQ(gn,1,gcd);
1947: divsn(NM(d0),gn,&qn); NTOQ(qn,1,a);
1948: for ( i = 0; i < n; i++ ) {
1949: mulq(a,q[i],&u);
1950: if ( r[i] ) {
1951: divsn(NM(r[i]),gn,&qn); NTOQ(qn,SGN(r[i]),u1);
1952: addq(u,u1,&q[i]);
1953: } else
1954: q[i] = u;
1955: }
1956: }
1957: for ( i = 0; i < n; i++ )
1958: c[i] = q[i];
1.16 noro 1959: }
1960:
1.19 noro 1961: void nd_mul_c(int mod,ND p,int mul)
1.1 noro 1962: {
1963: NM m;
1964: int c,c1;
1965:
1966: if ( !p )
1967: return;
1968: for ( m = BDY(p); m; m = NEXT(m) ) {
1.14 noro 1969: c1 = CM(m);
1.19 noro 1970: DMAR(c1,mul,0,mod,c);
1.14 noro 1971: CM(m) = c;
1.1 noro 1972: }
1973: }
1974:
1.16 noro 1975: void nd_mul_c_q(ND p,Q mul)
1976: {
1977: NM m;
1978: Q c;
1979:
1980: if ( !p )
1981: return;
1982: for ( m = BDY(p); m; m = NEXT(m) ) {
1983: mulq(CQ(m),mul,&c); CQ(m) = c;
1984: }
1985: }
1986:
1.1 noro 1987: void nd_free(ND p)
1988: {
1989: NM t,s;
1990:
1991: if ( !p )
1992: return;
1993: t = BDY(p);
1994: while ( t ) {
1995: s = NEXT(t);
1996: FREENM(t);
1997: t = s;
1998: }
1999: FREEND(p);
2000: }
2001:
2002: void nd_append_red(unsigned int *d,int td,int i)
2003: {
1.13 noro 2004: RHist m,m0;
1.1 noro 2005: int h;
2006:
1.13 noro 2007: NEWRHist(m);
1.1 noro 2008: h = ndl_hash_value(td,d);
1.13 noro 2009: m->index = i;
1.14 noro 2010: TD(m) = td;
2011: ndl_copy(d,DL(m));
1.1 noro 2012: NEXT(m) = nd_red[h];
2013: nd_red[h] = m;
2014: }
2015:
1.5 noro 2016: unsigned int *dp_compute_bound(DP p)
2017: {
2018: unsigned int *d,*d1,*d2,*t;
2019: MP m;
1.7 noro 2020: int i,l;
1.5 noro 2021:
2022: if ( !p )
2023: return 0;
2024: d1 = (unsigned int *)ALLOCA(nd_nvar*sizeof(unsigned int));
2025: d2 = (unsigned int *)ALLOCA(nd_nvar*sizeof(unsigned int));
2026: m = BDY(p);
1.14 noro 2027: bcopy(DL(m)->d,d1,nd_nvar*sizeof(unsigned int));
1.5 noro 2028: for ( m = NEXT(BDY(p)); m; m = NEXT(m) ) {
1.14 noro 2029: d = DL(m)->d;
1.5 noro 2030: for ( i = 0; i < nd_nvar; i++ )
2031: d2[i] = d[i] > d1[i] ? d[i] : d1[i];
2032: t = d1; d1 = d2; d2 = t;
2033: }
1.13 noro 2034: l = (nd_nvar+31);
1.7 noro 2035: t = (unsigned int *)MALLOC_ATOMIC(l*sizeof(unsigned int));
2036: bzero(t,l*sizeof(unsigned int));
1.5 noro 2037: bcopy(d1,t,nd_nvar*sizeof(unsigned int));
2038: return t;
2039: }
2040:
1.1 noro 2041: unsigned int *nd_compute_bound(ND p)
2042: {
2043: unsigned int *d1,*d2,*t;
1.9 noro 2044: int i,l;
1.1 noro 2045: NM m;
2046:
2047: if ( !p )
2048: return 0;
2049: d1 = (unsigned int *)ALLOCA(nd_wpd*sizeof(unsigned int));
2050: d2 = (unsigned int *)ALLOCA(nd_wpd*sizeof(unsigned int));
2051: bcopy(HDL(p),d1,nd_wpd*sizeof(unsigned int));
2052: for ( m = NEXT(BDY(p)); m; m = NEXT(m) ) {
1.14 noro 2053: ndl_lcm(DL(m),d1,d2);
1.1 noro 2054: t = d1; d1 = d2; d2 = t;
2055: }
1.12 noro 2056: l = nd_nvar+31;
1.9 noro 2057: t = (unsigned int *)MALLOC_ATOMIC(l*sizeof(unsigned int));
2058: bzero(t,l*sizeof(unsigned int));
1.5 noro 2059: for ( i = 0; i < nd_nvar; i++ )
2060: t[i] = (d1[i/nd_epw]>>((nd_epw-(i%nd_epw)-1)*nd_bpe))&nd_mask0;
1.1 noro 2061: return t;
2062: }
2063:
2064: void nd_setup_parameters() {
2065: int i;
2066:
2067: nd_epw = (sizeof(unsigned int)*8)/nd_bpe;
2068: nd_wpd = nd_nvar/nd_epw+(nd_nvar%nd_epw?1:0);
2069: if ( nd_bpe < 32 ) {
2070: nd_mask0 = (1<<nd_bpe)-1;
2071: } else {
2072: nd_mask0 = 0xffffffff;
2073: }
2074: bzero(nd_mask,sizeof(nd_mask));
2075: nd_mask1 = 0;
2076: for ( i = 0; i < nd_epw; i++ ) {
2077: nd_mask[nd_epw-i-1] = (nd_mask0<<(i*nd_bpe));
2078: nd_mask1 |= (1<<(nd_bpe-1))<<(i*nd_bpe);
2079: }
1.13 noro 2080: nm_adv = sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int);
1.3 noro 2081: nmv_adv = sizeof(struct oNMV)+(nd_wpd-1)*sizeof(unsigned int);
1.1 noro 2082: }
2083:
1.20 noro 2084: /* mod < 0 => realloc nd_ps and pd_psq */
2085:
2086: ND_pairs nd_reconstruct(int mod,int trace,ND_pairs d)
1.1 noro 2087: {
1.11 noro 2088: int i,obpe,oadv;
1.13 noro 2089: NM prev_nm_free_list;
2090: RHist mr0,mr;
2091: RHist r;
1.1 noro 2092: ND_pairs s0,s,t,prev_ndp_free_list;
1.15 noro 2093:
1.1 noro 2094: obpe = nd_bpe;
1.11 noro 2095: oadv = nmv_adv;
1.5 noro 2096: if ( obpe < 4 )
2097: nd_bpe = 4;
2098: else if ( obpe < 6 )
2099: nd_bpe = 6;
2100: else if ( obpe < 8 )
2101: nd_bpe = 8;
2102: else if ( obpe < 16 )
2103: nd_bpe = 16;
2104: else if ( obpe < 32 )
2105: nd_bpe = 32;
2106: else
2107: error("nd_reconstruct : exponent too large");
2108:
1.1 noro 2109: nd_setup_parameters();
2110: prev_nm_free_list = _nm_free_list;
2111: prev_ndp_free_list = _ndp_free_list;
2112: _nm_free_list = 0;
2113: _ndp_free_list = 0;
1.20 noro 2114: if ( mod != 0 )
2115: for ( i = nd_psn-1; i >= 0; i-- )
2116: ndv_realloc(nd_ps[i],obpe,oadv);
2117: if ( !mod || trace )
2118: for ( i = nd_psn-1; i >= 0; i-- )
2119: ndv_realloc(nd_psq[i],obpe,oadv);
1.1 noro 2120: s0 = 0;
2121: for ( t = d; t; t = NEXT(t) ) {
2122: NEXTND_pairs(s0,s);
2123: s->i1 = t->i1;
2124: s->i2 = t->i2;
1.14 noro 2125: TD(s) = TD(t);
2126: SG(s) = SG(t);
1.1 noro 2127: ndl_dup(obpe,t->lcm,s->lcm);
2128: }
1.6 noro 2129: for ( i = 0; i < REDTAB_LEN; i++ ) {
1.13 noro 2130: for ( mr0 = 0, r = nd_red[i]; r; r = NEXT(r) ) {
1.16 noro 2131: NEXTRHist(mr0,mr);
1.13 noro 2132: mr->index = r->index;
1.20 noro 2133: SG(mr) = SG(r);
1.14 noro 2134: TD(mr) = TD(r);
2135: ndl_dup(obpe,DL(r),DL(mr));
1.6 noro 2136: }
2137: if ( mr0 ) NEXT(mr) = 0;
2138: nd_red[i] = mr0;
2139: }
1.11 noro 2140: for ( i = 0; i < nd_psn; i++ ) {
1.20 noro 2141: NEWRHist(r); SG(r) = SG(nd_psh[i]);
2142: TD(r) = TD(nd_psh[i]); ndl_dup(obpe,DL(nd_psh[i]),DL(r));
1.13 noro 2143: nd_psh[i] = r;
1.11 noro 2144: }
1.1 noro 2145: if ( s0 ) NEXT(s) = 0;
2146: prev_nm_free_list = 0;
2147: prev_ndp_free_list = 0;
1.3 noro 2148: BDY(ndv_red) = (NMV)REALLOC(BDY(ndv_red),nmv_len*nmv_adv);
1.1 noro 2149: GC_gcollect();
2150: return s0;
2151: }
2152:
2153: void ndl_dup(int obpe,unsigned int *d,unsigned int *r)
2154: {
2155: int n,i,ei,oepw,cepw,cbpe;
2156:
2157: n = nd_nvar;
2158: oepw = (sizeof(unsigned int)*8)/obpe;
2159: cepw = nd_epw;
2160: cbpe = nd_bpe;
1.15 noro 2161: for ( i = 0; i < nd_wpd; i++ )
2162: r[i] = 0;
1.1 noro 2163: if ( is_rlex )
2164: for ( i = 0; i < n; i++ ) {
2165: ei = (d[(n-1-i)/oepw]>>((oepw-((n-1-i)%oepw)-1)*obpe))
2166: &((1<<obpe)-1);
2167: r[(n-1-i)/cepw] |= (ei<<((cepw-((n-1-i)%cepw)-1)*cbpe));
2168: }
2169: else
2170: for ( i = 0; i < n; i++ ) {
2171: ei = (d[i/oepw]>>((oepw-(i%oepw)-1)*obpe))
2172: &((1<<obpe)-1);
2173: r[i/cepw] |= (ei<<((cepw-(i%cepw)-1)*cbpe));
2174: }
2175: }
2176:
1.11 noro 2177: void nd_realloc(ND p,int obpe)
1.1 noro 2178: {
2179: NM m,mr,mr0;
2180:
1.11 noro 2181: if ( p ) {
2182: m = BDY(p);
1.1 noro 2183: for ( mr0 = 0; m; m = NEXT(m) ) {
2184: NEXTNM(mr0,mr);
1.14 noro 2185: CM(mr) = CM(m);
2186: TD(mr) = TD(m);
2187: ndl_dup(obpe,DL(m),DL(mr));
1.1 noro 2188: }
2189: NEXT(mr) = 0;
1.11 noro 2190: BDY(p) = mr0;
1.1 noro 2191: }
2192: }
1.3 noro 2193:
1.6 noro 2194: ND nd_copy(ND p)
2195: {
2196: NM m,mr,mr0;
2197: int c,n,s;
2198: ND r;
2199:
2200: if ( !p )
2201: return 0;
2202: else {
2203: s = sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int);
2204: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
2205: NEXTNM(mr0,mr);
1.14 noro 2206: CM(mr) = CM(m);
2207: TD(mr) = TD(m);
2208: ndl_copy(DL(m),DL(mr));
1.6 noro 2209: }
2210: NEXT(mr) = 0;
2211: MKND(NV(p),mr0,r);
1.14 noro 2212: SG(r) = SG(p);
1.6 noro 2213: return r;
2214: }
2215: }
2216:
1.16 noro 2217: int nd_sp(int mod,ND_pairs p,ND *rp)
1.11 noro 2218: {
2219: NM m;
2220: NDV p1,p2;
2221: ND t1,t2;
2222: unsigned int *lcm;
2223: int td;
2224:
1.20 noro 2225: if ( mod ) {
2226: p1 = nd_ps[p->i1]; p2 = nd_ps[p->i2];
2227: } else {
2228: p1 = nd_psq[p->i1]; p2 = nd_psq[p->i2];
2229: }
1.11 noro 2230: lcm = p->lcm;
1.14 noro 2231: td = TD(p);
1.11 noro 2232: NEWNM(m);
1.20 noro 2233: CQ(m) = HCQ(p2);
1.16 noro 2234: TD(m) = td-HTD(p1); ndl_sub(lcm,HDL(p1),DL(m));
1.14 noro 2235: if ( ndl_check_bound2(p->i1,DL(m)) )
1.11 noro 2236: return 0;
1.16 noro 2237: t1 = ndv_mul_nm_create(mod,p1,m);
2238: if ( mod )
2239: CM(m) = mod-HCM(p1);
2240: else
2241: chsgnq(HCQ(p1),&CQ(m));
2242: TD(m) = td-HTD(p2); ndl_sub(lcm,HDL(p2),DL(m));
1.14 noro 2243: if ( ndl_check_bound2(p->i2,DL(m)) ) {
1.11 noro 2244: nd_free(t1);
2245: return 0;
2246: }
1.16 noro 2247: ndv_mul_nm(mod,p2,m,ndv_red);
1.11 noro 2248: FREENM(m);
1.20 noro 2249: *rp = ndv_add(mod,t1,ndv_red);
1.11 noro 2250: return 1;
2251: }
2252:
1.19 noro 2253: void ndv_mul_c(int mod,NDV p,int mul)
1.11 noro 2254: {
2255: NMV m;
2256: int c,c1,len,i;
2257:
2258: if ( !p )
2259: return;
1.14 noro 2260: len = LEN(p);
1.11 noro 2261: for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {
1.14 noro 2262: c1 = CM(m);
1.19 noro 2263: DMAR(c1,mul,0,mod,c);
1.14 noro 2264: CM(m) = c;
1.11 noro 2265: }
2266: }
2267:
1.16 noro 2268: void ndv_mul_c_q(NDV p,Q mul)
2269: {
2270: NMV m;
2271: Q c;
2272: int len,i;
2273:
2274: if ( !p )
2275: return;
2276: len = LEN(p);
2277: for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {
2278: mulq(CQ(m),mul,&c); CQ(m) = c;
2279: }
2280: }
2281:
2282: void ndv_mul_nm(int mod,NDV p,NM m0,NDV r)
1.4 noro 2283: {
2284: NMV m,mr,mr0;
2285: unsigned int *d,*dt,*dm;
2286: int c,n,td,i,c1,c2,len;
1.16 noro 2287: Q q;
1.4 noro 2288:
2289: if ( !p )
2290: /* XXX */
1.14 noro 2291: LEN(r) = 0;
1.4 noro 2292: else {
1.14 noro 2293: n = NV(p); m = BDY(p); len = LEN(p);
1.16 noro 2294: d = DL(m0); td = TD(m0);
1.4 noro 2295: mr = BDY(r);
1.16 noro 2296: if ( mod ) {
2297: c = CM(m0);
2298: for ( ; len > 0; len--, NMV_ADV(m), NMV_ADV(mr) ) {
1.19 noro 2299: c1 = CM(m); DMAR(c1,c,0,mod,c2); CM(mr) = c2;
1.16 noro 2300: TD(mr) = TD(m)+td; ndl_add(DL(m),d,DL(mr));
2301: }
2302: } else {
2303: q = CQ(m0);
2304: for ( ; len > 0; len--, NMV_ADV(m), NMV_ADV(mr) ) {
2305: mulq(CQ(m),q,&CQ(mr));
2306: TD(mr) = TD(m)+td; ndl_add(DL(m),d,DL(mr));
2307: }
1.9 noro 2308: }
2309: NV(r) = NV(p);
1.14 noro 2310: LEN(r) = LEN(p);
2311: SG(r) = SG(p) + td;
1.9 noro 2312: }
2313: }
2314:
1.16 noro 2315: ND ndv_mul_nm_create(int mod,NDV p,NM m0)
1.9 noro 2316: {
2317: NM mr,mr0;
2318: NMV m;
2319: unsigned int *d,*dt,*dm;
2320: int c,n,td,i,c1,c2,len;
1.16 noro 2321: Q q;
1.9 noro 2322: ND r;
2323:
2324: if ( !p )
2325: return 0;
2326: else {
2327: n = NV(p); m = BDY(p);
1.16 noro 2328: d = DL(m0); td = TD(m0);
1.14 noro 2329: len = LEN(p);
1.9 noro 2330: mr0 = 0;
1.16 noro 2331: if ( mod ) {
2332: c = CM(m0);
2333: for ( i = 0; i < len; i++, NMV_ADV(m) ) {
2334: NEXTNM(mr0,mr);
2335: c1 = CM(m);
1.19 noro 2336: DMAR(c1,c,0,mod,c2);
1.16 noro 2337: CM(mr) = c2;
2338: TD(mr) = TD(m)+td;
2339: ndl_add(DL(m),d,DL(mr));
2340: }
2341: } else {
2342: q = CQ(m0);
2343: for ( i = 0; i < len; i++, NMV_ADV(m) ) {
2344: NEXTNM(mr0,mr);
2345: mulq(CQ(m),q,&CQ(mr));
2346: TD(mr) = TD(m)+td;
2347: ndl_add(DL(m),d,DL(mr));
2348: }
1.4 noro 2349: }
1.9 noro 2350: NEXT(mr) = 0;
2351: MKND(NV(p),mr0,r);
1.14 noro 2352: SG(r) = SG(p) + td;
1.9 noro 2353: return r;
1.4 noro 2354: }
2355: }
2356:
1.19 noro 2357: ND ndv_add(int mod,ND p1,NDV p2)
1.4 noro 2358: {
1.9 noro 2359: register NM prev,cur,new;
1.4 noro 2360: int c,c1,c2,t,td,td2,mul,len,i;
1.9 noro 2361: NM head;
1.4 noro 2362: unsigned int *d;
2363: NMV m2;
1.16 noro 2364: Q q;
1.4 noro 2365:
2366: if ( !p1 )
2367: return 0;
1.20 noro 2368: else if ( !mod )
2369: return ndv_add_q(p1,p2);
1.4 noro 2370: else {
2371: prev = 0; head = cur = BDY(p1);
1.14 noro 2372: NEWNM(new); len = LEN(p2);
1.9 noro 2373: for ( m2 = BDY(p2), i = 0; cur && i < len; ) {
1.14 noro 2374: td2 = TD(new) = TD(m2);
2375: if ( TD(cur) > td2 ) {
1.13 noro 2376: prev = cur; cur = NEXT(cur);
2377: continue;
1.14 noro 2378: } else if ( TD(cur) < td2 ) c = -1;
1.13 noro 2379: else if ( nd_wpd == 1 ) {
1.14 noro 2380: if ( DL(cur)[0] > DL(m2)[0] ) c = is_rlex ? -1 : 1;
2381: else if ( DL(cur)[0] < DL(m2)[0] ) c = is_rlex ? 1 : -1;
1.13 noro 2382: else c = 0;
2383: }
1.14 noro 2384: else c = ndl_compare(DL(cur),DL(m2));
1.4 noro 2385: switch ( c ) {
2386: case 0:
1.19 noro 2387: t = CM(m2)+CM(cur)-mod;
2388: if ( t < 0 ) t += mod;
1.14 noro 2389: if ( t ) CM(cur) = t;
1.4 noro 2390: else if ( !prev ) {
1.9 noro 2391: head = NEXT(cur); FREENM(cur); cur = head;
1.4 noro 2392: } else {
1.9 noro 2393: NEXT(prev) = NEXT(cur); FREENM(cur); cur = NEXT(prev);
1.4 noro 2394: }
2395: NMV_ADV(m2); i++;
2396: break;
2397: case 1:
1.9 noro 2398: prev = cur; cur = NEXT(cur);
1.4 noro 2399: break;
2400: case -1:
1.14 noro 2401: ndl_copy(DL(m2),DL(new));
1.16 noro 2402: CQ(new) = CQ(m2);
2403: if ( !prev ) {
2404: /* cur = head */
2405: prev = new; NEXT(prev) = head; head = prev;
2406: } else {
2407: NEXT(prev) = new; NEXT(new) = cur; prev = new;
2408: }
2409: NEWNM(new); NMV_ADV(m2); i++;
2410: break;
2411: }
2412: }
2413: for ( ; i < len; i++, NMV_ADV(m2) ) {
2414: td2 = TD(new) = TD(m2); CQ(new) = CQ(m2); ndl_copy(DL(m2),DL(new));
2415: if ( !prev ) {
2416: prev = new; NEXT(prev) = 0; head = prev;
2417: } else {
2418: NEXT(prev) = new; NEXT(new) = 0; prev = new;
2419: }
2420: NEWNM(new);
2421: }
2422: FREENM(new);
2423: if ( head ) {
2424: BDY(p1) = head; SG(p1) = MAX(SG(p1),SG(p2));
2425: return p1;
2426: } else {
2427: FREEND(p1);
2428: return 0;
2429: }
2430:
2431: }
2432: }
2433:
2434: ND ndv_add_q(ND p1,NDV p2)
2435: {
2436: register NM prev,cur,new;
2437: int c,c1,c2,t,td,td2,mul,len,i;
2438: NM head;
2439: unsigned int *d;
2440: NMV m2;
2441: Q q;
2442:
2443: if ( !p1 )
2444: return 0;
2445: else {
2446: prev = 0; head = cur = BDY(p1);
2447: NEWNM(new); len = LEN(p2);
2448: for ( m2 = BDY(p2), i = 0; cur && i < len; ) {
2449: td2 = TD(new) = TD(m2);
2450: if ( TD(cur) > td2 ) {
2451: prev = cur; cur = NEXT(cur);
2452: continue;
2453: } else if ( TD(cur) < td2 ) c = -1;
2454: else if ( nd_wpd == 1 ) {
2455: if ( DL(cur)[0] > DL(m2)[0] ) c = is_rlex ? -1 : 1;
2456: else if ( DL(cur)[0] < DL(m2)[0] ) c = is_rlex ? 1 : -1;
2457: else c = 0;
2458: }
2459: else c = ndl_compare(DL(cur),DL(m2));
2460: switch ( c ) {
2461: case 0:
2462: addq(CQ(cur),CQ(m2),&q);
2463: if ( q )
2464: CQ(cur) = q;
2465: else if ( !prev ) {
2466: head = NEXT(cur); FREENM(cur); cur = head;
2467: } else {
2468: NEXT(prev) = NEXT(cur); FREENM(cur); cur = NEXT(prev);
2469: }
2470: NMV_ADV(m2); i++;
2471: break;
2472: case 1:
2473: prev = cur; cur = NEXT(cur);
2474: break;
2475: case -1:
2476: ndl_copy(DL(m2),DL(new));
2477: CQ(new) = CQ(m2);
1.4 noro 2478: if ( !prev ) {
2479: /* cur = head */
1.9 noro 2480: prev = new; NEXT(prev) = head; head = prev;
1.4 noro 2481: } else {
1.9 noro 2482: NEXT(prev) = new; NEXT(new) = cur; prev = new;
1.4 noro 2483: }
1.9 noro 2484: NEWNM(new); NMV_ADV(m2); i++;
1.4 noro 2485: break;
2486: }
2487: }
1.9 noro 2488: for ( ; i < len; i++, NMV_ADV(m2) ) {
1.16 noro 2489: td2 = TD(new) = TD(m2); CQ(new) = CQ(m2); ndl_copy(DL(m2),DL(new));
1.9 noro 2490: if ( !prev ) {
2491: prev = new; NEXT(prev) = 0; head = prev;
2492: } else {
2493: NEXT(prev) = new; NEXT(new) = 0; prev = new;
2494: }
2495: NEWNM(new);
2496: }
1.4 noro 2497: FREENM(new);
2498: if ( head ) {
1.14 noro 2499: BDY(p1) = head; SG(p1) = MAX(SG(p1),SG(p2));
1.4 noro 2500: return p1;
2501: } else {
2502: FREEND(p1);
2503: return 0;
2504: }
2505:
2506: }
2507: }
2508:
1.11 noro 2509: void ndv_realloc(NDV p,int obpe,int oadv)
2510: {
1.13 noro 2511: NMV m,mr,mr0,t;
2512: int len,i,k;
1.11 noro 2513:
1.13 noro 2514: #define NMV_OPREV(m) (m = (NMV)(((char *)m)-oadv))
2515: #define NMV_PREV(m) (m = (NMV)(((char *)m)-nmv_adv))
1.11 noro 2516:
2517: if ( p ) {
1.14 noro 2518: m = BDY(p); len = LEN(p);
1.15 noro 2519: if ( nmv_adv > oadv )
2520: mr0 = (NMV)REALLOC(BDY(p),len*nmv_adv);
2521: else
2522: mr0 = BDY(p);
1.13 noro 2523: m = (NMV)((char *)mr0+(len-1)*oadv);
2524: mr = (NMV)((char *)mr0+(len-1)*nmv_adv);
2525: t = (NMV)ALLOCA(nmv_adv);
2526: for ( i = 0; i < len; i++, NMV_OPREV(m), NMV_PREV(mr) ) {
1.16 noro 2527: CQ(t) = CQ(m);
1.14 noro 2528: TD(t) = TD(m);
2529: for ( k = 0; k < nd_wpd; k++ ) DL(t)[k] = 0;
2530: ndl_dup(obpe,DL(m),DL(t));
1.16 noro 2531: CQ(mr) = CQ(t);
1.14 noro 2532: TD(mr) = TD(t);
2533: ndl_copy(DL(t),DL(mr));
1.11 noro 2534: }
2535: BDY(p) = mr0;
2536: }
2537: }
2538:
1.16 noro 2539: NDV ndtondv(int mod,ND p)
1.3 noro 2540: {
2541: NDV d;
2542: NMV m,m0;
2543: NM t;
2544: int i,len;
2545:
2546: if ( !p )
2547: return 0;
2548: len = nd_length(p);
1.16 noro 2549: if ( mod )
2550: m0 = m = (NMV)MALLOC_ATOMIC(len*nmv_adv);
2551: else
2552: m0 = m = (NMV)MALLOC(len*nmv_adv);
1.3 noro 2553: for ( t = BDY(p), i = 0; t; t = NEXT(t), i++, NMV_ADV(m) ) {
1.14 noro 2554: TD(m) = TD(t);
2555: ndl_copy(DL(t),DL(m));
1.16 noro 2556: CQ(m) = CQ(t);
1.3 noro 2557: }
2558: MKNDV(NV(p),m0,len,d);
1.14 noro 2559: SG(d) = SG(p);
1.3 noro 2560: return d;
2561: }
2562:
1.16 noro 2563: NDV dptondv(int mod,DP p)
1.11 noro 2564: {
2565: NDV d;
2566: NMV m0,m;
2567: MP t;
1.20 noro 2568: DP q;
1.11 noro 2569: int l,i,n;
2570:
2571: if ( !p )
2572: return 0;
2573: for ( t = BDY(p), l = 0; t; t = NEXT(t), l++ );
1.20 noro 2574: if ( mod ) {
2575: _dp_mod(p,mod,0,&q); p = q;
1.16 noro 2576: m0 = m = (NMV)MALLOC_ATOMIC(l*nmv_adv);
1.20 noro 2577: } else
1.16 noro 2578: m0 = m = (NMV)MALLOC(l*nmv_adv);
1.11 noro 2579: n = NV(p);
1.17 noro 2580: for ( t = BDY(p), i = 0; i < l; i++, t = NEXT(t), NMV_ADV(m) ) {
2581: if ( mod )
1.16 noro 2582: CM(m) = ITOS(C(t));
1.17 noro 2583: else
1.16 noro 2584: CQ(m) = (Q)C(t);
1.17 noro 2585: TD(m) = TD(DL(t));
2586: dltondl(n,DL(t),DL(m));
1.11 noro 2587: }
2588: MKNDV(n,m0,l,d);
1.14 noro 2589: SG(d) = SG(p);
1.11 noro 2590: return d;
2591: }
2592:
1.16 noro 2593: DP ndvtodp(int mod,NDV p)
1.11 noro 2594: {
2595: DP d;
2596: MP m0,m;
2597: NMV t;
2598: int len,i,n;
2599:
2600: if ( !p )
2601: return 0;
2602: m0 = 0;
1.14 noro 2603: len = LEN(p);
1.11 noro 2604: n = NV(p);
1.17 noro 2605: for ( t = BDY(p), i = 0; i < len; i++, NMV_ADV(t) ) {
2606: NEXTMP(m0,m);
2607: if ( mod )
1.16 noro 2608: C(m) = STOI(CM(t));
1.17 noro 2609: else
1.16 noro 2610: C(m) = (P)CQ(t);
1.17 noro 2611: DL(m) = ndltodl(n,TD(t),DL(t));
1.11 noro 2612: }
2613: NEXT(m) = 0;
2614: MKDP(NV(p),m0,d);
1.14 noro 2615: SG(d) = SG(p);
1.11 noro 2616: return d;
2617: }
2618:
1.3 noro 2619: void ndv_print(NDV p)
2620: {
2621: NMV m;
2622: int i,len;
2623:
2624: if ( !p )
2625: printf("0\n");
2626: else {
1.14 noro 2627: len = LEN(p);
1.3 noro 2628: for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {
1.14 noro 2629: printf("+%d*",CM(m));
1.16 noro 2630: ndl_print(DL(m));
2631: }
2632: printf("\n");
2633: }
2634: }
2635:
2636: void ndv_print_q(NDV p)
2637: {
2638: NMV m;
2639: int i,len;
2640:
2641: if ( !p )
2642: printf("0\n");
2643: else {
2644: len = LEN(p);
2645: for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {
2646: printf("+");
2647: printexpr(CO,CQ(m));
2648: printf("*");
1.14 noro 2649: ndl_print(DL(m));
1.3 noro 2650: }
2651: printf("\n");
2652: }
1.11 noro 2653: }