Annotation of OpenXM_contrib2/asir2000/builtin/fctr.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/builtin/fctr.c,v 1.1.1.1 1999/11/10 08:12:25 noro Exp $ */
! 2: #include "ca.h"
! 3: #include "parse.h"
! 4:
! 5: void Pfctr(), Pgcd(), Pgcdz(), Plcm(), Psqfr(), Pufctrhint();
! 6: void Pptozp(), Pcont();
! 7: void Pafctr(), Pagcd();
! 8: void Pmodsqfr(),Pmodfctr(),Pddd(),Pnewddd(),Pddd_tab();
! 9: void Pirred_check(), Pnfctr_mod();
! 10:
! 11: struct ftab fctr_tab[] = {
! 12: {"fctr",Pfctr,1},
! 13: {"gcd",Pgcd,-3},
! 14: {"gcdz",Pgcdz,2},
! 15: {"lcm",Plcm,2},
! 16: {"sqfr",Psqfr,1},
! 17: {"ufctrhint",Pufctrhint,2},
! 18: {"ptozp",Pptozp,1},
! 19: {"cont",Pcont,-2},
! 20: {"afctr",Pafctr,2},
! 21: {"agcd",Pagcd,3},
! 22: {"modsqfr",Pmodsqfr,2},
! 23: {"modfctr",Pmodfctr,2},
! 24: #if 0
! 25: {"ddd",Pddd,2},
! 26: {"newddd",Pnewddd,2},
! 27: #endif
! 28: {"ddd_tab",Pddd_tab,2},
! 29: {"irred_check",Pirred_check,2},
! 30: {"nfctr_mod",Pnfctr_mod,2},
! 31: {0,0,0},
! 32: };
! 33:
! 34: void Pfctr(arg,rp)
! 35: NODE arg;
! 36: LIST *rp;
! 37: {
! 38: DCP dc;
! 39:
! 40: asir_assert(ARG0(arg),O_P,"fctr");
! 41: fctrp(CO,(P)ARG0(arg),&dc);
! 42: dcptolist(dc,rp);
! 43: }
! 44:
! 45: void Pgcd(arg,rp)
! 46: NODE arg;
! 47: P *rp;
! 48: {
! 49: P p1,p2,g1,g2,g;
! 50: Num m;
! 51: int mod;
! 52:
! 53: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 54: asir_assert(p1,O_P,"gcd");
! 55: asir_assert(p2,O_P,"gcd");
! 56: if ( !p1 )
! 57: *rp = p2;
! 58: else if ( !p2 )
! 59: *rp = p1;
! 60: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 61: error("gcd : invalid argument");
! 62: else if ( argc(arg) == 2 )
! 63: ezgcdp(CO,p1,p2,rp);
! 64: else {
! 65: m = (Num)ARG2(arg);
! 66: asir_assert(m,O_P,"gcd");
! 67: mod = QTOS((Q)m);
! 68: ptomp(mod,p1,&g1); ptomp(mod,p2,&g2);
! 69: gcdprsmp(CO,mod,g1,g2,&g);
! 70: mptop(g,rp);
! 71: }
! 72: }
! 73:
! 74: void Pgcdz(arg,rp)
! 75: NODE arg;
! 76: P *rp;
! 77: {
! 78: P p1,p2,t;
! 79: Q c1,c2;
! 80: N n;
! 81:
! 82: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 83: asir_assert(p1,O_P,"gcdz");
! 84: asir_assert(p2,O_P,"gcdz");
! 85: if ( !p1 )
! 86: *rp = p2;
! 87: else if ( !p2 )
! 88: *rp = p1;
! 89: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 90: error("gcdz : invalid argument");
! 91: else if ( NUM(p1) || NUM(p2) ) {
! 92: if ( NUM(p1) )
! 93: c1 = (Q)p1;
! 94: else
! 95: ptozp(p1,1,&c1,&t);
! 96: if ( NUM(p2) )
! 97: c2 = (Q)p2;
! 98: else
! 99: ptozp(p2,1,&c2,&t);
! 100: gcdn(NM(c1),NM(c2),&n); NTOQ(n,1,c1); *rp = (P)c1;
! 101: } else {
! 102: #if 0
! 103: w[0] = p1; w[1] = p2; nezgcdnpz(CO,w,2,rp);
! 104: #endif
! 105: ezgcdpz(CO,p1,p2,rp);
! 106: }
! 107: }
! 108:
! 109: void Plcm(arg,rp)
! 110: NODE arg;
! 111: P *rp;
! 112: {
! 113: P t1,t2,p1,p2,g,q;
! 114: Q c;
! 115:
! 116: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 117: asir_assert(p1,O_P,"lcm");
! 118: asir_assert(p2,O_P,"lcm");
! 119: if ( !p1 || !p2 )
! 120: *rp = 0;
! 121: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 122: error("lcm : invalid argument");
! 123: else {
! 124: ptozp(p1,1,&c,&t1); ptozp(p2,1,&c,&t2);
! 125: ezgcdp(CO,t1,t2,&g); divsp(CO,t1,g,&q); mulp(CO,q,t2,rp);
! 126: }
! 127: }
! 128:
! 129: void Psqfr(arg,rp)
! 130: NODE arg;
! 131: LIST *rp;
! 132: {
! 133: DCP dc;
! 134:
! 135: asir_assert(ARG0(arg),O_P,"sqfr");
! 136: sqfrp(CO,(P)ARG0(arg),&dc);
! 137: dcptolist(dc,rp);
! 138: }
! 139:
! 140: void Pufctrhint(arg,rp)
! 141: NODE arg;
! 142: LIST *rp;
! 143: {
! 144: DCP dc;
! 145:
! 146: asir_assert(ARG0(arg),O_P,"ufctrhint");
! 147: asir_assert(ARG1(arg),O_N,"ufctrhint");
! 148: ufctr((P)ARG0(arg),QTOS((Q)ARG1(arg)),&dc);
! 149: dcptolist(dc,rp);
! 150: }
! 151:
! 152: #if 0
! 153: Pmgcd(arg,rp)
! 154: NODE arg;
! 155: Obj *rp;
! 156: {
! 157: NODE node,tn;
! 158: int i,m;
! 159: P *l;
! 160:
! 161: node = BDY((LIST)ARG0(arg));
! 162: for ( i = 0, tn = node; tn; tn = NEXT(tn), i++ );
! 163: m = i; l = (P *)ALLOCA(m*sizeof(P));
! 164: for ( i = 0, tn = node; i < m; tn = NEXT(tn), i++ )
! 165: l[i] = (P)BDY(tn);
! 166: nezgcdnpz(CO,l,m,rp);
! 167: }
! 168: #endif
! 169:
! 170: void Pcont(arg,rp)
! 171: NODE arg;
! 172: P *rp;
! 173: {
! 174: DCP dc;
! 175: int m;
! 176: P p,p1;
! 177: P *l;
! 178: V v;
! 179:
! 180: asir_assert(ARG0(arg),O_P,"cont");
! 181: p = (P)ARG0(arg);
! 182: if ( NUM(p) )
! 183: *rp = p;
! 184: else {
! 185: if ( argc(arg) == 2 ) {
! 186: v = VR((P)ARG1(arg));
! 187: change_mvar(CO,p,v,&p1);
! 188: if ( VR(p1) != v ) {
! 189: *rp = p1; return;
! 190: } else
! 191: p = p1;
! 192: }
! 193: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ );
! 194: l = (P *)ALLOCA(m*sizeof(P));
! 195: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ )
! 196: l[m] = COEF(dc);
! 197: nezgcdnpz(CO,l,m,rp);
! 198: }
! 199: }
! 200:
! 201: void Pptozp(arg,rp)
! 202: NODE arg;
! 203: P *rp;
! 204: {
! 205: Q t;
! 206:
! 207: asir_assert(ARG0(arg),O_P,"ptozp");
! 208: ptozp((P)ARG0(arg),1,&t,rp);
! 209: }
! 210:
! 211: void Pafctr(arg,rp)
! 212: NODE arg;
! 213: LIST *rp;
! 214: {
! 215: DCP dc;
! 216:
! 217: asir_assert(ARG0(arg),O_P,"afctr");
! 218: asir_assert(ARG1(arg),O_P,"afctr");
! 219: afctr(CO,(P)ARG0(arg),(P)ARG1(arg),&dc);
! 220: dcptolist(dc,rp);
! 221: }
! 222:
! 223: void Pagcd(arg,rp)
! 224: NODE arg;
! 225: P *rp;
! 226: {
! 227: asir_assert(ARG0(arg),O_P,"agcd");
! 228: asir_assert(ARG1(arg),O_P,"agcd");
! 229: asir_assert(ARG2(arg),O_P,"agcd");
! 230: gcda(CO,(P)ARG0(arg),(P)ARG1(arg),(P)ARG2(arg),rp);
! 231: }
! 232:
! 233: #if 1
! 234: #define Mulum mulum
! 235: #define Divum divum
! 236: #define Mulsum mulsum
! 237: #define Gcdum gcdum
! 238: #endif
! 239:
! 240: void Mulum(), Mulsum(), Gcdum();
! 241: int Divum();
! 242:
! 243: #define FCTR 0 /* berlekamp */
! 244: #define SQFR 1
! 245: #define DDD 2 /* Cantor-Zassenhauss */
! 246: #define NEWDDD 3 /* berlekamp + root-finding by Cantor-Zassenhauss */
! 247:
! 248: UM *resberle();
! 249:
! 250: void Pmodfctr(arg,rp)
! 251: NODE arg;
! 252: LIST *rp;
! 253: {
! 254: DCP dc;
! 255: int mod;
! 256:
! 257: mod = QTOS((Q)ARG1(arg));
! 258: if ( mod < 0 )
! 259: error("modfctr : invalid modulus");
! 260: modfctrp(ARG0(arg),mod,NEWDDD,&dc);
! 261: if ( !dc ) {
! 262: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 263: }
! 264: dcptolist(dc,rp);
! 265: }
! 266:
! 267: void Pmodsqfr(arg,rp)
! 268: NODE arg;
! 269: LIST *rp;
! 270: {
! 271: DCP dc;
! 272:
! 273: if ( !dc ) {
! 274: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 275: }
! 276: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),SQFR,&dc);
! 277: dcptolist(dc,rp);
! 278: }
! 279:
! 280: void Pddd(arg,rp)
! 281: NODE arg;
! 282: LIST *rp;
! 283: {
! 284: DCP dc;
! 285:
! 286: if ( !dc ) {
! 287: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 288: }
! 289: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),DDD,&dc);
! 290: dcptolist(dc,rp);
! 291: }
! 292:
! 293: void Pnewddd(arg,rp)
! 294: NODE arg;
! 295: LIST *rp;
! 296: {
! 297: DCP dc;
! 298:
! 299: if ( !dc ) {
! 300: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 301: }
! 302: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),NEWDDD,&dc);
! 303: dcptolist(dc,rp);
! 304: }
! 305:
! 306: void Pirred_check(arg,rp)
! 307: NODE arg;
! 308: Q *rp;
! 309: {
! 310: P p;
! 311: UM mp;
! 312: int r,mod;
! 313:
! 314: p = (P)ARG0(arg);
! 315: if ( !p ) {
! 316: *rp = 0; return;
! 317: }
! 318: mp = W_UMALLOC(UDEG(p));
! 319: mod = QTOS((Q)ARG1(arg));
! 320: ptoum(mod,p,mp);
! 321: r = irred_check(mp,mod);
! 322: if ( r )
! 323: *rp = ONE;
! 324: else
! 325: *rp = 0;
! 326: }
! 327:
! 328: void Pnfctr_mod(arg,rp)
! 329: NODE arg;
! 330: Q *rp;
! 331: {
! 332: P p;
! 333: UM mp;
! 334: int r,mod;
! 335:
! 336: p = (P)ARG0(arg);
! 337: if ( !p ) {
! 338: *rp = 0; return;
! 339: }
! 340: mp = W_UMALLOC(UDEG(p));
! 341: mod = QTOS((Q)ARG1(arg));
! 342: ptoum(mod,p,mp);
! 343: r = nfctr_mod(mp,mod);
! 344: STOQ(r,*rp);
! 345: }
! 346:
! 347: void Pddd_tab(arg,rp)
! 348: NODE arg;
! 349: VECT *rp;
! 350: {
! 351: P p;
! 352: UM mp,t,q,r1,w,w1;
! 353: UM *r,*s;
! 354: int dr,mod,n,i;
! 355: VECT result;
! 356: V v;
! 357:
! 358: p = (P)ARG0(arg); mod = QTOS((Q)ARG1(arg));
! 359: v = VR(p);
! 360: n = UDEG(p); mp = W_UMALLOC(n);
! 361: ptoum(mod,p,mp);
! 362: r = (UM *)W_ALLOC(n); s = (UM *)W_ALLOC(n);
! 363: r[0] = UMALLOC(0); DEG(r[0]) = 0; COEF(r[0])[0] = 1;
! 364: t = W_UMALLOC(mod); bzero(COEF(t),sizeof(int)*(mod+1));
! 365: DEG(t) = mod; COEF(t)[mod] = 1;
! 366: q = W_UMALLOC(mod);
! 367: dr = divum(mod,t,mp,q);
! 368: DEG(t) = dr; r[1] = r1 = UMALLOC(dr); cpyum(t,r1);
! 369: s[0] = W_UMALLOC(dr); cpyum(t,s[0]);
! 370: w = W_UMALLOC(n); bzero(COEF(w),sizeof(int)*(n+1));
! 371: w1 = W_UMALLOC(2*n); bzero(COEF(w1),sizeof(int)*(2*n+1));
! 372: for ( i = 1; i < n; i++ ) {
! 373: DEG(w) = i; COEF(w)[i-1] = 0; COEF(w)[i] = 1;
! 374: mulum(mod,r1,w,w1);
! 375: dr = divum(mod,w1,mp,q); DEG(w1) = dr;
! 376: s[i] = W_UMALLOC(dr); cpyum(w1,s[i]);
! 377: }
! 378: for ( i = 2; i < n; i++ ) {
! 379: mult_mod_tab(r[i-1],mod,s,w,n);
! 380: r[i] = UMALLOC(DEG(w)); cpyum(w,r[i]);
! 381: }
! 382: MKVECT(result,n);
! 383: for ( i = 0; i < n; i++ )
! 384: umtop(v,r[i],(P *)&BDY(result)[i]);
! 385: *rp = result;
! 386: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>