Annotation of OpenXM_contrib2/asir2000/engine/distm.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/engine/distm.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
! 2: #include "ca.h"
! 3: #include "inline.h"
! 4:
! 5: #define NV(p) ((p)->nv)
! 6: #define C(p) ((p)->c)
! 7: #if 0
! 8: #define ITOS(p) (((unsigned int)(p))>>1)
! 9: #define STOI(i) ((P)((((unsigned int)(i))<<1)|1))
! 10: #else
! 11: #define ITOS(p) (((unsigned int)(p))&0x7fffffff)
! 12: #define STOI(i) ((P)((unsigned int)(i)|0x80000000))
! 13: #endif
! 14:
! 15: extern int (*cmpdl)();
! 16:
! 17: void ptomd(vl,mod,dvl,p,pr)
! 18: VL vl,dvl;
! 19: int mod;
! 20: P p;
! 21: DP *pr;
! 22: {
! 23: P t;
! 24:
! 25: ptomp(mod,p,&t);
! 26: mptomd(vl,mod,dvl,t,pr);
! 27: }
! 28:
! 29: void mptomd(vl,mod,dvl,p,pr)
! 30: VL vl,dvl;
! 31: int mod;
! 32: P p;
! 33: DP *pr;
! 34: {
! 35: int n,i;
! 36: VL tvl;
! 37: V v;
! 38: DL d;
! 39: MP m;
! 40: DCP dc;
! 41: DP r,s,t,u;
! 42: P x,c;
! 43:
! 44: if ( !p )
! 45: *pr = 0;
! 46: else {
! 47: for ( n = 0, tvl = dvl; tvl; tvl = NEXT(tvl), n++ );
! 48: if ( NUM(p) ) {
! 49: NEWDL(d,n);
! 50: NEWMP(m); m->dl = d; C(m) = p; NEXT(m) = 0; MKDP(n,m,*pr);
! 51: } else {
! 52: for ( i = 0, tvl = dvl, v = VR(p);
! 53: tvl && tvl->v != v; tvl = NEXT(tvl), i++ );
! 54: if ( !tvl ) {
! 55: for ( dc = DC(p), s = 0, MKV(v,x); dc; dc = NEXT(dc) ) {
! 56: mptomd(vl,mod,dvl,COEF(dc),&t); pwrmp(vl,mod,x,DEG(dc),&c);
! 57: mulmdc(vl,mod,t,c,&r); addmd(vl,mod,r,s,&t); s = t;
! 58: }
! 59: *pr = s;
! 60: } else {
! 61: for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) {
! 62: mptomd(vl,mod,dvl,COEF(dc),&t);
! 63: NEWDL(d,n); d->td = QTOS(DEG(dc)); d->d[i] = d->td;
! 64: NEWMP(m); m->dl = d; C(m) = (P)ONEM; NEXT(m) = 0; MKDP(n,m,u);
! 65: mulmd(vl,mod,t,u,&r); addmd(vl,mod,r,s,&t); s = t;
! 66: }
! 67: *pr = s;
! 68: }
! 69: }
! 70: }
! 71: }
! 72:
! 73: void mdtop(vl,mod,dvl,p,pr)
! 74: VL vl,dvl;
! 75: int mod;
! 76: DP p;
! 77: P *pr;
! 78: {
! 79: int n,i;
! 80: DL d;
! 81: MP m;
! 82: P r,s,t,u,w;
! 83: Q q;
! 84: VL tvl;
! 85:
! 86: if ( !p )
! 87: *pr = 0;
! 88: else {
! 89: for ( n = p->nv, m = BDY(p), s = 0; m; m = NEXT(m) ) {
! 90: for ( i = 0, t = C(m), d = m->dl, tvl = dvl;
! 91: i < n; tvl = NEXT(tvl), i++ ) {
! 92: MKMV(tvl->v,r); STOQ(d->d[i],q); pwrmp(vl,mod,r,q,&u);
! 93: mulmp(vl,mod,t,u,&w); t = w;
! 94: }
! 95: addmp(vl,mod,s,t,&u); s = u;
! 96: }
! 97: mptop(s,pr);
! 98: }
! 99: }
! 100:
! 101: void addmd(vl,mod,p1,p2,pr)
! 102: VL vl;
! 103: int mod;
! 104: DP p1,p2,*pr;
! 105: {
! 106: int n;
! 107: MP m1,m2,mr,mr0;
! 108: P t;
! 109: int tmp;
! 110: MQ q;
! 111:
! 112: if ( !p1 )
! 113: *pr = p2;
! 114: else if ( !p2 )
! 115: *pr = p1;
! 116: else {
! 117: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; )
! 118: switch ( (*cmpdl)(n,m1->dl,m2->dl) ) {
! 119: case 0:
! 120: if ( NUM(C(m1)) && NUM(C(m2)) ) {
! 121: tmp = (CONT((MQ)C(m1))+CONT((MQ)C(m2))) - mod;
! 122: if ( tmp < 0 )
! 123: tmp += mod;
! 124: if ( tmp ) {
! 125: STOMQ(tmp,q); t = (P)q;
! 126: } else
! 127: t = 0;
! 128: } else
! 129: addmp(vl,mod,C(m1),C(m2),&t);
! 130: if ( t ) {
! 131: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = t;
! 132: }
! 133: m1 = NEXT(m1); m2 = NEXT(m2); break;
! 134: case 1:
! 135: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = C(m1);
! 136: m1 = NEXT(m1); break;
! 137: case -1:
! 138: NEXTMP(mr0,mr); mr->dl = m2->dl; C(mr) = C(m2);
! 139: m2 = NEXT(m2); break;
! 140: }
! 141: if ( !mr0 )
! 142: if ( m1 )
! 143: mr0 = m1;
! 144: else if ( m2 )
! 145: mr0 = m2;
! 146: else {
! 147: *pr = 0;
! 148: return;
! 149: }
! 150: else if ( m1 )
! 151: NEXT(mr) = m1;
! 152: else if ( m2 )
! 153: NEXT(mr) = m2;
! 154: else
! 155: NEXT(mr) = 0;
! 156: MKDP(NV(p1),mr0,*pr);
! 157: if ( *pr )
! 158: (*pr)->sugar = MAX(p1->sugar,p2->sugar);
! 159: }
! 160: }
! 161:
! 162: void submd(vl,mod,p1,p2,pr)
! 163: VL vl;
! 164: int mod;
! 165: DP p1,p2,*pr;
! 166: {
! 167: DP t;
! 168:
! 169: if ( !p2 )
! 170: *pr = p1;
! 171: else {
! 172: chsgnmd(mod,p2,&t); addmd(vl,mod,p1,t,pr);
! 173: }
! 174: }
! 175:
! 176: void chsgnmd(mod,p,pr)
! 177: int mod;
! 178: DP p,*pr;
! 179: {
! 180: MP m,mr,mr0;
! 181:
! 182: if ( !p )
! 183: *pr = 0;
! 184: else {
! 185: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 186: NEXTMP(mr0,mr); chsgnmp(mod,C(m),&C(mr)); mr->dl = m->dl;
! 187: }
! 188: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
! 189: if ( *pr )
! 190: (*pr)->sugar = p->sugar;
! 191: }
! 192: }
! 193:
! 194: void mulmd(vl,mod,p1,p2,pr)
! 195: VL vl;
! 196: int mod;
! 197: DP p1,p2,*pr;
! 198: {
! 199: MP m;
! 200: DP s,t,u;
! 201:
! 202: if ( !p1 || !p2 )
! 203: *pr = 0;
! 204: else if ( OID(p1) <= O_P )
! 205: mulmdc(vl,mod,p2,(P)p1,pr);
! 206: else if ( OID(p2) <= O_P )
! 207: mulmdc(vl,mod,p1,(P)p2,pr);
! 208: else {
! 209: for ( m = BDY(p2), s = 0; m; m = NEXT(m) ) {
! 210: mulmdm(vl,mod,p1,m,&t); addmd(vl,mod,s,t,&u); s = u;
! 211: }
! 212: *pr = s;
! 213: }
! 214: }
! 215:
! 216: void mulmdm(vl,mod,p,m0,pr)
! 217: VL vl;
! 218: int mod;
! 219: DP p;
! 220: MP m0;
! 221: DP *pr;
! 222: {
! 223: MP m,mr,mr0;
! 224: P c;
! 225: DL d;
! 226: int n,t;
! 227: MQ q;
! 228:
! 229: if ( !p )
! 230: *pr = 0;
! 231: else {
! 232: for ( mr0 = 0, m = BDY(p), c = C(m0), d = m0->dl, n = NV(p);
! 233: m; m = NEXT(m) ) {
! 234: NEXTMP(mr0,mr);
! 235: if ( NUM(C(m)) && NUM(c) ) {
! 236: t = dmar(((MQ)(C(m)))->cont,((MQ)c)->cont,0,mod);
! 237: STOMQ(t,q); C(mr) = (P)q;
! 238: } else
! 239: mulmp(vl,mod,C(m),c,&C(mr));
! 240: adddl(n,m->dl,d,&mr->dl);
! 241: }
! 242: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
! 243: if ( *pr )
! 244: (*pr)->sugar = p->sugar + m0->dl->td;
! 245: }
! 246: }
! 247:
! 248: void mulmdc(vl,mod,p,c,pr)
! 249: VL vl;
! 250: int mod;
! 251: DP p;
! 252: P c;
! 253: DP *pr;
! 254: {
! 255: MP m,mr,mr0;
! 256: int t;
! 257: MQ q;
! 258:
! 259: if ( !p || !c )
! 260: *pr = 0;
! 261: else {
! 262: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 263: NEXTMP(mr0,mr);
! 264: if ( NUM(C(m)) && NUM(c) ) {
! 265: t = dmar(((MQ)(C(m)))->cont,((MQ)c)->cont,0,mod);
! 266: STOMQ(t,q); C(mr) = (P)q;
! 267: } else
! 268: mulmp(vl,mod,C(m),c,&C(mr));
! 269: mr->dl = m->dl;
! 270: }
! 271: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
! 272: if ( *pr )
! 273: (*pr)->sugar = p->sugar;
! 274: }
! 275: }
! 276:
! 277: void divsmdc(vl,mod,p,c,pr)
! 278: VL vl;
! 279: int mod;
! 280: DP p;
! 281: P c;
! 282: DP *pr;
! 283: {
! 284: MP m,mr,mr0;
! 285:
! 286: if ( !c )
! 287: error("disvsdc : division by 0");
! 288: else if ( !p )
! 289: *pr = 0;
! 290: else {
! 291: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 292: NEXTMP(mr0,mr); divsmp(vl,mod,C(m),c,&C(mr)); mr->dl = m->dl;
! 293: }
! 294: NEXT(mr) = 0; MKDP(NV(p),mr0,*pr);
! 295: if ( *pr )
! 296: (*pr)->sugar = p->sugar;
! 297: }
! 298: }
! 299:
! 300: #define MKDPM(n,m,d) (NEWDP(d),(d)->nv=(n),BDY(d)=(m))
! 301:
! 302: void _mdtop(vl,mod,dvl,p,pr)
! 303: VL vl,dvl;
! 304: int mod;
! 305: DP p;
! 306: P *pr;
! 307: {
! 308: int n,i;
! 309: DL d;
! 310: MP m;
! 311: P r,s,t,u,w;
! 312: Q q;
! 313: MQ tq;
! 314: VL tvl;
! 315:
! 316: if ( !p )
! 317: *pr = 0;
! 318: else {
! 319: for ( n = p->nv, m = BDY(p), s = 0; m; m = NEXT(m) ) {
! 320: STOMQ(ITOS(C(m)),tq); t = (P)tq;
! 321: for ( i = 0, d = m->dl, tvl = dvl;
! 322: i < n; tvl = NEXT(tvl), i++ ) {
! 323: MKV(tvl->v,r); STOQ(d->d[i],q); pwrmp(vl,mod,r,q,&u);
! 324: mulmp(vl,mod,t,u,&w); t = w;
! 325: }
! 326: addmp(vl,mod,s,t,&u); s = u;
! 327: }
! 328: mptop(s,pr);
! 329: }
! 330: }
! 331:
! 332: void _addmd(vl,mod,p1,p2,pr)
! 333: VL vl;
! 334: int mod;
! 335: DP p1,p2,*pr;
! 336: {
! 337: int n;
! 338: MP m1,m2,mr,mr0;
! 339: int t;
! 340:
! 341: if ( !p1 )
! 342: *pr = p2;
! 343: else if ( !p2 )
! 344: *pr = p1;
! 345: else {
! 346: for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; )
! 347: switch ( (*cmpdl)(n,m1->dl,m2->dl) ) {
! 348: case 0:
! 349: t = (ITOS(C(m1))+ITOS(C(m2))) - mod;
! 350: if ( t < 0 )
! 351: t += mod;
! 352: if ( t ) {
! 353: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = STOI(t);
! 354: }
! 355: m1 = NEXT(m1); m2 = NEXT(m2); break;
! 356: case 1:
! 357: NEXTMP(mr0,mr); mr->dl = m1->dl; C(mr) = C(m1);
! 358: m1 = NEXT(m1); break;
! 359: case -1:
! 360: NEXTMP(mr0,mr); mr->dl = m2->dl; C(mr) = C(m2);
! 361: m2 = NEXT(m2); break;
! 362: }
! 363: if ( !mr0 )
! 364: if ( m1 )
! 365: mr0 = m1;
! 366: else if ( m2 )
! 367: mr0 = m2;
! 368: else {
! 369: *pr = 0;
! 370: return;
! 371: }
! 372: else if ( m1 )
! 373: NEXT(mr) = m1;
! 374: else if ( m2 )
! 375: NEXT(mr) = m2;
! 376: else
! 377: NEXT(mr) = 0;
! 378: MKDPM(NV(p1),mr0,*pr);
! 379: if ( *pr )
! 380: (*pr)->sugar = MAX(p1->sugar,p2->sugar);
! 381: }
! 382: }
! 383:
! 384: void _submd(vl,mod,p1,p2,pr)
! 385: VL vl;
! 386: int mod;
! 387: DP p1,p2,*pr;
! 388: {
! 389: DP t;
! 390:
! 391: if ( !p2 )
! 392: *pr = p1;
! 393: else {
! 394: _chsgnmd(mod,p2,&t); _addmd(vl,mod,p1,t,pr);
! 395: }
! 396: }
! 397:
! 398: void _chsgnmd(mod,p,pr)
! 399: int mod;
! 400: DP p,*pr;
! 401: {
! 402: MP m,mr,mr0;
! 403:
! 404: if ( !p )
! 405: *pr = 0;
! 406: else {
! 407: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 408: NEXTMP(mr0,mr); C(mr) = STOI(mod - ITOS(C(m))); mr->dl = m->dl;
! 409: }
! 410: NEXT(mr) = 0; MKDPM(NV(p),mr0,*pr);
! 411: if ( *pr )
! 412: (*pr)->sugar = p->sugar;
! 413: }
! 414: }
! 415:
! 416: void _mulmd(vl,mod,p1,p2,pr)
! 417: VL vl;
! 418: int mod;
! 419: DP p1,p2,*pr;
! 420: {
! 421: MP m;
! 422: DP s,t,u;
! 423:
! 424: if ( !p1 || !p2 )
! 425: *pr = 0;
! 426: else {
! 427: for ( m = BDY(p2), s = 0; m; m = NEXT(m) ) {
! 428: _mulmdm(vl,mod,p1,m,&t); _addmd(vl,mod,s,t,&u); s = u;
! 429: }
! 430: *pr = s;
! 431: }
! 432: }
! 433:
! 434: void _mulmdm(vl,mod,p,m0,pr)
! 435: VL vl;
! 436: int mod;
! 437: DP p;
! 438: MP m0;
! 439: DP *pr;
! 440: {
! 441: MP m,mr,mr0;
! 442: DL d;
! 443: int c,n,r;
! 444:
! 445: if ( !p )
! 446: *pr = 0;
! 447: else {
! 448: for ( mr0 = 0, m = BDY(p), c = ITOS(C(m0)), d = m0->dl, n = NV(p);
! 449: m; m = NEXT(m) ) {
! 450: NEXTMP(mr0,mr);
! 451: C(mr) = STOI(dmar(ITOS(C(m)),c,0,mod));
! 452: adddl(n,m->dl,d,&mr->dl);
! 453: }
! 454: NEXT(mr) = 0; MKDPM(NV(p),mr0,*pr);
! 455: if ( *pr )
! 456: (*pr)->sugar = p->sugar + m0->dl->td;
! 457: }
! 458: }
! 459:
! 460: void _dtop_mod(vl,dvl,p,pr)
! 461: VL vl,dvl;
! 462: DP p;
! 463: P *pr;
! 464: {
! 465: int n,i;
! 466: DL d;
! 467: MP m;
! 468: P r,s,t,u,w;
! 469: Q q;
! 470: VL tvl;
! 471:
! 472: if ( !p )
! 473: *pr = 0;
! 474: else {
! 475: for ( n = p->nv, m = BDY(p), s = 0; m; m = NEXT(m) ) {
! 476: i = ITOS(m->c); STOQ(i,q); t = (P)q;
! 477: for ( i = 0, d = m->dl, tvl = dvl;
! 478: i < n; tvl = NEXT(tvl), i++ ) {
! 479: MKV(tvl->v,r); STOQ(d->d[i],q); pwrp(vl,r,q,&u);
! 480: mulp(vl,t,u,&w); t = w;
! 481: }
! 482: addp(vl,s,t,&u); s = u;
! 483: }
! 484: *pr = s;
! 485: }
! 486: }
! 487:
! 488: void _dp_red_mod(p1,p2,mod,rp)
! 489: DP p1,p2;
! 490: int mod;
! 491: DP *rp;
! 492: {
! 493: int i,n;
! 494: DL d1,d2,d;
! 495: MP m;
! 496: DP t,s;
! 497: int c,c1;
! 498:
! 499: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 500: NEWDL(d,n); d->td = d1->td - d2->td;
! 501: for ( i = 0; i < n; i++ )
! 502: d->d[i] = d1->d[i]-d2->d[i];
! 503: c = invm(ITOS(BDY(p2)->c),mod); c1 = dmar(c,ITOS(BDY(p1)->c),0,mod);
! 504: NEWMP(m); m->dl = d; m->c = STOI(mod-c1); NEXT(m) = 0;
! 505: MKDP(n,m,s); s->sugar = d->td;
! 506: _mulmd(CO,mod,p2,s,&t); _addmd(CO,mod,p1,t,rp);
! 507: }
! 508:
! 509: void _dp_mod(p,mod,subst,rp)
! 510: DP p;
! 511: int mod;
! 512: NODE subst;
! 513: DP *rp;
! 514: {
! 515: MP m,mr,mr0;
! 516: P t,s;
! 517: NODE tn;
! 518:
! 519: if ( !p )
! 520: *rp = 0;
! 521: else {
! 522: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 523: for ( tn = subst, s = m->c; tn; tn = NEXT(NEXT(tn)) ) {
! 524: substp(CO,s,BDY(tn),BDY(NEXT(tn)),&t);
! 525: s = t;
! 526: }
! 527: ptomp(mod,s,&t);
! 528: if ( t ) {
! 529: NEXTMP(mr0,mr); mr->c = STOI(CONT((MQ)t)); mr->dl = m->dl;
! 530: }
! 531: }
! 532: if ( mr0 ) {
! 533: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 534: } else
! 535: *rp = 0;
! 536: }
! 537: }
! 538:
! 539: void _dp_sp_mod(p1,p2,mod,rp)
! 540: DP p1,p2;
! 541: int mod;
! 542: DP *rp;
! 543: {
! 544: int i,n,td;
! 545: int *w;
! 546: DL d1,d2,d;
! 547: MP m;
! 548: DP t,s,u;
! 549:
! 550: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 551: w = (int *)ALLOCA(n*sizeof(int));
! 552: for ( i = 0, td = 0; i < n; i++ ) {
! 553: w[i] = MAX(d1->d[i],d2->d[i]); td += w[i];
! 554: }
! 555: NEWDL(d,n); d->td = td - d1->td;
! 556: for ( i = 0; i < n; i++ )
! 557: d->d[i] = w[i] - d1->d[i];
! 558: NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
! 559: MKDP(n,m,s); s->sugar = d->td; _mulmd(CO,mod,p1,s,&t);
! 560: NEWDL(d,n); d->td = td - d2->td;
! 561: for ( i = 0; i < n; i++ )
! 562: d->d[i] = w[i] - d2->d[i];
! 563: NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
! 564: MKDP(n,m,s); s->sugar = d->td; _mulmd(CO,mod,p2,s,&u);
! 565: _addmd(CO,mod,t,u,rp);
! 566: }
! 567:
! 568: void _dp_sp_component_mod(p1,p2,mod,f1,f2)
! 569: DP p1,p2;
! 570: int mod;
! 571: DP *f1,*f2;
! 572: {
! 573: int i,n,td;
! 574: int *w;
! 575: DL d1,d2,d;
! 576: MP m;
! 577: DP t,s,u;
! 578:
! 579: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 580: w = (int *)ALLOCA(n*sizeof(int));
! 581: for ( i = 0, td = 0; i < n; i++ ) {
! 582: w[i] = MAX(d1->d[i],d2->d[i]); td += w[i];
! 583: }
! 584: NEWDL(d,n); d->td = td - d1->td;
! 585: for ( i = 0; i < n; i++ )
! 586: d->d[i] = w[i] - d1->d[i];
! 587: NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
! 588: MKDP(n,m,s); s->sugar = d->td; _mulmd(CO,mod,p1,s,f1);
! 589: NEWDL(d,n); d->td = td - d2->td;
! 590: for ( i = 0; i < n; i++ )
! 591: d->d[i] = w[i] - d2->d[i];
! 592: NEWMP(m); m->dl = d; m->c = BDY(p1)->c; NEXT(m) = 0;
! 593: MKDP(n,m,s); s->sugar = d->td; _mulmd(CO,mod,p2,s,f2);
! 594: }
! 595:
! 596: void _printdp(d)
! 597: DP d;
! 598: {
! 599: int n,i;
! 600: MP m;
! 601: DL dl;
! 602:
! 603: if ( !d ) {
! 604: printf("0"); return;
! 605: }
! 606: for ( n = d->nv, m = BDY(d); m; m = NEXT(m) ) {
! 607: printf("%d*<<",ITOS(m->c));
! 608: for ( i = 0, dl = m->dl; i < n-1; i++ )
! 609: printf("%d,",dl->d[i]);
! 610: printf("%d",dl->d[i]);
! 611: printf(">>");
! 612: if ( NEXT(m) )
! 613: printf("+");
! 614: }
! 615: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>