Annotation of OpenXM_contrib2/asir2000/engine/dist.c, Revision 1.2
1.2 ! noro 1: /* $OpenXM: OpenXM_contrib2/asir2000/engine/dist.c,v 1.1.1.1 1999/12/03 07:39:08 noro Exp $ */
1.1 noro 2: #include "ca.h"
3:
4: #define NV(p) ((p)->nv)
5: #define C(p) ((p)->c)
6:
7: #define ORD_REVGRADLEX 0
8: #define ORD_GRADLEX 1
9: #define ORD_LEX 2
10: #define ORD_BREVGRADLEX 3
11: #define ORD_BGRADLEX 4
12: #define ORD_BLEX 5
13: #define ORD_BREVREV 6
14: #define ORD_BGRADREV 7
15: #define ORD_BLEXREV 8
16: #define ORD_ELIM 9
17:
18: int (*cmpdl)()=cmpdl_revgradlex;
19: int (*primitive_cmpdl[3])() = {cmpdl_revgradlex,cmpdl_gradlex,cmpdl_lex};
20:
1.2 ! noro 21: void comm_muld(VL,DP,DP,DP *);
! 22: void weyl_muld(VL,DP,DP,DP *);
! 23: void weyl_muldm(VL,DP,MP,DP *);
! 24: void weyl_mulmm(VL,MP,MP,int,DP *);
! 25: void mkwc(int,int,Q *);
! 26:
! 27: int do_weyl;
! 28:
1.1 noro 29: int dp_nelim,dp_fcoeffs;
30: struct order_spec dp_current_spec;
31: int *dp_dl_work;
32:
33: int has_fcoef(DP);
34: int has_fcoef_p(P);
35:
36: int has_fcoef(f)
37: DP f;
38: {
39: MP t;
40:
41: if ( !f )
42: return 0;
43: for ( t = BDY(f); t; t = NEXT(t) )
44: if ( has_fcoef_p(t->c) )
45: break;
46: return t ? 1 : 0;
47: }
48:
49: int has_fcoef_p(f)
50: P f;
51: {
52: DCP dc;
53:
54: if ( !f )
55: return 0;
56: else if ( NUM(f) )
57: return (NID((Num)f) == N_LM || NID((Num)f) == N_GF2N) ? 1 : 0;
58: else {
59: for ( dc = DC(f); dc; dc = NEXT(dc) )
60: if ( has_fcoef_p(COEF(dc)) )
61: return 1;
62: return 0;
63: }
64: }
65:
66: void initd(spec)
67: struct order_spec *spec;
68: {
69: switch ( spec->id ) {
70: case 2:
71: cmpdl = cmpdl_matrix;
72: dp_dl_work = (int *)MALLOC_ATOMIC(spec->nv*sizeof(int));
73: break;
74: case 1:
75: cmpdl = cmpdl_order_pair;
76: break;
77: default:
78: switch ( spec->ord.simple ) {
79: case ORD_REVGRADLEX:
80: cmpdl = cmpdl_revgradlex; break;
81: case ORD_GRADLEX:
82: cmpdl = cmpdl_gradlex; break;
83: case ORD_BREVGRADLEX:
84: cmpdl = cmpdl_brevgradlex; break;
85: case ORD_BGRADLEX:
86: cmpdl = cmpdl_bgradlex; break;
87: case ORD_BLEX:
88: cmpdl = cmpdl_blex; break;
89: case ORD_BREVREV:
90: cmpdl = cmpdl_brevrev; break;
91: case ORD_BGRADREV:
92: cmpdl = cmpdl_bgradrev; break;
93: case ORD_BLEXREV:
94: cmpdl = cmpdl_blexrev; break;
95: case ORD_ELIM:
96: cmpdl = cmpdl_elim; break;
97: case ORD_LEX: default:
98: cmpdl = cmpdl_lex; break;
99: }
100: break;
101: }
102: dp_current_spec = *spec;
103: }
104:
105: void ptod(vl,dvl,p,pr)
106: VL vl,dvl;
107: P p;
108: DP *pr;
109: {
110: int isconst = 0;
111: int n,i;
112: VL tvl;
113: V v;
114: DL d;
115: MP m;
116: DCP dc;
117: DP r,s,t,u;
118: P x,c;
119:
120: if ( !p )
121: *pr = 0;
122: else {
123: for ( n = 0, tvl = dvl; tvl; tvl = NEXT(tvl), n++ );
124: if ( NUM(p) ) {
125: NEWDL(d,n);
126: NEWMP(m); m->dl = d; C(m) = p; NEXT(m) = 0; MKDP(n,m,*pr); (*pr)->sugar = 0;
127: } else {
128: for ( i = 0, tvl = dvl, v = VR(p);
129: tvl && tvl->v != v; tvl = NEXT(tvl), i++ );
130: if ( !tvl ) {
131: for ( dc = DC(p), s = 0, MKV(v,x); dc; dc = NEXT(dc) ) {
132: ptod(vl,dvl,COEF(dc),&t); pwrp(vl,x,DEG(dc),&c);
133: muldc(vl,t,c,&r); addd(vl,r,s,&t); s = t;
134: }
135: *pr = s;
136: } else {
137: for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) {
138: ptod(vl,dvl,COEF(dc),&t);
139: NEWDL(d,n); d->td = QTOS(DEG(dc)); d->d[i] = d->td;
140: NEWMP(m); m->dl = d; C(m) = (P)ONE; NEXT(m) = 0; MKDP(n,m,u); u->sugar = d->td;
1.2 ! noro 141: comm_muld(vl,t,u,&r); addd(vl,r,s,&t); s = t;
1.1 noro 142: }
143: *pr = s;
144: }
145: }
146: }
147: if ( !dp_fcoeffs && has_fcoef(*pr) )
148: dp_fcoeffs = 1;
149: }
150:
151: void dtop(vl,dvl,p,pr)
152: VL vl,dvl;
153: DP p;
154: P *pr;
155: {
156: int n,i;
157: DL d;
158: MP m;
159: P r,s,t,u,w;
160: Q q;
161: VL tvl;
162:
163: if ( !p )
164: *pr = 0;
165: else {
166: for ( n = p->nv, m = BDY(p), s = 0; m; m = NEXT(m) ) {
167: t = C(m);
168: if ( NUM(t) && NID((Num)t) == N_M ) {
169: mptop(t,&u); t = u;
170: }
171: for ( i = 0, d = m->dl, tvl = dvl;
172: i < n; tvl = NEXT(tvl), i++ ) {
173: MKV(tvl->v,r); STOQ(d->d[i],q); pwrp(vl,r,q,&u);
174: mulp(vl,t,u,&w); t = w;
175: }
176: addp(vl,s,t,&u); s = u;
177: }
178: *pr = s;
179: }
180: }
181:
182: void nodetod(node,dp)
183: NODE node;
184: DP *dp;
185: {
186: NODE t;
187: int len,i,td;
188: Q e;
189: DL d;
190: MP m;
191: DP u;
192:
193: for ( t = node, len = 0; t; t = NEXT(t), len++ );
194: NEWDL(d,len);
195: for ( t = node, i = 0, td = 0; i < len; t = NEXT(t), i++ ) {
196: e = (Q)BDY(t);
197: if ( !e )
198: d->d[i] = 0;
199: else if ( !NUM(e) || !RATN(e) || !INT(e) )
200: error("nodetod : invalid input");
201: else {
202: d->d[i] = QTOS((Q)e); td += d->d[i];
203: }
204: }
205: d->td = td;
206: NEWMP(m); m->dl = d; C(m) = (P)ONE; NEXT(m) = 0;
207: MKDP(len,m,u); u->sugar = td; *dp = u;
208: }
209:
210: int sugard(m)
211: MP m;
212: {
213: int s;
214:
215: for ( s = 0; m; m = NEXT(m) )
216: s = MAX(s,m->dl->td);
217: return s;
218: }
219:
220: void addd(vl,p1,p2,pr)
221: VL vl;
222: DP p1,p2,*pr;
223: {
224: int n;
225: MP m1,m2,mr,mr0;
226: P t;
227:
228: if ( !p1 )
229: *pr = p2;
230: else if ( !p2 )
231: *pr = p1;
232: else {
233: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; )
234: switch ( (*cmpdl)(n,m1->dl,m2->dl) ) {
235: case 0:
236: addp(vl,C(m1),C(m2),&t);
237: if ( t ) {
238: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = t;
239: }
240: m1 = NEXT(m1); m2 = NEXT(m2); break;
241: case 1:
242: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = C(m1);
243: m1 = NEXT(m1); break;
244: case -1:
245: NEXTMP(mr0,mr); mr->dl = m2->dl; C(mr) = C(m2);
246: m2 = NEXT(m2); break;
247: }
248: if ( !mr0 )
249: if ( m1 )
250: mr0 = m1;
251: else if ( m2 )
252: mr0 = m2;
253: else {
254: *pr = 0;
255: return;
256: }
257: else if ( m1 )
258: NEXT(mr) = m1;
259: else if ( m2 )
260: NEXT(mr) = m2;
261: else
262: NEXT(mr) = 0;
263: MKDP(NV(p1),mr0,*pr);
264: if ( *pr )
265: (*pr)->sugar = MAX(p1->sugar,p2->sugar);
266: }
267: }
268:
269: /* for F4 symbolic reduction */
270:
271: void symb_addd(p1,p2,pr)
272: DP p1,p2,*pr;
273: {
274: int n;
275: MP m1,m2,mr,mr0;
276: P t;
277:
278: if ( !p1 )
279: *pr = p2;
280: else if ( !p2 )
281: *pr = p1;
282: else {
283: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) {
284: NEXTMP(mr0,mr); C(mr) = (P)ONE;
285: switch ( (*cmpdl)(n,m1->dl,m2->dl) ) {
286: case 0:
287: mr->dl = m1->dl;
288: m1 = NEXT(m1); m2 = NEXT(m2); break;
289: case 1:
290: mr->dl = m1->dl;
291: m1 = NEXT(m1); break;
292: case -1:
293: mr->dl = m2->dl;
294: m2 = NEXT(m2); break;
295: }
296: }
297: if ( !mr0 )
298: if ( m1 )
299: mr0 = m1;
300: else if ( m2 )
301: mr0 = m2;
302: else {
303: *pr = 0;
304: return;
305: }
306: else if ( m1 )
307: NEXT(mr) = m1;
308: else if ( m2 )
309: NEXT(mr) = m2;
310: else
311: NEXT(mr) = 0;
312: MKDP(NV(p1),mr0,*pr);
313: if ( *pr )
314: (*pr)->sugar = MAX(p1->sugar,p2->sugar);
315: }
316: }
317:
318: void subd(vl,p1,p2,pr)
319: VL vl;
320: DP p1,p2,*pr;
321: {
322: DP t;
323:
324: if ( !p2 )
325: *pr = p1;
326: else {
327: chsgnd(p2,&t); addd(vl,p1,t,pr);
328: }
329: }
330:
331: void chsgnd(p,pr)
332: DP p,*pr;
333: {
334: MP m,mr,mr0;
335:
336: if ( !p )
337: *pr = 0;
338: else {
339: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
340: NEXTMP(mr0,mr); chsgnp(C(m),&C(mr)); mr->dl = m->dl;
341: }
342: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
343: if ( *pr )
344: (*pr)->sugar = p->sugar;
345: }
346: }
347:
348: void muld(vl,p1,p2,pr)
349: VL vl;
350: DP p1,p2,*pr;
351: {
1.2 ! noro 352: if ( ! do_weyl )
! 353: comm_muld(vl,p1,p2,pr);
! 354: else
! 355: weyl_muld(vl,p1,p2,pr);
! 356: }
! 357:
! 358: void comm_muld(vl,p1,p2,pr)
! 359: VL vl;
! 360: DP p1,p2,*pr;
! 361: {
1.1 noro 362: MP m;
363: DP s,t,u;
364:
365: if ( !p1 || !p2 )
366: *pr = 0;
367: else if ( OID(p1) <= O_P )
368: muldc(vl,p2,(P)p1,pr);
369: else if ( OID(p2) <= O_P )
370: muldc(vl,p1,(P)p2,pr);
371: else {
372: for ( m = BDY(p2), s = 0; m; m = NEXT(m) ) {
373: muldm(vl,p1,m,&t); addd(vl,s,t,&u); s = u;
374: }
375: *pr = s;
376: }
377: }
378:
379: void muldm(vl,p,m0,pr)
380: VL vl;
381: DP p;
382: MP m0;
383: DP *pr;
384: {
385: MP m,mr,mr0;
386: P c;
387: DL d;
388: int n;
389:
390: if ( !p )
391: *pr = 0;
392: else {
393: for ( mr0 = 0, m = BDY(p), c = C(m0), d = m0->dl, n = NV(p);
394: m; m = NEXT(m) ) {
395: NEXTMP(mr0,mr);
396: if ( NUM(C(m)) && RATN(C(m)) && NUM(c) && RATN(c) )
397: mulq((Q)C(m),(Q)c,(Q *)&C(mr));
398: else
399: mulp(vl,C(m),c,&C(mr));
400: adddl(n,m->dl,d,&mr->dl);
401: }
402: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
403: if ( *pr )
404: (*pr)->sugar = p->sugar + m0->dl->td;
1.2 ! noro 405: }
! 406: }
! 407:
! 408: void weyl_muld(vl,p1,p2,pr)
! 409: VL vl;
! 410: DP p1,p2,*pr;
! 411: {
! 412: MP m;
! 413: DP s,t,u;
! 414:
! 415: if ( !p1 || !p2 )
! 416: *pr = 0;
! 417: else if ( OID(p1) <= O_P )
! 418: muldc(vl,p2,(P)p1,pr);
! 419: else if ( OID(p2) <= O_P )
! 420: muldc(vl,p1,(P)p2,pr);
! 421: else {
! 422: for ( m = BDY(p2), s = 0; m; m = NEXT(m) ) {
! 423: weyl_muldm(vl,p1,m,&t); addd(vl,s,t,&u); s = u;
! 424: }
! 425: *pr = s;
! 426: }
! 427: }
! 428:
! 429: void weyl_muldm(vl,p,m0,pr)
! 430: VL vl;
! 431: DP p;
! 432: MP m0;
! 433: DP *pr;
! 434: {
! 435: DP r,t,t1;
! 436: MP m;
! 437: int n;
! 438:
! 439: if ( !p )
! 440: *pr = 0;
! 441: else {
! 442: for ( r = 0, m = BDY(p), n = NV(p); m; m = NEXT(m) ) {
! 443: weyl_mulmm(vl,m,m0,n,&t);
! 444: addd(vl,r,t,&t1); r = t1;
! 445: }
! 446: if ( r )
! 447: r->sugar = p->sugar + m0->dl->td;
! 448: *pr = r;
! 449: }
! 450: }
! 451:
! 452: /* m0 = x0^d0*x1^d1*... * dx0^d(n/2)*dx1^d(n/2+1)*... */
! 453:
! 454: void weyl_mulmm(vl,m0,m1,n,pr)
! 455: VL vl;
! 456: MP m0,m1;
! 457: int n;
! 458: DP *pr;
! 459: {
! 460: MP m,mr,mr0;
! 461: DP r,t,t1;
! 462: P c,c0,c1,cc;
! 463: DL d,d0,d1;
! 464: int i,j,a,b,k,l,n2,s,min;
! 465: Q *tab;
! 466:
! 467: if ( !m0 || !m1 )
! 468: *pr = 0;
! 469: else {
! 470: c0 = C(m0); c1 = C(m1);
! 471: mulp(vl,c0,c1,&c);
! 472: d0 = m0->dl; d1 = m1->dl;
! 473: n2 = n>>1;
! 474: for ( i = 0, r = (DP)ONE; i < n2; i++ ) {
! 475: a = d0->d[i]; b = d1->d[n2+i];
! 476: k = d0->d[n2+i]; l = d1->d[i];
! 477: /* compute xi^a*(Di^k*xi^l)*Di^b */
! 478: min = MIN(k,l);
! 479: tab = (Q *)ALLOCA((min+1)*sizeof(Q));
! 480: mkwc(k,l,tab);
! 481: for ( mr0 = 0, s = 0, j = 0; j <= min; j++ ) {
! 482: NEXTMP(mr0,mr);
! 483: NEWDL(d,n);
! 484: d->d[i] = l-j+a; d->d[n2+i] = k-j+b;
! 485: d->td = d->d[i]+d->d[n2+i]; /* XXX */
! 486: s = MAX(s,d->td); /* XXX */
! 487: mr->c = (P)tab[j];
! 488: mr->dl = d;
! 489: }
! 490: if ( mr0 )
! 491: NEXT(mr) = 0;
! 492: MKDP(n,mr0,t);
! 493: if ( t )
! 494: t->sugar = s;
! 495: comm_muld(vl,r,t,&t1); r = t1;
! 496: }
! 497: muldc(vl,r,c,pr);
1.1 noro 498: }
499: }
500:
501: void muldc(vl,p,c,pr)
502: VL vl;
503: DP p;
504: P c;
505: DP *pr;
506: {
507: MP m,mr,mr0;
508:
509: if ( !p || !c )
510: *pr = 0;
511: else if ( NUM(c) && UNIQ((Q)c) )
512: *pr = p;
513: else if ( NUM(c) && MUNIQ((Q)c) )
514: chsgnd(p,pr);
515: else {
516: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
517: NEXTMP(mr0,mr);
518: if ( NUM(C(m)) && RATN(C(m)) && NUM(c) && RATN(c) )
519: mulq((Q)C(m),(Q)c,(Q *)&C(mr));
520: else
521: mulp(vl,C(m),c,&C(mr));
522: mr->dl = m->dl;
523: }
524: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
525: if ( *pr )
526: (*pr)->sugar = p->sugar;
527: }
528: }
529:
530: void divsdc(vl,p,c,pr)
531: VL vl;
532: DP p;
533: P c;
534: DP *pr;
535: {
536: MP m,mr,mr0;
537:
538: if ( !c )
539: error("disvsdc : division by 0");
540: else if ( !p )
541: *pr = 0;
542: else {
543: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
544: NEXTMP(mr0,mr); divsp(vl,C(m),c,&C(mr)); mr->dl = m->dl;
545: }
546: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
547: if ( *pr )
548: (*pr)->sugar = p->sugar;
549: }
550: }
551:
552: void adddl(n,d1,d2,dr)
553: int n;
554: DL d1,d2;
555: DL *dr;
556: {
557: DL dt;
558: int i;
559:
560: if ( !d1->td )
561: *dr = d2;
562: else if ( !d2->td )
563: *dr = d1;
564: else {
565: *dr = dt = (DL)MALLOC_ATOMIC((n+1)*sizeof(int));
566: dt->td = d1->td + d2->td;
567: for ( i = 0; i < n; i++ )
568: dt->d[i] = d1->d[i]+d2->d[i];
569: }
570: }
571:
572: int compd(vl,p1,p2)
573: VL vl;
574: DP p1,p2;
575: {
576: int n,t;
577: MP m1,m2;
578:
579: if ( !p1 )
580: return p2 ? -1 : 0;
581: else if ( !p2 )
582: return 1;
583: else {
584: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2);
585: m1 && m2; m1 = NEXT(m1), m2 = NEXT(m2) )
586: if ( (t = (*cmpdl)(n,m1->dl,m2->dl)) ||
587: (t = compp(vl,C(m1),C(m2)) ) )
588: return t;
589: if ( m1 )
590: return 1;
591: else if ( m2 )
592: return -1;
593: else
594: return 0;
595: }
596: }
597:
598: int cmpdl_lex(n,d1,d2)
599: int n;
600: DL d1,d2;
601: {
602: int i;
603:
604: for ( i = 0; i < n && d1->d[i] == d2->d[i]; i++ );
605: return i == n ? 0 : (d1->d[i] > d2->d[i] ? 1 : -1);
606: }
607:
608: int cmpdl_revlex(n,d1,d2)
609: int n;
610: DL d1,d2;
611: {
612: int i;
613:
614: for ( i = n - 1; i >= 0 && d1->d[i] == d2->d[i]; i-- );
615: return i < 0 ? 0 : (d1->d[i] < d2->d[i] ? 1 : -1);
616: }
617:
618: int cmpdl_gradlex(n,d1,d2)
619: int n;
620: DL d1,d2;
621: {
622: if ( d1->td > d2->td )
623: return 1;
624: else if ( d1->td < d2->td )
625: return -1;
626: else
627: return cmpdl_lex(n,d1,d2);
628: }
629:
630: int cmpdl_revgradlex(n,d1,d2)
631: int n;
632: DL d1,d2;
633: {
634: if ( d1->td > d2->td )
635: return 1;
636: else if ( d1->td < d2->td )
637: return -1;
638: else
639: return cmpdl_revlex(n,d1,d2);
640: }
641:
642: int cmpdl_blex(n,d1,d2)
643: int n;
644: DL d1,d2;
645: {
646: int c;
647:
648: if ( c = cmpdl_lex(n-1,d1,d2) )
649: return c;
650: else {
651: c = d1->d[n-1] - d2->d[n-1];
652: return c > 0 ? 1 : c < 0 ? -1 : 0;
653: }
654: }
655:
656: int cmpdl_bgradlex(n,d1,d2)
657: int n;
658: DL d1,d2;
659: {
660: int e1,e2,c;
661:
662: e1 = d1->td - d1->d[n-1]; e2 = d2->td - d2->d[n-1];
663: if ( e1 > e2 )
664: return 1;
665: else if ( e1 < e2 )
666: return -1;
667: else {
668: c = cmpdl_lex(n-1,d1,d2);
669: if ( c )
670: return c;
671: else
672: return d1->td > d2->td ? 1 : d1->td < d2->td ? -1 : 0;
673: }
674: }
675:
676: int cmpdl_brevgradlex(n,d1,d2)
677: int n;
678: DL d1,d2;
679: {
680: int e1,e2,c;
681:
682: e1 = d1->td - d1->d[n-1]; e2 = d2->td - d2->d[n-1];
683: if ( e1 > e2 )
684: return 1;
685: else if ( e1 < e2 )
686: return -1;
687: else {
688: c = cmpdl_revlex(n-1,d1,d2);
689: if ( c )
690: return c;
691: else
692: return d1->td > d2->td ? 1 : d1->td < d2->td ? -1 : 0;
693: }
694: }
695:
696: int cmpdl_brevrev(n,d1,d2)
697: int n;
698: DL d1,d2;
699: {
700: int e1,e2,f1,f2,c,i;
701:
702: for ( i = 0, e1 = 0, e2 = 0; i < dp_nelim; i++ ) {
703: e1 += d1->d[i]; e2 += d2->d[i];
704: }
705: f1 = d1->td - e1; f2 = d2->td - e2;
706: if ( e1 > e2 )
707: return 1;
708: else if ( e1 < e2 )
709: return -1;
710: else {
711: c = cmpdl_revlex(dp_nelim,d1,d2);
712: if ( c )
713: return c;
714: else if ( f1 > f2 )
715: return 1;
716: else if ( f1 < f2 )
717: return -1;
718: else {
719: for ( i = n - 1; i >= dp_nelim && d1->d[i] == d2->d[i]; i-- );
720: return i < dp_nelim ? 0 : (d1->d[i] < d2->d[i] ? 1 : -1);
721: }
722: }
723: }
724:
725: int cmpdl_bgradrev(n,d1,d2)
726: int n;
727: DL d1,d2;
728: {
729: int e1,e2,f1,f2,c,i;
730:
731: for ( i = 0, e1 = 0, e2 = 0; i < dp_nelim; i++ ) {
732: e1 += d1->d[i]; e2 += d2->d[i];
733: }
734: f1 = d1->td - e1; f2 = d2->td - e2;
735: if ( e1 > e2 )
736: return 1;
737: else if ( e1 < e2 )
738: return -1;
739: else {
740: c = cmpdl_lex(dp_nelim,d1,d2);
741: if ( c )
742: return c;
743: else if ( f1 > f2 )
744: return 1;
745: else if ( f1 < f2 )
746: return -1;
747: else {
748: for ( i = n - 1; i >= dp_nelim && d1->d[i] == d2->d[i]; i-- );
749: return i < dp_nelim ? 0 : (d1->d[i] < d2->d[i] ? 1 : -1);
750: }
751: }
752: }
753:
754: int cmpdl_blexrev(n,d1,d2)
755: int n;
756: DL d1,d2;
757: {
758: int e1,e2,f1,f2,c,i;
759:
760: for ( i = 0, e1 = 0, e2 = 0; i < dp_nelim; i++ ) {
761: e1 += d1->d[i]; e2 += d2->d[i];
762: }
763: f1 = d1->td - e1; f2 = d2->td - e2;
764: c = cmpdl_lex(dp_nelim,d1,d2);
765: if ( c )
766: return c;
767: else if ( f1 > f2 )
768: return 1;
769: else if ( f1 < f2 )
770: return -1;
771: else {
772: for ( i = n - 1; i >= dp_nelim && d1->d[i] == d2->d[i]; i-- );
773: return i < dp_nelim ? 0 : (d1->d[i] < d2->d[i] ? 1 : -1);
774: }
775: }
776:
777: int cmpdl_elim(n,d1,d2)
778: int n;
779: DL d1,d2;
780: {
781: int e1,e2,i;
782:
783: for ( i = 0, e1 = 0, e2 = 0; i < dp_nelim; i++ ) {
784: e1 += d1->d[i]; e2 += d2->d[i];
785: }
786: if ( e1 > e2 )
787: return 1;
788: else if ( e1 < e2 )
789: return -1;
790: else
791: return cmpdl_revgradlex(n,d1,d2);
792: }
793:
794: int cmpdl_order_pair(n,d1,d2)
795: int n;
796: DL d1,d2;
797: {
798: int e1,e2,i,j,l;
799: int *t1,*t2;
800: int len;
801: struct order_pair *pair;
802:
803: len = dp_current_spec.ord.block.length;
804: pair = dp_current_spec.ord.block.order_pair;
805:
806: for ( i = 0, t1 = d1->d, t2 = d2->d; i < len; i++ ) {
807: l = pair[i].length;
808: switch ( pair[i].order ) {
809: case 0:
810: for ( j = 0, e1 = e2 = 0; j < l; j++ ) {
811: e1 += t1[j]; e2 += t2[j];
812: }
813: if ( e1 > e2 )
814: return 1;
815: else if ( e1 < e2 )
816: return -1;
817: else {
818: for ( j = l - 1; j >= 0 && t1[j] == t2[j]; j-- );
819: if ( j >= 0 )
820: return t1[j] < t2[j] ? 1 : -1;
821: }
822: break;
823: case 1:
824: for ( j = 0, e1 = e2 = 0; j < l; j++ ) {
825: e1 += t1[j]; e2 += t2[j];
826: }
827: if ( e1 > e2 )
828: return 1;
829: else if ( e1 < e2 )
830: return -1;
831: else {
832: for ( j = 0; j < l && t1[j] == t2[j]; j++ );
833: if ( j < l )
834: return t1[j] > t2[j] ? 1 : -1;
835: }
836: break;
837: case 2:
838: for ( j = 0; j < l && t1[j] == t2[j]; j++ );
839: if ( j < l )
840: return t1[j] > t2[j] ? 1 : -1;
841: break;
842: default:
843: error("cmpdl_order_pair : invalid order"); break;
844: }
845: t1 += l; t2 += l;
846: }
847: return 0;
848: }
849:
850: int cmpdl_matrix(n,d1,d2)
851: int n;
852: DL d1,d2;
853: {
854: int *v,*w,*t1,*t2;
855: int s,i,j,len;
856: int **matrix;
857:
858: for ( i = 0, t1 = d1->d, t2 = d2->d, w = dp_dl_work; i < n; i++ )
859: w[i] = t1[i]-t2[i];
860: len = dp_current_spec.ord.matrix.row;
861: matrix = dp_current_spec.ord.matrix.matrix;
862: for ( j = 0; j < len; j++ ) {
863: v = matrix[j];
864: for ( i = 0, s = 0; i < n; i++ )
865: s += v[i]*w[i];
866: if ( s > 0 )
867: return 1;
868: else if ( s < 0 )
869: return -1;
870: }
871: return 0;
872: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>