Annotation of OpenXM_contrib2/asir2018/builtin/pdiv.c, Revision 1.1
1.1 ! noro 1: /*
! 2: * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
! 3: * All rights reserved.
! 4: *
! 5: * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
! 6: * non-exclusive and royalty-free license to use, copy, modify and
! 7: * redistribute, solely for non-commercial and non-profit purposes, the
! 8: * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
! 9: * conditions of this Agreement. For the avoidance of doubt, you acquire
! 10: * only a limited right to use the SOFTWARE hereunder, and FLL or any
! 11: * third party developer retains all rights, including but not limited to
! 12: * copyrights, in and to the SOFTWARE.
! 13: *
! 14: * (1) FLL does not grant you a license in any way for commercial
! 15: * purposes. You may use the SOFTWARE only for non-commercial and
! 16: * non-profit purposes only, such as academic, research and internal
! 17: * business use.
! 18: * (2) The SOFTWARE is protected by the Copyright Law of Japan and
! 19: * international copyright treaties. If you make copies of the SOFTWARE,
! 20: * with or without modification, as permitted hereunder, you shall affix
! 21: * to all such copies of the SOFTWARE the above copyright notice.
! 22: * (3) An explicit reference to this SOFTWARE and its copyright owner
! 23: * shall be made on your publication or presentation in any form of the
! 24: * results obtained by use of the SOFTWARE.
! 25: * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
! 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
! 27: * for such modification or the source code of the modified part of the
! 28: * SOFTWARE.
! 29: *
! 30: * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
! 31: * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
! 32: * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
! 33: * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
! 34: * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
! 35: * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
! 36: * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
! 37: * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
! 38: * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
! 39: * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
! 40: * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
! 41: * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
! 42: * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
! 43: * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
! 44: * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
! 45: * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
! 46: * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
! 47: *
! 48: * $OpenXM$
! 49: */
! 50: #include "ca.h"
! 51: #include "parse.h"
! 52:
! 53: void Psdiv(), Psrem(), Ptdiv(), Psqr(), Pinva_mod(), Pprem();
! 54: void Psdiv_gf2n(), Psrem_gf2n(), Pgcd_gf2n();
! 55: void Psdivm(), Psremm(), Psqrm();
! 56: void Psrem_mod();
! 57: void Pugcd();
! 58: void Purem();
! 59: void Pudiv();
! 60:
! 61: struct ftab pdiv_tab[] = {
! 62: {"sdiv",Psdiv,-3},
! 63: {"srem",Psrem,-3},
! 64: {"prem",Pprem,-3},
! 65: {"sdiv_gf2n",Psdiv_gf2n,2},
! 66: {"srem_gf2n",Psrem_gf2n,2},
! 67: {"gcd_gf2n",Pgcd_gf2n,2},
! 68: {"sqr",Psqr,-3},
! 69: {"tdiv",Ptdiv,-3},
! 70: {"udiv",Pudiv,2},
! 71: {"sdivm",Psdivm,-4},
! 72: {"sremm",Psremm,-4},
! 73: {"sqrm",Psqrm,-4},
! 74: {"inva_mod",Pinva_mod,3},
! 75: {"srem_mod",Psrem_mod,3},
! 76: {"ugcd",Pugcd,2},
! 77: {"urem",Purem,2},
! 78: {0,0,0},
! 79: };
! 80:
! 81: void Psdiv(arg,rp)
! 82: NODE arg;
! 83: Obj *rp;
! 84: {
! 85: P q,r,dnd,dnd1,dvr,dvr1;
! 86: V v;
! 87: VL vl;
! 88:
! 89: asir_assert(ARG0(arg),O_P,"sdiv");
! 90: asir_assert(ARG1(arg),O_P,"sdiv");
! 91: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
! 92: if ( argc(arg) == 3 ) {
! 93: v = VR((P)ARG2(arg));
! 94: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 95: reordvar(CO,v,&vl);
! 96: divsrp(vl,dnd1,dvr1,&q,&r);
! 97: restore_mvar(CO,q,v,(P *)rp);
! 98: } else
! 99: divsrp(CO,dnd,dvr,(P *)rp,&r);
! 100: }
! 101:
! 102: void Psrem(arg,rp)
! 103: NODE arg;
! 104: Obj *rp;
! 105: {
! 106: P q,r,dnd,dnd1,dvr,dvr1;
! 107: V v;
! 108: VL vl;
! 109:
! 110: asir_assert(ARG0(arg),O_P,"srem");
! 111: asir_assert(ARG1(arg),O_P,"srem");
! 112: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
! 113: if ( argc(arg) == 3 ) {
! 114: v = VR((P)ARG2(arg));
! 115: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 116: reordvar(CO,v,&vl);
! 117: divsrp(vl,dnd1,dvr1,&q,&r);
! 118: restore_mvar(CO,r,v,(P *)rp);
! 119: } else
! 120: divsrp(CO,dnd,dvr,&q,(P *)rp);
! 121: }
! 122:
! 123: void Pprem(arg,rp)
! 124: NODE arg;
! 125: P *rp;
! 126: {
! 127: P q,r,dnd,dnd1,dvr,dvr1;
! 128: V v;
! 129: VL vl;
! 130:
! 131: asir_assert(ARG0(arg),O_P,"prem");
! 132: asir_assert(ARG1(arg),O_P,"prem");
! 133: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
! 134: if ( !dvr ) error("prem : division by 0");
! 135: if ( !dnd ) {
! 136: *rp = 0; return;
! 137: }
! 138: if ( argc(arg) == 3 ) {
! 139: v = VR((P)ARG2(arg));
! 140: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 141: reordvar(CO,v,&vl);
! 142: premp(vl,dnd1,dvr1,&r);
! 143: restore_mvar(CO,r,v,rp);
! 144: } else
! 145: premp(CO,dnd,dvr,rp);
! 146: }
! 147:
! 148: void Psqr(arg,rp)
! 149: NODE arg;
! 150: LIST *rp;
! 151: {
! 152: P q,q1,r,r1,dnd,dnd1,dvr,dvr1;
! 153: NODE n,tn;
! 154: V v;
! 155: VL vl;
! 156:
! 157: asir_assert(ARG0(arg),O_P,"sqr");
! 158: asir_assert(ARG1(arg),O_P,"sqr");
! 159: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
! 160: if ( argc(arg) == 3 ) {
! 161: v = VR((P)ARG2(arg));
! 162: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 163: reordvar(CO,v,&vl);
! 164: divsrp(vl,dnd1,dvr1,&q1,&r1);
! 165: restore_mvar(CO,q1,v,&q); restore_mvar(CO,r1,v,&r);
! 166: } else
! 167: divsrp(CO,dnd,dvr,&q,&r);
! 168: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
! 169: }
! 170:
! 171: void Psdiv_gf2n(arg,rp)
! 172: NODE arg;
! 173: GF2N *rp;
! 174: {
! 175: GF2N dnd,dvr;
! 176: UP2 q,r;
! 177:
! 178: dnd = (GF2N)ARG0(arg); dvr = (GF2N)ARG1(arg);
! 179: if ( !dvr )
! 180: error("sdiv_gf2n : division by 0");
! 181: else if ( !dnd )
! 182: *rp = 0;
! 183: else {
! 184: qrup2(dnd->body,dvr->body,&q,&r);
! 185: MKGF2N(q,*rp);
! 186: }
! 187: }
! 188:
! 189: void Psrem_gf2n(arg,rp)
! 190: NODE arg;
! 191: GF2N *rp;
! 192: {
! 193: GF2N dnd,dvr;
! 194: UP2 q,r;
! 195:
! 196: dnd = (GF2N)ARG0(arg); dvr = (GF2N)ARG1(arg);
! 197: if ( !dvr )
! 198: error("srem_gf2n : division by 0");
! 199: else if ( !dnd )
! 200: *rp = 0;
! 201: else {
! 202: qrup2(dnd->body,dvr->body,&q,&r);
! 203: MKGF2N(r,*rp);
! 204: }
! 205: }
! 206:
! 207: void Pgcd_gf2n(arg,rp)
! 208: NODE arg;
! 209: GF2N *rp;
! 210: {
! 211: GF2N p1,p2;
! 212: UP2 gcd;
! 213:
! 214: p1 = (GF2N)ARG0(arg); p2 = (GF2N)ARG1(arg);
! 215: if ( !p1 )
! 216: *rp = p2;
! 217: else if ( !p2 )
! 218: *rp = p1;
! 219: else {
! 220: gcdup2(p1->body,p2->body,&gcd);
! 221: MKGF2N(gcd,*rp);
! 222: }
! 223: }
! 224:
! 225: void Ptdiv(arg,rp)
! 226: NODE arg;
! 227: P *rp;
! 228: {
! 229: P p1,p2,q1,q2,q,c1,c2,c;
! 230: int m;
! 231:
! 232: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 233: asir_assert(p1,O_P,"tdiv");
! 234: asir_assert(p2,O_P,"tdiv");
! 235: if ( !p1 || !p2 )
! 236: *rp = 0;
! 237: else if ( (OID(p1) > O_P) || (OID(p2) > O_P ) )
! 238: *rp = 0;
! 239: else if ( argc(arg) == 3 ) {
! 240: m = QTOS((Q)ARG2(arg));
! 241: ptomp(m,p1,&q1); ptomp(m,p2,&q2);
! 242: if ( divtmp(CO,m,q1,q2,&q) )
! 243: mptop(q,rp);
! 244: else
! 245: *rp = 0;
! 246: } else if ( qpcheck((Obj)p1) && qpcheck((Obj)p2) ) {
! 247: ptozp(p1,1,(Q *)&c1,&q1); ptozp(p2,1,(Q *)&c2,&q2);
! 248: if ( divtpz(CO,q1,q2,&q) ) {
! 249: divq((Q)c1,(Q)c2,(Q *)&c); mulp(CO,q,c,rp);
! 250: } else
! 251: *rp = 0;
! 252: } else {
! 253: if ( !divtp(CO,p1,p2,rp) )
! 254: *rp = 0;
! 255: }
! 256: }
! 257:
! 258: void Pudiv(arg,rp)
! 259: NODE arg;
! 260: LIST *rp;
! 261: {
! 262: P q,r,dnd,dvr;
! 263: NODE n,tn;
! 264:
! 265: asir_assert(ARG0(arg),O_P,"udiv");
! 266: asir_assert(ARG1(arg),O_P,"udiv");
! 267: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
! 268: udivpz(dnd,dvr,&q,&r);
! 269: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
! 270: }
! 271:
! 272: void Psdivm(arg,rp)
! 273: NODE arg;
! 274: Obj *rp;
! 275: {
! 276: P q,r,dnd,dnd1,dndm,dvr,dvr1,dvrm,t;
! 277: V v;
! 278: VL vl;
! 279: int m;
! 280:
! 281: asir_assert(ARG0(arg),O_P,"sdivm");
! 282: asir_assert(ARG1(arg),O_P,"sdivm");
! 283: asir_assert(ARG2(arg),O_N,"sdivm");
! 284: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
! 285: if ( argc(arg) == 4 ) {
! 286: v = VR((P)ARG3(arg));
! 287: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 288: reordvar(CO,v,&vl);
! 289: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
! 290: divsrmp(vl,m,dndm,dvrm,&t,&r); mptop(t,&q);
! 291: restore_mvar(CO,q,v,(P *)rp);
! 292: } else {
! 293: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
! 294: divsrmp(CO,m,dndm,dvrm,&t,&r); mptop(t,(P *)rp);
! 295: }
! 296: }
! 297:
! 298: void Psremm(arg,rp)
! 299: NODE arg;
! 300: Obj *rp;
! 301: {
! 302: P q,r,dnd,dnd1,dndm,dvr,dvr1,dvrm,t;
! 303: V v;
! 304: VL vl;
! 305: int m;
! 306:
! 307: asir_assert(ARG0(arg),O_P,"sremm");
! 308: asir_assert(ARG1(arg),O_P,"sremm");
! 309: asir_assert(ARG2(arg),O_N,"sremm");
! 310: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
! 311: if ( argc(arg) == 4 ) {
! 312: v = VR((P)ARG3(arg));
! 313: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 314: reordvar(CO,v,&vl);
! 315: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
! 316: divsrmp(vl,m,dndm,dvrm,&q,&t); mptop(t,&r);
! 317: restore_mvar(CO,r,v,(P *)rp);
! 318: } else {
! 319: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
! 320: divsrmp(CO,m,dndm,dvrm,&q,&t); mptop(t,(P *)rp);
! 321: }
! 322: }
! 323:
! 324: void Psqrm(arg,rp)
! 325: NODE arg;
! 326: LIST *rp;
! 327: {
! 328: P q,q1,r,r1,dnd,dnd1,dndm,dvr,dvr1,dvrm;
! 329: NODE n,tn;
! 330: V v;
! 331: VL vl;
! 332: int m;
! 333:
! 334: asir_assert(ARG0(arg),O_P,"sqrm");
! 335: asir_assert(ARG1(arg),O_P,"sqrm");
! 336: asir_assert(ARG2(arg),O_N,"sqrm");
! 337: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
! 338: if ( argc(arg) == 4 ) {
! 339: v = VR((P)ARG3(arg));
! 340: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
! 341: reordvar(CO,v,&vl);
! 342: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
! 343: divsrmp(vl,m,dndm,dvrm,&q,&r); mptop(q,&q1); mptop(r,&r1);
! 344: restore_mvar(CO,q1,v,&q); restore_mvar(CO,r1,v,&r);
! 345: } else {
! 346: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
! 347: divsrmp(CO,m,dndm,dvrm,&q1,&r1); mptop(q1,&q); mptop(r1,&r);
! 348: }
! 349: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
! 350: }
! 351:
! 352: void Pinva_mod(arg,rp)
! 353: NODE arg;
! 354: P *rp;
! 355: {
! 356: P dp,f;
! 357: Z q;
! 358: int n,i;
! 359: int mod;
! 360: V v;
! 361: UM wf,wdp,winv;
! 362:
! 363: asir_assert(ARG0(arg),O_P,"gcda_mod");
! 364: asir_assert(ARG1(arg),O_N,"gcda_mod");
! 365: asir_assert(ARG2(arg),O_P,"gcda_mod");
! 366: dp = (P)ARG0(arg);
! 367: mod = QTOS((Q)ARG1(arg));
! 368: f = (P)ARG2(arg);
! 369: if ( NUM(f) ) {
! 370: i = invm(remqi((Q)f,mod),mod);
! 371: STOQ(i,q); *rp = (P)q;
! 372: } else {
! 373: v = VR(dp);
! 374: n = MAX(UDEG(dp),UDEG(f));
! 375: wf = W_UMALLOC(n); wdp = W_UMALLOC(n);
! 376: winv = W_UMALLOC(n);
! 377: ptoum(mod,f,wf); ptoum(mod,dp,wdp);
! 378: invum(mod,wdp,wf,winv);
! 379: if ( DEG(winv) < 0 )
! 380: *rp = 0;
! 381: else {
! 382: umtop(v,winv,rp);
! 383: }
! 384: }
! 385: }
! 386:
! 387: void Psrem_mod(arg,rp)
! 388: NODE arg;
! 389: P *rp;
! 390: {
! 391: P p1,p2;
! 392: int n,dr;
! 393: int mod;
! 394: V v;
! 395: UM wp1,wp2,q;
! 396:
! 397: asir_assert(ARG0(arg),O_P,"srem_mod");
! 398: asir_assert(ARG1(arg),O_P,"srem_mod");
! 399: asir_assert(ARG2(arg),O_N,"srem_mod");
! 400: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg); mod = QTOS((Q)ARG2(arg));
! 401: if ( !p1 || NUM(p1) )
! 402: *rp = p1;
! 403: else {
! 404: v = VR(p1);
! 405: n = MAX(UDEG(p1),UDEG(p2));
! 406: wp1 = W_UMALLOC(n); wp2 = W_UMALLOC(n); q = W_UMALLOC(n);
! 407: ptoum(mod,p1,wp1); ptoum(mod,p2,wp2);
! 408: dr = divum(mod,wp1,wp2,q);
! 409: if ( ( DEG(wp1) = dr ) == -1 )
! 410: *rp = 0;
! 411: else
! 412: umtop(v,wp1,rp);
! 413: }
! 414: }
! 415:
! 416: void Purem(arg,rp)
! 417: NODE arg;
! 418: P *rp;
! 419: {
! 420: asir_assert(ARG0(arg),O_P,"urem");
! 421: asir_assert(ARG1(arg),O_P,"urem");
! 422: uremp((P)ARG0(arg),(P)ARG1(arg),rp);
! 423: }
! 424:
! 425: void Pugcd(arg,rp)
! 426: NODE arg;
! 427: P *rp;
! 428: {
! 429: asir_assert(ARG0(arg),O_P,"ugcd");
! 430: asir_assert(ARG1(arg),O_P,"ugcd");
! 431: ugcdp((P)ARG0(arg),(P)ARG1(arg),rp);
! 432: }
! 433:
! 434: void invum(mod,dp,f,inv)
! 435: int mod;
! 436: UM dp,f,inv;
! 437: {
! 438: UM g1,g2,a1,a2,a3,wm,q,tum;
! 439: int d,dr;
! 440:
! 441: d = DEG(dp)+DEG(f)+10;
! 442: g1 = W_UMALLOC(d); g2 = W_UMALLOC(d); a1 = W_UMALLOC(d);
! 443: a2 = W_UMALLOC(d); a3 = W_UMALLOC(d); wm = W_UMALLOC(d);
! 444: q = W_UMALLOC(d);
! 445: DEG(a1) = 0; COEF(a1)[0] = 1; DEG(a2) = -1;
! 446: cpyum(f,g1); cpyum(dp,g2);
! 447: while ( 1 ) {
! 448: dr = divum(mod,g1,g2,q); tum = g1; g1 = g2; g2 = tum;
! 449: if ( ( DEG(g2) = dr ) == -1 )
! 450: break;
! 451: mulum(mod,a2,q,wm); subum(mod,a1,wm,a3); dr = divum(mod,a3,dp,q);
! 452: tum = a1; a1 = a2; a2 = a3; a3 = tum; DEG(a3) = dr;
! 453: }
! 454: if ( DEG(g1) != 0 )
! 455: DEG(inv) = -1;
! 456: else if ( COEF(g1)[0] != 1 )
! 457: mulsum(mod,a2,invm(COEF(g1)[0],mod),inv);
! 458: else
! 459: cpyum(a2,inv);
! 460: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>