Annotation of OpenXM_contrib2/asir2000/builtin/fctr.c, Revision 1.25
1.2 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
1.3 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.2 noro 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: *
1.25 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/fctr.c,v 1.24 2017/02/27 05:14:54 noro Exp $
1.2 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include "parse.h"
52:
53: void Pfctr(), Pgcd(), Pgcdz(), Plcm(), Psqfr(), Pufctrhint();
1.17 noro 54: void Pptozp(), Pcont(), Psfcont();
1.1 noro 55: void Pafctr(), Pagcd();
56: void Pmodsqfr(),Pmodfctr(),Pddd(),Pnewddd(),Pddd_tab();
1.15 noro 57: void Psfsqfr(),Psffctr(),Psfbfctr(),Psfufctr(),Psfmintdeg(),Psfgcd();
1.1 noro 58: void Pirred_check(), Pnfctr_mod();
1.16 noro 59: void Pbivariate_hensel_special();
1.1 noro 60:
1.11 noro 61: void sfmintdeg(VL vl,P fx,int dy,int c,P *fr);
62:
1.1 noro 63: struct ftab fctr_tab[] = {
1.25 ! noro 64: {"bivariate_hensel_special",Pbivariate_hensel_special,6},
! 65: {"fctr",Pfctr,-2},
! 66: {"gcd",Pgcd,-3},
! 67: {"gcdz",Pgcdz,2},
! 68: {"lcm",Plcm,2},
! 69: {"sqfr",Psqfr,1},
! 70: {"ufctrhint",Pufctrhint,2},
! 71: {"ptozp",Pptozp,1},
! 72: {"cont",Pcont,-2},
! 73: {"sfcont",Psfcont,-2},
! 74: {"afctr",Pafctr,2},
! 75: {"agcd",Pagcd,3},
! 76: {"modsqfr",Pmodsqfr,2},
! 77: {"modfctr",Pmodfctr,2},
! 78: {"sfsqfr",Psfsqfr,1},
! 79: {"sffctr",Psffctr,1},
! 80: {"sfufctr",Psfufctr,1},
! 81: {"sfbfctr",Psfbfctr,-4},
! 82: {"sfmintdeg",Psfmintdeg,5},
! 83: {"sfgcd",Psfgcd,2},
1.1 noro 84: #if 0
1.25 ! noro 85: {"ddd",Pddd,2},
! 86: {"newddd",Pnewddd,2},
1.1 noro 87: #endif
1.25 ! noro 88: {"ddd_tab",Pddd_tab,2},
! 89: {"irred_check",Pirred_check,2},
! 90: {"nfctr_mod",Pnfctr_mod,2},
! 91: {0,0,0},
1.1 noro 92: };
1.16 noro 93:
94: /* bivariate_hensel_special(f(x,y):monic in x,g0(x),h0(y),x,y,d) */
95:
96: void Pbivariate_hensel_special(arg,rp)
97: NODE arg;
98: LIST *rp;
99: {
1.25 ! noro 100: DCP dc;
! 101: struct oVN vn[2];
! 102: P f,g0,h0,ak,bk,gk,hk;
! 103: V vx,vy;
! 104: VL nvl;
! 105: Q qk,cbd,bb;
! 106: int d;
! 107: NODE n;
! 108:
! 109: f = (P)ARG0(arg);
! 110: g0 = (P)ARG1(arg);
! 111: h0 = (P)ARG2(arg);
! 112: vx = VR((P)ARG3(arg));
! 113: vy = VR((P)ARG4(arg));
! 114: d = QTOS((Q)ARG5(arg));
! 115: NEWVL(nvl); nvl->v = vx;
! 116: NEWVL(NEXT(nvl)); NEXT(nvl)->v = vy;
! 117: NEXT(NEXT(nvl)) = 0;
! 118: vn[0].v = vy; vn[0].n = 0;
! 119: vn[1].v = 0; vn[1].n = 0;
! 120: cbound(nvl,f,&cbd);
! 121: addq(cbd,cbd,&bb);
! 122: henzq1(g0,h0,bb,&bk,&ak,&qk);
! 123: henmv(nvl,vn,f,g0,h0,ak,bk,(P)ONE,(P)ONE,(P)ONE,(P)ONE,qk,d,&gk,&hk);
! 124: n = mknode(2,gk,hk);
! 125: MKLIST(*rp,n);
1.16 noro 126: }
1.1 noro 127:
128: void Pfctr(arg,rp)
129: NODE arg;
130: LIST *rp;
131: {
1.25 ! noro 132: DCP dc;
1.1 noro 133:
1.25 ! noro 134: asir_assert(ARG0(arg),O_P,"fctr");
! 135: if ( argc(arg) == 1 )
! 136: fctrp(CO,(P)ARG0(arg),&dc);
! 137: else {
! 138: asir_assert(ARG1(arg),O_P,"fctr");
! 139: fctr_wrt_v_p(CO,(P)ARG0(arg),VR((P)ARG1(arg)),&dc);
! 140: }
! 141: dcptolist(dc,rp);
1.1 noro 142: }
143:
144: void Pgcd(arg,rp)
145: NODE arg;
146: P *rp;
147: {
1.25 ! noro 148: P p1,p2,g1,g2,g;
! 149: Num m;
! 150: int mod;
! 151:
! 152: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 153: asir_assert(p1,O_P,"gcd");
! 154: asir_assert(p2,O_P,"gcd");
! 155: if ( !p1 )
! 156: *rp = p2;
! 157: else if ( !p2 )
! 158: *rp = p1;
! 159: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 160: gcdprsp(CO,p1,p2,rp);
! 161: else if ( argc(arg) == 2 )
! 162: ezgcdp(CO,p1,p2,rp);
! 163: else {
! 164: m = (Num)ARG2(arg);
! 165: asir_assert(m,O_P,"gcd");
! 166: mod = QTOS((Q)m);
! 167: ptomp(mod,p1,&g1); ptomp(mod,p2,&g2);
! 168: gcdprsmp(CO,mod,g1,g2,&g);
! 169: mptop(g,rp);
! 170: }
1.1 noro 171: }
172:
173: void Pgcdz(arg,rp)
174: NODE arg;
175: P *rp;
176: {
1.25 ! noro 177: P p1,p2,t;
! 178: Q c1,c2;
! 179: N n;
! 180:
! 181: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 182: asir_assert(p1,O_P,"gcdz");
! 183: asir_assert(p2,O_P,"gcdz");
! 184: if ( !p1 )
! 185: *rp = p2;
! 186: else if ( !p2 )
! 187: *rp = p1;
! 188: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 189: error("gcdz : invalid argument");
! 190: else if ( NUM(p1) || NUM(p2) ) {
! 191: if ( NUM(p1) )
! 192: c1 = (Q)p1;
! 193: else
! 194: ptozp(p1,1,&c1,&t);
! 195: if ( NUM(p2) )
! 196: c2 = (Q)p2;
! 197: else
! 198: ptozp(p2,1,&c2,&t);
! 199: gcdn(NM(c1),NM(c2),&n); NTOQ(n,1,c1); *rp = (P)c1;
! 200: } else {
1.1 noro 201: #if 0
1.25 ! noro 202: w[0] = p1; w[1] = p2; nezgcdnpz(CO,w,2,rp);
1.1 noro 203: #endif
1.25 ! noro 204: ezgcdpz(CO,p1,p2,rp);
! 205: }
1.1 noro 206: }
207:
208: void Plcm(arg,rp)
209: NODE arg;
210: P *rp;
211: {
1.25 ! noro 212: P t1,t2,p1,p2,g,q;
! 213: Q c;
1.1 noro 214:
1.25 ! noro 215: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
! 216: asir_assert(p1,O_P,"lcm");
! 217: asir_assert(p2,O_P,"lcm");
! 218: if ( !p1 || !p2 )
! 219: *rp = 0;
! 220: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
! 221: error("lcm : invalid argument");
! 222: else {
! 223: ptozp(p1,1,&c,&t1); ptozp(p2,1,&c,&t2);
! 224: ezgcdp(CO,t1,t2,&g); divsp(CO,t1,g,&q); mulp(CO,q,t2,rp);
! 225: }
1.1 noro 226: }
227:
228: void Psqfr(arg,rp)
229: NODE arg;
230: LIST *rp;
231: {
1.25 ! noro 232: DCP dc;
1.1 noro 233:
1.25 ! noro 234: asir_assert(ARG0(arg),O_P,"sqfr");
! 235: sqfrp(CO,(P)ARG0(arg),&dc);
! 236: dcptolist(dc,rp);
1.1 noro 237: }
238:
239: void Pufctrhint(arg,rp)
240: NODE arg;
241: LIST *rp;
242: {
1.25 ! noro 243: DCP dc;
1.1 noro 244:
1.25 ! noro 245: asir_assert(ARG0(arg),O_P,"ufctrhint");
! 246: asir_assert(ARG1(arg),O_N,"ufctrhint");
! 247: ufctr((P)ARG0(arg),QTOS((Q)ARG1(arg)),&dc);
! 248: dcptolist(dc,rp);
1.1 noro 249: }
250:
251: #if 0
252: Pmgcd(arg,rp)
253: NODE arg;
254: Obj *rp;
255: {
1.25 ! noro 256: NODE node,tn;
! 257: int i,m;
! 258: P *l;
! 259:
! 260: node = BDY((LIST)ARG0(arg));
! 261: for ( i = 0, tn = node; tn; tn = NEXT(tn), i++ );
! 262: m = i; l = (P *)ALLOCA(m*sizeof(P));
! 263: for ( i = 0, tn = node; i < m; tn = NEXT(tn), i++ )
! 264: l[i] = (P)BDY(tn);
! 265: nezgcdnpz(CO,l,m,rp);
1.1 noro 266: }
267: #endif
268:
269: void Pcont(arg,rp)
270: NODE arg;
271: P *rp;
272: {
1.25 ! noro 273: DCP dc;
! 274: int m;
! 275: P p,p1;
! 276: P *l;
! 277: V v;
! 278:
! 279: asir_assert(ARG0(arg),O_P,"cont");
! 280: p = (P)ARG0(arg);
! 281: if ( NUM(p) )
! 282: *rp = p;
! 283: else {
! 284: if ( argc(arg) == 2 ) {
! 285: v = VR((P)ARG1(arg));
! 286: change_mvar(CO,p,v,&p1);
! 287: if ( VR(p1) != v ) {
! 288: *rp = p1; return;
! 289: } else
! 290: p = p1;
! 291: }
! 292: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ );
! 293: l = (P *)ALLOCA(m*sizeof(P));
! 294: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ )
! 295: l[m] = COEF(dc);
! 296: nezgcdnpz(CO,l,m,rp);
! 297: }
1.17 noro 298: }
299:
300: void Psfcont(arg,rp)
301: NODE arg;
302: P *rp;
303: {
1.25 ! noro 304: DCP dc;
! 305: MP mp;
! 306: int m;
! 307: Obj obj;
! 308: P p,p1;
! 309: P *l;
! 310: V v;
! 311:
! 312: obj = (Obj)ARG0(arg);
! 313: if ( !obj || NUM(obj) )
! 314: *rp = (P)obj;
! 315: else if ( OID(obj) == O_P ) {
! 316: p = (P)obj;
! 317: if ( argc(arg) == 2 ) {
! 318: v = VR((P)ARG1(arg));
! 319: change_mvar(CO,p,v,&p1);
! 320: if ( VR(p1) != v ) {
! 321: *rp = p1; return;
! 322: } else
! 323: p = p1;
! 324: }
! 325: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ );
! 326: l = (P *)ALLOCA(m*sizeof(P));
! 327: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ )
! 328: l[m] = COEF(dc);
! 329: gcdsf(CO,l,m,rp);
! 330: } else if ( OID(obj) == O_DP ) {
! 331: for ( m = 0, mp = BDY((DP)obj); mp; mp = NEXT(mp), m++ );
! 332: l = (P *)ALLOCA(m*sizeof(P));
! 333: for ( m = 0, mp = BDY((DP)obj); mp; mp = NEXT(mp), m++)
! 334: l[m] = mp->c;
! 335: gcdsf(CO,l,m,rp);
! 336: }
1.1 noro 337: }
338:
339: void Pptozp(arg,rp)
340: NODE arg;
1.21 noro 341: Obj *rp;
1.1 noro 342: {
1.25 ! noro 343: Q t;
1.22 noro 344: NODE tt,p;
1.20 takayama 345: NODE n,n0;
346: char *key;
1.25 ! noro 347: P pp;
! 348: LIST list;
1.20 takayama 349: int get_factor=0;
350:
1.25 ! noro 351: asir_assert(ARG0(arg),O_P,"ptozp");
1.20 takayama 352:
353: /* analyze the option */
1.22 noro 354: if ( current_option ) {
355: for ( tt = current_option; tt; tt = NEXT(tt) ) {
1.20 takayama 356: p = BDY((LIST)BDY(tt));
357: key = BDY((STRING)BDY(p));
358: /* value = (Obj)BDY(NEXT(p)); */
359: if ( !strcmp(key,"factor") ) get_factor=1;
360: else {
361: error("ptozp: unknown option.");
362: }
363: }
364: }
365:
1.25 ! noro 366: ptozp((P)ARG0(arg),1,&t,&pp);
1.20 takayama 367:
368: /* printexpr(NULL,t); */
1.25 ! noro 369: /* if the option factor is given, then it returns the answer
1.20 takayama 370: in the format [zpoly, num] where num*zpoly is equal to the argument.*/
371: if (get_factor) {
1.25 ! noro 372: n0 = mknode(2,pp,t);
1.21 noro 373: MKLIST(list,n0);
1.25 ! noro 374: *rp = (Obj)list;
1.21 noro 375: } else
376: *rp = (Obj)pp;
1.1 noro 377: }
378:
379: void Pafctr(arg,rp)
380: NODE arg;
381: LIST *rp;
382: {
1.25 ! noro 383: DCP dc;
! 384:
! 385: asir_assert(ARG0(arg),O_P,"afctr");
! 386: asir_assert(ARG1(arg),O_P,"afctr");
! 387: afctr(CO,(P)ARG0(arg),(P)ARG1(arg),&dc);
! 388: dcptolist(dc,rp);
1.1 noro 389: }
390:
391: void Pagcd(arg,rp)
392: NODE arg;
393: P *rp;
394: {
1.25 ! noro 395: asir_assert(ARG0(arg),O_P,"agcd");
! 396: asir_assert(ARG1(arg),O_P,"agcd");
! 397: asir_assert(ARG2(arg),O_P,"agcd");
! 398: gcda(CO,(P)ARG0(arg),(P)ARG1(arg),(P)ARG2(arg),rp);
1.1 noro 399: }
400:
401: #if 1
402: #define Mulum mulum
403: #define Divum divum
404: #define Mulsum mulsum
405: #define Gcdum gcdum
406: #endif
407:
408: void Mulum(), Mulsum(), Gcdum();
409: int Divum();
410:
411: #define FCTR 0 /* berlekamp */
412: #define SQFR 1
413: #define DDD 2 /* Cantor-Zassenhauss */
414: #define NEWDDD 3 /* berlekamp + root-finding by Cantor-Zassenhauss */
415:
416: UM *resberle();
417:
1.18 noro 418: void reduce_sfdc(DCP sfdc, DCP *dc);
419:
1.1 noro 420: void Pmodfctr(arg,rp)
421: NODE arg;
422: LIST *rp;
423: {
1.25 ! noro 424: DCP dc,dcu;
! 425: int mod,i,t;
! 426: P p;
! 427: Obj u;
! 428: VL vl;
! 429:
! 430: mod = QTOS((Q)ARG1(arg));
! 431: if ( mod < 0 )
! 432: error("modfctr : invalid modulus");
! 433: p = (P)ARG0(arg);
! 434: clctv(CO,p,&vl);
! 435: if ( !vl ) {
! 436: NEWDC(dc); COEF(dc) = p; DEG(dc) = ONE; NEXT(dc) = 0;
! 437: } else if ( !NEXT(vl) )
! 438: modfctrp(ARG0(arg),mod,NEWDDD,&dc);
! 439: else {
! 440: /* XXX 16384 should be replaced by a macro */
! 441: for ( i = 1, t = mod; t*mod < 16384; t *= mod, i++ );
! 442: current_ff = FF_GFS;
! 443: setmod_sf(mod,i);
! 444: simp_ff((Obj)p,&u);
! 445: mfctrsf(CO,(P)u,&dcu);
! 446: reduce_sfdc(dcu,&dc);
! 447: }
! 448: if ( !dc ) {
! 449: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 450: }
! 451: dcptolist(dc,rp);
1.13 noro 452: }
453:
454: void Psfgcd(arg,rp)
455: NODE arg;
456: LIST *rp;
457: {
1.25 ! noro 458: P ps[2];
1.13 noro 459:
1.25 ! noro 460: ps[0] = (P)ARG0(arg);
! 461: ps[1] = (P)ARG1(arg);
! 462: gcdsf(CO,ps,2,rp);
1.6 noro 463: }
464:
1.15 noro 465: void Psffctr(arg,rp)
466: NODE arg;
467: LIST *rp;
468: {
1.25 ! noro 469: DCP dc;
1.15 noro 470:
1.25 ! noro 471: mfctrsf(CO,ARG0(arg),&dc);
! 472: dcptolist(dc,rp);
1.15 noro 473: }
474:
1.10 noro 475: void Psfsqfr(arg,rp)
476: NODE arg;
477: LIST *rp;
478: {
1.25 ! noro 479: DCP dc;
1.10 noro 480:
1.25 ! noro 481: sqfrsf(CO,ARG0(arg),&dc);
! 482: dcptolist(dc,rp);
1.10 noro 483: }
484:
485: void Psfufctr(arg,rp)
1.6 noro 486: NODE arg;
487: LIST *rp;
488: {
1.25 ! noro 489: DCP dc;
1.6 noro 490:
1.25 ! noro 491: ufctrsf(ARG0(arg),&dc);
! 492: dcptolist(dc,rp);
1.7 noro 493: }
494:
495: void Psfbfctr(arg,rp)
496: NODE arg;
497: LIST *rp;
498: {
1.25 ! noro 499: V x,y;
! 500: DCP dc,dct;
! 501: P t;
! 502: struct oVL vl1,vl2;
! 503: VL vl;
! 504: int degbound;
! 505:
! 506: x = VR((P)ARG1(arg));
! 507: y = VR((P)ARG2(arg));
! 508: vl1.v = x; vl1.next = &vl2;
! 509: vl2.v = y; vl2.next = 0;
! 510: vl = &vl1;
! 511: if ( argc(arg) == 4 )
! 512: degbound = QTOS((Q)ARG3(arg));
! 513: else
! 514: degbound = -1;
! 515:
! 516: sfbfctr((P)ARG0(arg),x,y,degbound,&dc);
! 517: for ( dct = dc; dct; dct = NEXT(dct) ) {
! 518: reorderp(CO,vl,COEF(dct),&t); COEF(dct) = t;
! 519: }
! 520: dcptolist(dc,rp);
1.1 noro 521: }
522:
1.11 noro 523: void Psfmintdeg(arg,rp)
524: NODE arg;
525: P *rp;
526: {
1.25 ! noro 527: V x,y;
! 528: P r;
! 529: struct oVL vl1,vl2;
! 530: VL vl;
! 531: int dy,c;
! 532:
! 533: x = VR((P)ARG1(arg));
! 534: y = VR((P)ARG2(arg));
! 535: vl1.v = x; vl1.next = &vl2;
! 536: vl2.v = y; vl2.next = 0;
! 537: vl = &vl1;
! 538: dy = QTOS((Q)ARG3(arg));
! 539: c = QTOS((Q)ARG4(arg));
! 540: sfmintdeg(vl,(P)ARG0(arg),dy,c,&r);
! 541: reorderp(CO,vl,r,rp);
1.11 noro 542: }
543:
1.1 noro 544: void Pmodsqfr(arg,rp)
545: NODE arg;
546: LIST *rp;
547: {
1.25 ! noro 548: DCP dc;
1.1 noro 549:
1.25 ! noro 550: if ( !ARG0(arg) ) {
! 551: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 552: } else
! 553: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),SQFR,&dc);
! 554: dcptolist(dc,rp);
1.1 noro 555: }
556:
557: void Pddd(arg,rp)
558: NODE arg;
559: LIST *rp;
560: {
1.25 ! noro 561: DCP dc;
1.1 noro 562:
1.25 ! noro 563: if ( !ARG0(arg) ) {
! 564: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 565: } else
! 566: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),DDD,&dc);
! 567: dcptolist(dc,rp);
1.1 noro 568: }
569:
570: void Pnewddd(arg,rp)
571: NODE arg;
572: LIST *rp;
573: {
1.25 ! noro 574: DCP dc=0;
1.1 noro 575:
1.25 ! noro 576: if ( !ARG0(arg) ) {
! 577: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
! 578: } else
! 579: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),NEWDDD,&dc);
! 580: dcptolist(dc,rp);
1.1 noro 581: }
582:
583: void Pirred_check(arg,rp)
584: NODE arg;
585: Q *rp;
586: {
1.25 ! noro 587: P p;
! 588: UM mp;
! 589: int r,mod;
! 590:
! 591: p = (P)ARG0(arg);
! 592: if ( !p ) {
! 593: *rp = 0; return;
! 594: }
! 595: mp = W_UMALLOC(UDEG(p));
! 596: mod = QTOS((Q)ARG1(arg));
! 597: ptoum(mod,p,mp);
! 598: r = irred_check(mp,mod);
! 599: if ( r )
! 600: *rp = ONE;
! 601: else
! 602: *rp = 0;
1.1 noro 603: }
604:
605: void Pnfctr_mod(arg,rp)
606: NODE arg;
607: Q *rp;
608: {
1.25 ! noro 609: P p;
! 610: UM mp;
! 611: int r,mod;
! 612:
! 613: p = (P)ARG0(arg);
! 614: if ( !p ) {
! 615: *rp = 0; return;
! 616: }
! 617: mp = W_UMALLOC(UDEG(p));
! 618: mod = QTOS((Q)ARG1(arg));
! 619: ptoum(mod,p,mp);
! 620: r = nfctr_mod(mp,mod);
! 621: STOQ(r,*rp);
1.1 noro 622: }
623:
624: void Pddd_tab(arg,rp)
625: NODE arg;
626: VECT *rp;
627: {
1.25 ! noro 628: P p;
! 629: UM mp,t,q,r1,w,w1;
! 630: UM *r,*s;
! 631: int dr,mod,n,i;
! 632: VECT result;
! 633: V v;
! 634:
! 635: p = (P)ARG0(arg); mod = QTOS((Q)ARG1(arg));
! 636: v = VR(p);
! 637: n = UDEG(p); mp = W_UMALLOC(n);
! 638: ptoum(mod,p,mp);
! 639: r = (UM *)W_ALLOC(n); s = (UM *)W_ALLOC(n);
! 640: r[0] = UMALLOC(0); DEG(r[0]) = 0; COEF(r[0])[0] = 1;
! 641: t = W_UMALLOC(mod); bzero(COEF(t),sizeof(int)*(mod+1));
! 642: DEG(t) = mod; COEF(t)[mod] = 1;
! 643: q = W_UMALLOC(mod);
! 644: dr = divum(mod,t,mp,q);
! 645: DEG(t) = dr; r[1] = r1 = UMALLOC(dr); cpyum(t,r1);
! 646: s[0] = W_UMALLOC(dr); cpyum(t,s[0]);
! 647: w = W_UMALLOC(n); bzero(COEF(w),sizeof(int)*(n+1));
! 648: w1 = W_UMALLOC(2*n); bzero(COEF(w1),sizeof(int)*(2*n+1));
! 649: for ( i = 1; i < n; i++ ) {
! 650: DEG(w) = i; COEF(w)[i-1] = 0; COEF(w)[i] = 1;
! 651: mulum(mod,r1,w,w1);
! 652: dr = divum(mod,w1,mp,q); DEG(w1) = dr;
! 653: s[i] = W_UMALLOC(dr); cpyum(w1,s[i]);
! 654: }
! 655: for ( i = 2; i < n; i++ ) {
! 656: mult_mod_tab(r[i-1],mod,s,w,n);
! 657: r[i] = UMALLOC(DEG(w)); cpyum(w,r[i]);
! 658: }
! 659: MKVECT(result,n);
! 660: for ( i = 0; i < n; i++ )
! 661: umtop(v,r[i],(P *)&BDY(result)[i]);
! 662: *rp = result;
1.18 noro 663: }
664:
665: void reduce_sfdc(DCP sfdc,DCP *dcr)
666: {
1.25 ! noro 667: P c,t,s,u,f;
! 668: DCP dc0,dc,tdc;
! 669: DCP *a;
! 670: int i,j,n;
! 671:
! 672: if ( !current_gfs_ext ) {
! 673: /* we simply apply sfptop() */
! 674: for ( dc0 = 0; sfdc; sfdc = NEXT(sfdc) ) {
! 675: NEXTDC(dc0,dc);
! 676: DEG(dc) = DEG(sfdc);
! 677: sfptop(COEF(sfdc),&COEF(dc));
! 678: }
! 679: NEXT(dc) = 0;
! 680: *dcr = dc0;
! 681: return;
! 682: }
! 683:
! 684: if ( NUM(COEF(sfdc)) ) {
! 685: sfptop(COEF(sfdc),&c);
! 686: sfdc = NEXT(sfdc);
! 687: } else
! 688: c = (P)ONE;
! 689:
! 690: for ( n = 0, tdc = sfdc; tdc; tdc = NEXT(tdc), n++ );
! 691: a = (DCP *)ALLOCA(n*sizeof(DCP));
! 692: for ( i = 0, tdc = sfdc; i < n; tdc = NEXT(tdc), i++ )
! 693: a[i] = tdc;
! 694:
! 695: dc0 = 0; NEXTDC(dc0,dc); DEG(dc) = ONE; COEF(dc) = c;
! 696: for ( i = 0; i < n; i++ ) {
! 697: if ( !a[i] )
! 698: continue;
! 699: t = COEF(a[i]);
! 700: f = t;
! 701: while ( 1 ) {
! 702: sf_galois_action(t,ONE,&s);
! 703: for ( j = i; j < n; j++ )
! 704: if ( a[j] && !compp(CO,s,COEF(a[j])) )
! 705: break;
! 706: if ( j == n )
! 707: error("reduce_sfdc : cannot happen");
! 708: if ( j == i ) {
! 709: NEXTDC(dc0,dc); DEG(dc) = DEG(a[i]);
! 710: sfptop(f,&COEF(dc));
! 711: break;
! 712: } else {
! 713: mulp(CO,f,s,&u); f = u;
! 714: t = s;
! 715: a[j] = 0;
! 716: }
! 717: }
! 718: }
! 719: *dcr = dc0;
1.1 noro 720: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>