![[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.27 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/dp-supp.c,v 1.26 2003/11/27 07:53:53 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:
65: /*
66: * content reduction
67: *
68: */
69:
1.20 noro 70: void dp_ptozp(DP p,DP *rp)
1.7 noro 71: {
72: MP m,mr,mr0;
73: int i,n;
74: Q *w;
75: Q dvr;
76: P t;
77:
78: if ( !p )
79: *rp = 0;
80: else {
81: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
82: w = (Q *)ALLOCA(n*sizeof(Q));
83: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
84: if ( NUM(m->c) )
85: w[i] = (Q)m->c;
86: else
87: ptozp(m->c,1,&w[i],&t);
88: sortbynm(w,n);
89: qltozl(w,n,&dvr);
90: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
91: NEXTMP(mr0,mr); divsp(CO,m->c,(P)dvr,&mr->c); mr->dl = m->dl;
92: }
93: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
94: }
95: }
96:
1.20 noro 97: void dp_ptozp2(DP p0,DP p1,DP *hp,DP *rp)
1.7 noro 98: {
99: DP t,s,h,r;
100: MP m,mr,mr0,m0;
101:
102: addd(CO,p0,p1,&t); dp_ptozp(t,&s);
103: if ( !p0 ) {
104: h = 0; r = s;
105: } else if ( !p1 ) {
106: h = s; r = 0;
107: } else {
108: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
109: m = NEXT(m), m0 = NEXT(m0) ) {
110: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
111: }
112: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
113: }
114: if ( h )
115: h->sugar = p0->sugar;
116: if ( r )
117: r->sugar = p1->sugar;
118: *hp = h; *rp = r;
119: }
1.1 noro 120:
1.20 noro 121: void dp_idiv(DP p,Q c,DP *rp)
1.1 noro 122: {
123: Q t;
124: N nm,q;
125: int sgn,s;
126: MP mr0,m,mr;
127:
128: if ( !p )
129: *rp = 0;
130: else if ( MUNIQ((Q)c) )
131: *rp = p;
132: else if ( MUNIQ((Q)c) )
133: chsgnd(p,rp);
134: else {
135: nm = NM(c); sgn = SGN(c);
136: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
137: NEXTMP(mr0,mr);
138:
139: divsn(NM((Q)(m->c)),nm,&q);
140: s = sgn*SGN((Q)(m->c));
141: NTOQ(q,s,t);
142: mr->c = (P)t;
143: mr->dl = m->dl;
144: }
145: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
146: if ( *rp )
147: (*rp)->sugar = p->sugar;
148: }
149: }
150:
1.20 noro 151: void dp_mbase(NODE hlist,NODE *mbase)
1.1 noro 152: {
153: DL *dl;
154: DL d;
155: int i,j,n,nvar,td;
156:
157: n = length(hlist); nvar = ((DP)BDY(hlist))->nv;
158: dl = (DL *)MALLOC(n*sizeof(DL));
159: for ( i = 0; i < n; i++, hlist = NEXT(hlist) )
160: dl[i] = BDY((DP)BDY(hlist))->dl;
161: NEWDL(d,nvar); *mbase = 0;
162: while ( 1 ) {
163: insert_to_node(d,mbase,nvar);
164: for ( i = nvar-1; i >= 0; ) {
1.21 noro 165: d->d[i]++;
166: d->td += MUL_WEIGHT(1,i);
1.1 noro 167: for ( j = 0; j < n; j++ ) {
168: if ( _dl_redble(dl[j],d,nvar) )
169: break;
170: }
171: if ( j < n ) {
172: for ( j = nvar-1; j >= i; j-- )
173: d->d[j] = 0;
174: for ( j = 0, td = 0; j < i; j++ )
1.21 noro 175: td += MUL_WEIGHT(d->d[j],j);
1.1 noro 176: d->td = td;
177: i--;
178: } else
179: break;
180: }
181: if ( i < 0 )
182: break;
183: }
184: }
185:
1.20 noro 186: int _dl_redble(DL d1,DL d2,int nvar)
1.1 noro 187: {
188: int i;
189:
190: if ( d1->td > d2->td )
191: return 0;
192: for ( i = 0; i < nvar; i++ )
193: if ( d1->d[i] > d2->d[i] )
194: break;
195: if ( i < nvar )
196: return 0;
197: else
198: return 1;
199: }
200:
1.20 noro 201: void insert_to_node(DL d,NODE *n,int nvar)
1.1 noro 202: {
203: DL d1;
204: MP m;
205: DP dp;
206: NODE n0,n1,n2;
207:
208: NEWDL(d1,nvar); d1->td = d->td;
209: bcopy((char *)d->d,(char *)d1->d,nvar*sizeof(int));
210: NEWMP(m); m->dl = d1; m->c = (P)ONE; NEXT(m) = 0;
211: MKDP(nvar,m,dp); dp->sugar = d->td;
212: if ( !(*n) ) {
213: MKNODE(n1,dp,0); *n = n1;
214: } else {
215: for ( n1 = *n, n0 = 0; n1; n0 = n1, n1 = NEXT(n1) )
216: if ( (*cmpdl)(nvar,d,BDY((DP)BDY(n1))->dl) > 0 ) {
217: MKNODE(n2,dp,n1);
218: if ( !n0 )
219: *n = n2;
220: else
221: NEXT(n0) = n2;
222: break;
223: }
224: if ( !n1 ) {
225: MKNODE(n2,dp,0); NEXT(n0) = n2;
226: }
227: }
228: }
229:
1.20 noro 230: void dp_vtod(Q *c,DP p,DP *rp)
1.1 noro 231: {
232: MP mr0,m,mr;
233: int i;
234:
235: if ( !p )
236: *rp = 0;
237: else {
238: for ( mr0 = 0, m = BDY(p), i = 0; m; m = NEXT(m), i++ ) {
239: NEXTMP(mr0,mr); mr->c = (P)c[i]; mr->dl = m->dl;
240: }
241: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
242: (*rp)->sugar = p->sugar;
243: }
244: }
245:
1.8 noro 246: extern int mpi_mag;
247: extern int PCoeffs;
248:
1.20 noro 249: void dp_ptozp_d(DP p,DP *rp)
1.1 noro 250: {
251: int i,j,k,l,n,nsep;
252: MP m;
253: NODE tn,n0,n1,n2,n3;
254: struct oVECT v;
255: VECT c,cs;
256: VECT qi,ri;
257: LIST *qr;
258: Obj dmy;
259: Q d0,d1,gcd,a,u,u1;
260: Q *q,*r;
261: STRING iqr_v;
262: pointer *b;
263: N qn,gn;
264: double get_rtime();
265: int blen;
1.8 noro 266: NODE dist;
267: int ndist;
1.1 noro 268: double t0;
269: double t_e,t_d,t_d1,t_c;
1.8 noro 270: extern int DP_NFStat;
271: extern LIST Dist;
1.20 noro 272: void Pox_rpc();
273: void Pox_pop_local();
1.1 noro 274:
275: if ( !p )
276: *rp = 0;
277: else {
1.8 noro 278: if ( PCoeffs ) {
279: dp_ptozp(p,rp); return;
280: }
1.9 noro 281: if ( !Dist || p_mag(BDY(p)->c) <= mpi_mag ) {
1.8 noro 282: dist = 0; ndist = 0;
283: if ( DP_NFStat ) fprintf(asir_out,"L");
284: } else {
285: dist = BDY(Dist); ndist = length(dist);
286: if ( DP_NFStat ) fprintf(asir_out,"D");
287: }
1.1 noro 288: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
289: nsep = ndist + 1;
290: if ( n <= nsep ) {
291: dp_ptozp(p,rp); return;
292: }
293: t0 = get_rtime();
294: dp_dtov(p,&c);
295: igcdv_estimate(c,&d0);
296: t_e = get_rtime()-t0;
297: t0 = get_rtime();
298: dp_dtov(p,&c);
299: sepvect(c,nsep,&cs);
300: MKSTR(iqr_v,"iqr");
301: qr = (LIST *)CALLOC(nsep,sizeof(LIST));
302: q = (Q *)CALLOC(n,sizeof(Q));
303: r = (Q *)CALLOC(n,sizeof(Q));
304: for ( i = 0, tn = dist, b = BDY(cs); i < ndist; i++, tn = NEXT(tn) ) {
305: MKNODE(n3,d0,0); MKNODE(n2,b[i],n3);
306: MKNODE(n1,iqr_v,n2); MKNODE(n0,BDY(tn),n1);
307: Pox_rpc(n0,&dmy);
308: }
309: iqrv(b[i],d0,&qr[i]);
310: dp_dtov(p,&c);
311: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
312: Pox_pop_local(tn,&qr[i]);
313: if ( OID(qr[i]) == O_ERR ) {
314: printexpr(CO,(Obj)qr[i]);
315: error("dp_ptozp_d : aborted");
316: }
317: }
318: t_d = get_rtime()-t0;
319: t_d1 = t_d/n;
320: t0 = get_rtime();
321: for ( i = j = 0; i < nsep; i++ ) {
322: tn = BDY(qr[i]); qi = (VECT)BDY(tn); ri = (VECT)BDY(NEXT(tn));
323: for ( k = 0, l = qi->len; k < l; k++, j++ ) {
324: q[j] = (Q)BDY(qi)[k]; r[j] = (Q)BDY(ri)[k];
325: }
326: }
327: v.id = O_VECT; v.len = n; v.body = (pointer *)r; igcdv(&v,&d1);
328: if ( d1 ) {
329: gcdn(NM(d0),NM(d1),&gn); NTOQ(gn,1,gcd);
330: divsn(NM(d0),gn,&qn); NTOQ(qn,1,a);
331: for ( i = 0; i < n; i++ ) {
332: mulq(a,q[i],&u);
333: if ( r[i] ) {
334: divsn(NM(r[i]),gn,&qn); NTOQ(qn,SGN(r[i]),u1);
335: addq(u,u1,&q[i]);
336: } else
337: q[i] = u;
338: }
339: } else
340: gcd = d0;
341: dp_vtod(q,p,rp);
342: t_c = get_rtime()-t0;
343: blen=p_mag((P)gcd);
344: pz_t_e += t_e; pz_t_d += t_d; pz_t_d1 += t_d1; pz_t_c += t_c;
345: if ( 0 )
346: fprintf(stderr,"(%d,%d)",p_mag((P)d0)-blen,blen);
347: }
348: }
349:
1.20 noro 350: void dp_ptozp2_d(DP p0,DP p1,DP *hp,DP *rp)
1.1 noro 351: {
352: DP t,s,h,r;
353: MP m,mr,mr0,m0;
354:
1.8 noro 355: addd(CO,p0,p1,&t); dp_ptozp_d(t,&s);
1.1 noro 356: if ( !p0 ) {
357: h = 0; r = s;
358: } else if ( !p1 ) {
359: h = s; r = 0;
360: } else {
361: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
362: m = NEXT(m), m0 = NEXT(m0) ) {
363: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
364: }
365: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
366: }
367: if ( h )
368: h->sugar = p0->sugar;
369: if ( r )
370: r->sugar = p1->sugar;
371: *hp = h; *rp = r;
1.5 noro 372: }
373:
1.22 noro 374: int have_sf_coef(P p)
375: {
376: DCP dc;
377:
378: if ( !p )
379: return 0;
380: else if ( NUM(p) )
381: return NID((Num)p) == N_GFS ? 1 : 0;
382: else {
383: for ( dc = DC(p); dc; dc = NEXT(dc) )
384: if ( have_sf_coef(COEF(dc)) )
385: return 1;
386: return 0;
387: }
388: }
389:
1.25 noro 390: void head_coef(P p,Num *c)
391: {
392: if ( !p )
393: *c = 0;
394: else if ( NUM(p) )
395: *c = (Num)p;
396: else
397: head_coef(COEF(DC(p)),c);
398: }
399:
400: void dp_monic_sf(DP p,DP *rp)
401: {
402: Num c;
403:
404: if ( !p )
405: *rp = 0;
406: else {
407: head_coef(BDY(p)->c,&c);
408: divsdc(CO,p,(P)c,rp);
409: }
410: }
411:
1.20 noro 412: void dp_prim(DP p,DP *rp)
1.5 noro 413: {
1.7 noro 414: P t,g;
415: DP p1;
416: MP m,mr,mr0;
417: int i,n;
418: P *w;
419: Q *c;
420: Q dvr;
1.5 noro 421:
1.7 noro 422: if ( !p )
423: *rp = 0;
1.23 noro 424: else if ( dp_fcoeffs == N_GFS ) {
425: for ( m = BDY(p); m; m = NEXT(m) )
1.22 noro 426: if ( OID(m->c) == O_N ) {
427: /* GCD of coeffs = 1 */
1.25 noro 428: dp_monic_sf(p,rp);
1.22 noro 429: return;
1.23 noro 430: } else break;
431: /* compute GCD over the finite fieid */
432: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
433: w = (P *)ALLOCA(n*sizeof(P));
434: for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ )
435: w[i] = m->c;
436: gcdsf(CO,w,n,&g);
437: if ( NUM(g) )
1.25 noro 438: dp_monic_sf(p,rp);
1.23 noro 439: else {
440: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
441: NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl;
1.22 noro 442: }
1.25 noro 443: NEXT(mr) = 0; MKDP(p->nv,mr0,p1); p1->sugar = p->sugar;
444: dp_monic_sf(p1,rp);
1.22 noro 445: }
1.23 noro 446: return;
447: } else if ( dp_fcoeffs )
1.7 noro 448: *rp = p;
1.23 noro 449: else if ( NoGCD )
1.7 noro 450: dp_ptozp(p,rp);
451: else {
452: dp_ptozp(p,&p1); p = p1;
453: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
454: if ( n == 1 ) {
455: m = BDY(p);
456: NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0;
457: MKDP(p->nv,mr,*rp); (*rp)->sugar = p->sugar;
458: return;
459: }
460: w = (P *)ALLOCA(n*sizeof(P));
461: c = (Q *)ALLOCA(n*sizeof(Q));
462: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
463: if ( NUM(m->c) ) {
464: c[i] = (Q)m->c; w[i] = (P)ONE;
465: } else
466: ptozp(m->c,1,&c[i],&w[i]);
467: qltozl(c,n,&dvr); heu_nezgcdnpz(CO,w,n,&t); mulp(CO,t,(P)dvr,&g);
468: if ( NUM(g) )
469: *rp = p;
470: else {
471: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
472: NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl;
473: }
474: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1.5 noro 475: }
1.7 noro 476: }
1.5 noro 477: }
478:
1.20 noro 479: void heu_nezgcdnpz(VL vl,P *pl,int m,P *pr)
1.5 noro 480: {
481: int i,r;
482: P gcd,t,s1,s2,u;
483: Q rq;
484:
485: while ( 1 ) {
486: for ( i = 0, s1 = 0; i < m; i++ ) {
487: r = random(); UTOQ(r,rq);
488: mulp(vl,pl[i],(P)rq,&t); addp(vl,s1,t,&u); s1 = u;
489: }
490: for ( i = 0, s2 = 0; i < m; i++ ) {
491: r = random(); UTOQ(r,rq);
492: mulp(vl,pl[i],(P)rq,&t); addp(vl,s2,t,&u); s2 = u;
493: }
494: ezgcdp(vl,s1,s2,&gcd);
495: for ( i = 0; i < m; i++ ) {
496: if ( !divtpz(vl,pl[i],gcd,&t) )
497: break;
498: }
499: if ( i == m )
500: break;
501: }
502: *pr = gcd;
503: }
504:
1.20 noro 505: void dp_prim_mod(DP p,int mod,DP *rp)
1.5 noro 506: {
507: P t,g;
508: MP m,mr,mr0;
509:
510: if ( !p )
511: *rp = 0;
512: else if ( NoGCD )
513: *rp = p;
514: else {
515: for ( m = BDY(p), g = m->c, m = NEXT(m); m; m = NEXT(m) ) {
516: gcdprsmp(CO,mod,g,m->c,&t); g = t;
517: }
518: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
519: NEXTMP(mr0,mr); divsmp(CO,mod,m->c,g,&mr->c); mr->dl = m->dl;
520: }
521: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
522: }
523: }
524:
1.20 noro 525: void dp_cont(DP p,Q *rp)
1.5 noro 526: {
1.7 noro 527: VECT v;
1.5 noro 528:
1.7 noro 529: dp_dtov(p,&v); igcdv(v,rp);
1.5 noro 530: }
531:
1.20 noro 532: void dp_dtov(DP dp,VECT *rp)
1.5 noro 533: {
1.7 noro 534: MP m,t;
535: int i,n;
536: VECT v;
537: pointer *p;
1.5 noro 538:
1.7 noro 539: m = BDY(dp);
540: for ( t = m, n = 0; t; t = NEXT(t), n++ );
541: MKVECT(v,n);
542: for ( i = 0, p = BDY(v), t = m; i < n; t = NEXT(t), i++ )
543: p[i] = (pointer)(t->c);
544: *rp = v;
1.5 noro 545: }
546:
1.7 noro 547: /*
548: * s-poly computation
549: *
550: */
1.5 noro 551:
1.20 noro 552: void dp_sp(DP p1,DP p2,DP *rp)
1.5 noro 553: {
1.7 noro 554: int i,n,td;
555: int *w;
556: DL d1,d2,d;
557: MP m;
558: DP t,s1,s2,u;
559: Q c,c1,c2;
560: N gn,tn;
1.5 noro 561:
1.7 noro 562: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
563: w = (int *)ALLOCA(n*sizeof(int));
564: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 565: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.5 noro 566: }
1.7 noro 567:
568: NEWDL(d,n); d->td = td - d1->td;
569: for ( i = 0; i < n; i++ )
570: d->d[i] = w[i] - d1->d[i];
571: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
572: if ( INT(c1) && INT(c2) ) {
573: gcdn(NM(c1),NM(c2),&gn);
574: if ( !UNIN(gn) ) {
575: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
576: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
1.5 noro 577: }
578: }
1.7 noro 579:
580: NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0;
581: MKDP(n,m,s1); s1->sugar = d->td; muld(CO,s1,p1,&t);
582:
583: NEWDL(d,n); d->td = td - d2->td;
584: for ( i = 0; i < n; i++ )
585: d->d[i] = w[i] - d2->d[i];
586: NEWMP(m); m->dl = d; m->c = (P)c1; NEXT(m) = 0;
587: MKDP(n,m,s2); s2->sugar = d->td; muld(CO,s2,p2,&u);
588:
589: subd(CO,t,u,rp);
1.14 noro 590: if ( GenTrace ) {
591: LIST hist;
592: NODE node;
593:
594: node = mknode(4,ONE,0,s1,ONE);
595: MKLIST(hist,node);
596: MKNODE(TraceList,hist,0);
597:
598: node = mknode(4,ONE,0,0,ONE);
599: chsgnd(s2,(DP *)&ARG2(node));
600: MKLIST(hist,node);
601: MKNODE(node,hist,TraceList); TraceList = node;
602: }
603: }
604:
1.20 noro 605: void _dp_sp_dup(DP p1,DP p2,DP *rp)
1.14 noro 606: {
607: int i,n,td;
608: int *w;
609: DL d1,d2,d;
610: MP m;
611: DP t,s1,s2,u;
612: Q c,c1,c2;
613: N gn,tn;
614:
615: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
616: w = (int *)ALLOCA(n*sizeof(int));
617: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 618: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.14 noro 619: }
620:
621: _NEWDL(d,n); d->td = td - d1->td;
622: for ( i = 0; i < n; i++ )
623: d->d[i] = w[i] - d1->d[i];
624: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
625: if ( INT(c1) && INT(c2) ) {
626: gcdn(NM(c1),NM(c2),&gn);
627: if ( !UNIN(gn) ) {
628: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
629: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
630: }
631: }
632:
633: _NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0;
634: _MKDP(n,m,s1); s1->sugar = d->td; _muld_dup(CO,s1,p1,&t); _free_dp(s1);
635:
636: _NEWDL(d,n); d->td = td - d2->td;
637: for ( i = 0; i < n; i++ )
638: d->d[i] = w[i] - d2->d[i];
639: _NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0;
640: _MKDP(n,m,s2); s2->sugar = d->td; _muld_dup(CO,s2,p2,&u); _free_dp(s2);
641:
642: _addd_destructive(CO,t,u,rp);
1.7 noro 643: if ( GenTrace ) {
644: LIST hist;
645: NODE node;
646:
647: node = mknode(4,ONE,0,s1,ONE);
648: MKLIST(hist,node);
649: MKNODE(TraceList,hist,0);
650:
651: node = mknode(4,ONE,0,0,ONE);
652: chsgnd(s2,(DP *)&ARG2(node));
653: MKLIST(hist,node);
654: MKNODE(node,hist,TraceList); TraceList = node;
655: }
656: }
657:
1.20 noro 658: void dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
1.7 noro 659: {
660: int i,n,td;
661: int *w;
662: DL d1,d2,d;
663: MP m;
664: DP t,s,u;
665:
666: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
667: w = (int *)ALLOCA(n*sizeof(int));
668: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 669: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 670: }
1.18 noro 671: NEWDL_NOINIT(d,n); d->td = td - d1->td;
1.7 noro 672: for ( i = 0; i < n; i++ )
673: d->d[i] = w[i] - d1->d[i];
674: NEWMP(m); m->dl = d; m->c = (P)BDY(p2)->c; NEXT(m) = 0;
675: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p1,s,&t);
1.18 noro 676: NEWDL_NOINIT(d,n); d->td = td - d2->td;
1.7 noro 677: for ( i = 0; i < n; i++ )
678: d->d[i] = w[i] - d2->d[i];
679: NEWMP(m); m->dl = d; m->c = (P)BDY(p1)->c; NEXT(m) = 0;
680: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p2,s,&u);
681: submd(CO,mod,t,u,rp);
682: }
683:
1.20 noro 684: void _dp_sp_mod_dup(DP p1,DP p2,int mod,DP *rp)
1.7 noro 685: {
686: int i,n,td;
687: int *w;
688: DL d1,d2,d;
689: MP m;
690: DP t,s,u;
691:
692: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
693: w = (int *)ALLOCA(n*sizeof(int));
694: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 695: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 696: }
697: _NEWDL(d,n); d->td = td - d1->td;
698: for ( i = 0; i < n; i++ )
699: d->d[i] = w[i] - d1->d[i];
700: _NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
701: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p1,&t); _free_dp(s);
702: _NEWDL(d,n); d->td = td - d2->td;
703: for ( i = 0; i < n; i++ )
704: d->d[i] = w[i] - d2->d[i];
705: _NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
706: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&u); _free_dp(s);
707: _addmd_destructive(mod,t,u,rp);
708: }
709:
1.20 noro 710: void _dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
1.7 noro 711: {
712: int i,n,td;
713: int *w;
714: DL d1,d2,d;
715: MP m;
716: DP t,s,u;
717:
718: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
719: w = (int *)ALLOCA(n*sizeof(int));
720: for ( i = 0, td = 0; i < n; i++ ) {
1.21 noro 721: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
1.7 noro 722: }
723: NEWDL(d,n); d->td = td - d1->td;
724: for ( i = 0; i < n; i++ )
725: d->d[i] = w[i] - d1->d[i];
726: NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
727: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p1,&t);
728: NEWDL(d,n); d->td = td - d2->td;
729: for ( i = 0; i < n; i++ )
730: d->d[i] = w[i] - d2->d[i];
731: NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
732: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p2,&u);
733: addmd_destructive(mod,t,u,rp);
734: }
735:
736: /*
737: * m-reduction
1.13 noro 738: * do content reduction over Z or Q(x,...)
739: * do nothing over finite fields
1.7 noro 740: *
741: */
742:
1.20 noro 743: void dp_red(DP p0,DP p1,DP p2,DP *head,DP *rest,P *dnp,DP *multp)
1.7 noro 744: {
745: int i,n;
746: DL d1,d2,d;
747: MP m;
748: DP t,s,r,h;
749: Q c,c1,c2;
750: N gn,tn;
751: P g,a;
1.23 noro 752: P p[2];
1.7 noro 753:
754: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
755: NEWDL(d,n); d->td = d1->td - d2->td;
756: for ( i = 0; i < n; i++ )
757: d->d[i] = d1->d[i]-d2->d[i];
758: c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c;
1.23 noro 759: if ( dp_fcoeffs == N_GFS ) {
760: p[0] = (P)c1; p[1] = (P)c2;
761: gcdsf(CO,p,2,&g);
762: divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a;
763: } else if ( dp_fcoeffs ) {
1.7 noro 764: /* do nothing */
765: } else if ( INT(c1) && INT(c2) ) {
766: gcdn(NM(c1),NM(c2),&gn);
767: if ( !UNIN(gn) ) {
768: divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c;
769: divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c;
770: }
771: } else {
772: ezgcdpz(CO,(P)c1,(P)c2,&g);
773: divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a;
774: }
775: NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
776: *multp = s;
777: muld(CO,s,p2,&t); muldc(CO,p1,(P)c2,&s); addd(CO,s,t,&r);
778: muldc(CO,p0,(P)c2,&h);
779: *head = h; *rest = r; *dnp = (P)c2;
780: }
781:
1.13 noro 782: /* m-reduction over a field */
783:
1.20 noro 784: void dp_red_f(DP p1,DP p2,DP *rest)
1.13 noro 785: {
786: int i,n;
787: DL d1,d2,d;
788: MP m;
1.20 noro 789: DP t,s;
1.13 noro 790: Obj a,b;
791:
792: n = p1->nv;
793: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
794:
795: NEWDL(d,n); d->td = d1->td - d2->td;
796: for ( i = 0; i < n; i++ )
797: d->d[i] = d1->d[i]-d2->d[i];
798:
799: NEWMP(m); m->dl = d;
800: divr(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); chsgnr(a,&b);
801: C(m) = (P)b;
802: NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
803:
804: muld(CO,s,p2,&t); addd(CO,p1,t,rest);
805: }
806:
1.20 noro 807: void dp_red_mod(DP p0,DP p1,DP p2,int mod,DP *head,DP *rest,P *dnp)
1.7 noro 808: {
809: int i,n;
810: DL d1,d2,d;
811: MP m;
812: DP t,s,r,h;
813: P c1,c2,g,u;
814:
815: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
816: NEWDL(d,n); d->td = d1->td - d2->td;
817: for ( i = 0; i < n; i++ )
818: d->d[i] = d1->d[i]-d2->d[i];
819: c1 = (P)BDY(p1)->c; c2 = (P)BDY(p2)->c;
820: gcdprsmp(CO,mod,c1,c2,&g);
821: divsmp(CO,mod,c1,g,&u); c1 = u; divsmp(CO,mod,c2,g,&u); c2 = u;
822: if ( NUM(c2) ) {
823: divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM;
824: }
825: NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,&m->c); NEXT(m) = 0;
1.11 noro 826: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,s,p2,&t);
1.7 noro 827: if ( NUM(c2) ) {
828: addmd(CO,mod,p1,t,&r); h = p0;
829: } else {
830: mulmdc(CO,mod,p1,c2,&s); addmd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h);
831: }
832: *head = h; *rest = r; *dnp = c2;
833: }
834:
1.10 noro 835: struct oEGT eg_red_mod;
836:
1.20 noro 837: void _dp_red_mod_destructive(DP p1,DP p2,int mod,DP *rp)
1.7 noro 838: {
839: int i,n;
840: DL d1,d2,d;
841: MP m;
842: DP t,s;
1.16 noro 843: int c,c1,c2;
844: extern int do_weyl;
1.7 noro 845:
846: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
847: _NEWDL(d,n); d->td = d1->td - d2->td;
848: for ( i = 0; i < n; i++ )
849: d->d[i] = d1->d[i]-d2->d[i];
1.16 noro 850: c = invm(ITOS(BDY(p2)->c),mod);
851: c2 = ITOS(BDY(p1)->c);
852: DMAR(c,c2,0,mod,c1);
1.7 noro 853: _NEWMP(m); m->dl = d; m->c = STOI(mod-c1); NEXT(m) = 0;
1.16 noro 854: #if 0
1.7 noro 855: _MKDP(n,m,s); s->sugar = d->td;
856: _mulmd_dup(mod,s,p2,&t); _free_dp(s);
1.16 noro 857: #else
858: if ( do_weyl ) {
1.19 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: _mulmdm_dup(mod,p2,m,&t); _FREEMP(m);
863: }
864: #endif
1.10 noro 865: /* get_eg(&t0); */
1.7 noro 866: _addmd_destructive(mod,p1,t,rp);
1.10 noro 867: /* get_eg(&t1); add_eg(&eg_red_mod,&t0,&t1); */
1.7 noro 868: }
869:
870: /*
871: * normal form computation
872: *
873: */
1.5 noro 874:
1.20 noro 875: void dp_true_nf(NODE b,DP g,DP *ps,int full,DP *rp,P *dnp)
1.5 noro 876: {
877: DP u,p,d,s,t,dmy;
878: NODE l;
879: MP m,mr;
880: int i,n;
881: int *wb;
882: int sugar,psugar;
883: P dn,tdn,tdn1;
884:
885: dn = (P)ONE;
886: if ( !g ) {
887: *rp = 0; *dnp = dn; return;
888: }
889: for ( n = 0, l = b; l; l = NEXT(l), n++ );
890: wb = (int *)ALLOCA(n*sizeof(int));
891: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
892: wb[i] = QTOS((Q)BDY(l));
893: sugar = g->sugar;
894: for ( d = 0; g; ) {
895: for ( u = 0, i = 0; i < n; i++ ) {
896: if ( dp_redble(g,p = ps[wb[i]]) ) {
897: dp_red(d,g,p,&t,&u,&tdn,&dmy);
898: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
899: sugar = MAX(sugar,psugar);
900: if ( !u ) {
901: if ( d )
902: d->sugar = sugar;
903: *rp = d; *dnp = dn; return;
904: } else {
905: d = t;
906: mulp(CO,dn,tdn,&tdn1); dn = tdn1;
907: }
908: break;
909: }
910: }
911: if ( u )
912: g = u;
913: else if ( !full ) {
914: if ( g ) {
915: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
916: }
917: *rp = g; *dnp = dn; return;
918: } else {
919: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
920: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
921: addd(CO,d,t,&s); d = s;
922: dp_rest(g,&t); g = t;
923: }
924: }
925: if ( d )
926: d->sugar = sugar;
927: *rp = d; *dnp = dn;
928: }
929:
1.13 noro 930: /* nf computation over Z */
931:
1.20 noro 932: void dp_nf_z(NODE b,DP g,DP *ps,int full,int multiple,DP *rp)
1.5 noro 933: {
934: DP u,p,d,s,t,dmy1;
935: P dmy;
936: NODE l;
937: MP m,mr;
938: int i,n;
939: int *wb;
940: int hmag;
941: int sugar,psugar;
942:
943: if ( !g ) {
944: *rp = 0; return;
945: }
946: for ( n = 0, l = b; l; l = NEXT(l), n++ );
947: wb = (int *)ALLOCA(n*sizeof(int));
948: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
949: wb[i] = QTOS((Q)BDY(l));
1.12 noro 950:
1.13 noro 951: hmag = multiple*HMAG(g);
1.5 noro 952: sugar = g->sugar;
1.12 noro 953:
1.5 noro 954: for ( d = 0; g; ) {
955: for ( u = 0, i = 0; i < n; i++ ) {
956: if ( dp_redble(g,p = ps[wb[i]]) ) {
957: dp_red(d,g,p,&t,&u,&dmy,&dmy1);
958: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
959: sugar = MAX(sugar,psugar);
960: if ( !u ) {
961: if ( d )
962: d->sugar = sugar;
963: *rp = d; return;
964: }
965: d = t;
966: break;
967: }
968: }
969: if ( u ) {
970: g = u;
971: if ( d ) {
1.13 noro 972: if ( multiple && HMAG(d) > hmag ) {
1.5 noro 973: dp_ptozp2(d,g,&t,&u); d = t; g = u;
974: hmag = multiple*HMAG(d);
975: }
976: } else {
1.13 noro 977: if ( multiple && HMAG(g) > hmag ) {
1.5 noro 978: dp_ptozp(g,&t); g = t;
979: hmag = multiple*HMAG(g);
980: }
981: }
982: }
983: else if ( !full ) {
984: if ( g ) {
985: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
986: }
987: *rp = g; return;
988: } else {
989: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
990: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
991: addd(CO,d,t,&s); d = s;
992: dp_rest(g,&t); g = t;
993:
994: }
995: }
996: if ( d )
997: d->sugar = sugar;
998: *rp = d;
999: }
1000:
1.13 noro 1001: /* nf computation over a field */
1002:
1.20 noro 1003: void dp_nf_f(NODE b,DP g,DP *ps,int full,DP *rp)
1.13 noro 1004: {
1005: DP u,p,d,s,t;
1006: NODE l;
1007: MP m,mr;
1008: int i,n;
1009: int *wb;
1010: int sugar,psugar;
1011:
1012: if ( !g ) {
1013: *rp = 0; return;
1014: }
1015: for ( n = 0, l = b; l; l = NEXT(l), n++ );
1016: wb = (int *)ALLOCA(n*sizeof(int));
1017: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
1018: wb[i] = QTOS((Q)BDY(l));
1019:
1020: sugar = g->sugar;
1021: for ( d = 0; g; ) {
1022: for ( u = 0, i = 0; i < n; i++ ) {
1023: if ( dp_redble(g,p = ps[wb[i]]) ) {
1024: dp_red_f(g,p,&u);
1025: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1026: sugar = MAX(sugar,psugar);
1027: if ( !u ) {
1028: if ( d )
1029: d->sugar = sugar;
1030: *rp = d; return;
1031: }
1032: break;
1033: }
1034: }
1035: if ( u )
1036: g = u;
1037: else if ( !full ) {
1038: if ( g ) {
1039: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1040: }
1041: *rp = g; return;
1042: } else {
1043: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1044: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1045: addd(CO,d,t,&s); d = s;
1046: dp_rest(g,&t); g = t;
1047: }
1048: }
1049: if ( d )
1050: d->sugar = sugar;
1051: *rp = d;
1052: }
1053:
1054: /* nf computation over GF(mod) (only for internal use) */
1055:
1.20 noro 1056: void dp_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
1.5 noro 1057: {
1058: DP u,p,d,s,t;
1059: P dmy;
1060: NODE l;
1061: MP m,mr;
1062: int sugar,psugar;
1063:
1064: if ( !g ) {
1065: *rp = 0; return;
1066: }
1067: sugar = g->sugar;
1068: for ( d = 0; g; ) {
1069: for ( u = 0, l = b; l; l = NEXT(l) ) {
1070: if ( dp_redble(g,p = ps[(int)BDY(l)]) ) {
1071: dp_red_mod(d,g,p,mod,&t,&u,&dmy);
1072: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1073: sugar = MAX(sugar,psugar);
1074: if ( !u ) {
1075: if ( d )
1076: d->sugar = sugar;
1077: *rp = d; return;
1078: }
1079: d = t;
1080: break;
1081: }
1082: }
1083: if ( u )
1084: g = u;
1085: else if ( !full ) {
1086: if ( g ) {
1087: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1088: }
1089: *rp = g; return;
1090: } else {
1091: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1092: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1093: addmd(CO,mod,d,t,&s); d = s;
1094: dp_rest(g,&t); g = t;
1095: }
1096: }
1097: if ( d )
1098: d->sugar = sugar;
1099: *rp = d;
1100: }
1101:
1.20 noro 1102: void dp_true_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp,P *dnp)
1.5 noro 1103: {
1104: DP u,p,d,s,t;
1105: NODE l;
1106: MP m,mr;
1107: int i,n;
1108: int *wb;
1109: int sugar,psugar;
1110: P dn,tdn,tdn1;
1111:
1112: dn = (P)ONEM;
1113: if ( !g ) {
1114: *rp = 0; *dnp = dn; return;
1115: }
1116: for ( n = 0, l = b; l; l = NEXT(l), n++ );
1117: wb = (int *)ALLOCA(n*sizeof(int));
1118: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
1119: wb[i] = QTOS((Q)BDY(l));
1120: sugar = g->sugar;
1121: for ( d = 0; g; ) {
1122: for ( u = 0, i = 0; i < n; i++ ) {
1123: if ( dp_redble(g,p = ps[wb[i]]) ) {
1124: dp_red_mod(d,g,p,mod,&t,&u,&tdn);
1125: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1126: sugar = MAX(sugar,psugar);
1127: if ( !u ) {
1128: if ( d )
1129: d->sugar = sugar;
1130: *rp = d; *dnp = dn; return;
1131: } else {
1132: d = t;
1133: mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1;
1134: }
1135: break;
1136: }
1137: }
1138: if ( u )
1139: g = u;
1140: else if ( !full ) {
1141: if ( g ) {
1142: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
1143: }
1144: *rp = g; *dnp = dn; return;
1145: } else {
1146: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
1147: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
1148: addmd(CO,mod,d,t,&s); d = s;
1149: dp_rest(g,&t); g = t;
1150: }
1151: }
1152: if ( d )
1153: d->sugar = sugar;
1154: *rp = d; *dnp = dn;
1155: }
1156:
1.20 noro 1157: void _dp_nf_mod_destructive(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
1.5 noro 1158: {
1.20 noro 1159: DP u,p,d;
1.7 noro 1160: NODE l;
1.20 noro 1161: MP m,mrd;
1162: int sugar,psugar,n,h_reducible;
1.5 noro 1163:
1.7 noro 1164: if ( !g ) {
1165: *rp = 0; return;
1.5 noro 1166: }
1.7 noro 1167: sugar = g->sugar;
1168: n = g->nv;
1169: for ( d = 0; g; ) {
1170: for ( h_reducible = 0, l = b; l; l = NEXT(l) ) {
1171: if ( dp_redble(g,p = ps[(int)BDY(l)]) ) {
1172: h_reducible = 1;
1173: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
1174: _dp_red_mod_destructive(g,p,mod,&u); g = u;
1175: sugar = MAX(sugar,psugar);
1176: if ( !g ) {
1177: if ( d )
1178: d->sugar = sugar;
1179: _dptodp(d,rp); _free_dp(d); return;
1180: }
1181: break;
1182: }
1183: }
1184: if ( !h_reducible ) {
1185: /* head term is not reducible */
1186: if ( !full ) {
1187: if ( g )
1188: g->sugar = sugar;
1189: _dptodp(g,rp); _free_dp(g); return;
1190: } else {
1191: m = BDY(g);
1192: if ( NEXT(m) ) {
1193: BDY(g) = NEXT(m); NEXT(m) = 0;
1194: } else {
1195: _FREEDP(g); g = 0;
1196: }
1197: if ( d ) {
1198: for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) );
1199: NEXT(mrd) = m;
1200: } else {
1201: _MKDP(n,m,d);
1202: }
1203: }
1204: }
1.5 noro 1205: }
1.7 noro 1206: if ( d )
1207: d->sugar = sugar;
1208: _dptodp(d,rp); _free_dp(d);
1.5 noro 1209: }
1.13 noro 1210:
1211: /* reduction by linear base over a field */
1212:
1.20 noro 1213: void dp_lnf_f(DP p1,DP p2,NODE g,DP *r1p,DP *r2p)
1.13 noro 1214: {
1215: DP r1,r2,b1,b2,t,s;
1216: Obj c,c1,c2;
1217: NODE l,b;
1218: int n;
1219:
1220: if ( !p1 ) {
1221: *r1p = p1; *r2p = p2; return;
1222: }
1223: n = p1->nv;
1224: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
1225: if ( !r1 ) {
1226: *r1p = r1; *r2p = r2; return;
1227: }
1228: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
1229: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
1230: b2 = (DP)BDY(NEXT(b));
1231: divr(CO,(Obj)ONE,(Obj)BDY(b1)->c,&c1);
1232: mulr(CO,c1,(Obj)BDY(r1)->c,&c2); chsgnr(c2,&c);
1233: muldc(CO,b1,(P)c,&t); addd(CO,r1,t,&s); r1 = s;
1234: muldc(CO,b2,(P)c,&t); addd(CO,r2,t,&s); r2 = s;
1235: }
1236: }
1237: *r1p = r1; *r2p = r2;
1238: }
1239:
1240: /* reduction by linear base over GF(mod) */
1.5 noro 1241:
1.20 noro 1242: void dp_lnf_mod(DP p1,DP p2,NODE g,int mod,DP *r1p,DP *r2p)
1.5 noro 1243: {
1.7 noro 1244: DP r1,r2,b1,b2,t,s;
1245: P c;
1246: MQ 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: invmq(mod,(MQ)BDY(b1)->c,&c1);
1262: mulmq(mod,c1,(MQ)BDY(r1)->c,&c2); chsgnmp(mod,(P)c2,&c);
1263: mulmdc(CO,mod,b1,c,&t); addmd(CO,mod,r1,t,&s); r1 = s;
1264: mulmdc(CO,mod,b2,c,&t); addmd(CO,mod,r2,t,&s); r2 = s;
1265: }
1266: }
1267: *r1p = r1; *r2p = r2;
1.5 noro 1268: }
1269:
1.20 noro 1270: void dp_nf_tab_mod(DP p,LIST *tab,int mod,DP *rp)
1.5 noro 1271: {
1.7 noro 1272: DP s,t,u;
1273: MP m;
1274: DL h;
1275: int i,n;
1276:
1277: if ( !p ) {
1278: *rp = p; return;
1279: }
1280: n = p->nv;
1281: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
1282: h = m->dl;
1283: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
1284: i++;
1285: mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
1286: addmd(CO,mod,s,t,&u); s = u;
1.24 noro 1287: }
1288: *rp = s;
1289: }
1290:
1291: void dp_nf_tab_f(DP p,LIST *tab,DP *rp)
1292: {
1293: DP s,t,u;
1294: MP m;
1295: DL h;
1296: int i,n;
1297:
1298: if ( !p ) {
1299: *rp = p; return;
1300: }
1301: n = p->nv;
1302: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
1303: h = m->dl;
1304: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
1305: i++;
1306: muldc(CO,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
1307: addd(CO,s,t,&u); s = u;
1.7 noro 1308: }
1309: *rp = s;
1.5 noro 1310: }
1311:
1.7 noro 1312: /*
1313: * setting flags
1314: *
1315: */
1316:
1.27 ! noro 1317: int create_order_spec(VL vl,Obj obj,struct order_spec **specp)
1.5 noro 1318: {
1.7 noro 1319: int i,j,n,s,row,col;
1.27 ! noro 1320: struct order_spec *spec;
1.7 noro 1321: struct order_pair *l;
1322: NODE node,t,tn;
1323: MAT m;
1324: pointer **b;
1325: int **w;
1.5 noro 1326:
1.27 ! noro 1327: if ( vl && obj && OID(obj) == O_LIST )
! 1328: return create_composite_order_spec(vl,(LIST)obj,specp);
! 1329:
! 1330: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
1.7 noro 1331: if ( !obj || NUM(obj) ) {
1332: spec->id = 0; spec->obj = obj;
1333: spec->ord.simple = QTOS((Q)obj);
1334: return 1;
1335: } else if ( OID(obj) == O_LIST ) {
1336: node = BDY((LIST)obj);
1337: for ( n = 0, t = node; t; t = NEXT(t), n++ );
1338: l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair));
1339: for ( i = 0, t = node, s = 0; i < n; t = NEXT(t), i++ ) {
1340: tn = BDY((LIST)BDY(t)); l[i].order = QTOS((Q)BDY(tn));
1341: tn = NEXT(tn); l[i].length = QTOS((Q)BDY(tn));
1342: s += l[i].length;
1343: }
1344: spec->id = 1; spec->obj = obj;
1345: spec->ord.block.order_pair = l;
1346: spec->ord.block.length = n; spec->nv = s;
1347: return 1;
1348: } else if ( OID(obj) == O_MAT ) {
1349: m = (MAT)obj; row = m->row; col = m->col; b = BDY(m);
1350: w = almat(row,col);
1351: for ( i = 0; i < row; i++ )
1352: for ( j = 0; j < col; j++ )
1353: w[i][j] = QTOS((Q)b[i][j]);
1354: spec->id = 2; spec->obj = obj;
1355: spec->nv = col; spec->ord.matrix.row = row;
1356: spec->ord.matrix.matrix = w;
1357: return 1;
1358: } else
1.5 noro 1359: return 0;
1360: }
1361:
1.27 ! noro 1362: /* order = [w_or_b, w_or_b, ... ] */
! 1363: /* w_or_b = w or b */
! 1364: /* w = [1,2,...] or [x,1,y,2,...] */
! 1365: /* b = [@lex,x,y,...,z] etc */
! 1366:
! 1367: int create_composite_order_spec(VL vl,LIST order,struct order_spec **specp)
! 1368: {
! 1369: NODE wb,t,p;
! 1370: struct order_spec *spec;
! 1371: VL tvl;
! 1372: int n,i,j,k,l,len;
! 1373: int *dw;
! 1374: struct sparse_weight *sw;
! 1375: struct weight_or_block *w_or_b;
! 1376: Obj a0;
! 1377: NODE a;
! 1378: V v;
! 1379:
! 1380: /* l = number of vars in vl */
! 1381: for ( l = 0, tvl = vl; tvl; tvl = NEXT(tvl), l++ );
! 1382: /* n = number of primitives in order */
! 1383: wb = BDY(order);
! 1384: n = length(wb);
! 1385: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 1386: spec->id = 3;
! 1387: spec->obj = (Obj)order;
! 1388: spec->nv = l;
! 1389: spec->ord.composite.length = n;
! 1390: spec->ord.composite.w_or_b = (struct weight_or_block *)
! 1391: MALLOC(sizeof(struct weight_or_block)*n);
! 1392: for ( t = wb, i = 0; t; t = NEXT(t) ) {
! 1393: a = BDY((LIST)BDY(wb));
! 1394: len = length(a);
! 1395: a0 = (Obj)BDY(a);
! 1396: if ( !a0 || OID(a0) == O_N ) {
! 1397: /* a is dense weight */
! 1398: dw = (int *)MALLOC(sizeof(int)*len);
! 1399: for ( j = 0, p = a; j < len; p = NEXT(p), j++ )
! 1400: dw[j] = QTOS((Q)BDY(p));
! 1401: w_or_b[i].type = IS_DENSE_WEIGHT;
! 1402: w_or_b[i].length = len;
! 1403: w_or_b[i].body.dense_weight = dw;
! 1404: } else if ( OID(a0) == O_P ) {
! 1405: sw = (struct sparse_weight *)
! 1406: MALLOC(sizeof(struct sparse_weight)*len);
! 1407: for ( j = 0, p = a; j < len; j++ ) {
! 1408: v = VR((P)BDY(a)); a = NEXT(a);
! 1409: for ( tvl = vl, k = 0; tvl && tvl->v != v;
! 1410: k++, tvl = NEXT(tvl) );
! 1411: if ( !tvl )
! 1412: error("invalid variable name");
! 1413: sw[j].pos = k;
! 1414: sw[j].value = QTOS((Q)BDY(a));
! 1415: }
! 1416: w_or_b[i].type = IS_SPARSE_WEIGHT;
! 1417: w_or_b[i].length = len;
! 1418: w_or_b[i].body.sparse_weight = sw;
! 1419: } else {
! 1420: error("not implemented yet");
! 1421: }
! 1422: }
! 1423: }
! 1424:
1.7 noro 1425: /*
1426: * converters
1427: *
1428: */
1429:
1.20 noro 1430: void dp_homo(DP p,DP *rp)
1.5 noro 1431: {
1.7 noro 1432: MP m,mr,mr0;
1433: int i,n,nv,td;
1434: DL dl,dlh;
1.5 noro 1435:
1.7 noro 1436: if ( !p )
1437: *rp = 0;
1438: else {
1439: n = p->nv; nv = n + 1;
1440: m = BDY(p); td = sugard(m);
1441: for ( mr0 = 0; m; m = NEXT(m) ) {
1442: NEXTMP(mr0,mr); mr->c = m->c;
1443: dl = m->dl;
1444: mr->dl = dlh = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
1445: dlh->td = td;
1446: for ( i = 0; i < n; i++ )
1447: dlh->d[i] = dl->d[i];
1448: dlh->d[n] = td - dl->td;
1449: }
1450: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
1.5 noro 1451: }
1452: }
1453:
1.20 noro 1454: void dp_dehomo(DP p,DP *rp)
1.5 noro 1455: {
1.7 noro 1456: MP m,mr,mr0;
1457: int i,n,nv;
1458: DL dl,dlh;
1.5 noro 1459:
1.7 noro 1460: if ( !p )
1461: *rp = 0;
1462: else {
1463: n = p->nv; nv = n - 1;
1464: m = BDY(p);
1465: for ( mr0 = 0; m; m = NEXT(m) ) {
1466: NEXTMP(mr0,mr); mr->c = m->c;
1467: dlh = m->dl;
1468: mr->dl = dl = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
1469: dl->td = dlh->td - dlh->d[nv];
1470: for ( i = 0; i < nv; i++ )
1471: dl->d[i] = dlh->d[i];
1472: }
1473: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
1474: }
1.5 noro 1475: }
1476:
1.20 noro 1477: void dp_mod(DP p,int mod,NODE subst,DP *rp)
1.5 noro 1478: {
1.7 noro 1479: MP m,mr,mr0;
1480: P t,s,s1;
1481: V v;
1482: NODE tn;
1.5 noro 1483:
1.7 noro 1484: if ( !p )
1485: *rp = 0;
1486: else {
1487: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
1488: for ( tn = subst, s = m->c; tn; tn = NEXT(tn) ) {
1489: v = VR((P)BDY(tn)); tn = NEXT(tn);
1490: substp(CO,s,v,(P)BDY(tn),&s1); s = s1;
1491: }
1492: ptomp(mod,s,&t);
1493: if ( t ) {
1494: NEXTMP(mr0,mr); mr->c = t; mr->dl = m->dl;
1495: }
1496: }
1497: if ( mr0 ) {
1498: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1499: } else
1500: *rp = 0;
1501: }
1.5 noro 1502: }
1503:
1.20 noro 1504: void dp_rat(DP p,DP *rp)
1.5 noro 1505: {
1.7 noro 1506: MP m,mr,mr0;
1.5 noro 1507:
1.7 noro 1508: if ( !p )
1509: *rp = 0;
1510: else {
1511: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
1512: NEXTMP(mr0,mr); mptop(m->c,&mr->c); mr->dl = m->dl;
1.5 noro 1513: }
1.7 noro 1514: if ( mr0 ) {
1515: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
1516: } else
1517: *rp = 0;
1.5 noro 1518: }
1519: }
1520:
1521:
1.27 ! noro 1522: void homogenize_order(struct order_spec *old,int n,struct order_spec **newp)
1.5 noro 1523: {
1.7 noro 1524: struct order_pair *l;
1525: int length,nv,row,i,j;
1526: int **newm,**oldm;
1.27 ! noro 1527: struct order_spec *new;
1.5 noro 1528:
1.27 ! noro 1529: *newp = new = (struct order_spec *)MALLOC(sizeof(struct order_spec));
1.7 noro 1530: switch ( old->id ) {
1531: case 0:
1532: switch ( old->ord.simple ) {
1533: case 0:
1534: new->id = 0; new->ord.simple = 0; break;
1535: case 1:
1536: l = (struct order_pair *)
1537: MALLOC_ATOMIC(2*sizeof(struct order_pair));
1538: l[0].length = n; l[0].order = 1;
1539: l[1].length = 1; l[1].order = 2;
1540: new->id = 1;
1541: new->ord.block.order_pair = l;
1542: new->ord.block.length = 2; new->nv = n+1;
1543: break;
1544: case 2:
1545: new->id = 0; new->ord.simple = 1; break;
1546: case 3: case 4: case 5:
1547: new->id = 0; new->ord.simple = old->ord.simple+3;
1548: dp_nelim = n-1; break;
1549: case 6: case 7: case 8: case 9:
1550: new->id = 0; new->ord.simple = old->ord.simple; break;
1551: default:
1552: error("homogenize_order : invalid input");
1553: }
1554: break;
1555: case 1:
1556: length = old->ord.block.length;
1557: l = (struct order_pair *)
1558: MALLOC_ATOMIC((length+1)*sizeof(struct order_pair));
1559: bcopy((char *)old->ord.block.order_pair,(char *)l,length*sizeof(struct order_pair));
1560: l[length].order = 2; l[length].length = 1;
1561: new->id = 1; new->nv = n+1;
1562: new->ord.block.order_pair = l;
1563: new->ord.block.length = length+1;
1564: break;
1565: case 2:
1566: nv = old->nv; row = old->ord.matrix.row;
1567: oldm = old->ord.matrix.matrix; newm = almat(row+1,nv+1);
1568: for ( i = 0; i <= nv; i++ )
1569: newm[0][i] = 1;
1570: for ( i = 0; i < row; i++ ) {
1571: for ( j = 0; j < nv; j++ )
1572: newm[i+1][j] = oldm[i][j];
1573: newm[i+1][j] = 0;
1574: }
1575: new->id = 2; new->nv = nv+1;
1576: new->ord.matrix.row = row+1; new->ord.matrix.matrix = newm;
1577: break;
1578: default:
1579: error("homogenize_order : invalid input");
1.5 noro 1580: }
1.7 noro 1581: }
1582:
1.20 noro 1583: void qltozl(Q *w,int n,Q *dvr)
1.7 noro 1584: {
1585: N nm,dn;
1586: N g,l1,l2,l3;
1587: Q c,d;
1588: int i;
1589: struct oVECT v;
1.5 noro 1590:
1591: for ( i = 0; i < n; i++ )
1.7 noro 1592: if ( w[i] && !INT(w[i]) )
1593: break;
1594: if ( i == n ) {
1595: v.id = O_VECT; v.len = n; v.body = (pointer *)w;
1596: igcdv(&v,dvr); return;
1597: }
1598: c = w[0]; nm = NM(c); dn = INT(c) ? ONEN : DN(c);
1599: for ( i = 1; i < n; i++ ) {
1600: c = w[i]; l1 = INT(c) ? ONEN : DN(c);
1601: gcdn(nm,NM(c),&g); nm = g;
1602: gcdn(dn,l1,&l2); muln(dn,l1,&l3); divsn(l3,l2,&dn);
1.5 noro 1603: }
1.7 noro 1604: if ( UNIN(dn) )
1605: NTOQ(nm,1,d);
1606: else
1607: NDTOQ(nm,dn,1,d);
1608: *dvr = d;
1609: }
1.5 noro 1610:
1.20 noro 1611: int comp_nm(Q *a,Q *b)
1.7 noro 1612: {
1613: return cmpn((*a)?NM(*a):0,(*b)?NM(*b):0);
1614: }
1615:
1.20 noro 1616: void sortbynm(Q *w,int n)
1.7 noro 1617: {
1618: qsort(w,n,sizeof(Q),(int (*)(const void *,const void *))comp_nm);
1619: }
1.5 noro 1620:
1621:
1.7 noro 1622: /*
1623: * simple operations
1624: *
1625: */
1.5 noro 1626:
1.20 noro 1627: int dp_redble(DP p1,DP p2)
1.7 noro 1628: {
1629: int i,n;
1630: DL d1,d2;
1.5 noro 1631:
1.7 noro 1632: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
1633: if ( d1->td < d2->td )
1634: return 0;
1635: else {
1636: for ( i = 0, n = p1->nv; i < n; i++ )
1637: if ( d1->d[i] < d2->d[i] )
1638: return 0;
1639: return 1;
1.5 noro 1640: }
1641: }
1642:
1.20 noro 1643: void dp_subd(DP p1,DP p2,DP *rp)
1.5 noro 1644: {
1.7 noro 1645: int i,n;
1.5 noro 1646: DL d1,d2,d;
1647: MP m;
1.7 noro 1648: DP s;
1.5 noro 1649:
1650: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
1.7 noro 1651: NEWDL(d,n); d->td = d1->td - d2->td;
1.5 noro 1652: for ( i = 0; i < n; i++ )
1.7 noro 1653: d->d[i] = d1->d[i]-d2->d[i];
1654: NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0;
1655: MKDP(n,m,s); s->sugar = d->td;
1656: *rp = s;
1657: }
1658:
1.20 noro 1659: void dltod(DL d,int n,DP *rp)
1.7 noro 1660: {
1661: MP m;
1662: DP s;
1663:
1664: NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0;
1665: MKDP(n,m,s); s->sugar = d->td;
1666: *rp = s;
1.5 noro 1667: }
1668:
1.20 noro 1669: void dp_hm(DP p,DP *rp)
1.5 noro 1670: {
1671: MP m,mr;
1672:
1673: if ( !p )
1674: *rp = 0;
1675: else {
1676: m = BDY(p);
1677: NEWMP(mr); mr->dl = m->dl; mr->c = m->c; NEXT(mr) = 0;
1678: MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
1679: }
1680: }
1681:
1.20 noro 1682: void dp_rest(DP p,DP *rp)
1.5 noro 1683: {
1684: MP m;
1685:
1686: m = BDY(p);
1687: if ( !NEXT(m) )
1688: *rp = 0;
1689: else {
1690: MKDP(p->nv,NEXT(m),*rp);
1691: if ( *rp )
1692: (*rp)->sugar = p->sugar;
1693: }
1694: }
1695:
1.20 noro 1696: DL lcm_of_DL(int nv,DL dl1,DL dl2,DL dl)
1.5 noro 1697: {
1.21 noro 1698: register int i, *d1, *d2, *d, td;
1.5 noro 1699:
1700: if ( !dl ) NEWDL(dl,nv);
1701: d = dl->d, d1 = dl1->d, d2 = dl2->d;
1.21 noro 1702: for ( td = 0, i = 0; i < nv; d1++, d2++, d++, i++ ) {
1703: *d = *d1 > *d2 ? *d1 : *d2;
1704: td += MUL_WEIGHT(*d,i);
1705: }
1.5 noro 1706: dl->td = td;
1707: return dl;
1708: }
1709:
1.20 noro 1710: int dl_equal(int nv,DL dl1,DL dl2)
1.5 noro 1711: {
1712: register int *d1, *d2, n;
1713:
1714: if ( dl1->td != dl2->td ) return 0;
1715: for ( d1 = dl1->d, d2 = dl2->d, n = nv; --n >= 0; d1++, d2++ )
1716: if ( *d1 != *d2 ) return 0;
1717: return 1;
1718: }
1719:
1.20 noro 1720: int dp_nt(DP p)
1.5 noro 1721: {
1722: int i;
1723: MP m;
1724:
1725: if ( !p )
1726: return 0;
1727: else {
1728: for ( i = 0, m = BDY(p); m; m = NEXT(m), i++ );
1729: return i;
1730: }
1731: }
1732:
1.20 noro 1733: int dp_homogeneous(DP p)
1.15 noro 1734: {
1735: MP m;
1736: int d;
1737:
1738: if ( !p )
1739: return 1;
1740: else {
1741: m = BDY(p);
1742: d = m->dl->td;
1743: m = NEXT(m);
1744: for ( ; m; m = NEXT(m) ) {
1745: if ( m->dl->td != d )
1746: return 0;
1747: }
1748: return 1;
1749: }
1.16 noro 1750: }
1751:
1.20 noro 1752: void _print_mp(int nv,MP m)
1.16 noro 1753: {
1754: int i;
1755:
1.17 noro 1756: if ( !m )
1.16 noro 1757: return;
1758: for ( ; m; m = NEXT(m) ) {
1759: fprintf(stderr,"%d<",ITOS(C(m)));
1760: for ( i = 0; i < nv; i++ ) {
1761: fprintf(stderr,"%d",m->dl->d[i]);
1762: if ( i != nv-1 )
1763: fprintf(stderr," ");
1764: }
1765: fprintf(stderr,">",C(m));
1766: }
1767: fprintf(stderr,"\n");
1.15 noro 1768: }
1.26 noro 1769:
1770: static int cmp_mp_nvar;
1771:
1772: int comp_mp(MP *a,MP *b)
1773: {
1774: return -(*cmpdl)(cmp_mp_nvar,(*a)->dl,(*b)->dl);
1775: }
1776:
1777: void dp_sort(DP p,DP *rp)
1778: {
1779: MP t,mp,mp0;
1780: int i,n;
1781: DP r;
1782: MP *w;
1783:
1784: if ( !p ) {
1785: *rp = 0;
1786: return;
1787: }
1788: for ( t = BDY(p), n = 0; t; t = NEXT(t), n++ );
1789: w = (MP *)ALLOCA(n*sizeof(MP));
1790: for ( t = BDY(p), i = 0; i < n; t = NEXT(t), i++ )
1791: w[i] = t;
1792: cmp_mp_nvar = NV(p);
1793: qsort(w,n,sizeof(MP),(int (*)(const void *,const void *))comp_mp);
1794: mp0 = 0;
1795: for ( i = n-1; i >= 0; i-- ) {
1796: NEWMP(mp); mp->dl = w[i]->dl; C(mp) = C(w[i]);
1797: NEXT(mp) = mp0; mp0 = mp;
1798: }
1799: MKDP(p->nv,mp0,r);
1800: r->sugar = p->sugar;
1801: *rp = r;
1802: }
1803: