Annotation of OpenXM_contrib2/asir2000/builtin/dp-supp.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/builtin/dp-supp.c,v 1.2 1999/11/18 05:42:01 noro Exp $ */
! 2: #include "ca.h"
! 3: #include "base.h"
! 4: #include "parse.h"
! 5: #include "ox.h"
! 6:
! 7: extern int (*cmpdl)();
! 8: double pz_t_e,pz_t_d,pz_t_d1,pz_t_c;
! 9:
! 10: void Pdp_mbase(),Pdp_lnf_mod(),Pdp_nf_tab_mod(),Pdp_mdtod();
! 11: void Pdp_vtoe(), Pdp_etov(), Pdp_dtov(), Pdp_idiv(), Pdp_sep();
! 12: void Pdp_cont(), Pdp_igcdv_hist(), Pdp_idivv_hist();
! 13:
! 14: struct ftab dp_supp_tab[] = {
! 15: {"dp_mbase",Pdp_mbase,1},
! 16: {"dp_lnf_mod",Pdp_lnf_mod,3},
! 17: {"dp_nf_tab_mod",Pdp_nf_tab_mod,3},
! 18: {"dp_etov",Pdp_etov,1},
! 19: {"dp_vtoe",Pdp_vtoe,1},
! 20: {"dp_dtov",Pdp_dtov,1},
! 21: {"dp_idiv",Pdp_idiv,2},
! 22: {"dp_cont",Pdp_cont,1},
! 23: {"dp_sep",Pdp_sep,2},
! 24: {"dp_mdtod",Pdp_mdtod,1},
! 25: {"dp_igcdv_hist",Pdp_igcdv_hist,1},
! 26: {"dp_idivv_hist",Pdp_idivv_hist,1},
! 27: {0,0,0}
! 28: };
! 29:
! 30: void Pdp_mdtod(arg,rp)
! 31: NODE arg;
! 32: DP *rp;
! 33: {
! 34: MP m,mr,mr0;
! 35: DP p;
! 36: P t;
! 37:
! 38: p = (DP)ARG0(arg);
! 39: if ( !p )
! 40: *rp = 0;
! 41: else {
! 42: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 43: mptop(m->c,&t); NEXTMP(mr0,mr); mr->c = t; mr->dl = m->dl;
! 44: }
! 45: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 46: }
! 47: }
! 48:
! 49: void Pdp_sep(arg,rp)
! 50: NODE arg;
! 51: VECT *rp;
! 52: {
! 53: DP p,r;
! 54: MP m,t;
! 55: MP *w0,*w;
! 56: int i,n,d,nv,sugar;
! 57: VECT v;
! 58: pointer *pv;
! 59:
! 60: p = (DP)ARG0(arg); m = BDY(p);
! 61: d = QTOS((Q)ARG1(arg));
! 62: for ( t = m, n = 0; t; t = NEXT(t), n++ );
! 63: if ( d > n )
! 64: d = n;
! 65: MKVECT(v,d); *rp = v;
! 66: pv = BDY(v); nv = p->nv; sugar = p->sugar;
! 67: w0 = (MP *)MALLOC(d*sizeof(MP)); bzero(w0,d*sizeof(MP));
! 68: w = (MP *)MALLOC(d*sizeof(MP)); bzero(w,d*sizeof(MP));
! 69: for ( t = BDY(p), i = 0; t; t = NEXT(t), i++, i %= d ) {
! 70: NEXTMP(w0[i],w[i]); w[i]->c = t->c; w[i]->dl = t->dl;
! 71: }
! 72: for ( i = 0; i < d; i++ ) {
! 73: NEXT(w[i]) = 0; MKDP(nv,w0[i],r); r->sugar = sugar;
! 74: pv[i] = (pointer)r;
! 75: }
! 76: }
! 77:
! 78: void Pdp_idiv(arg,rp)
! 79: NODE arg;
! 80: DP *rp;
! 81: {
! 82: dp_idiv((DP)ARG0(arg),(Q)ARG1(arg),rp);
! 83: }
! 84:
! 85: void Pdp_cont(arg,rp)
! 86: NODE arg;
! 87: Q *rp;
! 88: {
! 89: dp_cont((DP)ARG0(arg),rp);
! 90: }
! 91:
! 92: void Pdp_dtov(arg,rp)
! 93: NODE arg;
! 94: VECT *rp;
! 95: {
! 96: dp_dtov((DP)ARG0(arg),rp);
! 97: }
! 98:
! 99: void Pdp_mbase(arg,rp)
! 100: NODE arg;
! 101: LIST *rp;
! 102: {
! 103: NODE mb;
! 104:
! 105: asir_assert(ARG0(arg),O_LIST,"dp_mbase");
! 106: dp_mbase(BDY((LIST)ARG0(arg)),&mb);
! 107: MKLIST(*rp,mb);
! 108: }
! 109:
! 110: void Pdp_etov(arg,rp)
! 111: NODE arg;
! 112: VECT *rp;
! 113: {
! 114: DP dp;
! 115: int n,i;
! 116: int *d;
! 117: VECT v;
! 118: Q t;
! 119:
! 120: dp = (DP)ARG0(arg);
! 121: asir_assert(dp,O_DP,"dp_etov");
! 122: n = dp->nv; d = BDY(dp)->dl->d;
! 123: MKVECT(v,n);
! 124: for ( i = 0; i < n; i++ ) {
! 125: STOQ(d[i],t); v->body[i] = (pointer)t;
! 126: }
! 127: *rp = v;
! 128: }
! 129:
! 130: void Pdp_vtoe(arg,rp)
! 131: NODE arg;
! 132: DP *rp;
! 133: {
! 134: DP dp;
! 135: DL dl;
! 136: MP m;
! 137: int n,i,td;
! 138: int *d;
! 139: VECT v;
! 140:
! 141: v = (VECT)ARG0(arg);
! 142: asir_assert(v,O_VECT,"dp_vtoe");
! 143: n = v->len;
! 144: NEWDL(dl,n); d = dl->d;
! 145: for ( i = 0, td = 0; i < n; i++ ) {
! 146: d[i] = QTOS((Q)(v->body[i])); td += d[i];
! 147: }
! 148: dl->td = td;
! 149: NEWMP(m); m->dl = dl; m->c = (P)ONE; NEXT(m) = 0;
! 150: MKDP(n,m,dp); dp->sugar = td;
! 151: *rp = dp;
! 152: }
! 153:
! 154: void Pdp_lnf_mod(arg,rp)
! 155: NODE arg;
! 156: LIST *rp;
! 157: {
! 158: DP r1,r2;
! 159: NODE b,g,n;
! 160: int mod;
! 161:
! 162: asir_assert(ARG0(arg),O_LIST,"dp_lnf_mod");
! 163: asir_assert(ARG1(arg),O_LIST,"dp_lnf_mod");
! 164: asir_assert(ARG2(arg),O_N,"dp_lnf_mod");
! 165: b = BDY((LIST)ARG0(arg)); g = BDY((LIST)ARG1(arg));
! 166: mod = QTOS((Q)ARG2(arg));
! 167: dp_lnf_mod((DP)BDY(b),(DP)BDY(NEXT(b)),g,mod,&r1,&r2);
! 168: NEWNODE(n); BDY(n) = (pointer)r1;
! 169: NEWNODE(NEXT(n)); BDY(NEXT(n)) = (pointer)r2;
! 170: NEXT(NEXT(n)) = 0; MKLIST(*rp,n);
! 171: }
! 172:
! 173: void Pdp_nf_tab_mod(arg,rp)
! 174: NODE arg;
! 175: DP *rp;
! 176: {
! 177: asir_assert(ARG0(arg),O_DP,"dp_nf_tab_mod");
! 178: asir_assert(ARG1(arg),O_VECT,"dp_nf_tab_mod");
! 179: asir_assert(ARG2(arg),O_N,"dp_nf_tab_mod");
! 180: dp_nf_tab_mod((DP)ARG0(arg),(LIST *)BDY((VECT)ARG1(arg)),
! 181: QTOS((Q)ARG2(arg)),rp);
! 182: }
! 183:
! 184: void Pdp_igcdv_hist(arg,rp)
! 185: NODE arg;
! 186: Q *rp;
! 187: {
! 188: dp_igcdv_hist((DP)ARG0(arg),rp);
! 189: }
! 190:
! 191: void Pdp_idivv_hist(arg,rp)
! 192: NODE arg;
! 193: DP *rp;
! 194: {
! 195: dp_idivv_hist((Q)ARG0(arg),rp);
! 196: }
! 197:
! 198: void dp_idiv(p,c,rp)
! 199: DP p;
! 200: Q c;
! 201: DP *rp;
! 202: {
! 203: Q t;
! 204: N nm,q;
! 205: int sgn,s;
! 206: MP mr0,m,mr;
! 207:
! 208: if ( !p )
! 209: *rp = 0;
! 210: else if ( MUNIQ((Q)c) )
! 211: *rp = p;
! 212: else if ( MUNIQ((Q)c) )
! 213: chsgnd(p,rp);
! 214: else {
! 215: nm = NM(c); sgn = SGN(c);
! 216: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 217: NEXTMP(mr0,mr);
! 218:
! 219: divsn(NM((Q)(m->c)),nm,&q);
! 220: s = sgn*SGN((Q)(m->c));
! 221: NTOQ(q,s,t);
! 222: mr->c = (P)t;
! 223: mr->dl = m->dl;
! 224: }
! 225: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
! 226: if ( *rp )
! 227: (*rp)->sugar = p->sugar;
! 228: }
! 229: }
! 230:
! 231: void dp_cont(p,rp)
! 232: DP p;
! 233: Q *rp;
! 234: {
! 235: VECT v;
! 236:
! 237: dp_dtov(p,&v); igcdv(v,rp);
! 238: }
! 239:
! 240: void dp_dtov(dp,rp)
! 241: DP dp;
! 242: VECT *rp;
! 243: {
! 244: MP m,t;
! 245: int i,n;
! 246: VECT v;
! 247: pointer *p;
! 248:
! 249: m = BDY(dp);
! 250: for ( t = m, n = 0; t; t = NEXT(t), n++ );
! 251: MKVECT(v,n);
! 252: for ( i = 0, p = BDY(v), t = m; i < n; t = NEXT(t), i++ )
! 253: p[i] = (pointer)(t->c);
! 254: *rp = v;
! 255: }
! 256:
! 257: void dp_mbase(hlist,mbase)
! 258: NODE hlist;
! 259: NODE *mbase;
! 260: {
! 261: DL *dl;
! 262: DL d;
! 263: int i,j,n,nvar,td;
! 264:
! 265: n = length(hlist); nvar = ((DP)BDY(hlist))->nv;
! 266: dl = (DL *)MALLOC(n*sizeof(DL));
! 267: for ( i = 0; i < n; i++, hlist = NEXT(hlist) )
! 268: dl[i] = BDY((DP)BDY(hlist))->dl;
! 269: NEWDL(d,nvar); *mbase = 0;
! 270: while ( 1 ) {
! 271: insert_to_node(d,mbase,nvar);
! 272: for ( i = nvar-1; i >= 0; ) {
! 273: d->d[i]++; d->td++;
! 274: for ( j = 0; j < n; j++ ) {
! 275: if ( _dl_redble(dl[j],d,nvar) )
! 276: break;
! 277: }
! 278: if ( j < n ) {
! 279: for ( j = nvar-1; j >= i; j-- )
! 280: d->d[j] = 0;
! 281: for ( j = 0, td = 0; j < i; j++ )
! 282: td += d->d[j];
! 283: d->td = td;
! 284: i--;
! 285: } else
! 286: break;
! 287: }
! 288: if ( i < 0 )
! 289: break;
! 290: }
! 291: }
! 292:
! 293: int _dl_redble(d1,d2,nvar)
! 294: DL d1,d2;
! 295: int nvar;
! 296: {
! 297: int i;
! 298:
! 299: if ( d1->td > d2->td )
! 300: return 0;
! 301: for ( i = 0; i < nvar; i++ )
! 302: if ( d1->d[i] > d2->d[i] )
! 303: break;
! 304: if ( i < nvar )
! 305: return 0;
! 306: else
! 307: return 1;
! 308: }
! 309:
! 310: void insert_to_node(d,n,nvar)
! 311: DL d;
! 312: NODE *n;
! 313: int nvar;
! 314: {
! 315: DL d1;
! 316: MP m;
! 317: DP dp;
! 318: NODE n0,n1,n2;
! 319:
! 320: NEWDL(d1,nvar); d1->td = d->td;
! 321: bcopy((char *)d->d,(char *)d1->d,nvar*sizeof(int));
! 322: NEWMP(m); m->dl = d1; m->c = (P)ONE; NEXT(m) = 0;
! 323: MKDP(nvar,m,dp); dp->sugar = d->td;
! 324: if ( !(*n) ) {
! 325: MKNODE(n1,dp,0); *n = n1;
! 326: } else {
! 327: for ( n1 = *n, n0 = 0; n1; n0 = n1, n1 = NEXT(n1) )
! 328: if ( (*cmpdl)(nvar,d,BDY((DP)BDY(n1))->dl) > 0 ) {
! 329: MKNODE(n2,dp,n1);
! 330: if ( !n0 )
! 331: *n = n2;
! 332: else
! 333: NEXT(n0) = n2;
! 334: break;
! 335: }
! 336: if ( !n1 ) {
! 337: MKNODE(n2,dp,0); NEXT(n0) = n2;
! 338: }
! 339: }
! 340: }
! 341:
! 342: void dp_lnf_mod(p1,p2,g,mod,r1p,r2p)
! 343: DP p1,p2;
! 344: NODE g;
! 345: int mod;
! 346: DP *r1p,*r2p;
! 347: {
! 348: DP r1,r2,b1,b2,t,s;
! 349: P c;
! 350: MQ c1,c2;
! 351: NODE l,b;
! 352: int n;
! 353:
! 354: if ( !p1 ) {
! 355: *r1p = p1; *r2p = p2; return;
! 356: }
! 357: n = p1->nv;
! 358: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
! 359: if ( !r1 ) {
! 360: *r1p = r1; *r2p = r2; return;
! 361: }
! 362: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
! 363: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
! 364: b2 = (DP)BDY(NEXT(b));
! 365: invmq(mod,(MQ)BDY(b1)->c,&c1);
! 366: mulmq(mod,c1,(MQ)BDY(r1)->c,&c2); chsgnmp(mod,(P)c2,&c);
! 367: mulmdc(CO,mod,b1,c,&t); addmd(CO,mod,r1,t,&s); r1 = s;
! 368: mulmdc(CO,mod,b2,c,&t); addmd(CO,mod,r2,t,&s); r2 = s;
! 369: }
! 370: }
! 371: *r1p = r1; *r2p = r2;
! 372: }
! 373:
! 374: void dp_nf_tab_mod(p,tab,mod,rp)
! 375: DP p;
! 376: LIST *tab;
! 377: int mod;
! 378: DP *rp;
! 379: {
! 380: DP s,t,u;
! 381: MP m;
! 382: DL h;
! 383: int i,n;
! 384:
! 385: if ( !p ) {
! 386: *rp = p; return;
! 387: }
! 388: n = p->nv;
! 389: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
! 390: h = m->dl;
! 391: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
! 392: i++;
! 393: mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
! 394: addmd(CO,mod,s,t,&u); s = u;
! 395: }
! 396: *rp = s;
! 397: }
! 398:
! 399: void dp_ptozp(p,rp)
! 400: DP p,*rp;
! 401: {
! 402: MP m,mr,mr0;
! 403: int i,n;
! 404: Q *w;
! 405: Q dvr;
! 406: P t;
! 407:
! 408: if ( !p )
! 409: *rp = 0;
! 410: else {
! 411: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
! 412: w = (Q *)ALLOCA(n*sizeof(Q));
! 413: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 414: if ( NUM(m->c) )
! 415: w[i] = (Q)m->c;
! 416: else
! 417: ptozp(m->c,1,&w[i],&t);
! 418: sortbynm(w,n);
! 419: qltozl(w,n,&dvr);
! 420: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 421: NEXTMP(mr0,mr); divsp(CO,m->c,(P)dvr,&mr->c); mr->dl = m->dl;
! 422: }
! 423: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 424: }
! 425: }
! 426:
! 427: void dp_ptozp2(p0,p1,hp,rp)
! 428: DP p0,p1;
! 429: DP *hp,*rp;
! 430: {
! 431: DP t,s,h,r;
! 432: MP m,mr,mr0,m0;
! 433:
! 434: addd(CO,p0,p1,&t); dp_ptozp(t,&s);
! 435: if ( !p0 ) {
! 436: h = 0; r = s;
! 437: } else if ( !p1 ) {
! 438: h = s; r = 0;
! 439: } else {
! 440: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
! 441: m = NEXT(m), m0 = NEXT(m0) ) {
! 442: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
! 443: }
! 444: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
! 445: }
! 446: if ( h )
! 447: h->sugar = p0->sugar;
! 448: if ( r )
! 449: r->sugar = p1->sugar;
! 450: *hp = h; *rp = r;
! 451: }
! 452:
! 453: struct igcdv_hist {
! 454: int len;
! 455: DP dp;
! 456: Q d,d1;
! 457: Q *q,*q1;
! 458: };
! 459:
! 460: static struct igcdv_hist current_igcdv_hist;
! 461:
! 462: void dp_igcdv_hist(p,rp)
! 463: DP p;
! 464: Q *rp;
! 465: {
! 466: VECT c;
! 467: Q *a,*q,*r,*q1,*r1;
! 468: int i,n,nz;
! 469: Q d,d1;
! 470: N qn,rn,gn;
! 471: struct oVECT v;
! 472:
! 473: dp_dtov(p,&c); a = (Q *)c->body;
! 474: n = c->len;
! 475: q = (Q *)MALLOC(n*sizeof(Q)); r = (Q *)MALLOC(n*sizeof(Q));
! 476: igcdv_estimate(c,&d);
! 477: for ( i = 0, nz = 0; i < n; i++ ) {
! 478: divn(NM(a[i]),NM(d),&qn,&rn);
! 479: NTOQ(qn,SGN(a[i]),q[i]); NTOQ(rn,SGN(a[i]),r[i]);
! 480: if ( r[i] )
! 481: nz = 1;
! 482: }
! 483: current_igcdv_hist.len = n;
! 484: current_igcdv_hist.dp = p;
! 485: current_igcdv_hist.d = d;
! 486: current_igcdv_hist.q = q;
! 487: if ( nz ) {
! 488: v.id = O_VECT; v.len = n; v.body = (pointer *)r;
! 489: igcdv(&v,&d1);
! 490: q1 = (Q *)MALLOC(n*sizeof(Q)); r1 = (Q *)MALLOC(n*sizeof(Q));
! 491: for ( i = 0; i < n; i++ ) {
! 492: if ( r[i] )
! 493: divn(NM(r[i]),NM(d1),&qn,&rn);
! 494: else {
! 495: qn = rn = 0;
! 496: }
! 497: NTOQ(qn,SGN(r[i]),q1[i]); NTOQ(rn,SGN(r[i]),r1[i]);
! 498: }
! 499: gcdn(NM(d),NM(d1),&gn); NTOQ(gn,1,*rp);
! 500: } else {
! 501: d1 = 0; q1 = 0; *rp = d;
! 502: }
! 503: current_igcdv_hist.d1 = d1;
! 504: current_igcdv_hist.q1 = q1;
! 505: }
! 506:
! 507: void dp_idivv_hist(g,rp)
! 508: Q g;
! 509: DP *rp;
! 510: {
! 511: int i,n;
! 512: Q *q,*q1,*s;
! 513: Q t,u,m;
! 514: N qn;
! 515:
! 516: n = current_igcdv_hist.len;
! 517: q = current_igcdv_hist.q; q1 = current_igcdv_hist.q1;
! 518: divsn(NM(current_igcdv_hist.d),NM(g),&qn); NTOQ(qn,1,m);
! 519: s = (Q *)MALLOC(n*sizeof(Q));
! 520: for ( i = 0; i < n; i++ )
! 521: mulq(m,q[i],&s[i]);
! 522: if ( q1 ) {
! 523: divsn(NM(current_igcdv_hist.d1),NM(g),&qn); NTOQ(qn,1,m);
! 524: for ( i = 0; i < n; i++ ) {
! 525: mulq(m,q1[i],&t); addq(t,s[i],&u); s[i] = u;
! 526: }
! 527: }
! 528: dp_vtod(s,current_igcdv_hist.dp,rp);
! 529: }
! 530:
! 531: void dp_vtod(c,p,rp)
! 532: Q *c;
! 533: DP p;
! 534: DP *rp;
! 535: {
! 536: MP mr0,m,mr;
! 537: int i;
! 538:
! 539: if ( !p )
! 540: *rp = 0;
! 541: else {
! 542: for ( mr0 = 0, m = BDY(p), i = 0; m; m = NEXT(m), i++ ) {
! 543: NEXTMP(mr0,mr); mr->c = (P)c[i]; mr->dl = m->dl;
! 544: }
! 545: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
! 546: (*rp)->sugar = p->sugar;
! 547: }
! 548: }
! 549:
! 550: #if INET
! 551: void dp_ptozp_d_old(dist,ndist,p,rp)
! 552: NODE dist;
! 553: int ndist;
! 554: DP p,*rp;
! 555: {
! 556: DP d,s,u,su;
! 557: MP m,mr,mr0;
! 558: MP *w0,*w;
! 559: Q *g;
! 560: int i,n,nv,nsep;
! 561: P t;
! 562: NODE tn,n0,n1,n2;
! 563: struct oVECT v;
! 564: int id,sindex;
! 565: Obj dmy;
! 566: Q gcd;
! 567: STRING cont,div;
! 568:
! 569: if ( !p )
! 570: *rp = 0;
! 571: else {
! 572: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 573: nsep = ndist + 1;
! 574: if ( n <= nsep ) {
! 575: dp_ptozp(p,rp); return;
! 576: }
! 577: nv = p->nv;
! 578: w0 = (MP *)MALLOC(nsep*sizeof(MP)); bzero(w0,nsep*sizeof(MP));
! 579: w = (MP *)MALLOC(nsep*sizeof(MP)); bzero(w,nsep*sizeof(MP));
! 580: g = (Q *)MALLOC(nsep*sizeof(Q)); bzero(g,nsep*sizeof(Q));
! 581: for ( m = BDY(p), i = 0; m; m = NEXT(m), i++, i %= nsep ) {
! 582: NEXTMP(w0[i],w[i]); w[i]->c = m->c; w[i]->dl = m->dl;
! 583: }
! 584: MKSTR(cont,"reg_dp_cont"); MKSTR(div,"reg_dp_idiv");
! 585:
! 586: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 587: NEXT(w[i]) = 0; MKDP(nv,w0[i],d); d->sugar = p->sugar;
! 588: MKNODE(n2,d,0); MKNODE(n1,cont,n2); MKNODE(n0,BDY(tn),n1);
! 589: Pox_rpc(n0,&dmy);
! 590: }
! 591: NEXT(w[i]) = 0; MKDP(nv,w0[i],d); d->sugar = p->sugar;
! 592: dp_cont(d,&g[i]);
! 593:
! 594: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) )
! 595: Pox_pop_local(tn,&g[i]);
! 596: v.id = O_VECT; v.len = nsep; v.body = (pointer *)g; igcdv(&v,&gcd);
! 597:
! 598: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 599: MKNODE(n2,gcd,0); MKNODE(n1,div,n2); MKNODE(n0,BDY(tn),n1);
! 600: Pox_rpc(n0,&dmy);
! 601: }
! 602: dp_idiv(d,gcd,&s);
! 603:
! 604: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 605: Pox_pop_local(tn,&u);
! 606: addd(CO,s,u,&su); s = su;
! 607: }
! 608: *rp = s;
! 609: }
! 610: }
! 611:
! 612: void dp_ptozp_d(dist,ndist,p,rp)
! 613: NODE dist;
! 614: int ndist;
! 615: DP p,*rp;
! 616: {
! 617: int i,j,k,l,n,nsep;
! 618: MP m;
! 619: NODE tn,n0,n1,n2,n3;
! 620: struct oVECT v;
! 621: VECT c,cs;
! 622: VECT qi,ri;
! 623: LIST *qr;
! 624: int s,id;
! 625: Obj dmy;
! 626: Q d0,d1,gcd,a,u,u1;
! 627: Q *q,*r;
! 628: STRING iqr_v;
! 629: pointer *b;
! 630: N qn,gn;
! 631: double get_rtime();
! 632: int blen;
! 633: double t0;
! 634: double t_e,t_d,t_d1,t_c;
! 635:
! 636: if ( !p )
! 637: *rp = 0;
! 638: else {
! 639: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 640: nsep = ndist + 1;
! 641: if ( n <= nsep ) {
! 642: dp_ptozp(p,rp); return;
! 643: }
! 644: t0 = get_rtime();
! 645: dp_dtov(p,&c);
! 646: igcdv_estimate(c,&d0);
! 647: t_e = get_rtime()-t0;
! 648: t0 = get_rtime();
! 649: dp_dtov(p,&c);
! 650: sepvect(c,nsep,&cs);
! 651: MKSTR(iqr_v,"iqr");
! 652: qr = (LIST *)CALLOC(nsep,sizeof(LIST));
! 653: q = (Q *)CALLOC(n,sizeof(Q));
! 654: r = (Q *)CALLOC(n,sizeof(Q));
! 655: for ( i = 0, tn = dist, b = BDY(cs); i < ndist; i++, tn = NEXT(tn) ) {
! 656: MKNODE(n3,d0,0); MKNODE(n2,b[i],n3);
! 657: MKNODE(n1,iqr_v,n2); MKNODE(n0,BDY(tn),n1);
! 658: Pox_rpc(n0,&dmy);
! 659: }
! 660: iqrv(b[i],d0,&qr[i]);
! 661: dp_dtov(p,&c);
! 662: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 663: Pox_pop_local(tn,&qr[i]);
! 664: if ( OID(qr[i]) == O_ERR ) {
! 665: printexpr(CO,(Obj)qr[i]);
! 666: error("dp_ptozp_d : aborted");
! 667: }
! 668: }
! 669: t_d = get_rtime()-t0;
! 670: t_d1 = t_d/n;
! 671: t0 = get_rtime();
! 672: for ( i = j = 0; i < nsep; i++ ) {
! 673: tn = BDY(qr[i]); qi = (VECT)BDY(tn); ri = (VECT)BDY(NEXT(tn));
! 674: for ( k = 0, l = qi->len; k < l; k++, j++ ) {
! 675: q[j] = (Q)BDY(qi)[k]; r[j] = (Q)BDY(ri)[k];
! 676: }
! 677: }
! 678: v.id = O_VECT; v.len = n; v.body = (pointer *)r; igcdv(&v,&d1);
! 679: if ( d1 ) {
! 680: gcdn(NM(d0),NM(d1),&gn); NTOQ(gn,1,gcd);
! 681: divsn(NM(d0),gn,&qn); NTOQ(qn,1,a);
! 682: for ( i = 0; i < n; i++ ) {
! 683: mulq(a,q[i],&u);
! 684: if ( r[i] ) {
! 685: divsn(NM(r[i]),gn,&qn); NTOQ(qn,SGN(r[i]),u1);
! 686: addq(u,u1,&q[i]);
! 687: } else
! 688: q[i] = u;
! 689: }
! 690: } else
! 691: gcd = d0;
! 692: dp_vtod(q,p,rp);
! 693: t_c = get_rtime()-t0;
! 694: blen=p_mag((P)gcd);
! 695: pz_t_e += t_e; pz_t_d += t_d; pz_t_d1 += t_d1; pz_t_c += t_c;
! 696: if ( 0 )
! 697: fprintf(stderr,"(%d,%d)",p_mag((P)d0)-blen,blen);
! 698: }
! 699: }
! 700:
! 701: void dp_ptozp2_d(dist,ndist,p0,p1,hp,rp)
! 702: NODE dist;
! 703: int ndist;
! 704: DP p0,p1;
! 705: DP *hp,*rp;
! 706: {
! 707: DP t,s,h,r;
! 708: MP m,mr,mr0,m0;
! 709:
! 710: addd(CO,p0,p1,&t); dp_ptozp_d(dist,ndist,t,&s);
! 711: if ( !p0 ) {
! 712: h = 0; r = s;
! 713: } else if ( !p1 ) {
! 714: h = s; r = 0;
! 715: } else {
! 716: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
! 717: m = NEXT(m), m0 = NEXT(m0) ) {
! 718: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
! 719: }
! 720: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
! 721: }
! 722: if ( h )
! 723: h->sugar = p0->sugar;
! 724: if ( r )
! 725: r->sugar = p1->sugar;
! 726: *hp = h; *rp = r;
! 727: }
! 728: #endif
! 729:
! 730: /*
! 731: * Old codes
! 732: */
! 733:
! 734: #if 0
! 735: void Pdp_sep(arg,rp)
! 736: NODE arg;
! 737: VECT *rp;
! 738: {
! 739: DP p,s;
! 740: MP m,t,mr,mr0;
! 741: int i,j,n,d,r,q,q1,nv;
! 742: VECT w;
! 743: pointer *pw;
! 744:
! 745: p = (DP)ARG0(arg); m = BDY(p);
! 746: d = QTOS((Q)ARG1(arg));
! 747: for ( t = m, n = 0; t; t = NEXT(t), n++ );
! 748: if ( d > n )
! 749: d = n;
! 750: q = n/d; r = n%d; q1 = q+1;
! 751: MKVECT(w,d); *rp = w;
! 752: t = m; pw = BDY(w); nv = p->nv;
! 753: for ( i = 0; i < r; i++ ) {
! 754: for ( mr0 = 0, j = 0; j < q1; j++, t = NEXT(t) ) {
! 755: NEXTMP(mr0,mr);
! 756: mr->c = t->c;
! 757: mr->dl = t->dl;
! 758: }
! 759: NEXT(mr) = 0; MKDP(nv,mr0,s); s->sugar = p->sugar;
! 760: pw[i] = (pointer)s;
! 761: }
! 762: for ( ; i < d; i++ ) {
! 763: for ( mr0 = 0, j = 0; j < q; j++, t = NEXT(t) ) {
! 764: NEXTMP(mr0,mr);
! 765: mr->c = t->c;
! 766: mr->dl = t->dl;
! 767: }
! 768: NEXT(mr) = 0; MKDP(nv,mr0,s); s->sugar = p->sugar;
! 769: pw[i] = (pointer)s;
! 770: }
! 771: }
! 772:
! 773: void dp_ptozp_mpi(dist,ndist,p,rp)
! 774: NODE dist;
! 775: int ndist;
! 776: DP p,*rp;
! 777: {
! 778: int i,j,k,l,n,nsep;
! 779: MP m;
! 780: NODE tn,n0,n1,n2,n3;
! 781: struct oVECT v;
! 782: VECT c,cs;
! 783: VECT qi,ri;
! 784: LIST *qr;
! 785: int id;
! 786: Obj dmy;
! 787: Q d0,d1,gcd,a,u,u1;
! 788: Q *q,*r;
! 789: STRING iqr_v;
! 790: pointer *b;
! 791: N qn,gn;
! 792: double get_rtime();
! 793: int blen;
! 794: double t0;
! 795: double t_e,t_d,t_d1,t_c;
! 796:
! 797: if ( !p )
! 798: *rp = 0;
! 799: else {
! 800: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 801: nsep = ndist + 1;
! 802: if ( n <= nsep ) {
! 803: dp_ptozp(p,rp); return;
! 804: }
! 805: t0 = get_rtime();
! 806: dp_dtov(p,&c);
! 807: igcdv_estimate(c,&d0);
! 808: t_e = get_rtime()-t0;
! 809: t0 = get_rtime();
! 810: dp_dtov(p,&c);
! 811: sepvect(c,nsep,&cs);
! 812: MKSTR(iqr_v,"iqr");
! 813: qr = (LIST *)CALLOC(nsep,sizeof(LIST));
! 814: q = (Q *)CALLOC(n,sizeof(Q));
! 815: r = (Q *)CALLOC(n,sizeof(Q));
! 816: for ( i = 0, tn = dist, b = BDY(cs); i < ndist;
! 817: i++, tn = NEXT(tn) ) {
! 818: div_vect_mpi(QTOS((Q)BDY(tn)),b[i],d0);
! 819: }
! 820: iqrv(b[i],d0,&qr[i]);
! 821: dp_dtov(p,&c);
! 822: for ( i = j = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 823: get_div_vect_mpi(QTOS((Q)BDY(tn)),&qr[i]);
! 824: }
! 825: t_d = get_rtime()-t0;
! 826: t_d1 = t_d/n;
! 827: t0 = get_rtime();
! 828: for ( i = j = 0; i < nsep; i++ ) {
! 829: tn = BDY(qr[i]); qi = (VECT)BDY(tn); ri = (VECT)BDY(NEXT(tn));
! 830: for ( k = 0, l = qi->len; k < l; k++, j++ ) {
! 831: q[j] = (Q)BDY(qi)[k]; r[j] = (Q)BDY(ri)[k];
! 832: }
! 833: }
! 834: v.id = O_VECT; v.len = n; v.body = (pointer *)r; igcdv(&v,&d1);
! 835: if ( d1 ) {
! 836: gcdn(NM(d0),NM(d1),&gn); NTOQ(gn,1,gcd);
! 837: divsn(NM(d0),gn,&qn); NTOQ(qn,1,a);
! 838: for ( i = 0; i < n; i++ ) {
! 839: mulq(a,q[i],&u);
! 840: if ( r[i] ) {
! 841: divsn(NM(r[i]),gn,&qn); NTOQ(qn,SGN(r[i]),u1);
! 842: addq(u,u1,&q[i]);
! 843: } else
! 844: q[i] = u;
! 845: }
! 846: } else
! 847: gcd = d0;
! 848: dp_vtod(q,p,rp);
! 849: t_c = get_rtime()-t0;
! 850: blen=p_mag(gcd);
! 851: pz_t_e += t_e; pz_t_d += t_d; pz_t_d1 += t_d1; pz_t_c += t_c;
! 852: if ( 0 )
! 853: fprintf(stderr,"(%d,%d)",p_mag(d0)-blen,blen);
! 854: }
! 855: }
! 856:
! 857: void dp_ptozp_d(dist,ndist,p,rp)
! 858: NODE dist;
! 859: int ndist;
! 860: DP p,*rp;
! 861: {
! 862: DP d,s,u,su;
! 863: MP m,mr,mr0;
! 864: MP *w0,*w;
! 865: Q *g;
! 866: int i,n,nv,nsep;
! 867: P t;
! 868: NODE tn,n0,n1,n2;
! 869: struct oVECT v;
! 870: int sindex,id;
! 871: Obj dmy;
! 872: Q gcd;
! 873: STRING gcd_h,div_h;
! 874:
! 875: if ( !p )
! 876: *rp = 0;
! 877: else {
! 878: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 879: nsep = ndist + 1;
! 880: if ( n <= nsep ) {
! 881: dp_ptozp(p,rp); return;
! 882: }
! 883: nv = p->nv;
! 884: w0 = (MP *)MALLOC(nsep*sizeof(MP)); bzero(w0,nsep*sizeof(MP));
! 885: w = (MP *)MALLOC(nsep*sizeof(MP)); bzero(w,nsep*sizeof(MP));
! 886: g = (Q *)MALLOC(nsep*sizeof(Q)); bzero(g,nsep*sizeof(Q));
! 887: for ( m = BDY(p), i = 0; m; m = NEXT(m), i++, i %= nsep ) {
! 888: NEXTMP(w0[i],w[i]); w[i]->c = m->c; w[i]->dl = m->dl;
! 889: }
! 890: MKSTR(gcd_h,"dp_igcdv_hist"); MKSTR(div_h,"dp_idivv_hist");
! 891:
! 892: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 893: NEXT(w[i]) = 0; MKDP(nv,w0[i],d); d->sugar = p->sugar;
! 894: MKNODE(n2,d,0); MKNODE(n1,gcd_h,n2); MKNODE(n0,BDY(tn),n1);
! 895: Pox_rpc(n0,&dmy);
! 896: }
! 897: NEXT(w[i]) = 0; MKDP(nv,w0[i],d); d->sugar = p->sugar;
! 898: dp_igcdv_hist(d,&g[i]);
! 899:
! 900: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) )
! 901: Pox_pop_local(tn,&g[i]);
! 902: v.id = O_VECT; v.len = nsep; v.body = (pointer *)g; igcdv(&v,&gcd);
! 903:
! 904: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 905: MKNODE(n2,gcd,0); MKNODE(n1,div_h,n2); MKNODE(n0,BDY(tn),n1);
! 906: Pox_rpc(n0,&dmy);
! 907: }
! 908: dp_idivv_hist(gcd,&s);
! 909:
! 910: for ( i = 0, tn = dist; i < ndist; i++, tn = NEXT(tn) ) {
! 911: Pox_pop_local(tn,&u);
! 912: addd(CO,s,u,&su); s = su;
! 913: }
! 914: *rp = s;
! 915: }
! 916: }
! 917: #endif
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>