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