Annotation of OpenXM_contrib2/asir2000/builtin/fctr.c, Revision 1.9
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.9 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/fctr.c,v 1.8 2001/06/26 03:00:40 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();
54: void Pptozp(), Pcont();
55: void Pafctr(), Pagcd();
56: void Pmodsqfr(),Pmodfctr(),Pddd(),Pnewddd(),Pddd_tab();
1.7 noro 57: void Psffctr(),Psfbfctr();
1.1 noro 58: void Pirred_check(), Pnfctr_mod();
59:
60: struct ftab fctr_tab[] = {
1.5 noro 61: {"fctr",Pfctr,-2},
1.1 noro 62: {"gcd",Pgcd,-3},
63: {"gcdz",Pgcdz,2},
64: {"lcm",Plcm,2},
65: {"sqfr",Psqfr,1},
66: {"ufctrhint",Pufctrhint,2},
67: {"ptozp",Pptozp,1},
68: {"cont",Pcont,-2},
69: {"afctr",Pafctr,2},
70: {"agcd",Pagcd,3},
71: {"modsqfr",Pmodsqfr,2},
72: {"modfctr",Pmodfctr,2},
1.6 noro 73: {"sffctr",Psffctr,1},
1.7 noro 74: {"sfbfctr",Psfbfctr,3},
1.1 noro 75: #if 0
76: {"ddd",Pddd,2},
77: {"newddd",Pnewddd,2},
78: #endif
79: {"ddd_tab",Pddd_tab,2},
80: {"irred_check",Pirred_check,2},
81: {"nfctr_mod",Pnfctr_mod,2},
82: {0,0,0},
83: };
84:
85: void Pfctr(arg,rp)
86: NODE arg;
87: LIST *rp;
88: {
89: DCP dc;
90:
91: asir_assert(ARG0(arg),O_P,"fctr");
1.5 noro 92: if ( argc(arg) == 1 )
93: fctrp(CO,(P)ARG0(arg),&dc);
94: else {
95: asir_assert(ARG1(arg),O_P,"fctr");
96: fctr_wrt_v_p(CO,(P)ARG0(arg),VR((P)ARG1(arg)),&dc);
97: }
1.1 noro 98: dcptolist(dc,rp);
99: }
100:
101: void Pgcd(arg,rp)
102: NODE arg;
103: P *rp;
104: {
105: P p1,p2,g1,g2,g;
106: Num m;
107: int mod;
108:
109: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
110: asir_assert(p1,O_P,"gcd");
111: asir_assert(p2,O_P,"gcd");
112: if ( !p1 )
113: *rp = p2;
114: else if ( !p2 )
115: *rp = p1;
116: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
1.4 noro 117: gcdprsp(CO,p1,p2,rp);
1.1 noro 118: else if ( argc(arg) == 2 )
119: ezgcdp(CO,p1,p2,rp);
120: else {
121: m = (Num)ARG2(arg);
122: asir_assert(m,O_P,"gcd");
123: mod = QTOS((Q)m);
124: ptomp(mod,p1,&g1); ptomp(mod,p2,&g2);
125: gcdprsmp(CO,mod,g1,g2,&g);
126: mptop(g,rp);
127: }
128: }
129:
130: void Pgcdz(arg,rp)
131: NODE arg;
132: P *rp;
133: {
134: P p1,p2,t;
135: Q c1,c2;
136: N n;
137:
138: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
139: asir_assert(p1,O_P,"gcdz");
140: asir_assert(p2,O_P,"gcdz");
141: if ( !p1 )
142: *rp = p2;
143: else if ( !p2 )
144: *rp = p1;
145: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
146: error("gcdz : invalid argument");
147: else if ( NUM(p1) || NUM(p2) ) {
148: if ( NUM(p1) )
149: c1 = (Q)p1;
150: else
151: ptozp(p1,1,&c1,&t);
152: if ( NUM(p2) )
153: c2 = (Q)p2;
154: else
155: ptozp(p2,1,&c2,&t);
156: gcdn(NM(c1),NM(c2),&n); NTOQ(n,1,c1); *rp = (P)c1;
157: } else {
158: #if 0
159: w[0] = p1; w[1] = p2; nezgcdnpz(CO,w,2,rp);
160: #endif
161: ezgcdpz(CO,p1,p2,rp);
162: }
163: }
164:
165: void Plcm(arg,rp)
166: NODE arg;
167: P *rp;
168: {
169: P t1,t2,p1,p2,g,q;
170: Q c;
171:
172: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
173: asir_assert(p1,O_P,"lcm");
174: asir_assert(p2,O_P,"lcm");
175: if ( !p1 || !p2 )
176: *rp = 0;
177: else if ( !qpcheck((Obj)p1) || !qpcheck((Obj)p2) )
178: error("lcm : invalid argument");
179: else {
180: ptozp(p1,1,&c,&t1); ptozp(p2,1,&c,&t2);
181: ezgcdp(CO,t1,t2,&g); divsp(CO,t1,g,&q); mulp(CO,q,t2,rp);
182: }
183: }
184:
185: void Psqfr(arg,rp)
186: NODE arg;
187: LIST *rp;
188: {
189: DCP dc;
190:
191: asir_assert(ARG0(arg),O_P,"sqfr");
192: sqfrp(CO,(P)ARG0(arg),&dc);
193: dcptolist(dc,rp);
194: }
195:
196: void Pufctrhint(arg,rp)
197: NODE arg;
198: LIST *rp;
199: {
200: DCP dc;
201:
202: asir_assert(ARG0(arg),O_P,"ufctrhint");
203: asir_assert(ARG1(arg),O_N,"ufctrhint");
204: ufctr((P)ARG0(arg),QTOS((Q)ARG1(arg)),&dc);
205: dcptolist(dc,rp);
206: }
207:
208: #if 0
209: Pmgcd(arg,rp)
210: NODE arg;
211: Obj *rp;
212: {
213: NODE node,tn;
214: int i,m;
215: P *l;
216:
217: node = BDY((LIST)ARG0(arg));
218: for ( i = 0, tn = node; tn; tn = NEXT(tn), i++ );
219: m = i; l = (P *)ALLOCA(m*sizeof(P));
220: for ( i = 0, tn = node; i < m; tn = NEXT(tn), i++ )
221: l[i] = (P)BDY(tn);
222: nezgcdnpz(CO,l,m,rp);
223: }
224: #endif
225:
226: void Pcont(arg,rp)
227: NODE arg;
228: P *rp;
229: {
230: DCP dc;
231: int m;
232: P p,p1;
233: P *l;
234: V v;
235:
236: asir_assert(ARG0(arg),O_P,"cont");
237: p = (P)ARG0(arg);
238: if ( NUM(p) )
239: *rp = p;
240: else {
241: if ( argc(arg) == 2 ) {
242: v = VR((P)ARG1(arg));
243: change_mvar(CO,p,v,&p1);
244: if ( VR(p1) != v ) {
245: *rp = p1; return;
246: } else
247: p = p1;
248: }
249: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ );
250: l = (P *)ALLOCA(m*sizeof(P));
251: for ( m = 0, dc = DC(p); dc; dc = NEXT(dc), m++ )
252: l[m] = COEF(dc);
253: nezgcdnpz(CO,l,m,rp);
254: }
255: }
256:
257: void Pptozp(arg,rp)
258: NODE arg;
259: P *rp;
260: {
261: Q t;
262:
263: asir_assert(ARG0(arg),O_P,"ptozp");
264: ptozp((P)ARG0(arg),1,&t,rp);
265: }
266:
267: void Pafctr(arg,rp)
268: NODE arg;
269: LIST *rp;
270: {
271: DCP dc;
272:
273: asir_assert(ARG0(arg),O_P,"afctr");
274: asir_assert(ARG1(arg),O_P,"afctr");
275: afctr(CO,(P)ARG0(arg),(P)ARG1(arg),&dc);
276: dcptolist(dc,rp);
277: }
278:
279: void Pagcd(arg,rp)
280: NODE arg;
281: P *rp;
282: {
283: asir_assert(ARG0(arg),O_P,"agcd");
284: asir_assert(ARG1(arg),O_P,"agcd");
285: asir_assert(ARG2(arg),O_P,"agcd");
286: gcda(CO,(P)ARG0(arg),(P)ARG1(arg),(P)ARG2(arg),rp);
287: }
288:
289: #if 1
290: #define Mulum mulum
291: #define Divum divum
292: #define Mulsum mulsum
293: #define Gcdum gcdum
294: #endif
295:
296: void Mulum(), Mulsum(), Gcdum();
297: int Divum();
298:
299: #define FCTR 0 /* berlekamp */
300: #define SQFR 1
301: #define DDD 2 /* Cantor-Zassenhauss */
302: #define NEWDDD 3 /* berlekamp + root-finding by Cantor-Zassenhauss */
303:
304: UM *resberle();
305:
306: void Pmodfctr(arg,rp)
307: NODE arg;
308: LIST *rp;
309: {
310: DCP dc;
311: int mod;
312:
313: mod = QTOS((Q)ARG1(arg));
314: if ( mod < 0 )
315: error("modfctr : invalid modulus");
316: modfctrp(ARG0(arg),mod,NEWDDD,&dc);
1.6 noro 317: if ( !dc ) {
318: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
319: }
320: dcptolist(dc,rp);
321: }
322:
323: void Psffctr(arg,rp)
324: NODE arg;
325: LIST *rp;
326: {
327: DCP dc;
328:
329: fctrsf(ARG0(arg),&dc);
1.1 noro 330: if ( !dc ) {
331: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
332: }
333: dcptolist(dc,rp);
1.7 noro 334: }
335:
336: void Psfbfctr(arg,rp)
337: NODE arg;
338: LIST *rp;
339: {
340: V x,y;
1.8 noro 341: DCP dc,dct;
342: P t;
343: struct oVL vl1,vl2;
344: VL vl;
1.7 noro 345:
346: x = VR((P)ARG1(arg));
347: y = VR((P)ARG2(arg));
1.8 noro 348: vl1.v = x; vl1.next = &vl2;
349: vl2.v = y; vl2.next = 0;
350: vl = &vl1;
351:
352: sfbfctr((P)ARG0(arg),x,y,&dc);
353: for ( dct = dc; dct; dct = NEXT(dct) ) {
354: reorderp(CO,vl,COEF(dct),&t); COEF(dct) = t;
1.7 noro 355: }
1.8 noro 356: dcptolist(dc,rp);
1.1 noro 357: }
358:
359: void Pmodsqfr(arg,rp)
360: NODE arg;
361: LIST *rp;
362: {
363: DCP dc;
364:
1.9 ! noro 365: if ( !ARG0(arg) ) {
1.1 noro 366: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
1.9 ! noro 367: } else
! 368: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),SQFR,&dc);
1.1 noro 369: dcptolist(dc,rp);
370: }
371:
372: void Pddd(arg,rp)
373: NODE arg;
374: LIST *rp;
375: {
376: DCP dc;
377:
1.9 ! noro 378: if ( !ARG0(arg) ) {
1.1 noro 379: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
1.9 ! noro 380: } else
! 381: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),DDD,&dc);
1.1 noro 382: dcptolist(dc,rp);
383: }
384:
385: void Pnewddd(arg,rp)
386: NODE arg;
387: LIST *rp;
388: {
1.9 ! noro 389: DCP dc=0;
1.1 noro 390:
1.9 ! noro 391: if ( !ARG0(arg) ) {
1.1 noro 392: NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; NEXT(dc) = 0;
1.9 ! noro 393: } else
! 394: modfctrp(ARG0(arg),QTOS((Q)ARG1(arg)),NEWDDD,&dc);
1.1 noro 395: dcptolist(dc,rp);
396: }
397:
398: void Pirred_check(arg,rp)
399: NODE arg;
400: Q *rp;
401: {
402: P p;
403: UM mp;
404: int r,mod;
405:
406: p = (P)ARG0(arg);
407: if ( !p ) {
408: *rp = 0; return;
409: }
410: mp = W_UMALLOC(UDEG(p));
411: mod = QTOS((Q)ARG1(arg));
412: ptoum(mod,p,mp);
413: r = irred_check(mp,mod);
414: if ( r )
415: *rp = ONE;
416: else
417: *rp = 0;
418: }
419:
420: void Pnfctr_mod(arg,rp)
421: NODE arg;
422: Q *rp;
423: {
424: P p;
425: UM mp;
426: int r,mod;
427:
428: p = (P)ARG0(arg);
429: if ( !p ) {
430: *rp = 0; return;
431: }
432: mp = W_UMALLOC(UDEG(p));
433: mod = QTOS((Q)ARG1(arg));
434: ptoum(mod,p,mp);
435: r = nfctr_mod(mp,mod);
436: STOQ(r,*rp);
437: }
438:
439: void Pddd_tab(arg,rp)
440: NODE arg;
441: VECT *rp;
442: {
443: P p;
444: UM mp,t,q,r1,w,w1;
445: UM *r,*s;
446: int dr,mod,n,i;
447: VECT result;
448: V v;
449:
450: p = (P)ARG0(arg); mod = QTOS((Q)ARG1(arg));
451: v = VR(p);
452: n = UDEG(p); mp = W_UMALLOC(n);
453: ptoum(mod,p,mp);
454: r = (UM *)W_ALLOC(n); s = (UM *)W_ALLOC(n);
455: r[0] = UMALLOC(0); DEG(r[0]) = 0; COEF(r[0])[0] = 1;
456: t = W_UMALLOC(mod); bzero(COEF(t),sizeof(int)*(mod+1));
457: DEG(t) = mod; COEF(t)[mod] = 1;
458: q = W_UMALLOC(mod);
459: dr = divum(mod,t,mp,q);
460: DEG(t) = dr; r[1] = r1 = UMALLOC(dr); cpyum(t,r1);
461: s[0] = W_UMALLOC(dr); cpyum(t,s[0]);
462: w = W_UMALLOC(n); bzero(COEF(w),sizeof(int)*(n+1));
463: w1 = W_UMALLOC(2*n); bzero(COEF(w1),sizeof(int)*(2*n+1));
464: for ( i = 1; i < n; i++ ) {
465: DEG(w) = i; COEF(w)[i-1] = 0; COEF(w)[i] = 1;
466: mulum(mod,r1,w,w1);
467: dr = divum(mod,w1,mp,q); DEG(w1) = dr;
468: s[i] = W_UMALLOC(dr); cpyum(w1,s[i]);
469: }
470: for ( i = 2; i < n; i++ ) {
471: mult_mod_tab(r[i-1],mod,s,w,n);
472: r[i] = UMALLOC(DEG(w)); cpyum(w,r[i]);
473: }
474: MKVECT(result,n);
475: for ( i = 0; i < n; i++ )
476: umtop(v,r[i],(P *)&BDY(result)[i]);
477: *rp = result;
478: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to fctr.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / builtin