![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to dp-supp.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / builtin|
Return to dp-supp.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / builtin |
1.2 noro 1: /*
2: * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
3: * All rights reserved.
4: *
5: * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
6: * non-exclusive and royalty-free license to use, copy, modify and
7: * redistribute, solely for non-commercial and non-profit purposes, the
8: * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
9: * conditions of this Agreement. For the avoidance of doubt, you acquire
10: * only a limited right to use the SOFTWARE hereunder, and FLL or any
11: * third party developer retains all rights, including but not limited to
12: * copyrights, in and to the SOFTWARE.
13: *
14: * (1) FLL does not grant you a license in any way for commercial
15: * purposes. You may use the SOFTWARE only for non-commercial and
16: * non-profit purposes only, such as academic, research and internal
17: * business use.
18: * (2) The SOFTWARE is protected by the Copyright Law of Japan and
19: * international copyright treaties. If you make copies of the SOFTWARE,
20: * with or without modification, as permitted hereunder, you shall affix
21: * to all such copies of the SOFTWARE the above copyright notice.
22: * (3) An explicit reference to this SOFTWARE and its copyright owner
23: * shall be made on your publication or presentation in any form of the
24: * results obtained by use of the SOFTWARE.
25: * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
1.3 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.2 noro 27: * for such modification or the source code of the modified part of the
28: * SOFTWARE.
29: *
30: * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
31: * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
32: * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
33: * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
34: * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
35: * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
36: * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
37: * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
38: * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
39: * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
40: * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
41: * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
42: * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
43: * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
44: * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
45: * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
46: * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
47: *
1.38 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/dp-supp.c,v 1.37 2004/09/15 06:06:42 noro Exp $
1.2 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include "base.h"
1.16 noro 52: #include "inline.h"
1.1 noro 53: #include "parse.h"
54: #include "ox.h"
55:
1.5 noro 56: #define HMAG(p) (p_mag(BDY(p)->c))
57:
1.1 noro 58: extern int (*cmpdl)();
1.5 noro 59: extern double pz_t_e,pz_t_d,pz_t_d1,pz_t_c;
60: extern int dp_nelim,dp_fcoeffs;
1.7 noro 61: extern int NoGCD;
62: extern int GenTrace;
63: extern NODE TraceList;
64:
1.37 noro 65: int show_orderspec;
66:
67: void print_composite_order_spec(struct order_spec *spec);
68:
1.7 noro 69: /*
70: * content reduction
71: *
72: */
73:
1.20 noro 74: void dp_ptozp(DP p,DP *rp)
1.7 noro 75: {
76: MP m,mr,mr0;
77: int i,n;
78: Q *w;
79: Q dvr;
80: P t;
81:
82: if ( !p )
83: *rp = 0;
84: else {
85: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
86: w = (Q *)ALLOCA(n*sizeof(Q));
87: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
88: if ( NUM(m->c) )
89: w[i] = (Q)m->c;
90: else
91: ptozp(m->c,1,&w[i],&t);
92: sortbynm(w,n);
93: qltozl(w,n,&dvr);
94: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
95: NEXTMP(mr0,mr); divsp(CO,m->c,(P)dvr,&mr->c); mr->dl = m->dl;
96: }
97: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
98: }
99: }
100:
1.20 noro 101: void dp_ptozp2(DP p0,DP p1,DP *hp,DP *rp)
1.7 noro 102: {
103: DP t,s,h,r;
104: MP m,mr,mr0,m0;
105:
106: addd(CO,p0,p1,&t); dp_ptozp(t,&s);
107: if ( !p0 ) {
108: h = 0; r = s;
109: } else if ( !p1 ) {
110: h = s; r = 0;
111: } else {
112: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
113: m = NEXT(m), m0 = NEXT(m0) ) {
114: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
115: }
116: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
117: }
118: if ( h )
119: h->sugar = p0->sugar;
120: if ( r )
121: r->sugar = p1->sugar;
122: *hp = h; *rp = r;
123: }
1.1 noro 124:
1.20 noro 125: void dp_idiv(DP p,Q c,DP *rp)
1.1 noro 126: {
127: Q t;
128: N nm,q;
129: int sgn,s;
130: MP mr0,m,mr;
131:
132: if ( !p )
133: *rp = 0;
134: else if ( MUNIQ((Q)c) )
135: *rp = p;
136: else if ( MUNIQ((Q)c) )
137: chsgnd(p,rp);
138: else {
139: nm = NM(c); sgn = SGN(c);
140: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
141: NEXTMP(mr0,mr);
142:
143: divsn(NM((Q)(m->c)),nm,&q);
144: s = sgn*SGN((Q)(m->c));
145: NTOQ(q,s,t);
146: mr->c = (P)t;
147: mr->dl = m->dl;
148: }
149: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
150: if ( *rp )
151: (*rp)->sugar = p->sugar;
152: }
153: }
154:
1.20 noro 155: void dp_mbase(NODE hlist,NODE *mbase)
1.1 noro 156: {
157: DL *dl;
158: DL d;
159: int i,j,n,nvar,td;
160:
161: n = length(hlist); nvar = ((DP)BDY(hlist))->nv;
162: dl = (DL *)MALLOC(n*sizeof(DL));
163: for ( i = 0; i < n; i++, hlist = NEXT(hlist) )
164: dl[i] = BDY((DP)BDY(hlist))->dl;
165: NEWDL(d,nvar); *mbase = 0;
166: while ( 1 ) {
167: insert_to_node(d,mbase,nvar);
168: for ( i = nvar-1; i >= 0; ) {
1.21 noro 169: d->d[i]++;
170: d->td += MUL_WEIGHT(1,i);
1.1 noro 171: for ( j = 0; j < n; j++ ) {
172: if ( _dl_redble(dl[j],d,nvar) )
173: break;
174: }
175: if ( j < n ) {
176: for ( j = nvar-1; j >= i; j-- )
177: d->d[j] = 0;
178: for ( j = 0, td = 0; j < i; j++ )
1.21 noro 179: td += MUL_WEIGHT(d->d[j],j);
1.1 noro 180: d->td = td;
181: i--;
182: } else
183: break;
184: }
185: if ( i < 0 )
186: break;
187: }
188: }
189:
1.20 noro 190: int _dl_redble(DL d1,DL d2,int nvar)
1.1 noro 191: {
192: int i;
193:
194: if ( d1->td > d2->td )
195: return 0;
196: for ( i = 0; i < nvar; i++ )
197: if ( d1->d[i] > d2->d[i] )
198: break;
199: if ( i < nvar )
200: return 0;
201: else
202: return 1;
203: }
204:
1.20 noro 205: void insert_to_node(DL d,NODE *n,int nvar)
1.1 noro 206: {
207: DL d1;
208: MP m;
209: DP dp;
210: NODE n0,n1,n2;
211:
212: NEWDL(d1,nvar); d1->td = d->td;
213: bcopy((char *)d->d,(char *)d1->d,nvar*sizeof(int));
214: NEWMP(m); m->dl = d1; m->c = (P)ONE; NEXT(m) = 0;
215: MKDP(nvar,m,dp); dp->sugar = d->td;
216: if ( !(*n) ) {
217: MKNODE(n1,dp,0); *n = n1;
218: } else {
219: for ( n1 = *n, n0 = 0; n1; n0 = n1, n1 = NEXT(n1) )
220: if ( (*cmpdl)(nvar,d,BDY((DP)BDY(n1))->dl) > 0 ) {
221: MKNODE(n2,dp,n1);
222: if ( !n0 )
223: *n = n2;
224: else
225: NEXT(n0) = n2;
226: break;
227: }
228: if ( !n1 ) {
229: MKNODE(n2,dp,0); NEXT(n0) = n2;
230: }
231: }
232: }
233:
1.20 noro 234: void dp_vtod(Q *c,DP p,DP *rp)
1.1 noro 235: {
236: MP mr0,m,mr;
237: int i;
238:
239: if ( !p )
240: *rp = 0;
241: else {
242: for ( mr0 = 0, m = BDY(p), i = 0; m; m = NEXT(m), i++ ) {
243: NEXTMP(mr0,mr); mr->c = (P)c[i]; mr->dl = m->dl;
244: }
245: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
246: (*rp)->sugar = p->sugar;
247: }
248: }
249:
1.8 noro 250: extern int mpi_mag;
251: extern int PCoeffs;
252:
1.20 noro 253: void dp_ptozp_d(DP p,DP *rp)
1.1 noro 254: {
255: int i,j,k,l,n,nsep;
256: MP m;
257: NODE tn,n0,n1,n2,n3;
258: struct oVECT v;
259: VECT c,cs;
260: VECT qi,ri;
261: LIST *qr;
262: Obj dmy;
263: Q d0,d1,gcd,a,u,u1;
264: Q *q,*r;
265: STRING iqr_v;
266: pointer *b;
267: N qn,gn;
268: double get_rtime();
269: int blen;
1.8 noro 270: NODE dist;
271: int ndist;
1.1 noro 272: double t0;
273: double t_e,t_d,t_d1,t_c;
1.8 noro 274: extern int DP_NFStat;
275: extern LIST Dist;
1.20 noro 276: void Pox_rpc();
277: void Pox_pop_local();
1.1 noro 278:
279: if ( !p )
280: *rp = 0;
281: else {
1.8 noro 282: if ( PCoeffs ) {
283: dp_ptozp(p,rp); return;
284: }
1.9 noro 285: if ( !Dist || p_mag(BDY(p)->c) <= mpi_mag ) {
1.8 noro 286: dist = 0; ndist = 0;
287: if ( DP_NFStat ) fprintf(asir_out,"L");
288: } else {
289: dist = BDY(Dist); ndist = length(dist);
290: if ( DP_NFStat ) fprintf(asir_out,"D");
291: }
1.1 noro 292: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
293: nsep = ndist + 1;
294: if ( n <= nsep ) {
295: dp_ptozp(p,rp); return;
296: }
297: t0 = get_rtime();
298: dp_dtov(p,&c);
299: igcdv_estimate(c,&d0);
300: t_e = get_rtime()-t0;
301: t0 = get_rtime();
302: dp_dtov(p,&c);
303: sepvect(c,nsep,&cs);
304: MKSTR(iqr_v,"iqr");
305: qr = (LIST *)CALLOC(nsep,sizeof(LIST));
306: q = (Q *)CALLOC(n,sizeof(Q));
307: r = (Q *)CALLOC(n,sizeof(Q));
308: for ( i = 0, tn = dist, b = BDY(cs); i < ndist; i++, tn = NEXT(tn) ) {
309: MKNODE(n3,d0,0); MKNODE(n2,b[i],n3);
310: MKNODE(n1,iqr_v,n2); MKNODE(n0,BDY(tn),n1);
311: Pox_rpc(n0,&dmy);
312: }
313: iqrv(b[i],d0,&qr[i]);
314: dp_dtov(p,&c);
315: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
316: Pox_pop_local(tn,&qr[i]);
317: if ( OID(qr[i]) == O_ERR ) {
318: printexpr(CO,(Obj)qr[i]);
319: error("dp_ptozp_d : aborted");
320: }
321: }
322: t_d = get_rtime()-t0;
323: t_d1 = t_d/n;
324: t0 = get_rtime();
325: for ( i = j = 0; i < nsep; i++ ) {
326: tn = BDY(qr[i]); qi = (VECT)BDY(tn); ri = (VECT)BDY(NEXT(tn));
327: for ( k = 0, l = qi->len; k < l; k++, j++ ) {
328: q[j] = (Q)BDY(qi)[k]; r[j] = (Q)BDY(ri)[k];
329: }
330: }
331: v.id = O_VECT; v.len = n; v.body = (pointer *)r; igcdv(&v,&d1);
332: if ( d1 ) {
333: gcdn(NM(d0),NM(d1),&gn); NTOQ(gn,1,gcd);
334: divsn(NM(d0),gn,&qn); NTOQ(qn,1,a);
335: for ( i = 0; i < n; i++ ) {
336: mulq(a,q[i],&u);
337: if ( r[i] ) {
338: divsn(NM(r[i]),gn,&qn); NTOQ(qn,SGN(r[i]),u1);
339: addq(u,u1,&q[i]);
340: } else
341: q[i] = u;
342: }
343: } else
344: gcd = d0;
345: dp_vtod(q,p,rp);
346: t_c = get_rtime()-t0;
347: blen=p_mag((P)gcd);
348: pz_t_e += t_e; pz_t_d += t_d; pz_t_d1 += t_d1; pz_t_c += t_c;
349: if ( 0 )
350: fprintf(stderr,"(%d,%d)",p_mag((P)d0)-blen,blen);
351: }
352: }
353:
1.20 noro 354: void dp_ptozp2_d(DP p0,DP p1,DP *hp,DP *rp)
1.1 noro 355: {
356: DP t,s,h,r;
357: MP m,mr,mr0,m0;
358:
1.8 noro 359: addd(CO,p0,p1,&t); dp_ptozp_d(t,&s);
1.1 noro 360: if ( !p0 ) {
361: h = 0; r = s;
362: } else if ( !p1 ) {
363: h = s; r = 0;
364: } else {
365: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
366: m = NEXT(m), m0 = NEXT(m0) ) {
367: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
368: }
369: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
370: }
371: if ( h )
372: h->sugar = p0->sugar;
373: if ( r )
374: r->sugar = p1->sugar;
375: *hp = h; *rp = r;
1.5 noro 376: }
377:
1.22 noro 378: int have_sf_coef(P p)
379: {
380: DCP dc;
381:
382: if ( !p )
383: return 0;
384: else if ( NUM(p) )
385: return NID((Num)p) == N_GFS ? 1 : 0;
386: else {
387: for ( dc = DC(p); dc; dc = NEXT(dc) )
388: if ( have_sf_coef(COEF(dc)) )
389: return 1;
390: return 0;
391: }
392: }
393:
1.25 noro 394: void head_coef(P p,Num *c)
395: {
396: if ( !p )
397: *c = 0;
398: else if ( NUM(p) )
399: *c = (Num)p;
400: else
401: head_coef(COEF(DC(p)),c);
402: }
403:
404: void dp_monic_sf(DP p,DP *rp)
405: {
406: Num c;
407:
408: if ( !p )
409: *rp = 0;
410: else {
411: head_coef(BDY(p)->c,&c);
412: divsdc(CO,p,(P)c,rp);
413: }
414: }
415:
1.20 noro 416: void dp_prim(DP p,DP *rp)
1.5 noro 417: {
1.7 noro 418: P t,g;
419: DP p1;
420: MP m,mr,mr0;
421: int i,n;
422: P *w;
423: Q *c;
424: Q dvr;
1.5 noro 425:
1.7 noro 426: if ( !p )
427: *rp = 0;
1.23 noro 428: else if ( dp_fcoeffs == N_GFS ) {
429: for ( m = BDY(p); m; m = NEXT(m) )
1.22 noro 430: if ( OID(m->c) == O_N ) {
431: /* GCD of coeffs = 1 */
1.25 noro 432: dp_monic_sf(p,rp);
1.22 noro 433: return;
1.23 noro 434: } else break;
435: /* compute GCD over the finite fieid */
436: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
437: w = (P *)ALLOCA(n*sizeof(P));
438: for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ )
439: w[i] = m->c;
440: gcdsf(CO,w,n,&g);
441: if ( NUM(g) )
1.25 noro 442: dp_monic_sf(p,rp);
1.23 noro 443: else {
444: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
445: NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl;
1.22 noro 446: }
1.25 noro 447: NEXT(mr) = 0; MKDP(p->nv,mr0,p1); p1->sugar = p->sugar;
448: dp_monic_sf(p1,rp);
1.22 noro 449: }
1.23 noro 450: return;
451: } else if ( dp_fcoeffs )
1.7 noro 452: *rp = p;
1.23 noro 453: else if ( NoGCD )
1.7 noro 454: dp_ptozp(p,rp);
455: else {
456: dp_ptozp(p,&p1); p = p1;
457: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
458: if ( n == 1 ) {
459: m = BDY(p);
460: NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0;
461: MKDP(p->nv,mr,*rp); (*rp)->sugar = p->sugar;
462: return;
463: }
464: w = (P *)ALLOCA(n*sizeof(P));
465: c = (Q *)ALLOCA(n*sizeof(Q));
466: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
467: if ( NUM(m->c) ) {
468: c[i] = (Q)m->c; w[i] = (P)ONE;
469: } else
470: ptozp(m->c,1,&c[i],&w[i]);
471: qltozl(c,n,&dvr); heu_nezgcdnpz(CO,w,n,&t); mulp(CO,t,(P)dvr,&g);
472: if ( NUM(g) )
473: *rp = p;
474: else {
475: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
476: NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl;
477: }
478: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1.5 noro 479: }
1.7 noro 480: }
1.5 noro 481: }
482:
1.20 noro 483: void heu_nezgcdnpz(VL vl,P *pl,int m,P *pr)
1.5 noro 484: {
485: int i,r;
486: P gcd,t,s1,s2,u;
487: Q rq;
488:
489: while ( 1 ) {
490: for ( i = 0, s1 = 0; i < m; i++ ) {
491: r = random(); UTOQ(r,rq);
492: mulp(vl,pl[i],(P)rq,&t); addp(vl,s1,t,&u); s1 = u;
493: }
494: for ( i = 0, s2 = 0; i < m; i++ ) {
495: r = random(); UTOQ(r,rq);
496: mulp(vl,pl[i],(P)rq,&t); addp(vl,s2,t,&u); s2 = u;
497: }
498: ezgcdp(vl,s1,s2,&gcd);
499: for ( i = 0; i < m; i++ ) {
500: if ( !divtpz(vl,pl[i],gcd,&t) )
501: break;
502: }
503: if ( i == m )
504: break;
505: }
506: *pr = gcd;
507: }
508:
1.20 noro 509: void dp_prim_mod(DP p,int mod,DP *rp)
1.5 noro 510: {
511: P t,g;
512: MP m,mr,mr0;
513:
514: if ( !p )
515: *rp = 0;
516: else if ( NoGCD )
517: *rp = p;
518: else {
519: for ( m = BDY(p), g = m->c, m = NEXT(m); m; m = NEXT(m) ) {
520: gcdprsmp(CO,mod,g,m->c,&t); g = t;
521: }
522: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
523: NEXTMP(mr0,mr); divsmp(CO,mod,m->c,g,&mr->c); mr->dl = m->dl;
524: }
525: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
526: }
527: }
528:
1.20 noro 529: void dp_cont(DP p,Q *rp)
1.5 noro 530: {
1.7 noro 531: VECT v;
1.5 noro 532:
1.7 noro 533: dp_dtov(p,&v); igcdv(v,rp);
1.5 noro 534: }
535:
1.20 noro 536: void dp_dtov(DP dp,VECT *rp)
1.5 noro 537: {
1.7 noro 538: MP m,t;
539: int i,n;
540: VECT v;
541: pointer *p;
1.5 noro 542:
1.7 noro 543: m = BDY(dp);
544: for ( t = m, n = 0; t; t = NEXT(t), n++ );
545: MKVECT(v,n);
546: for ( i = 0, p = BDY(v), t = m; i < n; t = NEXT(t), i++ )
547: p[i] = (pointer)(t->c);
548: *rp = v;
1.5 noro 549: }
550:
1.7 noro 551: /*
552: * s-poly computation
553: *
554: */
1.5 noro 555:
1.20 noro 556: void dp_sp(DP p1,DP p2,DP *rp)
1.5 noro 557: {
1.7 noro 558: int i,n,td;
559: int *w;
560: DL d1,d2,d;
561: MP m;
562: DP t,s1,s2,u;
563: Q c,c1,c2;
564: N gn,tn;
1.5 noro 565:
1.7 noro 566: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
567: w = (int *)ALLOCA(n*sizeof(int));
568: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 569: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.5 noro 570: }
1.7 noro 571:
572: NEWDL(d,n); d->td = td - d1->td;
573: for ( i = 0; i < n; i++ )
574: d->d[i] = w[i] - d1->d[i];
575: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
576: if ( INT(c1) && INT(c2) ) {
577: gcdn(NM(c1),NM(c2),&gn);
578: if ( !UNIN(gn) ) {
579: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
580: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
1.5 noro 581: }
582: }
1.7 noro 583:
584: NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0;
585: MKDP(n,m,s1); s1->sugar = d->td; muld(CO,s1,p1,&t);
586:
587: NEWDL(d,n); d->td = td - d2->td;
588: for ( i = 0; i < n; i++ )
589: d->d[i] = w[i] - d2->d[i];
590: NEWMP(m); m->dl = d; m->c = (P)c1; NEXT(m) = 0;
591: MKDP(n,m,s2); s2->sugar = d->td; muld(CO,s2,p2,&u);
592:
593: subd(CO,t,u,rp);
1.14 noro 594: if ( GenTrace ) {
595: LIST hist;
596: NODE node;
597:
598: node = mknode(4,ONE,0,s1,ONE);
599: MKLIST(hist,node);
600: MKNODE(TraceList,hist,0);
601:
602: node = mknode(4,ONE,0,0,ONE);
603: chsgnd(s2,(DP *)&ARG2(node));
604: MKLIST(hist,node);
605: MKNODE(node,hist,TraceList); TraceList = node;
606: }
607: }
608:
1.20 noro 609: void _dp_sp_dup(DP p1,DP p2,DP *rp)
1.14 noro 610: {
611: int i,n,td;
612: int *w;
613: DL d1,d2,d;
614: MP m;
615: DP t,s1,s2,u;
616: Q c,c1,c2;
617: N gn,tn;
618:
619: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
620: w = (int *)ALLOCA(n*sizeof(int));
621: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 622: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.14 noro 623: }
624:
625: _NEWDL(d,n); d->td = td - d1->td;
626: for ( i = 0; i < n; i++ )
627: d->d[i] = w[i] - d1->d[i];
628: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
629: if ( INT(c1) && INT(c2) ) {
630: gcdn(NM(c1),NM(c2),&gn);
631: if ( !UNIN(gn) ) {
632: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
633: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
634: }
635: }
636:
637: _NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0;
638: _MKDP(n,m,s1); s1->sugar = d->td; _muld_dup(CO,s1,p1,&t); _free_dp(s1);
639:
640: _NEWDL(d,n); d->td = td - d2->td;
641: for ( i = 0; i < n; i++ )
642: d->d[i] = w[i] - d2->d[i];
643: _NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0;
644: _MKDP(n,m,s2); s2->sugar = d->td; _muld_dup(CO,s2,p2,&u); _free_dp(s2);
645:
646: _addd_destructive(CO,t,u,rp);
1.7 noro 647: if ( GenTrace ) {
648: LIST hist;
649: NODE node;
650:
651: node = mknode(4,ONE,0,s1,ONE);
652: MKLIST(hist,node);
653: MKNODE(TraceList,hist,0);
654:
655: node = mknode(4,ONE,0,0,ONE);
656: chsgnd(s2,(DP *)&ARG2(node));
657: MKLIST(hist,node);
658: MKNODE(node,hist,TraceList); TraceList = node;
659: }
660: }
661:
1.20 noro 662: void dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
1.7 noro 663: {
664: int i,n,td;
665: int *w;
666: DL d1,d2,d;
667: MP m;
668: DP t,s,u;
669:
670: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
671: w = (int *)ALLOCA(n*sizeof(int));
672: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 673: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 674: }
1.18 noro 675: NEWDL_NOINIT(d,n); d->td = td - d1->td;
1.7 noro 676: for ( i = 0; i < n; i++ )
677: d->d[i] = w[i] - d1->d[i];
678: NEWMP(m); m->dl = d; m->c = (P)BDY(p2)->c; NEXT(m) = 0;
679: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p1,s,&t);
1.18 noro 680: NEWDL_NOINIT(d,n); d->td = td - d2->td;
1.7 noro 681: for ( i = 0; i < n; i++ )
682: d->d[i] = w[i] - d2->d[i];
683: NEWMP(m); m->dl = d; m->c = (P)BDY(p1)->c; NEXT(m) = 0;
684: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p2,s,&u);
685: submd(CO,mod,t,u,rp);
686: }
687:
1.20 noro 688: void _dp_sp_mod_dup(DP p1,DP p2,int mod,DP *rp)
1.7 noro 689: {
690: int i,n,td;
691: int *w;
692: DL d1,d2,d;
693: MP m;
694: DP t,s,u;
695:
696: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
697: w = (int *)ALLOCA(n*sizeof(int));
698: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 699: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 700: }
701: _NEWDL(d,n); d->td = td - d1->td;
702: for ( i = 0; i < n; i++ )
703: d->d[i] = w[i] - d1->d[i];
704: _NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
705: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p1,&t); _free_dp(s);
706: _NEWDL(d,n); d->td = td - d2->td;
707: for ( i = 0; i < n; i++ )
708: d->d[i] = w[i] - d2->d[i];
709: _NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
710: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&u); _free_dp(s);
711: _addmd_destructive(mod,t,u,rp);
712: }
713:
1.20 noro 714: void _dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
1.7 noro 715: {
716: int i,n,td;
717: int *w;
718: DL d1,d2,d;
719: MP m;
720: DP t,s,u;
721:
722: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
723: w = (int *)ALLOCA(n*sizeof(int));
724: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 725: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 726: }
727: NEWDL(d,n); d->td = td - d1->td;
728: for ( i = 0; i < n; i++ )
729: d->d[i] = w[i] - d1->d[i];
730: NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
731: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p1,&t);
732: NEWDL(d,n); d->td = td - d2->td;
733: for ( i = 0; i < n; i++ )
734: d->d[i] = w[i] - d2->d[i];
735: NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
736: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p2,&u);
737: addmd_destructive(mod,t,u,rp);
738: }
739:
740: /*
741: * m-reduction
1.13 noro 742: * do content reduction over Z or Q(x,...)
743: * do nothing over finite fields
1.7 noro 744: *
745: */
746:
1.20 noro 747: void dp_red(DP p0,DP p1,DP p2,DP *head,DP *rest,P *dnp,DP *multp)
1.7 noro 748: {
749: int i,n;
750: DL d1,d2,d;
751: MP m;
752: DP t,s,r,h;
753: Q c,c1,c2;
754: N gn,tn;
755: P g,a;
1.23 noro 756: P p[2];
1.7 noro 757:
758: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
759: NEWDL(d,n); d->td = d1->td - d2->td;
760: for ( i = 0; i < n; i++ )
761: d->d[i] = d1->d[i]-d2->d[i];
762: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
1.23 noro 763: if ( dp_fcoeffs == N_GFS ) {
764: p[0] = (P)c1; p[1] = (P)c2;
765: gcdsf(CO,p,2,&g);
766: divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a;
767: } else if ( dp_fcoeffs ) {
1.7 noro 768: /* do nothing */
769: } else if ( INT(c1) && INT(c2) ) {
770: gcdn(NM(c1),NM(c2),&gn);
771: if ( !UNIN(gn) ) {
772: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
773: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
774: }
775: } else {
776: ezgcdpz(CO,(P)c1,(P)c2,&g);
777: divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a;
778: }
779: NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
780: *multp = s;
781: muld(CO,s,p2,&t); muldc(CO,p1,(P)c2,&s); addd(CO,s,t,&r);
782: muldc(CO,p0,(P)c2,&h);
783: *head = h; *rest = r; *dnp = (P)c2;
784: }
785:
1.13 noro 786: /* m-reduction over a field */
787:
1.20 noro 788: void dp_red_f(DP p1,DP p2,DP *rest)
1.13 noro 789: {
790: int i,n;
791: DL d1,d2,d;
792: MP m;
1.20 noro 793: DP t,s;
1.13 noro 794: Obj a,b;
795:
796: n = p1->nv;
797: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
798:
799: NEWDL(d,n); d->td = d1->td - d2->td;
800: for ( i = 0; i < n; i++ )
801: d->d[i] = d1->d[i]-d2->d[i];
802:
803: NEWMP(m); m->dl = d;
804: divr(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); chsgnr(a,&b);
805: C(m) = (P)b;
806: NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
807:
808: muld(CO,s,p2,&t); addd(CO,p1,t,rest);
809: }
810:
1.20 noro 811: void dp_red_mod(DP p0,DP p1,DP p2,int mod,DP *head,DP *rest,P *dnp)
1.7 noro 812: {
813: int i,n;
814: DL d1,d2,d;
815: MP m;
816: DP t,s,r,h;
817: P c1,c2,g,u;
818:
819: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
820: NEWDL(d,n); d->td = d1->td - d2->td;
821: for ( i = 0; i < n; i++ )
822: d->d[i] = d1->d[i]-d2->d[i];
823: c1 = (P)BDY(p1)->c; c2 = (P)BDY(p2)->c;
824: gcdprsmp(CO,mod,c1,c2,&g);
825: divsmp(CO,mod,c1,g,&u); c1 = u; divsmp(CO,mod,c2,g,&u); c2 = u;
826: if ( NUM(c2) ) {
827: divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM;
828: }
829: NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,&m->c); NEXT(m) = 0;
1.11 noro 830: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,s,p2,&t);
1.7 noro 831: if ( NUM(c2) ) {
832: addmd(CO,mod,p1,t,&r); h = p0;
833: } else {
834: mulmdc(CO,mod,p1,c2,&s); addmd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h);
835: }
836: *head = h; *rest = r; *dnp = c2;
837: }
838:
1.10 noro 839: struct oEGT eg_red_mod;
840:
1.20 noro 841: void _dp_red_mod_destructive(DP p1,DP p2,int mod,DP *rp)
1.7 noro 842: {
843: int i,n;
844: DL d1,d2,d;
845: MP m;
846: DP t,s;
1.16 noro 847: int c,c1,c2;
848: extern int do_weyl;
1.7 noro 849:
850: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
851: _NEWDL(d,n); d->td = d1->td - d2->td;
852: for ( i = 0; i < n; i++ )
853: d->d[i] = d1->d[i]-d2->d[i];
1.16 noro 854: c = invm(ITOS(BDY(p2)->c),mod);
855: c2 = ITOS(BDY(p1)->c);
856: DMAR(c,c2,0,mod,c1);
1.7 noro 857: _NEWMP(m); m->dl = d; m->c = STOI(mod-c1); NEXT(m) = 0;
1.16 noro 858: #if 0
1.7 noro 859: _MKDP(n,m,s); s->sugar = d->td;
860: _mulmd_dup(mod,s,p2,&t); _free_dp(s);
1.16 noro 861: #else
862: if ( do_weyl ) {
1.19 noro 863: _MKDP(n,m,s); s->sugar = d->td;
864: _mulmd_dup(mod,s,p2,&t); _free_dp(s);
1.16 noro 865: } else {
866: _mulmdm_dup(mod,p2,m,&t); _FREEMP(m);
867: }
868: #endif
1.10 noro 869: /* get_eg(&t0); */
1.7 noro 870: _addmd_destructive(mod,p1,t,rp);
1.10 noro 871: /* get_eg(&t1); add_eg(&eg_red_mod,&t0,&t1); */
1.7 noro 872: }
873:
874: /*
875: * normal form computation
876: *
877: */
1.5 noro 878:
1.20 noro 879: void dp_true_nf(NODE b,DP g,DP *ps,int full,DP *rp,P *dnp)
1.5 noro 880: {
881: DP u,p,d,s,t,dmy;
882: NODE l;
883: MP m,mr;
884: int i,n;
885: int *wb;
886: int sugar,psugar;
887: P dn,tdn,tdn1;
888:
889: dn = (P)ONE;
890: if ( !g ) {
891: *rp = 0; *dnp = dn; return;
892: }
893: for ( n = 0, l = b; l; l = NEXT(l), n++ );
894: wb = (int *)ALLOCA(n*sizeof(int));
895: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
896: wb[i] = QTOS((Q)BDY(l));
897: sugar = g->sugar;
898: for ( d = 0; g; ) {
899: for ( u = 0, i = 0; i < n; i++ ) {
900: if ( dp_redble(g,p = ps[wb[i]]) ) {
901: dp_red(d,g,p,&t,&u,&tdn,&dmy);
902: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
903: sugar = MAX(sugar,psugar);
904: if ( !u ) {
905: if ( d )
906: d->sugar = sugar;
907: *rp = d; *dnp = dn; return;
908: } else {
909: d = t;
910: mulp(CO,dn,tdn,&tdn1); dn = tdn1;
911: }
912: break;
913: }
914: }
915: if ( u )
916: g = u;
917: else if ( !full ) {
918: if ( g ) {
919: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
920: }
921: *rp = g; *dnp = dn; return;
922: } else {
923: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
924: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
925: addd(CO,d,t,&s); d = s;
926: dp_rest(g,&t); g = t;
927: }
928: }
929: if ( d )
930: d->sugar = sugar;
931: *rp = d; *dnp = dn;
932: }
933:
1.13 noro 934: /* nf computation over Z */
935:
1.20 noro 936: void dp_nf_z(NODE b,DP g,DP *ps,int full,int multiple,DP *rp)
1.5 noro 937: {
938: DP u,p,d,s,t,dmy1;
939: P dmy;
940: NODE l;
941: MP m,mr;
942: int i,n;
943: int *wb;
944: int hmag;
945: int sugar,psugar;
946:
947: if ( !g ) {
948: *rp = 0; return;
949: }
950: for ( n = 0, l = b; l; l = NEXT(l), n++ );
951: wb = (int *)ALLOCA(n*sizeof(int));
952: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
953: wb[i] = QTOS((Q)BDY(l));
1.12 noro 954:
1.13 noro 955: hmag = multiple*HMAG(g);
1.5 noro 956: sugar = g->sugar;
1.12 noro 957:
1.5 noro 958: for ( d = 0; g; ) {
959: for ( u = 0, i = 0; i < n; i++ ) {
960: if ( dp_redble(g,p = ps[wb[i]]) ) {
961: dp_red(d,g,p,&t,&u,&dmy,&dmy1);
962: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
963: sugar = MAX(sugar,psugar);
964: if ( !u ) {
965: if ( d )
966: d->sugar = sugar;
967: *rp = d; return;
968: }
969: d = t;
970: break;
971: }
972: }
973: if ( u ) {
974: g = u;
975: if ( d ) {
1.13 noro 976: if ( multiple && HMAG(d) > hmag ) {
1.5 noro 977: dp_ptozp2(d,g,&t,&u); d = t; g = u;
978: hmag = multiple*HMAG(d);
979: }
980: } else {
1.13 noro 981: if ( multiple && HMAG(g) > hmag ) {
1.5 noro 982: dp_ptozp(g,&t); g = t;
983: hmag = multiple*HMAG(g);
984: }
985: }
986: }
987: else if ( !full ) {
988: if ( g ) {
989: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
990: }
991: *rp = g; return;
992: } else {
993: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
994: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
995: addd(CO,d,t,&s); d = s;
996: dp_rest(g,&t); g = t;
997:
998: }
999: }
1000: if ( d )
1001: d->sugar = sugar;
1002: *rp = d;
1003: }
1004:
1.13 noro 1005: /* nf computation over a field */
1006:
1.20 noro 1007: void dp_nf_f(NODE b,DP g,DP *ps,int full,DP *rp)
1.13 noro 1008: {
1009: DP u,p,d,s,t;
1010: NODE l;
1011: MP m,mr;
1012: int i,n;
1013: int *wb;
1014: int sugar,psugar;
1015:
1016: if ( !g ) {
1017: *rp = 0; return;
1018: }
1019: for ( n = 0, l = b; l; l = NEXT(l), n++ );
1020: wb = (int *)ALLOCA(n*sizeof(int));
1021: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
1022: wb[i] = QTOS((Q)BDY(l));
1023:
1024: sugar = g->sugar;
1025: for ( d = 0; g; ) {
1026: for ( u = 0, i = 0; i < n; i++ ) {
1027: if ( dp_redble(g,p = ps[wb[i]]) ) {
1028: dp_red_f(g,p,&u);
1029: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1030: sugar = MAX(sugar,psugar);
1031: if ( !u ) {
1032: if ( d )
1033: d->sugar = sugar;
1034: *rp = d; return;
1035: }
1036: break;
1037: }
1038: }
1039: if ( u )
1040: g = u;
1041: else if ( !full ) {
1042: if ( g ) {
1043: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1044: }
1045: *rp = g; return;
1046: } else {
1047: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1048: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1049: addd(CO,d,t,&s); d = s;
1050: dp_rest(g,&t); g = t;
1051: }
1052: }
1053: if ( d )
1054: d->sugar = sugar;
1055: *rp = d;
1056: }
1057:
1058: /* nf computation over GF(mod) (only for internal use) */
1059:
1.20 noro 1060: void dp_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
1.5 noro 1061: {
1062: DP u,p,d,s,t;
1063: P dmy;
1064: NODE l;
1065: MP m,mr;
1066: int sugar,psugar;
1067:
1068: if ( !g ) {
1069: *rp = 0; return;
1070: }
1071: sugar = g->sugar;
1072: for ( d = 0; g; ) {
1073: for ( u = 0, l = b; l; l = NEXT(l) ) {
1074: if ( dp_redble(g,p = ps[(int)BDY(l)]) ) {
1075: dp_red_mod(d,g,p,mod,&t,&u,&dmy);
1076: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1077: sugar = MAX(sugar,psugar);
1078: if ( !u ) {
1079: if ( d )
1080: d->sugar = sugar;
1081: *rp = d; return;
1082: }
1083: d = t;
1084: break;
1085: }
1086: }
1087: if ( u )
1088: g = u;
1089: else if ( !full ) {
1090: if ( g ) {
1091: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1092: }
1093: *rp = g; return;
1094: } else {
1095: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1096: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1097: addmd(CO,mod,d,t,&s); d = s;
1098: dp_rest(g,&t); g = t;
1099: }
1100: }
1101: if ( d )
1102: d->sugar = sugar;
1103: *rp = d;
1104: }
1105:
1.20 noro 1106: void dp_true_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp,P *dnp)
1.5 noro 1107: {
1108: DP u,p,d,s,t;
1109: NODE l;
1110: MP m,mr;
1111: int i,n;
1112: int *wb;
1113: int sugar,psugar;
1114: P dn,tdn,tdn1;
1115:
1116: dn = (P)ONEM;
1117: if ( !g ) {
1118: *rp = 0; *dnp = dn; return;
1119: }
1120: for ( n = 0, l = b; l; l = NEXT(l), n++ );
1121: wb = (int *)ALLOCA(n*sizeof(int));
1122: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
1123: wb[i] = QTOS((Q)BDY(l));
1124: sugar = g->sugar;
1125: for ( d = 0; g; ) {
1126: for ( u = 0, i = 0; i < n; i++ ) {
1127: if ( dp_redble(g,p = ps[wb[i]]) ) {
1128: dp_red_mod(d,g,p,mod,&t,&u,&tdn);
1129: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1130: sugar = MAX(sugar,psugar);
1131: if ( !u ) {
1132: if ( d )
1133: d->sugar = sugar;
1134: *rp = d; *dnp = dn; return;
1135: } else {
1136: d = t;
1137: mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1;
1138: }
1139: break;
1140: }
1141: }
1142: if ( u )
1143: g = u;
1144: else if ( !full ) {
1145: if ( g ) {
1146: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1147: }
1148: *rp = g; *dnp = dn; return;
1149: } else {
1150: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1151: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1152: addmd(CO,mod,d,t,&s); d = s;
1153: dp_rest(g,&t); g = t;
1154: }
1155: }
1156: if ( d )
1157: d->sugar = sugar;
1158: *rp = d; *dnp = dn;
1159: }
1160:
1.20 noro 1161: void _dp_nf_mod_destructive(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
1.5 noro 1162: {
1.20 noro 1163: DP u,p,d;
1.7 noro 1164: NODE l;
1.20 noro 1165: MP m,mrd;
1166: int sugar,psugar,n,h_reducible;
1.5 noro 1167:
1.7 noro 1168: if ( !g ) {
1169: *rp = 0; return;
1.5 noro 1170: }
1.7 noro 1171: sugar = g->sugar;
1172: n = g->nv;
1173: for ( d = 0; g; ) {
1174: for ( h_reducible = 0, l = b; l; l = NEXT(l) ) {
1175: if ( dp_redble(g,p = ps[(int)BDY(l)]) ) {
1176: h_reducible = 1;
1177: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1178: _dp_red_mod_destructive(g,p,mod,&u); g = u;
1179: sugar = MAX(sugar,psugar);
1180: if ( !g ) {
1181: if ( d )
1182: d->sugar = sugar;
1183: _dptodp(d,rp); _free_dp(d); return;
1184: }
1185: break;
1186: }
1187: }
1188: if ( !h_reducible ) {
1189: /* head term is not reducible */
1190: if ( !full ) {
1191: if ( g )
1192: g->sugar = sugar;
1193: _dptodp(g,rp); _free_dp(g); return;
1194: } else {
1195: m = BDY(g);
1196: if ( NEXT(m) ) {
1197: BDY(g) = NEXT(m); NEXT(m) = 0;
1198: } else {
1199: _FREEDP(g); g = 0;
1200: }
1201: if ( d ) {
1202: for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) );
1203: NEXT(mrd) = m;
1204: } else {
1205: _MKDP(n,m,d);
1206: }
1207: }
1208: }
1.5 noro 1209: }
1.7 noro 1210: if ( d )
1211: d->sugar = sugar;
1212: _dptodp(d,rp); _free_dp(d);
1.5 noro 1213: }
1.13 noro 1214:
1215: /* reduction by linear base over a field */
1216:
1.20 noro 1217: void dp_lnf_f(DP p1,DP p2,NODE g,DP *r1p,DP *r2p)
1.13 noro 1218: {
1219: DP r1,r2,b1,b2,t,s;
1220: Obj c,c1,c2;
1221: NODE l,b;
1222: int n;
1223:
1224: if ( !p1 ) {
1225: *r1p = p1; *r2p = p2; return;
1226: }
1227: n = p1->nv;
1228: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
1229: if ( !r1 ) {
1230: *r1p = r1; *r2p = r2; return;
1231: }
1232: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
1233: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
1234: b2 = (DP)BDY(NEXT(b));
1235: divr(CO,(Obj)ONE,(Obj)BDY(b1)->c,&c1);
1236: mulr(CO,c1,(Obj)BDY(r1)->c,&c2); chsgnr(c2,&c);
1237: muldc(CO,b1,(P)c,&t); addd(CO,r1,t,&s); r1 = s;
1238: muldc(CO,b2,(P)c,&t); addd(CO,r2,t,&s); r2 = s;
1239: }
1240: }
1241: *r1p = r1; *r2p = r2;
1242: }
1243:
1244: /* reduction by linear base over GF(mod) */
1.5 noro 1245:
1.20 noro 1246: void dp_lnf_mod(DP p1,DP p2,NODE g,int mod,DP *r1p,DP *r2p)
1.5 noro 1247: {
1.7 noro 1248: DP r1,r2,b1,b2,t,s;
1249: P c;
1250: MQ c1,c2;
1251: NODE l,b;
1252: int n;
1253:
1254: if ( !p1 ) {
1255: *r1p = p1; *r2p = p2; return;
1256: }
1257: n = p1->nv;
1258: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
1259: if ( !r1 ) {
1260: *r1p = r1; *r2p = r2; return;
1261: }
1262: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
1263: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
1264: b2 = (DP)BDY(NEXT(b));
1265: invmq(mod,(MQ)BDY(b1)->c,&c1);
1266: mulmq(mod,c1,(MQ)BDY(r1)->c,&c2); chsgnmp(mod,(P)c2,&c);
1267: mulmdc(CO,mod,b1,c,&t); addmd(CO,mod,r1,t,&s); r1 = s;
1268: mulmdc(CO,mod,b2,c,&t); addmd(CO,mod,r2,t,&s); r2 = s;
1269: }
1270: }
1271: *r1p = r1; *r2p = r2;
1.5 noro 1272: }
1273:
1.20 noro 1274: void dp_nf_tab_mod(DP p,LIST *tab,int mod,DP *rp)
1.5 noro 1275: {
1.7 noro 1276: DP s,t,u;
1277: MP m;
1278: DL h;
1279: int i,n;
1280:
1281: if ( !p ) {
1282: *rp = p; return;
1283: }
1284: n = p->nv;
1285: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
1286: h = m->dl;
1287: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
1288: i++;
1289: mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
1290: addmd(CO,mod,s,t,&u); s = u;
1.24 noro 1291: }
1292: *rp = s;
1293: }
1294:
1295: void dp_nf_tab_f(DP p,LIST *tab,DP *rp)
1296: {
1297: DP s,t,u;
1298: MP m;
1299: DL h;
1300: int i,n;
1301:
1302: if ( !p ) {
1303: *rp = p; return;
1304: }
1305: n = p->nv;
1306: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
1307: h = m->dl;
1308: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
1309: i++;
1310: muldc(CO,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
1311: addd(CO,s,t,&u); s = u;
1.7 noro 1312: }
1313: *rp = s;
1.5 noro 1314: }
1315:
1.7 noro 1316: /*
1317: * setting flags
1.30 noro 1318: * call create_order_spec with vl=0 to set old type order.
1.7 noro 1319: *
1320: */
1321:
1.27 noro 1322: int create_order_spec(VL vl,Obj obj,struct order_spec **specp)
1.5 noro 1323: {
1.37 noro 1324: int i,j,n,s,row,col,ret;
1.27 noro 1325: struct order_spec *spec;
1.7 noro 1326: struct order_pair *l;
1327: NODE node,t,tn;
1328: MAT m;
1329: pointer **b;
1330: int **w;
1.5 noro 1331:
1.37 noro 1332: if ( vl && obj && OID(obj) == O_LIST ) {
1333: ret = create_composite_order_spec(vl,(LIST)obj,specp);
1334: if ( show_orderspec )
1335: print_composite_order_spec(*specp);
1336: return ret;
1337: }
1.27 noro 1338:
1339: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
1.7 noro 1340: if ( !obj || NUM(obj) ) {
1341: spec->id = 0; spec->obj = obj;
1342: spec->ord.simple = QTOS((Q)obj);
1343: return 1;
1344: } else if ( OID(obj) == O_LIST ) {
1345: node = BDY((LIST)obj);
1346: for ( n = 0, t = node; t; t = NEXT(t), n++ );
1347: l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair));
1348: for ( i = 0, t = node, s = 0; i < n; t = NEXT(t), i++ ) {
1349: tn = BDY((LIST)BDY(t)); l[i].order = QTOS((Q)BDY(tn));
1350: tn = NEXT(tn); l[i].length = QTOS((Q)BDY(tn));
1351: s += l[i].length;
1352: }
1353: spec->id = 1; spec->obj = obj;
1354: spec->ord.block.order_pair = l;
1355: spec->ord.block.length = n; spec->nv = s;
1356: return 1;
1357: } else if ( OID(obj) == O_MAT ) {
1358: m = (MAT)obj; row = m->row; col = m->col; b = BDY(m);
1359: w = almat(row,col);
1360: for ( i = 0; i < row; i++ )
1361: for ( j = 0; j < col; j++ )
1362: w[i][j] = QTOS((Q)b[i][j]);
1363: spec->id = 2; spec->obj = obj;
1364: spec->nv = col; spec->ord.matrix.row = row;
1365: spec->ord.matrix.matrix = w;
1366: return 1;
1367: } else
1.5 noro 1368: return 0;
1369: }
1370:
1.28 noro 1371: void print_composite_order_spec(struct order_spec *spec)
1372: {
1373: int nv,n,len,i,j,k,start;
1374: struct weight_or_block *worb;
1375:
1376: nv = spec->nv;
1377: n = spec->ord.composite.length;
1378: worb = spec->ord.composite.w_or_b;
1379: for ( i = 0; i < n; i++, worb++ ) {
1380: len = worb->length;
1381: printf("[ ");
1382: switch ( worb->type ) {
1383: case IS_DENSE_WEIGHT:
1384: for ( j = 0; j < len; j++ )
1385: printf("%d ",worb->body.dense_weight[j]);
1386: for ( ; j < nv; j++ )
1387: printf("0 ");
1388: break;
1389: case IS_SPARSE_WEIGHT:
1390: for ( j = 0, k = 0; j < nv; j++ )
1391: if ( j == worb->body.sparse_weight[k].pos )
1392: printf("%d ",worb->body.sparse_weight[k++].value);
1393: else
1394: printf("0 ");
1395: break;
1396: case IS_BLOCK:
1397: start = worb->body.block.start;
1398: for ( j = 0; j < start; j++ ) printf("0 ");
1399: switch ( worb->body.block.order ) {
1400: case 0:
1401: for ( k = 0; k < len; k++, j++ ) printf("R ");
1402: break;
1403: case 1:
1404: for ( k = 0; k < len; k++, j++ ) printf("G ");
1405: break;
1406: case 2:
1407: for ( k = 0; k < len; k++, j++ ) printf("L ");
1408: break;
1409: }
1410: for ( ; j < nv; j++ ) printf("0 ");
1411: break;
1412: }
1413: printf("]\n");
1414: }
1.38 ! noro 1415: }
! 1416:
! 1417: struct order_spec *append_block(struct order_spec *spec,
! 1418: int nv,int nalg,int ord)
! 1419: {
! 1420: MAT m,mat;
! 1421: int i,j,row,col,n;
! 1422: Q **b,**wp;
! 1423: int **w;
! 1424: NODE t,s,s0;
! 1425: struct order_pair *l,*l0;
! 1426: int n0,nv0;
! 1427: LIST list0,list1,list;
! 1428: Q oq,nq;
! 1429: struct order_spec *r;
! 1430:
! 1431: r = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 1432: switch ( spec->id ) {
! 1433: case 0:
! 1434: STOQ(spec->ord.simple,oq); STOQ(nv,nq);
! 1435: t = mknode(2,oq,nq); MKLIST(list0,t);
! 1436: STOQ(ord,oq); STOQ(nalg,nq);
! 1437: t = mknode(2,oq,nq); MKLIST(list1,t);
! 1438: t = mknode(2,list0,list1); MKLIST(list,t);
! 1439: l = (struct order_pair *)MALLOC_ATOMIC(2*sizeof(struct order_pair));
! 1440: l[0].order = spec->ord.simple; l[0].length = nv;
! 1441: l[1].order = ord; l[1].length = nalg;
! 1442: r->id = 1; r->obj = (Obj)list;
! 1443: r->ord.block.order_pair = l;
! 1444: r->ord.block.length = 2;
! 1445: r->nv = nv+nalg;
! 1446: break;
! 1447: case 1:
! 1448: if ( spec->nv != nv )
! 1449: error("append_block : number of variables mismatch");
! 1450: l0 = spec->ord.block.order_pair;
! 1451: n0 = spec->ord.block.length;
! 1452: nv0 = spec->nv;
! 1453: list0 = (LIST)spec->obj;
! 1454: n = n0+1;
! 1455: l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair));
! 1456: for ( i = 0; i < n0; i++ )
! 1457: l[i] = l0[i];
! 1458: l[i].order = ord; l[i].length = nalg;
! 1459: for ( t = BDY(list0), s0 = 0; t; t = NEXT(t) ) {
! 1460: NEXTNODE(s0,s); BDY(s) = BDY(t);
! 1461: }
! 1462: STOQ(ord,oq); STOQ(nalg,nq);
! 1463: t = mknode(2,oq,nq); MKLIST(list,t);
! 1464: NEXTNODE(s0,s); BDY(s) = (pointer)list; NEXT(s) = 0;
! 1465: MKLIST(list,s0);
! 1466: r->id = 1; r->obj = (Obj)list;
! 1467: r->ord.block.order_pair = l;
! 1468: r->ord.block.length = n;
! 1469: r->nv = nv+nalg;
! 1470: break;
! 1471: case 2:
! 1472: if ( spec->nv != nv )
! 1473: error("append_block : number of variables mismatch");
! 1474: m = (MAT)spec->obj;
! 1475: row = m->row; col = m->col; b = (Q **)BDY(m);
! 1476: w = almat(row+nalg,col+nalg);
! 1477: MKMAT(mat,row+nalg,col+nalg); wp = (Q **)BDY(mat);
! 1478: for ( i = 0; i < row; i++ )
! 1479: for ( j = 0; j < col; j++ ) {
! 1480: w[i][j] = QTOS(b[i][j]);
! 1481: wp[i][j] = b[i][j];
! 1482: }
! 1483: for ( i = 0; i < nalg; i++ ) {
! 1484: w[i+row][i+col] = 1;
! 1485: wp[i+row][i+col] = ONE;
! 1486: }
! 1487: r->id = 2; r->obj = (Obj)mat;
! 1488: r->nv = col+nalg; r->ord.matrix.row = row+nalg;
! 1489: r->ord.matrix.matrix = w;
! 1490: break;
! 1491: case 3:
! 1492: default:
! 1493: /* XXX */
! 1494: error("append_block : not implemented yet");
! 1495: }
! 1496: return r;
1.28 noro 1497: }
1498:
1.37 noro 1499: int comp_sw(struct sparse_weight *a, struct sparse_weight *b)
1500: {
1501: if ( a->pos > b->pos ) return 1;
1502: else if ( a->pos < b->pos ) return -1;
1503: else return 0;
1504: }
1505:
1.27 noro 1506: /* order = [w_or_b, w_or_b, ... ] */
1507: /* w_or_b = w or b */
1508: /* w = [1,2,...] or [x,1,y,2,...] */
1509: /* b = [@lex,x,y,...,z] etc */
1510:
1511: int create_composite_order_spec(VL vl,LIST order,struct order_spec **specp)
1512: {
1513: NODE wb,t,p;
1514: struct order_spec *spec;
1515: VL tvl;
1.29 noro 1516: int n,i,j,k,l,start,end,len,w;
1.27 noro 1517: int *dw;
1518: struct sparse_weight *sw;
1519: struct weight_or_block *w_or_b;
1520: Obj a0;
1521: NODE a;
1.29 noro 1522: V v,sv,ev;
1523: SYMBOL sym;
1524: int *top;
1.27 noro 1525:
1526: /* l = number of vars in vl */
1527: for ( l = 0, tvl = vl; tvl; tvl = NEXT(tvl), l++ );
1528: /* n = number of primitives in order */
1529: wb = BDY(order);
1530: n = length(wb);
1531: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
1532: spec->id = 3;
1533: spec->obj = (Obj)order;
1534: spec->nv = l;
1535: spec->ord.composite.length = n;
1.28 noro 1536: w_or_b = spec->ord.composite.w_or_b = (struct weight_or_block *)
1.29 noro 1537: MALLOC(sizeof(struct weight_or_block)*(n+1));
1538:
1539: /* top : register the top variable in each w_or_b specification */
1540: top = (int *)ALLOCA(l*sizeof(int));
1541: for ( i = 0; i < l; i++ ) top[i] = 0;
1542:
1.28 noro 1543: for ( t = wb, i = 0; t; t = NEXT(t), i++ ) {
1.30 noro 1544: if ( !BDY(t) || OID((Obj)BDY(t)) != O_LIST )
1545: error("a list of lists must be specified for the key \"order\"");
1.28 noro 1546: a = BDY((LIST)BDY(t));
1.27 noro 1547: len = length(a);
1548: a0 = (Obj)BDY(a);
1549: if ( !a0 || OID(a0) == O_N ) {
1.28 noro 1550: /* a is a dense weight vector */
1.27 noro 1551: dw = (int *)MALLOC(sizeof(int)*len);
1.30 noro 1552: for ( j = 0, p = a; j < len; p = NEXT(p), j++ ) {
1553: if ( !INT((Q)BDY(p)) )
1554: error("a dense weight vector must be specified as a list of integers");
1.27 noro 1555: dw[j] = QTOS((Q)BDY(p));
1.30 noro 1556: }
1.27 noro 1557: w_or_b[i].type = IS_DENSE_WEIGHT;
1558: w_or_b[i].length = len;
1559: w_or_b[i].body.dense_weight = dw;
1.29 noro 1560:
1561: /* find the top */
1562: for ( k = 0; k < len && !dw[k]; k++ );
1563: if ( k < len ) top[k] = 1;
1564:
1.27 noro 1565: } else if ( OID(a0) == O_P ) {
1.28 noro 1566: /* a is a sparse weight vector */
1567: len >>= 1;
1.27 noro 1568: sw = (struct sparse_weight *)
1569: MALLOC(sizeof(struct sparse_weight)*len);
1570: for ( j = 0, p = a; j < len; j++ ) {
1.30 noro 1571: if ( !BDY(p) || OID((P)BDY(p)) != O_P )
1572: error("a sparse weight vector must be specified as [var1,weight1,...]");
1.28 noro 1573: v = VR((P)BDY(p)); p = NEXT(p);
1.27 noro 1574: for ( tvl = vl, k = 0; tvl && tvl->v != v;
1575: k++, tvl = NEXT(tvl) );
1576: if ( !tvl )
1.30 noro 1577: error("invalid variable name in a sparse weight vector");
1.27 noro 1578: sw[j].pos = k;
1.30 noro 1579: if ( !INT((Q)BDY(p)) )
1580: error("a sparse weight vector must be specified as [var1,weight1,...]");
1.28 noro 1581: sw[j].value = QTOS((Q)BDY(p)); p = NEXT(p);
1.27 noro 1582: }
1.37 noro 1583: qsort(sw,len,sizeof(struct sparse_weight),
1584: (int (*)(const void *,const void *))comp_sw);
1.27 noro 1585: w_or_b[i].type = IS_SPARSE_WEIGHT;
1586: w_or_b[i].length = len;
1587: w_or_b[i].body.sparse_weight = sw;
1.29 noro 1588:
1589: /* find the top */
1590: for ( k = 0; k < len && !sw[k].value; k++ );
1591: if ( k < len ) top[sw[k].pos] = 1;
1592: } else if ( OID(a0) == O_RANGE ) {
1593: /* [range(v1,v2),w] */
1594: sv = VR((P)(((RANGE)a0)->start));
1595: ev = VR((P)(((RANGE)a0)->end));
1596: for ( tvl = vl, start = 0; tvl && tvl->v != sv; start++, tvl = NEXT(tvl) );
1597: if ( !tvl )
1598: error("invalid range");
1599: for ( end = start; tvl && tvl->v != ev; end++, tvl = NEXT(tvl) );
1600: if ( !tvl )
1601: error("invalid range");
1602: len = end-start+1;
1603: sw = (struct sparse_weight *)
1604: MALLOC(sizeof(struct sparse_weight)*len);
1605: w = QTOS((Q)BDY(NEXT(a)));
1606: for ( tvl = vl, k = 0; k < start; k++, tvl = NEXT(tvl) );
1607: for ( j = 0 ; k <= end; k++, tvl = NEXT(tvl), j++ ) {
1608: sw[j].pos = k;
1609: sw[j].value = w;
1610: }
1611: w_or_b[i].type = IS_SPARSE_WEIGHT;
1612: w_or_b[i].length = len;
1613: w_or_b[i].body.sparse_weight = sw;
1614:
1615: /* register the top */
1616: if ( w ) top[start] = 1;
1.28 noro 1617: } else if ( OID(a0) == O_SYMBOL ) {
1618: /* a is a block */
1.29 noro 1619: sym = (SYMBOL)a0; a = NEXT(a); len--;
1620: if ( OID((Obj)BDY(a)) == O_RANGE ) {
1621: sv = VR((P)(((RANGE)BDY(a))->start));
1622: ev = VR((P)(((RANGE)BDY(a))->end));
1623: for ( tvl = vl, start = 0; tvl && tvl->v != sv; start++, tvl = NEXT(tvl) );
1624: if ( !tvl )
1625: error("invalid range");
1626: for ( end = start; tvl && tvl->v != ev; end++, tvl = NEXT(tvl) );
1627: if ( !tvl )
1628: error("invalid range");
1629: len = end-start+1;
1630: } else {
1631: for ( start = 0, tvl = vl; tvl->v != VR((P)BDY(a));
1.28 noro 1632: tvl = NEXT(tvl), start++ );
1.29 noro 1633: for ( p = NEXT(a), tvl = NEXT(tvl); p;
1.30 noro 1634: p = NEXT(p), tvl = NEXT(tvl) ) {
1635: if ( !BDY(p) || OID((P)BDY(p)) != O_P )
1636: error("a block must be specified as [ordsymbol,var1,var2,...]");
1.29 noro 1637: if ( tvl->v != VR((P)BDY(p)) ) break;
1.30 noro 1638: }
1.29 noro 1639: if ( p )
1.30 noro 1640: error("a block must be contiguous in the variable list");
1.29 noro 1641: }
1.28 noro 1642: w_or_b[i].type = IS_BLOCK;
1643: w_or_b[i].length = len;
1644: w_or_b[i].body.block.start = start;
1645: if ( !strcmp(sym->name,"@grlex") )
1646: w_or_b[i].body.block.order = 0;
1647: else if ( !strcmp(sym->name,"@glex") )
1648: w_or_b[i].body.block.order = 1;
1649: else if ( !strcmp(sym->name,"@lex") )
1650: w_or_b[i].body.block.order = 2;
1651: else
1.29 noro 1652: error("invalid ordername");
1653: /* register the tops */
1654: for ( j = 0, k = start; j < len; j++, k++ )
1655: top[k] = 1;
1.28 noro 1656: }
1.29 noro 1657: }
1658: for ( k = 0; k < l && top[k]; k++ );
1659: if ( k < l ) {
1660: /* incomplete order specification; add @grlex */
1661: w_or_b[n].type = IS_BLOCK;
1662: w_or_b[n].length = l;
1663: w_or_b[n].body.block.start = 0;
1664: w_or_b[n].body.block.order = 0;
1665: spec->ord.composite.length = n+1;
1.27 noro 1666: }
1667: }
1668:
1.35 noro 1669: /* module order spec */
1670:
1671: void create_modorder_spec(int id,LIST shift,struct modorder_spec **s)
1672: {
1673: struct modorder_spec *spec;
1674: NODE n,t;
1675: LIST list;
1676: int *ds;
1677: int i,l;
1678: Q q;
1679:
1680: *s = spec = (struct modorder_spec *)MALLOC(sizeof(struct modorder_spec));
1681: spec->id = id;
1682: if ( shift ) {
1683: n = BDY(shift);
1684: spec->len = l = length(n);
1685: spec->degree_shift = ds = (int *)MALLOC_ATOMIC(l*sizeof(int));
1686: for ( t = n, i = 0; t; t = NEXT(t), i++ )
1687: ds[i] = QTOS((Q)BDY(t));
1688: } else {
1689: spec->len = 0;
1690: spec->degree_shift = 0;
1691: }
1692: STOQ(id,q);
1693: n = mknode(2,q,shift);
1694: MKLIST(list,n);
1695: spec->obj = (Obj)list;
1696: }
1697:
1.7 noro 1698: /*
1699: * converters
1700: *
1701: */
1702:
1.20 noro 1703: void dp_homo(DP p,DP *rp)
1.5 noro 1704: {
1.7 noro 1705: MP m,mr,mr0;
1706: int i,n,nv,td;
1707: DL dl,dlh;
1.5 noro 1708:
1.7 noro 1709: if ( !p )
1710: *rp = 0;
1711: else {
1712: n = p->nv; nv = n + 1;
1713: m = BDY(p); td = sugard(m);
1714: for ( mr0 = 0; m; m = NEXT(m) ) {
1715: NEXTMP(mr0,mr); mr->c = m->c;
1716: dl = m->dl;
1717: mr->dl = dlh = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
1718: dlh->td = td;
1719: for ( i = 0; i < n; i++ )
1720: dlh->d[i] = dl->d[i];
1721: dlh->d[n] = td - dl->td;
1722: }
1723: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
1.5 noro 1724: }
1725: }
1726:
1.20 noro 1727: void dp_dehomo(DP p,DP *rp)
1.5 noro 1728: {
1.7 noro 1729: MP m,mr,mr0;
1730: int i,n,nv;
1731: DL dl,dlh;
1.5 noro 1732:
1.7 noro 1733: if ( !p )
1734: *rp = 0;
1735: else {
1736: n = p->nv; nv = n - 1;
1737: m = BDY(p);
1738: for ( mr0 = 0; m; m = NEXT(m) ) {
1739: NEXTMP(mr0,mr); mr->c = m->c;
1740: dlh = m->dl;
1741: mr->dl = dl = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
1742: dl->td = dlh->td - dlh->d[nv];
1743: for ( i = 0; i < nv; i++ )
1744: dl->d[i] = dlh->d[i];
1745: }
1746: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
1747: }
1.5 noro 1748: }
1749:
1.20 noro 1750: void dp_mod(DP p,int mod,NODE subst,DP *rp)
1.5 noro 1751: {
1.7 noro 1752: MP m,mr,mr0;
1753: P t,s,s1;
1754: V v;
1755: NODE tn;
1.5 noro 1756:
1.7 noro 1757: if ( !p )
1758: *rp = 0;
1759: else {
1760: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
1761: for ( tn = subst, s = m->c; tn; tn = NEXT(tn) ) {
1762: v = VR((P)BDY(tn)); tn = NEXT(tn);
1763: substp(CO,s,v,(P)BDY(tn),&s1); s = s1;
1764: }
1765: ptomp(mod,s,&t);
1766: if ( t ) {
1767: NEXTMP(mr0,mr); mr->c = t; mr->dl = m->dl;
1768: }
1769: }
1770: if ( mr0 ) {
1771: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1772: } else
1773: *rp = 0;
1774: }
1.5 noro 1775: }
1776:
1.20 noro 1777: void dp_rat(DP p,DP *rp)
1.5 noro 1778: {
1.7 noro 1779: MP m,mr,mr0;
1.5 noro 1780:
1.7 noro 1781: if ( !p )
1782: *rp = 0;
1783: else {
1784: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
1785: NEXTMP(mr0,mr); mptop(m->c,&mr->c); mr->dl = m->dl;
1.5 noro 1786: }
1.7 noro 1787: if ( mr0 ) {
1788: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1789: } else
1790: *rp = 0;
1.5 noro 1791: }
1792: }
1793:
1794:
1.27 noro 1795: void homogenize_order(struct order_spec *old,int n,struct order_spec **newp)
1.5 noro 1796: {
1.7 noro 1797: struct order_pair *l;
1798: int length,nv,row,i,j;
1799: int **newm,**oldm;
1.27 noro 1800: struct order_spec *new;
1.31 noro 1801: int onv,nnv,nlen,olen,owlen;
1802: struct weight_or_block *owb,*nwb;
1.5 noro 1803:
1.27 noro 1804: *newp = new = (struct order_spec *)MALLOC(sizeof(struct order_spec));
1.7 noro 1805: switch ( old->id ) {
1806: case 0:
1807: switch ( old->ord.simple ) {
1808: case 0:
1809: new->id = 0; new->ord.simple = 0; break;
1810: case 1:
1811: l = (struct order_pair *)
1812: MALLOC_ATOMIC(2*sizeof(struct order_pair));
1813: l[0].length = n; l[0].order = 1;
1814: l[1].length = 1; l[1].order = 2;
1815: new->id = 1;
1816: new->ord.block.order_pair = l;
1817: new->ord.block.length = 2; new->nv = n+1;
1818: break;
1819: case 2:
1820: new->id = 0; new->ord.simple = 1; break;
1821: case 3: case 4: case 5:
1822: new->id = 0; new->ord.simple = old->ord.simple+3;
1823: dp_nelim = n-1; break;
1824: case 6: case 7: case 8: case 9:
1825: new->id = 0; new->ord.simple = old->ord.simple; break;
1826: default:
1827: error("homogenize_order : invalid input");
1828: }
1829: break;
1830: case 1:
1831: length = old->ord.block.length;
1832: l = (struct order_pair *)
1833: MALLOC_ATOMIC((length+1)*sizeof(struct order_pair));
1834: bcopy((char *)old->ord.block.order_pair,(char *)l,length*sizeof(struct order_pair));
1835: l[length].order = 2; l[length].length = 1;
1836: new->id = 1; new->nv = n+1;
1837: new->ord.block.order_pair = l;
1838: new->ord.block.length = length+1;
1839: break;
1840: case 2:
1841: nv = old->nv; row = old->ord.matrix.row;
1842: oldm = old->ord.matrix.matrix; newm = almat(row+1,nv+1);
1843: for ( i = 0; i <= nv; i++ )
1844: newm[0][i] = 1;
1845: for ( i = 0; i < row; i++ ) {
1846: for ( j = 0; j < nv; j++ )
1847: newm[i+1][j] = oldm[i][j];
1848: newm[i+1][j] = 0;
1849: }
1850: new->id = 2; new->nv = nv+1;
1851: new->ord.matrix.row = row+1; new->ord.matrix.matrix = newm;
1.31 noro 1852: break;
1853: case 3:
1854: onv = old->nv;
1855: nnv = onv+1;
1856: olen = old->ord.composite.length;
1857: nlen = olen+1;
1858: owb = old->ord.composite.w_or_b;
1859: nwb = (struct weight_or_block *)
1860: MALLOC(nlen*sizeof(struct weight_or_block));
1861: for ( i = 0; i < olen; i++ ) {
1862: nwb[i].type = owb[i].type;
1863: switch ( owb[i].type ) {
1864: case IS_DENSE_WEIGHT:
1865: owlen = owb[i].length;
1866: nwb[i].length = owlen+1;
1867: nwb[i].body.dense_weight = (int *)MALLOC((owlen+1)*sizeof(int));
1868: for ( j = 0; j < owlen; j++ )
1869: nwb[i].body.dense_weight[j] = owb[i].body.dense_weight[j];
1870: nwb[i].body.dense_weight[owlen] = 0;
1871: break;
1872: case IS_SPARSE_WEIGHT:
1873: nwb[i].length = owb[i].length;
1874: nwb[i].body.sparse_weight = owb[i].body.sparse_weight;
1875: break;
1876: case IS_BLOCK:
1877: nwb[i].length = owb[i].length;
1878: nwb[i].body.block = owb[i].body.block;
1879: break;
1880: }
1881: }
1882: nwb[i].type = IS_SPARSE_WEIGHT;
1883: nwb[i].body.sparse_weight =
1884: (struct sparse_weight *)MALLOC(sizeof(struct sparse_weight));
1885: nwb[i].body.sparse_weight[0].pos = onv;
1886: nwb[i].body.sparse_weight[0].value = 1;
1887: new->id = 3;
1888: new->nv = nnv;
1889: new->ord.composite.length = nlen;
1890: new->ord.composite.w_or_b = nwb;
1891: print_composite_order_spec(new);
1.7 noro 1892: break;
1893: default:
1894: error("homogenize_order : invalid input");
1.5 noro 1895: }
1.7 noro 1896: }
1897:
1.20 noro 1898: void qltozl(Q *w,int n,Q *dvr)
1.7 noro 1899: {
1900: N nm,dn;
1901: N g,l1,l2,l3;
1902: Q c,d;
1903: int i;
1904: struct oVECT v;
1.5 noro 1905:
1906: for ( i = 0; i < n; i++ )
1.7 noro 1907: if ( w[i] && !INT(w[i]) )
1908: break;
1909: if ( i == n ) {
1910: v.id = O_VECT; v.len = n; v.body = (pointer *)w;
1911: igcdv(&v,dvr); return;
1912: }
1913: c = w[0]; nm = NM(c); dn = INT(c) ? ONEN : DN(c);
1914: for ( i = 1; i < n; i++ ) {
1915: c = w[i]; l1 = INT(c) ? ONEN : DN(c);
1916: gcdn(nm,NM(c),&g); nm = g;
1917: gcdn(dn,l1,&l2); muln(dn,l1,&l3); divsn(l3,l2,&dn);
1.5 noro 1918: }
1.7 noro 1919: if ( UNIN(dn) )
1920: NTOQ(nm,1,d);
1921: else
1922: NDTOQ(nm,dn,1,d);
1923: *dvr = d;
1924: }
1.5 noro 1925:
1.20 noro 1926: int comp_nm(Q *a,Q *b)
1.7 noro 1927: {
1928: return cmpn((*a)?NM(*a):0,(*b)?NM(*b):0);
1929: }
1930:
1.20 noro 1931: void sortbynm(Q *w,int n)
1.7 noro 1932: {
1933: qsort(w,n,sizeof(Q),(int (*)(const void *,const void *))comp_nm);
1934: }
1.5 noro 1935:
1936:
1.7 noro 1937: /*
1938: * simple operations
1939: *
1940: */
1.5 noro 1941:
1.20 noro 1942: int dp_redble(DP p1,DP p2)
1.7 noro 1943: {
1944: int i,n;
1945: DL d1,d2;
1.5 noro 1946:
1.7 noro 1947: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
1948: if ( d1->td < d2->td )
1949: return 0;
1950: else {
1951: for ( i = 0, n = p1->nv; i < n; i++ )
1952: if ( d1->d[i] < d2->d[i] )
1953: return 0;
1954: return 1;
1.5 noro 1955: }
1956: }
1957:
1.20 noro 1958: void dp_subd(DP p1,DP p2,DP *rp)
1.5 noro 1959: {
1.7 noro 1960: int i,n;
1.5 noro 1961: DL d1,d2,d;
1962: MP m;
1.7 noro 1963: DP s;
1.5 noro 1964:
1965: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
1.7 noro 1966: NEWDL(d,n); d->td = d1->td - d2->td;
1.5 noro 1967: for ( i = 0; i < n; i++ )
1.7 noro 1968: d->d[i] = d1->d[i]-d2->d[i];
1969: NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0;
1970: MKDP(n,m,s); s->sugar = d->td;
1971: *rp = s;
1972: }
1973:
1.20 noro 1974: void dltod(DL d,int n,DP *rp)
1.7 noro 1975: {
1976: MP m;
1977: DP s;
1978:
1979: NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0;
1980: MKDP(n,m,s); s->sugar = d->td;
1981: *rp = s;
1.5 noro 1982: }
1983:
1.20 noro 1984: void dp_hm(DP p,DP *rp)
1.5 noro 1985: {
1986: MP m,mr;
1987:
1988: if ( !p )
1989: *rp = 0;
1990: else {
1991: m = BDY(p);
1992: NEWMP(mr); mr->dl = m->dl; mr->c = m->c; NEXT(mr) = 0;
1993: MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
1994: }
1995: }
1996:
1.35 noro 1997: void dp_ht(DP p,DP *rp)
1998: {
1999: MP m,mr;
2000:
2001: if ( !p )
2002: *rp = 0;
2003: else {
2004: m = BDY(p);
2005: NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0;
2006: MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
2007: }
2008: }
2009:
1.20 noro 2010: void dp_rest(DP p,DP *rp)
1.5 noro 2011: {
2012: MP m;
2013:
2014: m = BDY(p);
2015: if ( !NEXT(m) )
2016: *rp = 0;
2017: else {
2018: MKDP(p->nv,NEXT(m),*rp);
2019: if ( *rp )
2020: (*rp)->sugar = p->sugar;
2021: }
2022: }
2023:
1.20 noro 2024: DL lcm_of_DL(int nv,DL dl1,DL dl2,DL dl)
1.5 noro 2025: {
1.21 noro 2026: register int i, *d1, *d2, *d, td;
1.5 noro 2027:
2028: if ( !dl ) NEWDL(dl,nv);
2029: d = dl->d, d1 = dl1->d, d2 = dl2->d;
1.21 noro 2030: for ( td = 0, i = 0; i < nv; d1++, d2++, d++, i++ ) {
2031: *d = *d1 > *d2 ? *d1 : *d2;
2032: td += MUL_WEIGHT(*d,i);
2033: }
1.5 noro 2034: dl->td = td;
2035: return dl;
2036: }
2037:
1.20 noro 2038: int dl_equal(int nv,DL dl1,DL dl2)
1.5 noro 2039: {
2040: register int *d1, *d2, n;
2041:
2042: if ( dl1->td != dl2->td ) return 0;
2043: for ( d1 = dl1->d, d2 = dl2->d, n = nv; --n >= 0; d1++, d2++ )
2044: if ( *d1 != *d2 ) return 0;
2045: return 1;
2046: }
2047:
1.20 noro 2048: int dp_nt(DP p)
1.5 noro 2049: {
2050: int i;
2051: MP m;
2052:
2053: if ( !p )
2054: return 0;
2055: else {
2056: for ( i = 0, m = BDY(p); m; m = NEXT(m), i++ );
2057: return i;
2058: }
2059: }
2060:
1.20 noro 2061: int dp_homogeneous(DP p)
1.15 noro 2062: {
2063: MP m;
2064: int d;
2065:
2066: if ( !p )
2067: return 1;
2068: else {
2069: m = BDY(p);
2070: d = m->dl->td;
2071: m = NEXT(m);
2072: for ( ; m; m = NEXT(m) ) {
2073: if ( m->dl->td != d )
2074: return 0;
2075: }
2076: return 1;
2077: }
1.16 noro 2078: }
2079:
1.20 noro 2080: void _print_mp(int nv,MP m)
1.16 noro 2081: {
2082: int i;
2083:
1.17 noro 2084: if ( !m )
1.16 noro 2085: return;
2086: for ( ; m; m = NEXT(m) ) {
2087: fprintf(stderr,"%d<",ITOS(C(m)));
2088: for ( i = 0; i < nv; i++ ) {
2089: fprintf(stderr,"%d",m->dl->d[i]);
2090: if ( i != nv-1 )
2091: fprintf(stderr," ");
2092: }
2093: fprintf(stderr,">",C(m));
2094: }
2095: fprintf(stderr,"\n");
1.15 noro 2096: }
1.26 noro 2097:
2098: static int cmp_mp_nvar;
2099:
2100: int comp_mp(MP *a,MP *b)
2101: {
2102: return -(*cmpdl)(cmp_mp_nvar,(*a)->dl,(*b)->dl);
2103: }
2104:
2105: void dp_sort(DP p,DP *rp)
2106: {
2107: MP t,mp,mp0;
2108: int i,n;
2109: DP r;
2110: MP *w;
2111:
2112: if ( !p ) {
2113: *rp = 0;
2114: return;
2115: }
2116: for ( t = BDY(p), n = 0; t; t = NEXT(t), n++ );
2117: w = (MP *)ALLOCA(n*sizeof(MP));
2118: for ( t = BDY(p), i = 0; i < n; t = NEXT(t), i++ )
2119: w[i] = t;
2120: cmp_mp_nvar = NV(p);
2121: qsort(w,n,sizeof(MP),(int (*)(const void *,const void *))comp_mp);
2122: mp0 = 0;
2123: for ( i = n-1; i >= 0; i-- ) {
2124: NEWMP(mp); mp->dl = w[i]->dl; C(mp) = C(w[i]);
2125: NEXT(mp) = mp0; mp0 = mp;
2126: }
2127: MKDP(p->nv,mp0,r);
2128: r->sugar = p->sugar;
2129: *rp = r;
2130: }
2131:
1.32 noro 2132: DP extract_initial_term_from_dp(DP p,int *weight,int n);
2133: LIST extract_initial_term(LIST f,int *weight,int n);
2134:
2135: DP extract_initial_term_from_dp(DP p,int *weight,int n)
2136: {
1.34 noro 2137: int w,t,i,top;
1.32 noro 2138: MP m,r0,r;
2139: DP dp;
2140:
2141: if ( !p ) return 0;
1.34 noro 2142: top = 1;
1.32 noro 2143: for ( m = BDY(p); m; m = NEXT(m) ) {
2144: for ( i = 0, t = 0; i < n; i++ )
2145: t += weight[i]*m->dl->d[i];
1.34 noro 2146: if ( top || t > w ) {
1.32 noro 2147: r0 = 0;
2148: w = t;
1.34 noro 2149: top = 0;
1.32 noro 2150: }
2151: if ( t == w ) {
2152: NEXTMP(r0,r);
2153: r->dl = m->dl;
2154: r->c = m->c;
2155: }
2156: }
2157: NEXT(r) = 0;
2158: MKDP(p->nv,r0,dp);
2159: return dp;
2160: }
2161:
2162: LIST extract_initial_term(LIST f,int *weight,int n)
2163: {
2164: NODE nd,r0,r;
2165: Obj p;
2166: LIST l;
2167:
2168: nd = BDY(f);
2169: for ( r0 = 0; nd; nd = NEXT(nd) ) {
2170: NEXTNODE(r0,r);
2171: p = (Obj)BDY(nd);
2172: BDY(r) = (pointer)extract_initial_term_from_dp((DP)p,weight,n);
2173: }
2174: if ( r0 ) NEXT(r) = 0;
2175: MKLIST(l,r0);
2176: return l;
2177: }
2178:
2179: LIST dp_initial_term(LIST f,struct order_spec *ord)
2180: {
2181: int n,l,i;
2182: struct weight_or_block *worb;
2183: int *weight;
2184:
2185: switch ( ord->id ) {
2186: case 2: /* matrix order */
2187: /* extract the first row */
2188: n = ord->nv;
2189: weight = ord->ord.matrix.matrix[0];
2190: return extract_initial_term(f,weight,n);
2191: case 3: /* composite order */
2192: /* the first w_or_b */
2193: worb = ord->ord.composite.w_or_b;
2194: switch ( worb->type ) {
2195: case IS_DENSE_WEIGHT:
2196: n = worb->length;
2197: weight = worb->body.dense_weight;
2198: return extract_initial_term(f,weight,n);
2199: case IS_SPARSE_WEIGHT:
2200: n = ord->nv;
2201: weight = (int *)ALLOCA(n*sizeof(int));
1.33 noro 2202: for ( i = 0; i < n; i++ ) weight[i] = 0;
1.32 noro 2203: l = worb->length;
2204: for ( i = 0; i < l; i++ )
2205: weight[worb->body.sparse_weight[i].pos]
2206: = worb->body.sparse_weight[i].value;
2207: return extract_initial_term(f,weight,n);
2208: default:
2209: error("dp_initial_term : unsupported order");
2210: }
2211: default:
2212: error("dp_initial_term : unsupported order");
2213: }
2214: }
2215:
2216: int highest_order_dp(DP p,int *weight,int n);
2217: LIST highest_order(LIST f,int *weight,int n);
2218:
2219: int highest_order_dp(DP p,int *weight,int n)
2220: {
1.34 noro 2221: int w,t,i,top;
1.32 noro 2222: MP m;
2223:
2224: if ( !p ) return -1;
1.34 noro 2225: top = 1;
1.32 noro 2226: for ( m = BDY(p); m; m = NEXT(m) ) {
2227: for ( i = 0, t = 0; i < n; i++ )
2228: t += weight[i]*m->dl->d[i];
1.34 noro 2229: if ( top || t > w ) {
1.32 noro 2230: w = t;
1.34 noro 2231: top = 0;
2232: }
1.32 noro 2233: }
2234: return w;
2235: }
2236:
2237: LIST highest_order(LIST f,int *weight,int n)
2238: {
2239: int h;
2240: NODE nd,r0,r;
2241: Obj p;
2242: LIST l;
2243: Q q;
2244:
2245: nd = BDY(f);
2246: for ( r0 = 0; nd; nd = NEXT(nd) ) {
2247: NEXTNODE(r0,r);
2248: p = (Obj)BDY(nd);
2249: h = highest_order_dp((DP)p,weight,n);
2250: STOQ(h,q);
2251: BDY(r) = (pointer)q;
2252: }
2253: if ( r0 ) NEXT(r) = 0;
2254: MKLIST(l,r0);
2255: return l;
2256: }
2257:
2258: LIST dp_order(LIST f,struct order_spec *ord)
2259: {
2260: int n,l,i;
2261: struct weight_or_block *worb;
2262: int *weight;
2263:
2264: switch ( ord->id ) {
2265: case 2: /* matrix order */
2266: /* extract the first row */
2267: n = ord->nv;
2268: weight = ord->ord.matrix.matrix[0];
2269: return highest_order(f,weight,n);
2270: case 3: /* composite order */
2271: /* the first w_or_b */
2272: worb = ord->ord.composite.w_or_b;
2273: switch ( worb->type ) {
2274: case IS_DENSE_WEIGHT:
2275: n = worb->length;
2276: weight = worb->body.dense_weight;
2277: return highest_order(f,weight,n);
2278: case IS_SPARSE_WEIGHT:
2279: n = ord->nv;
2280: weight = (int *)ALLOCA(n*sizeof(int));
1.33 noro 2281: for ( i = 0; i < n; i++ ) weight[i] = 0;
1.32 noro 2282: l = worb->length;
2283: for ( i = 0; i < l; i++ )
2284: weight[worb->body.sparse_weight[i].pos]
2285: = worb->body.sparse_weight[i].value;
2286: return highest_order(f,weight,n);
2287: default:
2288: error("dp_initial_term : unsupported order");
2289: }
2290: default:
2291: error("dp_initial_term : unsupported order");
1.35 noro 2292: }
2293: }
2294:
2295: int dpv_ht(DPV p,DP *h)
2296: {
2297: int len,max,maxi,i,t;
2298: DP *e;
2299: MP m,mr;
2300:
2301: len = p->len;
2302: e = p->body;
2303: max = -1;
2304: maxi = -1;
2305: for ( i = 0; i < len; i++ )
2306: if ( e[i] && (t = BDY(e[i])->dl->td) > max ) {
2307: max = t;
2308: maxi = i;
2309: }
2310: if ( max < 0 ) {
2311: *h = 0;
2312: return -1;
2313: } else {
2314: m = BDY(e[maxi]);
2315: NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0;
2316: MKDP(e[maxi]->nv,mr,*h); (*h)->sugar = mr->dl->td; /* XXX */
2317: return maxi;
1.32 noro 2318: }
2319: }