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