Annotation of OpenXM_contrib2/asir2000/engine/P.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/engine/P.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
! 2: #ifndef FBASE
! 3: #define FBASE
! 4: #endif
! 5:
! 6: #include "b.h"
! 7: #include "ca.h"
! 8:
! 9: void D_ADDP(vl,p1,p2,pr)
! 10: VL vl;
! 11: P p1,p2,*pr;
! 12: {
! 13: register DCP dc1,dc2,dcr0,dcr;
! 14: V v1,v2;
! 15: P t;
! 16:
! 17: if ( !p1 )
! 18: *pr = p2;
! 19: else if ( !p2 )
! 20: *pr = p1;
! 21: else if ( NUM(p1) )
! 22: if ( NUM(p2) )
! 23: ADDNUM(p1,p2,pr);
! 24: else
! 25: ADDPQ(p2,p1,pr);
! 26: else if ( NUM(p2) )
! 27: ADDPQ(p1,p2,pr);
! 28: else if ( ( v1 = VR(p1) ) == ( v2 = VR(p2) ) ) {
! 29: for ( dc1 = DC(p1), dc2 = DC(p2), dcr0 = 0; dc1 && dc2; )
! 30: switch ( cmpq(DEG(dc1),DEG(dc2)) ) {
! 31: case 0:
! 32: ADDP(vl,COEF(dc1),COEF(dc2),&t);
! 33: if ( t ) {
! 34: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc1); COEF(dcr) = t;
! 35: }
! 36: dc1 = NEXT(dc1); dc2 = NEXT(dc2); break;
! 37: case 1:
! 38: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc1); COEF(dcr) = COEF(dc1);
! 39: dc1 = NEXT(dc1); break;
! 40: case -1:
! 41: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc2); COEF(dcr) = COEF(dc2);
! 42: dc2 = NEXT(dc2); break;
! 43: }
! 44: if ( !dcr0 )
! 45: if ( dc1 )
! 46: dcr0 = dc1;
! 47: else if ( dc2 )
! 48: dcr0 = dc2;
! 49: else {
! 50: *pr = 0;
! 51: return;
! 52: }
! 53: else
! 54: if ( dc1 )
! 55: NEXT(dcr) = dc1;
! 56: else if ( dc2 )
! 57: NEXT(dcr) = dc2;
! 58: else
! 59: NEXT(dcr) = 0;
! 60: MKP(v1,dcr0,*pr);
! 61: } else {
! 62: while ( v1 != VR(vl) && v2 != VR(vl) )
! 63: vl = NEXT(vl);
! 64: if ( v1 == VR(vl) )
! 65: ADDPTOC(vl,p1,p2,pr);
! 66: else
! 67: ADDPTOC(vl,p2,p1,pr);
! 68: }
! 69: }
! 70:
! 71: void D_ADDPQ(p,q,pr)
! 72: P p,q,*pr;
! 73: {
! 74: DCP dc,dcr,dcr0;
! 75: P t;
! 76:
! 77: if ( NUM(p) )
! 78: ADDNUM(p,q,pr);
! 79: else {
! 80: for ( dcr0 = 0, dc = DC(p); dc && DEG(dc); dc = NEXT(dc) ) {
! 81: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = COEF(dc);
! 82: }
! 83: if ( !dc ) {
! 84: NEWDC(NEXT(dcr)); dcr = NEXT(dcr); DEG(dcr) = 0; COEF(dcr) = q;
! 85: } else {
! 86: ADDPQ(COEF(dc),q,&t);
! 87: if ( t ) {
! 88: NEWDC(NEXT(dcr)); dcr = NEXT(dcr); DEG(dcr) = 0; COEF(dcr) = t;
! 89: }
! 90: }
! 91: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
! 92: }
! 93: }
! 94:
! 95: void D_ADDPTOC(vl,p,c,pr)
! 96: VL vl;
! 97: P p,c,*pr;
! 98: {
! 99: DCP dc,dcr,dcr0;
! 100: P t;
! 101:
! 102: for ( dcr0 = 0, dc = DC(p); dc && DEG(dc); dc = NEXT(dc) ) {
! 103: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = COEF(dc);
! 104: }
! 105: if ( !dc ) {
! 106: NEWDC(NEXT(dcr)); dcr = NEXT(dcr); DEG(dcr) = 0; COEF(dcr) = c;
! 107: } else {
! 108: ADDP(vl,COEF(dc),c,&t);
! 109: if ( t ) {
! 110: NEWDC(NEXT(dcr)); dcr = NEXT(dcr); DEG(dcr) = 0; COEF(dcr) = t;
! 111: }
! 112: }
! 113: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
! 114: }
! 115:
! 116: void D_SUBP(vl,p1,p2,pr)
! 117: VL vl;
! 118: P p1,p2,*pr;
! 119: {
! 120: P t;
! 121:
! 122: if ( !p2 )
! 123: *pr = p1;
! 124: else {
! 125: CHSGNP(p2,&t); ADDP(vl,p1,t,pr);
! 126: }
! 127: }
! 128:
! 129: void D_MULP(vl,p1,p2,pr)
! 130: VL vl;
! 131: P p1,p2,*pr;
! 132: {
! 133: register DCP dc,dct,dcr,dcr0;
! 134: V v1,v2;
! 135: P t,s,u;
! 136: int n1,n2;
! 137:
! 138: if ( !p1 || !p2 ) *pr = 0;
! 139: else if ( NUM(p1) )
! 140: MULPQ(p2,p1,pr);
! 141: else if ( NUM(p2) )
! 142: MULPQ(p1,p2,pr);
! 143: else if ( ( v1 = VR(p1) ) == ( v2 = VR(p2) ) ) {
! 144: for ( dc = DC(p1), n1 = 0; dc; dc = NEXT(dc), n1++ );
! 145: for ( dc = DC(p2), n2 = 0; dc; dc = NEXT(dc), n2++ );
! 146: if ( n1 > n2 )
! 147: for ( dc = DC(p2), s = 0; dc; dc = NEXT(dc) ) {
! 148: for ( dcr0 = 0, dct = DC(p1); dct; dct = NEXT(dct) ) {
! 149: NEXTDC(dcr0,dcr); MULP(vl,COEF(dct),COEF(dc),&COEF(dcr));
! 150: addq(DEG(dct),DEG(dc),&DEG(dcr));
! 151: }
! 152: NEXT(dcr) = 0; MKP(v1,dcr0,t);
! 153: ADDP(vl,s,t,&u); s = u; t = u = 0;
! 154: }
! 155: else
! 156: for ( dc = DC(p1), s = 0; dc; dc = NEXT(dc) ) {
! 157: for ( dcr0 = 0, dct = DC(p2); dct; dct = NEXT(dct) ) {
! 158: NEXTDC(dcr0,dcr); MULP(vl,COEF(dct),COEF(dc),&COEF(dcr));
! 159: addq(DEG(dct),DEG(dc),&DEG(dcr));
! 160: }
! 161: NEXT(dcr) = 0; MKP(v1,dcr0,t);
! 162: ADDP(vl,s,t,&u); s = u; t = u = 0;
! 163: }
! 164: *pr = s;
! 165: } else {
! 166: while ( v1 != VR(vl) && v2 != VR(vl) )
! 167: vl = NEXT(vl);
! 168: if ( v1 == VR(vl) )
! 169: MULPC(vl,p1,p2,pr);
! 170: else
! 171: MULPC(vl,p2,p1,pr);
! 172: }
! 173: }
! 174:
! 175: void D_MULPQ(p,q,pr)
! 176: P p,q,*pr;
! 177: {
! 178: DCP dc,dcr,dcr0;
! 179: P t;
! 180:
! 181: if (!p || !q)
! 182: *pr = 0;
! 183: else if ( Uniq(q) )
! 184: *pr = p;
! 185: else if ( NUM(p) )
! 186: MULNUM(p,q,pr);
! 187: else {
! 188: for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
! 189: MULPQ(COEF(dc),q,&t);
! 190: if ( t ) {
! 191: NEXTDC(dcr0,dcr); COEF(dcr) = t; DEG(dcr) = DEG(dc);
! 192: }
! 193: }
! 194: if ( dcr0 ) {
! 195: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
! 196: } else
! 197: *pr = 0;
! 198: }
! 199: }
! 200:
! 201: void D_MULPC(vl,p,c,pr)
! 202: VL vl;
! 203: P p,c,*pr;
! 204: {
! 205: DCP dc,dcr,dcr0;
! 206: P t;
! 207:
! 208: if ( NUM(c) )
! 209: MULPQ(p,c,pr);
! 210: else {
! 211: for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
! 212: MULP(vl,COEF(dc),c,&t);
! 213: if ( t ) {
! 214: NEXTDC(dcr0,dcr); COEF(dcr) = t; DEG(dcr) = DEG(dc);
! 215: }
! 216: }
! 217: if ( dcr0 ) {
! 218: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
! 219: } else
! 220: *pr = 0;
! 221: }
! 222: }
! 223:
! 224: void D_PWRP(vl,p,q,pr)
! 225: VL vl;
! 226: P p,*pr;
! 227: Q q;
! 228: {
! 229: DCP dc,dcr;
! 230: int n,i;
! 231: P *x,*y;
! 232: P t,s,u;
! 233: DCP dct;
! 234: P *pt;
! 235:
! 236: if ( !q ) {
! 237: *pr = (P)One;
! 238: } else if ( !p )
! 239: *pr = 0;
! 240: else if ( UNIQ(q) )
! 241: *pr = p;
! 242: else if ( NUM(p) )
! 243: PWRNUM(p,q,pr);
! 244: else {
! 245: dc = DC(p);
! 246: if ( !NEXT(dc) ) {
! 247: NEWDC(dcr);
! 248: PWRP(vl,COEF(dc),q,&COEF(dcr)); mulq(DEG(dc),q,&DEG(dcr));
! 249: NEXT(dcr) = 0; MKP(VR(p),dcr,*pr);
! 250: } else if ( !INT(q) ) {
! 251: error("pwrp: can't calculate fractional power."); *pr = 0;
! 252: } else if ( PL(NM(q)) == 1 ) {
! 253: n = QTOS(q); x = (P *)ALLOCA((n+1)*sizeof(pointer));
! 254: NEWDC(dct); DEG(dct) = DEG(dc); COEF(dct) = COEF(dc);
! 255: NEXT(dct) = 0; MKP(VR(p),dct,t);
! 256: for ( i = 0, u = (P)One; i < n; i++ ) {
! 257: x[i] = u; MULP(vl,u,t,&s); u = s;
! 258: }
! 259: x[n] = u; y = (P *)ALLOCA((n+1)*sizeof(pointer));
! 260: if ( DEG(NEXT(dc)) ) {
! 261: dct = NEXT(dc); MKP(VR(p),dct,t);
! 262: } else
! 263: t = COEF(NEXT(dc));
! 264: for ( i = 0, u = (P)One; i < n; i++ ) {
! 265: y[i] = u; MULP(vl,u,t,&s); u = s;
! 266: }
! 267: y[n] = u;
! 268: pt = (P *)ALLOCA((n+1)*sizeof(pointer)); MKBC(n,pt);
! 269: for ( i = 0, u = 0; i <= n; i++ ) {
! 270: MULP(vl,x[i],y[n-i],&t); MULP(vl,t,pt[i],&s);
! 271: ADDP(vl,u,s,&t); u = t;
! 272: }
! 273: *pr = u;
! 274: } else {
! 275: error("exponent too big");
! 276: *pr = 0;
! 277: }
! 278: }
! 279: }
! 280:
! 281: void D_CHSGNP(p,pr)
! 282: P p,*pr;
! 283: {
! 284: register DCP dc,dcr,dcr0;
! 285: P t;
! 286:
! 287: if ( !p )
! 288: *pr = NULL;
! 289: else if ( NUM(p) ) {
! 290: #if defined(THINK_C) || defined(_PA_RISC1_1) || defined(__alpha) || defined(mips)
! 291: #ifdef FBASE
! 292: chsgnnum((Num)p,(Num *)pr);
! 293: #else
! 294: MQ mq;
! 295:
! 296: NEWMQ(mq); CONT(mq)=mod-CONT((MQ)p); *pr = (P)mq;
! 297: #endif
! 298: #else
! 299: CHSGNNUM(p,t); *pr = t;
! 300: #endif
! 301: } else {
! 302: for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
! 303: NEXTDC(dcr0,dcr); CHSGNP(COEF(dc),&COEF(dcr)); DEG(dcr) = DEG(dc);
! 304: }
! 305: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
! 306: }
! 307: }
! 308:
! 309:
! 310: #ifdef FBASE
! 311: void ptozp(p,sgn,c,pr)
! 312: P p;
! 313: int sgn;
! 314: Q *c;
! 315: P *pr;
! 316: {
! 317: N nm,dn;
! 318:
! 319: if ( !p ) {
! 320: *c = 0; *pr = 0;
! 321: } else {
! 322: lgp(p,&nm,&dn);
! 323: if ( UNIN(dn) )
! 324: NTOQ(nm,sgn,*c);
! 325: else
! 326: NDTOQ(nm,dn,sgn,*c);
! 327: divcp(p,*c,pr);
! 328: }
! 329: }
! 330:
! 331: void lgp(p,g,l)
! 332: P p;
! 333: N *g,*l;
! 334: {
! 335: N g1,g2,l1,l2,l3,l4,tmp;
! 336: DCP dc;
! 337:
! 338: if ( NUM(p) ) {
! 339: *g = NM((Q)p);
! 340: if ( INT((Q)p) )
! 341: *l = ONEN;
! 342: else
! 343: *l = DN((Q)p);
! 344: } else {
! 345: dc = DC(p); lgp(COEF(dc),g,l);
! 346: for ( dc = NEXT(dc); dc; dc = NEXT(dc) ) {
! 347: lgp(COEF(dc),&g1,&l1); gcdn(*g,g1,&g2); *g = g2;
! 348: gcdn(*l,l1,&l2); kmuln(*l,l1,&l3); divn(l3,l2,&l4,&tmp); *l = l4;
! 349: }
! 350: }
! 351: }
! 352:
! 353: void divcp(f,q,rp)
! 354: P f;
! 355: Q q;
! 356: P *rp;
! 357: {
! 358: DCP dc,dcr,dcr0;
! 359:
! 360: if ( !f )
! 361: *rp = 0;
! 362: else if ( NUM(f) )
! 363: DIVNUM(f,q,rp);
! 364: else {
! 365: for ( dc = DC(f), dcr0 = 0; dc; dc = NEXT(dc) ) {
! 366: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc);
! 367: divcp(COEF(dc),q,&COEF(dcr));
! 368: }
! 369: NEXT(dcr) = 0; MKP(VR(f),dcr0,*rp);
! 370: }
! 371: }
! 372:
! 373: void diffp(vl,p,v,r)
! 374: VL vl;
! 375: P p;
! 376: V v;
! 377: P *r;
! 378: {
! 379: P t;
! 380: DCP dc,dcr,dcr0;
! 381:
! 382: if ( !p || NUM(p) )
! 383: *r = 0;
! 384: else {
! 385: if ( v == VR(p) ) {
! 386: for ( dc = DC(p), dcr0 = 0;
! 387: dc && DEG(dc); dc = NEXT(dc) ) {
! 388: MULPQ(COEF(dc),(P)DEG(dc),&t);
! 389: if ( t ) {
! 390: NEXTDC(dcr0,dcr); SUBQ(DEG(dc),ONE,&DEG(dcr));
! 391: COEF(dcr) = t;
! 392: }
! 393: }
! 394: if ( !dcr0 )
! 395: *r = 0;
! 396: else {
! 397: NEXT(dcr) = 0; MKP(v,dcr0,*r);
! 398: }
! 399: } else {
! 400: for ( dc = DC(p), dcr0 = 0; dc; dc = NEXT(dc) ) {
! 401: diffp(vl,COEF(dc),v,&t);
! 402: if ( t ) {
! 403: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = t;
! 404: }
! 405: }
! 406: if ( !dcr0 )
! 407: *r = 0;
! 408: else {
! 409: NEXT(dcr) = 0; MKP(VR(p),dcr0,*r);
! 410: }
! 411: }
! 412: }
! 413: }
! 414:
! 415: void coefp(p,d,pr)
! 416: P p;
! 417: int d;
! 418: P *pr;
! 419: {
! 420: DCP dc;
! 421: int sgn;
! 422: Q dq;
! 423:
! 424: if ( NUM(p) )
! 425: if ( d == 0 )
! 426: *pr = p;
! 427: else
! 428: *pr = 0;
! 429: else {
! 430: for ( STOQ(d,dq), dc = DC(p); dc; dc = NEXT(dc) )
! 431: if ( (sgn = cmpq(DEG(dc),dq)) > 0 )
! 432: continue;
! 433: else if ( sgn == 0 ) {
! 434: *pr = COEF(dc);
! 435: return;
! 436: } else {
! 437: *pr = 0;
! 438: break;
! 439: }
! 440: *pr = 0;
! 441: }
! 442: }
! 443:
! 444: int compp(vl,p1,p2)
! 445: VL vl;
! 446: P p1,p2;
! 447: {
! 448: DCP dc1,dc2;
! 449: V v1,v2;
! 450:
! 451: if ( !p1 )
! 452: return p2 ? -1 : 0;
! 453: else if ( !p2 )
! 454: return 1;
! 455: else if ( NUM(p1) )
! 456: return NUM(p2) ? (*cmpnumt[MAX(NID(p1),NID(p2))])(p1,p2) : -1;
! 457: else if ( NUM(p2) )
! 458: return 1;
! 459: if ( (v1 = VR(p1)) == (v2 = VR(p2)) ) {
! 460: for ( dc1 = DC(p1), dc2 = DC(p2);
! 461: dc1 && dc2; dc1 = NEXT(dc1), dc2 = NEXT(dc2) )
! 462: switch ( cmpq(DEG(dc1),DEG(dc2)) ) {
! 463: case 1:
! 464: return 1; break;
! 465: case -1:
! 466: return -1; break;
! 467: default:
! 468: switch ( compp(vl,COEF(dc1),COEF(dc2)) ) {
! 469: case 1:
! 470: return 1; break;
! 471: case -1:
! 472: return -1; break;
! 473: default:
! 474: break;
! 475: }
! 476: break;
! 477: }
! 478: return dc1 ? 1 : (dc2 ? -1 : 0);
! 479: } else {
! 480: for ( ; v1 != VR(vl) && v2 != VR(vl); vl = NEXT(vl) );
! 481: return v1 == VR(vl) ? 1 : -1;
! 482: }
! 483: }
! 484:
! 485: void csump(vl,p,c)
! 486: VL vl;
! 487: P p;
! 488: Q *c;
! 489: {
! 490: DCP dc;
! 491: Q s,s1,s2;
! 492:
! 493: if ( !p || NUM(p) )
! 494: *c = (Q)p;
! 495: else {
! 496: for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) {
! 497: csump(vl,COEF(dc),&s1); addq(s,s1,&s2); s = s2;
! 498: }
! 499: *c = s;
! 500: }
! 501: }
! 502:
! 503: void degp(v,p,d)
! 504: V v;
! 505: P p;
! 506: Q *d;
! 507: {
! 508: DCP dc;
! 509: Q m,m1;
! 510:
! 511: if ( !p || NUM(p) )
! 512: *d = 0;
! 513: else if ( v == VR(p) )
! 514: *d = DEG(DC(p));
! 515: else {
! 516: for ( dc = DC(p), m = 0; dc; dc = NEXT(dc) ) {
! 517: degp(v,COEF(dc),&m1);
! 518: if ( cmpq(m,m1) < 0 )
! 519: m = m1;
! 520: }
! 521: *d = m;
! 522: }
! 523: }
! 524: #endif
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>