![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to poly.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / builtin|
Return to poly.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / builtin |
1.4 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.5 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.4 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.16 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/poly.c,v 1.15 2001/09/03 07:01:05 noro Exp $
1.4 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include "parse.h"
52: #include "base.h"
53:
54: void Pranp();
55:
56: void Pumul(),Pumul_ff(),Pusquare(),Pusquare_ff(),Putmul(),Putmul_ff();
57: void Pkmul(),Pksquare(),Pktmul();
58: void Pord(), Pcoef0(), Pcoef(), Pdeg(), Pmindeg(), Psetmod();
59: void Pcoef_gf2n();
60: void getcoef(), getdeglist(), mergedeglist(), change_mvar(), restore_mvar();
61:
1.2 noro 62: void Pp_mag(),Pmaxblen();
1.1 noro 63: void Pmergelist(), Pch_mv(), Pre_mv(), Pdeglist();
64: void Pptomp(),Pmptop();
65: void Pptolmp(),Plmptop();
1.12 noro 66: void Pptosfp(),Psfptop(),Psf_galois_action(),Psf_embed(),Psf_find_root();
1.13 noro 67: void Psf_minipoly(),Psf_log();
1.1 noro 68: void Pptogf2n(),Pgf2ntop(),Pgf2ntovect();
69: void Pptogfpn(),Pgfpntop();
70: void Pfind_root_gf2n();
1.2 noro 71: void Pumul_specialmod(),Pusquare_specialmod(),Putmul_specialmod();
1.1 noro 72:
73: void Pureverse(),Putrunc(),Pudecomp(),Purembymul(),Purembymul_precomp();
74: void Puinvmod(),Purevinvmod();
75: void Ppwrmod_ff(),Ppwrtab_ff(),Pgeneric_pwrmod_ff();
76: void Pkpwrmod_lm(),Pkpwrtab_lm(),Pkgeneric_pwrmod_lm();
77:
78: void Pkmulum();
79: void Pksquareum();
80:
81: void Pfmultest();
82: void Phfmul_lm();
83: void Plazy_lm();
84:
85: void Psetmod_ff();
86: void Psimp_ff();
87: void Pextdeg_ff();
88: void Pcharacteristic_ff();
89: void Pfield_type_ff();
90: void Pfield_order_ff();
91: void Prandom_ff();
92:
93: void Ptracemod_gf2n();
94: void Psparsemod_gf2n();
95: void Pmultest_gf2n();
96: void Psquaretest_gf2n();
97: void Pinvtest_gf2n();
98: void Pbininv_gf2n();
99: void Prinvtest_gf2n();
100: void Pis_irred_gf2();
101: void Pis_irred_ddd_gf2();
1.2 noro 102: void Pget_next_fft_prime();
1.3 noro 103: void Puadj_coef();
1.9 noro 104: void Preorder();
105: void Phomogeneous_part();
106: void Phomogeneous_deg();
1.1 noro 107: void simp_ff(Obj,Obj *);
108: void ranp(int,UP *);
109: void field_order_ff(N *);
110:
111: int current_ff;
112:
113: struct ftab poly_tab[] = {
1.9 noro 114: {"homogeneous_deg",Phomogeneous_deg,-2},
115: {"homogeneous_part",Phomogeneous_part,-3},
116: {"reorder",Preorder,3},
1.3 noro 117: {"uadj_coef",Puadj_coef,3},
1.1 noro 118: {"ranp",Pranp,2},
119: {"p_mag",Pp_mag,1},
1.2 noro 120: {"maxblen",Pmaxblen,1},
1.1 noro 121: {"ord",Pord,-1},
122: {"coef0",Pcoef0,-3},
123: {"coef",Pcoef,-3},
124: {"coef_gf2n",Pcoef_gf2n,2},
125: {"deg",Pdeg,2},
126: {"mindeg",Pmindeg,2},
127: {"setmod",Psetmod,-1},
128:
129: {"sparsemod_gf2n",Psparsemod_gf2n,-1},
130:
1.15 noro 131: {"setmod_ff",Psetmod_ff,-3},
1.1 noro 132: {"simp_ff",Psimp_ff,1},
133: {"extdeg_ff",Pextdeg_ff,0},
134: {"characteristic_ff",Pcharacteristic_ff,0},
135: {"field_type_ff",Pfield_type_ff,0},
136: {"field_order_ff",Pfield_order_ff,0},
137: {"random_ff",Prandom_ff,0},
138:
139: {"deglist",Pdeglist,2},
140: {"mergelist",Pmergelist,2},
141: {"ch_mv",Pch_mv,2},
142: {"re_mv",Pre_mv,2},
143:
144: {"ptomp",Pptomp,2},
145: {"mptop",Pmptop,1},
146:
147: {"ptolmp",Pptolmp,1},
148: {"lmptop",Plmptop,1},
149:
1.11 noro 150: {"sf_galois_action",Psf_galois_action,2},
1.12 noro 151: {"sf_find_root",Psf_find_root,1},
152: {"sf_minipoly",Psf_minipoly,2},
153: {"sf_embed",Psf_embed,3},
1.13 noro 154: {"sf_log",Psf_log,1},
155:
1.11 noro 156: {"ptosfp",Pptosfp,1},
1.6 noro 157: {"sfptop",Psfptop,1},
1.1 noro 158: {"ptogf2n",Pptogf2n,1},
159: {"gf2ntop",Pgf2ntop,-2},
160: {"gf2ntovect",Pgf2ntovect,1},
161:
162: {"ptogfpn",Pptogfpn,1},
163: {"gfpntop",Pgfpntop,-2},
164:
165: {"kmul",Pkmul,2},
166: {"ksquare",Pksquare,1},
167: {"ktmul",Pktmul,3},
168:
169: {"umul",Pumul,2},
170: {"usquare",Pusquare,1},
171: {"ureverse_inv_as_power_series",Purevinvmod,2},
172: {"uinv_as_power_series",Puinvmod,2},
173:
1.2 noro 174: {"umul_specialmod",Pumul_specialmod,3},
175: {"usquare_specialmod",Pusquare_specialmod,2},
176: {"utmul_specialmod",Putmul_specialmod,4},
177:
1.1 noro 178: {"utmul",Putmul,3},
179: {"umul_ff",Pumul_ff,2},
180: {"usquare_ff",Pusquare_ff,1},
181: {"utmul_ff",Putmul_ff,3},
182:
183: /* for historical reason */
184: {"trunc",Putrunc,2},
185: {"decomp",Pudecomp,2},
186:
187: {"utrunc",Putrunc,2},
188: {"udecomp",Pudecomp,2},
1.3 noro 189: {"ureverse",Pureverse,-2},
1.1 noro 190: {"urembymul",Purembymul,2},
191: {"urembymul_precomp",Purembymul_precomp,3},
192:
193: {"lazy_lm",Plazy_lm,1},
194: {"lazy_ff",Plazy_lm,1},
195:
196: {"pwrmod_ff",Ppwrmod_ff,1},
197: {"generic_pwrmod_ff",Pgeneric_pwrmod_ff,3},
198: {"pwrtab_ff",Ppwrtab_ff,2},
199:
200: {"tracemod_gf2n",Ptracemod_gf2n,3},
201: {"b_find_root_gf2n",Pfind_root_gf2n,1},
202:
203: {"is_irred_gf2",Pis_irred_gf2,1},
204: {"is_irred_ddd_gf2",Pis_irred_ddd_gf2,1},
205:
206: {"kpwrmod_lm",Pkpwrmod_lm,1},
207: {"kgeneric_pwrmod_lm",Pkgeneric_pwrmod_lm,3},
208: {"kpwrtab_lm",Pkpwrtab_lm,2},
209:
210: {"kmulum",Pkmulum,3},
211: {"ksquareum",Pksquareum,2},
212:
213: {"fmultest",Pfmultest,3},
214: {"hfmul_lm",Phfmul_lm,2},
215:
216: {"multest_gf2n",Pmultest_gf2n,2},
217: {"squaretest_gf2n",Psquaretest_gf2n,1},
218: {"bininv_gf2n",Pbininv_gf2n,2},
219: {"invtest_gf2n",Pinvtest_gf2n,1},
220: {"rinvtest_gf2n",Prinvtest_gf2n,0},
1.2 noro 221: {"get_next_fft_prime",Pget_next_fft_prime,2},
1.1 noro 222: {0,0,0},
223: };
224:
1.16 ! noro 225: void Phomogeneous_part(NODE arg,P *rp)
1.9 noro 226: {
227: if ( argc(arg) == 2 )
228: exthp(CO,(P)ARG0(arg),QTOS((Q)ARG1(arg)),rp);
229: else
230: exthpc_generic(CO,(P)ARG0(arg),QTOS((Q)ARG2(arg)),
231: VR((P)ARG1(arg)),rp);
232: }
233:
1.16 ! noro 234: void Phomogeneous_deg(NODE arg,Q *rp)
1.9 noro 235: {
236: int d;
237:
238: if ( argc(arg) == 1 )
239: d = homdeg((P)ARG0(arg));
240: else
241: d = getchomdeg(VR((P)ARG1(arg)),(P)ARG0(arg));
242: STOQ(d,*rp);
243: }
244:
245: /*
246: p1 = reorder(p,ovl,nvl) => p1 is 'sorted accoding to nvl.
247: */
248:
1.16 ! noro 249: void Preorder(NODE arg,P *rp)
1.9 noro 250: {
251: VL ovl,nvl,tvl;
252: NODE n;
253:
254: for ( ovl = 0, n = BDY((LIST)ARG1(arg)); n; n = NEXT(n) ) {
255: if ( !ovl ) {
256: NEWVL(ovl); tvl = ovl;
257: } else {
258: NEWVL(NEXT(tvl)); tvl = NEXT(tvl);
259: }
260: VR(tvl) = VR((P)BDY(n));
261: }
262: for ( nvl = 0, n = BDY((LIST)ARG2(arg)); n; n = NEXT(n) ) {
263: if ( !nvl ) {
264: NEWVL(nvl); tvl = nvl;
265: } else {
266: NEWVL(NEXT(tvl)); tvl = NEXT(tvl);
267: }
268: VR(tvl) = VR((P)BDY(n));
269: }
270: reorderp(nvl,ovl,(P)ARG0(arg),rp);
271: }
1.3 noro 272:
273: /*
274: uadj_coef(F,M,M2)
275: if ( F is a non-negative integer )
276: return F > M2 ? F-M : M2;
277: else
278: F = CN*V^N+...+C0
279: return uadj_coef(CN,M,M2)*V^N+...+uadj_coef(C0,M,M2);
280: */
281:
1.16 ! noro 282: void Puadj_coef(NODE arg,P *rp)
1.3 noro 283: {
284: UP f,r;
285: N m,m2;
286:
287: ptoup((P)ARG0(arg),&f);
288: m = NM((Q)ARG1(arg));
289: m2 = NM((Q)ARG2(arg));
290: adj_coefup(f,m,m2,&r);
291: uptop(r,rp);
292: }
293:
1.2 noro 294: /*
295: get_next_fft_prime(StartIndex,Bits)
296: tries to find smallest Index >= StartIndex s.t.
297: 2^(Bits-1)|FFTprime[Index]-1
298: return [Index,Mod] or 0 (not exist)
299: */
300:
1.16 ! noro 301: void Pget_next_fft_prime(NODE arg,LIST *rp)
1.2 noro 302: {
303: unsigned int mod,d;
304: int start,bits,i;
305: NODE n;
306: Q q,ind;
307:
308: start = QTOS((Q)ARG0(arg));
309: bits = QTOS((Q)ARG1(arg));
310: for ( i = start; ; i++ ) {
311: get_fft_prime(i,&mod,&d);
312: if ( !mod ) {
313: *rp = 0; return;
314: }
1.16 ! noro 315: if ( bits <= (int)d ) {
1.2 noro 316: UTOQ(mod,q);
317: UTOQ(i,ind);
318: n = mknode(2,ind,q);
319: MKLIST(*rp,n);
320: return;
321: }
322: }
323: }
1.1 noro 324:
1.16 ! noro 325: void Pranp(NODE arg,P *rp)
1.1 noro 326: {
327: int n;
328: UP c;
329:
330: n = QTOS((Q)ARG0(arg));
331: ranp(n,&c);
332: if ( c ) {
333: up_var = VR((P)ARG1(arg));
334: uptop(c,rp);
335: } else
336: *rp = 0;
337: }
338:
1.16 ! noro 339: void ranp(int n,UP *nr)
1.1 noro 340: {
341: int i;
342: unsigned int r;
343: Q q;
344: UP c;
345:
346: *nr = c = UPALLOC(n);
347: for ( i = 0; i <= n; i++ ) {
348: r = random();
349: UTOQ(r,q);
350: c->c[i] = (Num)q;
351: }
352: for ( i = n; i >= 0 && !c->c[i]; i-- );
353: if ( i >= 0 )
354: c->d = i;
355: else
356: *nr = 0;
357: }
358:
1.16 ! noro 359: void Pmaxblen(NODE arg,Q *rp)
1.2 noro 360: {
361: int l;
362: l = maxblenp(ARG0(arg));
363: STOQ(l,*rp);
364: }
365:
1.16 ! noro 366: void Pp_mag(NODE arg,Q *rp)
1.1 noro 367: {
368: int l;
369: l = p_mag(ARG0(arg));
370: STOQ(l,*rp);
371: }
372:
1.16 ! noro 373: void Pord(NODE arg,LIST *listp)
1.1 noro 374: {
375: NODE n,tn;
376: LIST l;
377: VL vl,tvl,svl;
378: P t;
379: int i,j;
380: V *va;
381: V v;
382:
383: if ( argc(arg) ) {
384: asir_assert(ARG0(arg),O_LIST,"ord");
385: for ( vl = 0, i = 0, n = BDY((LIST)ARG0(arg));
386: n; n = NEXT(n), i++ ) {
387: if ( !vl ) {
388: NEWVL(vl); tvl = vl;
389: } else {
390: NEWVL(NEXT(tvl)); tvl = NEXT(tvl);
391: }
392: if ( !(t = (P)BDY(n)) || (OID(t) != O_P) )
393: error("ord : invalid argument");
394: VR(tvl) = VR(t);
395: }
396: va = (V *)ALLOCA(i*sizeof(V));
397: for ( j = 0, svl = vl; j < i; j++, svl = NEXT(svl) )
398: va[j] = VR(svl);
399: for ( svl = CO; svl; svl = NEXT(svl) ) {
400: v = VR(svl);
401: for ( j = 0; j < i; j++ )
402: if ( v == va[j] )
403: break;
404: if ( j == i ) {
405: if ( !vl ) {
406: NEWVL(vl); tvl = vl;
407: } else {
408: NEWVL(NEXT(tvl)); tvl = NEXT(tvl);
409: }
410: VR(tvl) = v;
411: }
412: }
413: if ( vl )
414: NEXT(tvl) = 0;
415: CO = vl;
416: }
417: for ( n = 0, vl = CO; vl; vl = NEXT(vl) ) {
418: NEXTNODE(n,tn); MKV(VR(vl),t); BDY(tn) = (pointer)t;
419: }
420: NEXT(tn) = 0; MKLIST(l,n); *listp = l;
421: }
422:
1.16 ! noro 423: void Pcoef0(NODE arg,Obj *rp)
1.1 noro 424: {
425: Obj t,n;
426: P s;
427: DCP dc;
428: int id;
429: V v;
430: VL vl;
431:
432: if ( !(t = (Obj)ARG0(arg)) || ((id = OID(ARG0(arg))) > O_P) )
433: *rp = 0;
434: else if ( (n = (Obj)ARG1(arg)) && (OID(n) > O_N) )
435: *rp = 0;
436: else if ( id == O_N )
437: if ( !n )
438: *rp = t;
439: else
440: *rp = 0;
441: else {
442: if ( argc(arg) == 3 ) {
443: if ( (v = VR((P)ARG2(arg))) != VR((P)t) ) {
444: reordvar(CO,v,&vl); reorderp(vl,CO,(P)t,&s);
445: } else
446: s = (P)t;
447: if ( VR(s) != v ) {
448: if ( n )
449: *rp = 0;
450: else
451: *rp = t;
452: return;
453: }
454: } else
455: s = (P)t;
456: for ( dc = DC(s); dc && cmpq(DEG(dc),(Q)n); dc = NEXT(dc) );
457: if ( dc )
458: *rp = (Obj)COEF(dc);
459: else
460: *rp = 0;
461: }
462: }
463:
1.16 ! noro 464: void Pcoef(NODE arg,Obj *rp)
1.1 noro 465: {
466: Obj t,n;
467: P s;
468: DCP dc;
469: int id;
470: V v;
471:
472: if ( !(t = (Obj)ARG0(arg)) || ((id = OID(ARG0(arg))) > O_P) )
473: *rp = 0;
474: else if ( (n = (Obj)ARG1(arg)) && (OID(n) > O_N) )
475: *rp = 0;
476: else if ( id == O_N ) {
477: if ( !n )
478: *rp = t;
479: else
480: *rp = 0;
481: } else {
482: if ( argc(arg) == 3 ) {
483: if ( (v = VR((P)ARG2(arg))) != VR((P)t) ) {
484: getcoef(CO,(P)t,v,(Q)n,(P *)rp); return;
485: } else
486: s = (P)t;
487: if ( VR(s) != v ) {
488: if ( n )
489: *rp = 0;
490: else
491: *rp = t;
492: return;
493: }
494: } else
495: s = (P)t;
496: for ( dc = DC(s); dc && cmpq(DEG(dc),(Q)n); dc = NEXT(dc) );
497: if ( dc )
498: *rp = (Obj)COEF(dc);
499: else
500: *rp = 0;
501: }
502: }
503:
1.16 ! noro 504: void Pcoef_gf2n(NODE arg,Obj *rp)
1.1 noro 505: {
506: Obj t,n;
507: int id,d;
508: UP2 up2;
509:
510: if ( !(t = (Obj)ARG0(arg)) || ((id = OID(ARG0(arg))) > O_P) )
511: *rp = 0;
512: else if ( (n = (Obj)ARG1(arg)) && (OID(n) > O_N) )
513: *rp = 0;
514: else if ( id == O_N && NID((Num)t) == N_GF2N ) {
515: d = QTOS((Q)n);
516: up2 = ((GF2N)t)->body;
517: if ( d > degup2(up2) )
518: *rp = 0;
519: else
520: *rp = (Obj)(up2->b[d/BSH]&(((unsigned long)1)<<(d%BSH))?ONE:0);
521: } else
522: *rp = 0;
523: }
524:
1.16 ! noro 525: void Pdeg(NODE arg,Q *rp)
1.1 noro 526: {
527: Obj t,v;
528: int d;
529:
530: #if 0
531: if ( !(t = (Obj)ARG0(arg)) || (OID(t) != O_P) )
532: *rp = 0;
533: else if ( !(v = (Obj)ARG1(arg)) || (VR((P)v) != VR((P)t)) )
534: *rp = 0;
535: else
536: *rp = (Obj)DEG(DC((P)t));
537: #endif
538: if ( !(t = (Obj)ARG0(arg)) )
539: STOQ(-1,*rp);
540: else if ( OID(t) != O_P ) {
541: if ( OID(t) == O_N && NID(t) == N_GF2N
542: && (v=(Obj)ARG1(arg)) && OID(v)== O_N && NID(v) == N_GF2N ) {
543: d = degup2(((GF2N)t)->body);
544: STOQ(d,*rp);
545: } else
546: *rp = 0;
547: } else
548: degp(VR((P)ARG1(arg)),(P)ARG0(arg),rp);
549: }
550:
1.16 ! noro 551: void Pmindeg(NODE arg,Q *rp)
1.1 noro 552: {
553: Obj t;
554:
555: if ( !(t = (Obj)ARG0(arg)) || (OID(t) != O_P) )
556: *rp = 0;
557: else
558: getmindeg(VR((P)ARG1(arg)),(P)ARG0(arg),rp);
559: }
560:
1.16 ! noro 561: void Psetmod(NODE arg,Q *rp)
1.1 noro 562: {
563: if ( arg ) {
564: asir_assert(ARG0(arg),O_N,"setmod");
565: current_mod = QTOS((Q)ARG0(arg));
566: }
567: STOQ(current_mod,*rp);
568: }
569:
1.16 ! noro 570: void Psparsemod_gf2n(NODE arg,Q *rp)
1.1 noro 571: {
572: int id;
573:
574: if ( arg && current_mod_gf2n )
575: current_mod_gf2n->id = ARG0(arg)?1:0;
576: if ( !current_mod_gf2n )
577: id = -1;
578: else
579: id = current_mod_gf2n->id;
580: STOQ(id,*rp);
581: }
582:
1.16 ! noro 583: void Pmultest_gf2n(NODE arg,GF2N *rp)
1.1 noro 584: {
585: GF2N a,b,c;
586: int i;
587:
588: a = (GF2N)ARG0(arg); b = (GF2N)ARG0(arg);
589: for ( i = 0; i < 10000; i++ )
590: mulgf2n(a,b,&c);
591: *rp = c;
592: }
593:
1.16 ! noro 594: void Psquaretest_gf2n(NODE arg,GF2N *rp)
1.1 noro 595: {
596: GF2N a,c;
597: int i;
598:
599: a = (GF2N)ARG0(arg);
600: for ( i = 0; i < 10000; i++ )
601: squaregf2n(a,&c);
602: *rp = c;
603: }
604:
1.16 ! noro 605: void Pinvtest_gf2n(NODE arg,GF2N *rp)
1.1 noro 606: {
607: GF2N a,c;
608: int i;
609:
610: a = (GF2N)ARG0(arg);
611: for ( i = 0; i < 10000; i++ )
612: invgf2n(a,&c);
613: *rp = c;
614: }
615:
1.16 ! noro 616: void Pbininv_gf2n(NODE arg,GF2N *rp)
1.1 noro 617: {
618: UP2 a,inv;
619: int n;
620:
621: a = ((GF2N)ARG0(arg))->body;
622: n = QTOS((Q)ARG1(arg));
623: type1_bin_invup2(a,n,&inv);
624: MKGF2N(inv,*rp);
625: }
626:
1.16 ! noro 627: void Prinvtest_gf2n(Real *rp)
1.1 noro 628: {
629: GF2N *a;
630: GF2N c;
631: double t0,t1,r;
632: int i;
633: double get_clock();
634:
635: a = (GF2N *)ALLOCA(1000*sizeof(GF2N));
636: for ( i = 0; i < 1000; i++ ) {
637: randomgf2n(&a[i]);
638: }
639: t0 = get_clock();
640: for ( i = 0; i < 1000; i++ )
641: invgf2n(a[i],&c);
642: t1 = get_clock();
643: r = (t1-t0)/1000;
644: MKReal(r,*rp);
645: }
646:
1.16 ! noro 647: void Pfind_root_gf2n(NODE arg,GF2N *rp)
1.1 noro 648: {
649:
650: #if 0
651: UP p;
652:
653: ptoup((P)ARG0(arg),&p);
654: find_root_gf2n(p,rp);
655: #else
656: UP2 p;
657:
658: ptoup2((P)ARG0(arg),&p);
659: find_root_up2(p,rp);
660: #endif
661: }
662:
1.16 ! noro 663: void Pis_irred_gf2(NODE arg,Q *rp)
1.1 noro 664: {
665: UP2 t;
666:
667: ptoup2(ARG0(arg),&t);
668: *rp = irredcheckup2(t) ? ONE : 0;
669: }
670:
1.16 ! noro 671: void Pis_irred_ddd_gf2(NODE arg,Q *rp)
1.1 noro 672: {
673: UP2 t;
674: int ret;
675:
676: ptoup2(ARG0(arg),&t);
677: ret = irredcheck_dddup2(t);
678: STOQ(ret,*rp);
679: }
680:
1.16 ! noro 681: void Psetmod_ff(NODE arg,Obj *rp)
1.1 noro 682: {
683: int ac;
1.7 noro 684: int d;
1.1 noro 685: Obj mod,defpoly;
686: N n;
687: UP up;
688: UP2 up2;
1.14 noro 689: UM dp;
1.6 noro 690: Q q,r;
1.15 noro 691: P p,p1,y;
1.1 noro 692: NODE n0,n1;
693: LIST list;
694:
695: ac = argc(arg);
696: if ( ac == 1 ) {
697: mod = (Obj)ARG0(arg);
698: if ( !mod )
699: error("setmod_ff : invalid argument");
700: switch ( OID(mod) ) {
701: case O_N:
1.7 noro 702: current_ff = FF_GFP;
703: setmod_lm(NM((Q)mod));
1.6 noro 704: break;
1.1 noro 705: case O_P:
706: current_ff = FF_GF2N;
707: setmod_gf2n((P)mod); break;
708: default:
709: error("setmod_ff : invalid argument");
710: }
711: } else if ( ac == 2 ) {
1.7 noro 712: if ( OID(ARG0(arg)) == O_N ) {
713: /* small finite field; primitive root representation */
714: current_ff = FF_GFS;
715: setmod_sf(QTOS((Q)ARG0(arg)),QTOS((Q)ARG1(arg)));
1.6 noro 716: } else {
1.7 noro 717: mod = (Obj)ARG1(arg);
1.6 noro 718: current_ff = FF_GFPN;
719: defpoly = (Obj)ARG0(arg);
720: if ( !mod || !defpoly )
721: error("setmod_ff : invalid argument");
722: setmod_lm(NM((Q)mod));
723: setmod_gfpn((P)defpoly);
724: }
1.14 noro 725: } else if ( ac == 3 ) {
726: /* finite extension of a small finite field */
727: current_ff = FF_GFS;
728: setmod_sf(QTOS((Q)ARG0(arg)),QTOS((Q)ARG1(arg)));
729: d = QTOS((Q)ARG2(arg));
730: generate_defpoly_sfum(d,&dp);
1.15 noro 731: setmod_gfsn(dp);
732: current_ff = FF_GFSN;
1.1 noro 733: }
734: switch ( current_ff ) {
735: case FF_GFP:
736: getmod_lm(&n); NTOQ(n,1,q); *rp = (Obj)q; break;
737: case FF_GF2N:
738: getmod_gf2n(&up2); up2top(up2,&p); *rp = (Obj)p; break;
739: case FF_GFPN:
740: getmod_lm(&n); NTOQ(n,1,q);
741: getmod_gfpn(&up); uptop(up,&p);
742: MKNODE(n1,q,0); MKNODE(n0,p,n1);
743: MKLIST(list,n0);
744: *rp = (Obj)list; break;
1.6 noro 745: case FF_GFS:
1.15 noro 746: case FF_GFSN:
1.7 noro 747: STOQ(current_gfs_p,q);
1.15 noro 748: if ( current_gfs_ext )
1.7 noro 749: enc_to_p(current_gfs_p,current_gfs_iton[1],
750: VR(current_gfs_ext),&p);
1.15 noro 751: else {
1.11 noro 752: if ( current_gfs_p == 2 )
753: r = ONE;
754: else
755: STOQ(current_gfs_iton[1],r);
1.15 noro 756: p = (P)r;
757: }
758: switch ( current_ff ) {
759: case FF_GFS:
760: n0 = mknode(3,q,current_gfs_ext,p);
761: break;
762: case FF_GFSN:
763: getmod_gfsn(&dp);
764: makevar("y",&y);
765: sfumtop(VR(y),dp,&p1);
766: n0 = mknode(4,q,current_gfs_ext,p,p1);
767: break;
1.7 noro 768: }
1.6 noro 769: MKLIST(list,n0);
770: *rp = (Obj)list; break;
1.1 noro 771: default:
772: *rp = 0; break;
773: }
774: }
775:
1.16 ! noro 776: void Pextdeg_ff(Q *rp)
1.1 noro 777: {
778: int d;
779: UP2 up2;
780: UP up;
1.15 noro 781: UM dp;
1.1 noro 782:
783: switch ( current_ff ) {
784: case FF_GFP:
785: *rp = ONE; break;
786: case FF_GF2N:
787: getmod_gf2n(&up2); d = degup2(up2); STOQ(d,*rp); break;
788: case FF_GFPN:
789: getmod_gfpn(&up); STOQ(up->d,*rp); break;
1.10 noro 790: case FF_GFS:
791: if ( !current_gfs_ext )
792: *rp = ONE;
793: else
794: *rp = DEG(DC(current_gfs_ext));
795: break;
1.15 noro 796: case FF_GFSN:
797: getmod_gfsn(&dp);
798: STOQ(DEG(dp),*rp);
799: break;
1.1 noro 800: default:
801: error("extdeg_ff : current_ff is not set");
802: }
803: }
804:
1.16 ! noro 805: void Pcharacteristic_ff(Q *rp)
1.1 noro 806: {
807: N lm;
808:
809: switch ( current_ff ) {
810: case FF_GFP:
811: case FF_GFPN:
812: getmod_lm(&lm); NTOQ(lm,1,*rp); break;
813: case FF_GF2N:
814: STOQ(2,*rp); break;
1.10 noro 815: case FF_GFS:
1.15 noro 816: case FF_GFSN:
1.10 noro 817: STOQ(current_gfs_p,*rp); break;
1.1 noro 818: default:
819: error("characteristic_ff : current_ff is not set");
820: }
821: }
822:
1.16 ! noro 823: void Pfield_type_ff(Q *rp)
1.1 noro 824: {
825: STOQ(current_ff,*rp);
826: }
827:
1.16 ! noro 828: void Pfield_order_ff(Q *rp)
1.1 noro 829: {
830: N n;
831:
832: field_order_ff(&n);
833: NTOQ(n,1,*rp);
834: }
835:
1.16 ! noro 836: void Prandom_ff(Obj *rp)
1.1 noro 837: {
838: LM l;
839: GF2N g;
840: GFPN p;
1.7 noro 841: GFS s;
1.15 noro 842: GFSN spn;
1.1 noro 843:
844: switch ( current_ff ) {
845: case FF_GFP:
846: random_lm(&l); *rp = (Obj)l; break;
847: case FF_GF2N:
848: randomgf2n(&g); *rp = (Obj)g; break;
849: case FF_GFPN:
850: randomgfpn(&p); *rp = (Obj)p; break;
1.7 noro 851: case FF_GFS:
852: randomgfs(&s); *rp = (Obj)s; break;
1.15 noro 853: case FF_GFSN:
854: randomgfsn(&spn); *rp = (Obj)spn; break;
1.1 noro 855: default:
856: error("random_ff : current_ff is not set");
857: }
858: }
859:
1.16 ! noro 860: void Psimp_ff(NODE arg,Obj *rp)
! 861: {
1.1 noro 862: simp_ff((Obj)ARG0(arg),rp);
863: }
864:
1.16 ! noro 865: void getcoef(VL vl,P p,V v,Q d,P *r)
1.1 noro 866: {
867: P s,t,u,a,b,x;
868: DCP dc;
869: V w;
870:
871: if ( !p )
872: *r = 0;
873: else if ( NUM(p) )
874: *r = d ? 0 : p;
875: else if ( (w=VR(p)) == v ) {
876: for ( dc = DC(p); dc && cmpq(DEG(dc),d); dc = NEXT(dc) );
877: *r = dc ? COEF(dc) : 0;
878: } else {
879: MKV(w,x);
880: for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) {
881: getcoef(vl,COEF(dc),v,d,&t);
882: if ( t ) {
883: pwrp(vl,x,DEG(dc),&u); mulp(vl,t,u,&a);
884: addp(vl,s,a,&b); s = b;
885: }
886: }
887: *r = s;
888: }
889: }
890:
1.16 ! noro 891: void Pdeglist(NODE arg,LIST *rp)
1.1 noro 892: {
893: NODE d;
894:
895: getdeglist((P)ARG0(arg),VR((P)ARG1(arg)),&d);
896: MKLIST(*rp,d);
897: }
898:
1.16 ! noro 899: void Pch_mv(NODE arg,P *rp)
1.1 noro 900: {
901: change_mvar(CO,(P)ARG0(arg),VR((P)ARG1(arg)),rp);
902: }
903:
1.16 ! noro 904: void Pre_mv(NODE arg,P *rp)
1.1 noro 905: {
906: restore_mvar(CO,(P)ARG0(arg),VR((P)ARG1(arg)),rp);
907: }
908:
1.16 ! noro 909: void change_mvar(VL vl,P p,V v,P *r)
1.1 noro 910: {
911: Q d;
912: DCP dc,dc0;
913: NODE dl;
914:
915: if ( !p || NUM(p) || (VR(p) == v) )
916: *r = p;
917: else {
918: getdeglist(p,v,&dl);
919: for ( dc0 = 0; dl; dl = NEXT(dl) ) {
920: NEXTDC(dc0,dc); DEG(dc) = d = (Q)BDY(dl);
921: getcoef(vl,p,v,d,&COEF(dc));
922: }
923: NEXT(dc) = 0; MKP(v,dc0,*r);
924: }
925: }
926:
1.16 ! noro 927: void restore_mvar(VL vl,P p,V v,P *r)
1.1 noro 928: {
929: P s,u,a,b,x;
930: DCP dc;
931:
932: if ( !p || NUM(p) || (VR(p) != v) )
933: *r = p;
934: else {
935: MKV(v,x);
936: for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) {
937: pwrp(vl,x,DEG(dc),&u); mulp(vl,COEF(dc),u,&a);
938: addp(vl,s,a,&b); s = b;
939: }
940: *r = s;
941: }
942: }
943:
1.16 ! noro 944: void getdeglist(P p,V v,NODE *d)
1.1 noro 945: {
946: NODE n,n0,d0,d1,d2;
947: DCP dc;
948:
949: if ( !p || NUM(p) ) {
950: MKNODE(n,0,0); *d = n;
951: } else if ( VR(p) == v ) {
952: for ( n0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
953: NEXTNODE(n0,n); BDY(n) = (pointer)DEG(dc);
954: }
955: NEXT(n) = 0; *d = n0;
956: } else {
957: for ( dc = DC(p), d0 = 0; dc; dc = NEXT(dc) ) {
958: getdeglist(COEF(dc),v,&d1); mergedeglist(d0,d1,&d2); d0 = d2;
959: }
960: *d = d0;
961: }
962: }
1.16 ! noro 963:
! 964: void Pmergelist(NODE arg,LIST *rp)
1.1 noro 965: {
966: NODE n;
967:
968: asir_assert(ARG0(arg),O_LIST,"mergelist");
969: asir_assert(ARG1(arg),O_LIST,"mergelist");
970: mergedeglist(BDY((LIST)ARG0(arg)),BDY((LIST)ARG1(arg)),&n);
971: MKLIST(*rp,n);
972: }
973:
1.16 ! noro 974: void mergedeglist(NODE d0,NODE d1,NODE *dr)
1.1 noro 975: {
976: NODE t0,t,dt;
977: Q d;
978: int c;
979:
980: if ( !d0 )
981: *dr = d1;
982: else {
983: while ( d1 ) {
984: dt = d1; d1 = NEXT(d1); d = (Q)BDY(dt);
985: for ( t0 = 0, t = d0; t; t0 = t, t = NEXT(t) ) {
986: c = cmpq(d,(Q)BDY(t));
987: if ( !c )
988: break;
989: else if ( c > 0 ) {
990: if ( !t0 ) {
991: NEXT(dt) = d0; d0 = dt;
992: } else {
993: NEXT(t0) = dt; NEXT(dt) = t;
994: }
995: break;
996: }
997: }
998: if ( !t ) {
999: NEXT(t0) = dt; *dr = d0; return;
1000: }
1001: }
1002: *dr = d0;
1003: }
1004: }
1005:
1.16 ! noro 1006: void Pptomp(NODE arg,P *rp)
1.1 noro 1007: {
1008: ptomp(QTOS((Q)ARG1(arg)),(P)ARG0(arg),rp);
1009: }
1010:
1.16 ! noro 1011: void Pmptop(NODE arg,P *rp)
1.1 noro 1012: {
1013: mptop((P)ARG0(arg),rp);
1014: }
1015:
1.16 ! noro 1016: void Pptolmp(NODE arg,P *rp)
1.1 noro 1017: {
1018: ptolmp((P)ARG0(arg),rp);
1019: }
1020:
1.16 ! noro 1021: void Plmptop(NODE arg,P *rp)
1.1 noro 1022: {
1023: lmptop((P)ARG0(arg),rp);
1.11 noro 1024: }
1025:
1.16 ! noro 1026: void Psf_galois_action(NODE arg,P *rp)
1.11 noro 1027: {
1028: sf_galois_action(ARG0(arg),ARG1(arg),rp);
1.12 noro 1029: }
1030:
1031: /*
1032: sf_embed(F,B,PM)
1033: F : an element of GF(pn)
1034: B : the image of the primitive root of GF(pn)
1035: PM : order of GF(pm)
1036: */
1037:
1.16 ! noro 1038: void Psf_embed(NODE arg,P *rp)
1.12 noro 1039: {
1040: int k,pm;
1041:
1042: /* GF(pn)={0,1,a,a^2,...}->GF(pm)={0,1,b,b^2,...}; a->b^k */
1043: k = CONT((GFS)ARG1(arg));
1044: pm = QTOS((Q)ARG2(arg));
1045: sf_embed((P)ARG0(arg),k,pm,rp);
1.13 noro 1046: }
1047:
1.16 ! noro 1048: void Psf_log(NODE arg,Q *rp)
1.13 noro 1049: {
1050: int k;
1051:
1052: if ( !ARG0(arg) )
1053: error("sf_log : invalid armument");
1054: k = CONT((GFS)ARG0(arg));
1055: STOQ(k,*rp);
1.12 noro 1056: }
1057:
1.16 ! noro 1058: void Psf_find_root(NODE arg,GFS *rp)
1.12 noro 1059: {
1060: P p;
1061: Obj t;
1062: int d;
1063: UM u;
1064: int *root;
1065:
1066: p = (P)ARG0(arg);
1067: simp_ff((Obj)p,&t); p = (P)t;
1068: d = getdeg(VR(p),p);
1069: u = W_UMALLOC(d);
1070: ptosfum(p,u);
1071: root = (int *)ALLOCA(d*sizeof(int));
1072: find_rootsf(u,root);
1073: MKGFS(IFTOF(root[0]),*rp);
1074: }
1075:
1.16 ! noro 1076: void Psf_minipoly(NODE arg,P *rp)
1.12 noro 1077: {
1078: Obj t;
1079: P p1,p2;
1080: int d1,d2;
1081: UM up1,up2,m;
1082:
1083: p1 = (P)ARG0(arg); simp_ff((Obj)p1,&t); p1 = (P)t;
1084: p2 = (P)ARG1(arg); simp_ff((Obj)p2,&t); p2 = (P)t;
1085: d1 = getdeg(VR(p1),p1); up1 = W_UMALLOC(d1); ptosfum(p1,up1);
1086: d2 = getdeg(VR(p2),p2); up2 = W_UMALLOC(d2); ptosfum(p2,up2);
1087: m = W_UMALLOC(d2);
1088: minipolysf(up1,up2,m);
1089: sfumtop(VR(p2),m,&p1);
1090: sfptop(p1,rp);
1.11 noro 1091: }
1092:
1.16 ! noro 1093: void Pptosfp(NODE arg,P *rp)
1.11 noro 1094: {
1095: ptosfp(ARG0(arg),rp);
1.6 noro 1096: }
1097:
1.16 ! noro 1098: void Psfptop(NODE arg,P *rp)
1.6 noro 1099: {
1100: sfptop((P)ARG0(arg),rp);
1.1 noro 1101: }
1102:
1.16 ! noro 1103: void Pptogf2n(NODE arg,GF2N *rp)
1.1 noro 1104: {
1105: ptogf2n((Obj)ARG0(arg),rp);
1106: }
1107:
1.16 ! noro 1108: void Pgf2ntop(NODE arg,P *rp)
1.1 noro 1109: {
1110: if ( argc(arg) == 2 )
1111: up2_var = VR((P)ARG1(arg));
1112: gf2ntop((GF2N)ARG0(arg),rp);
1113: }
1114:
1.16 ! noro 1115: void Pgf2ntovect(NODE arg,VECT *rp)
1.1 noro 1116: {
1117: gf2ntovect((GF2N)ARG0(arg),rp);
1118: }
1119:
1.16 ! noro 1120: void Pptogfpn(NODE arg,GFPN *rp)
1.1 noro 1121: {
1122: ptogfpn((Obj)ARG0(arg),rp);
1123: }
1124:
1.16 ! noro 1125: void Pgfpntop(NODE arg,P *rp)
1.1 noro 1126: {
1127: if ( argc(arg) == 2 )
1128: up_var = VR((P)ARG1(arg));
1129: gfpntop((GFPN)ARG0(arg),rp);
1130: }
1131:
1.16 ! noro 1132: void Pureverse(NODE arg,P *rp)
1.1 noro 1133: {
1134: UP p,r;
1135:
1136: ptoup((P)ARG0(arg),&p);
1.3 noro 1137: if ( argc(arg) == 1 )
1138: reverseup(p,p->d,&r);
1139: else
1140: reverseup(p,QTOS((Q)ARG1(arg)),&r);
1.1 noro 1141: uptop(r,rp);
1142: }
1143:
1.16 ! noro 1144: void Putrunc(NODE arg,P *rp)
1.1 noro 1145: {
1146: UP p,r;
1147:
1148: ptoup((P)ARG0(arg),&p);
1149: truncup(p,QTOS((Q)ARG1(arg))+1,&r);
1150: uptop(r,rp);
1151: }
1152:
1.16 ! noro 1153: void Pudecomp(NODE arg,LIST *rp)
1.1 noro 1154: {
1155: P u,l;
1156: UP p,up,low;
1157: NODE n0,n1;
1158:
1159: ptoup((P)ARG0(arg),&p);
1160: decompup(p,QTOS((Q)ARG1(arg))+1,&low,&up);
1161: uptop(low,&l);
1162: uptop(up,&u);
1163: MKNODE(n1,u,0); MKNODE(n0,l,n1);
1164: MKLIST(*rp,n0);
1165: }
1166:
1.16 ! noro 1167: void Purembymul(NODE arg,P *rp)
1.1 noro 1168: {
1169: UP p1,p2,r;
1170:
1171: if ( !ARG0(arg) || !ARG1(arg) )
1172: *rp = 0;
1173: else {
1174: ptoup((P)ARG0(arg),&p1);
1175: ptoup((P)ARG1(arg),&p2);
1176: rembymulup(p1,p2,&r);
1177: uptop(r,rp);
1178: }
1179: }
1180:
1181: /*
1182: * d1 = deg(p1), d2 = deg(p2)
1183: * d1 <= 2*d2-1
1184: * p2*inv = 1 mod x^d2
1185: */
1186:
1.16 ! noro 1187: void Purembymul_precomp(NODE arg,P *rp)
1.1 noro 1188: {
1189: UP p1,p2,inv,r;
1190:
1191: if ( !ARG0(arg) || !ARG1(arg) )
1192: *rp = 0;
1193: else {
1194: ptoup((P)ARG0(arg),&p1);
1195: ptoup((P)ARG1(arg),&p2);
1196: ptoup((P)ARG2(arg),&inv);
1197: if ( p1->d >= 2*p2->d ) {
1198: error("urembymul_precomp : degree of 1st arg is too large");
1199: /* fprintf(stderr,"urembymul_precomp : degree of 1st arg is too large"); */
1200: remup(p1,p2,&r);
1201: } else
1202: hybrid_rembymulup_special(current_ff,p1,p2,inv,&r);
1203: uptop(r,rp);
1204: }
1205: }
1206:
1.16 ! noro 1207: void Puinvmod(NODE arg,P *rp)
1.1 noro 1208: {
1209: UP p,r;
1210:
1211: ptoup((P)ARG0(arg),&p);
1212: invmodup(p,QTOS((Q)ARG1(arg)),&r);
1213: uptop(r,rp);
1214: }
1215:
1.16 ! noro 1216: void Purevinvmod(NODE arg,P *rp)
1.1 noro 1217: {
1218: UP p,pr,r;
1219:
1220: ptoup((P)ARG0(arg),&p);
1221: reverseup(p,p->d,&pr);
1222: invmodup(pr,QTOS((Q)ARG1(arg)),&r);
1223: uptop(r,rp);
1224: }
1225:
1.16 ! noro 1226: void Ppwrmod_ff(NODE arg,P *rp)
1.1 noro 1227: {
1228: UP p1,p2;
1229:
1230: ptoup((P)ARG0(arg),&p1);
1231: switch ( current_ff ) {
1232: case FF_GFP:
1233: hybrid_powermodup(p1,&p2); break;
1234: case FF_GF2N:
1235: powermodup_gf2n(p1,&p2); break;
1236: case FF_GFPN:
1.10 noro 1237: case FF_GFS:
1.15 noro 1238: case FF_GFSN:
1.1 noro 1239: powermodup(p1,&p2); break;
1240: default:
1241: error("pwrmod_ff : current_ff is not set");
1242: }
1243: uptop(p2,rp);
1244: }
1245:
1.16 ! noro 1246: void Pgeneric_pwrmod_ff(NODE arg,P *rp)
1.1 noro 1247: {
1248: UP g,f,r;
1249:
1250: ptoup((P)ARG0(arg),&g);
1251: ptoup((P)ARG1(arg),&f);
1252: switch ( current_ff ) {
1253: case FF_GFP:
1254: hybrid_generic_powermodup(g,f,(Q)ARG2(arg),&r); break;
1255: case FF_GF2N:
1256: generic_powermodup_gf2n(g,f,(Q)ARG2(arg),&r); break;
1257: case FF_GFPN:
1.10 noro 1258: case FF_GFS:
1.15 noro 1259: case FF_GFSN:
1.1 noro 1260: generic_powermodup(g,f,(Q)ARG2(arg),&r); break;
1261: default:
1262: error("generic_pwrmod_ff : current_ff is not set");
1263: }
1264: uptop(r,rp);
1265: }
1266:
1.16 ! noro 1267: void Ppwrtab_ff(NODE arg,VECT *rp)
1.1 noro 1268: {
1269: UP f,xp;
1270: UP *tab;
1271: VECT r;
1272: int i,d;
1273:
1274: ptoup((P)ARG0(arg),&f);
1275: ptoup((P)ARG1(arg),&xp);
1276: d = f->d;
1277:
1278: tab = (UP *)ALLOCA(d*sizeof(UP));
1279: switch ( current_ff ) {
1280: case FF_GFP:
1281: hybrid_powertabup(f,xp,tab); break;
1282: case FF_GF2N:
1283: powertabup_gf2n(f,xp,tab); break;
1284: case FF_GFPN:
1.10 noro 1285: case FF_GFS:
1.15 noro 1286: case FF_GFSN:
1.1 noro 1287: powertabup(f,xp,tab); break;
1288: default:
1289: error("pwrtab_ff : current_ff is not set");
1290: }
1291: MKVECT(r,d); *rp = r;
1292: for ( i = 0; i < d; i++ )
1293: uptop(tab[i],(P *)&BDY(r)[i]);
1294: }
1295:
1.16 ! noro 1296: void Pkpwrmod_lm(NODE arg,P *rp)
1.1 noro 1297: {
1298: UP p1,p2;
1299:
1300: ptoup((P)ARG0(arg),&p1);
1301: powermodup(p1,&p2);
1302: uptop(p2,rp);
1303: }
1304:
1.16 ! noro 1305: void Pkgeneric_pwrmod_lm(NODE arg,P *rp)
1.1 noro 1306: {
1307: UP g,f,r;
1308:
1309: ptoup((P)ARG0(arg),&g);
1310: ptoup((P)ARG1(arg),&f);
1311: generic_powermodup(g,f,(Q)ARG2(arg),&r);
1312: uptop(r,rp);
1313: }
1314:
1.16 ! noro 1315: void Pkpwrtab_lm(NODE arg,VECT *rp)
1.1 noro 1316: {
1317: UP f,xp;
1318: UP *tab;
1319: VECT r;
1320: int i,d;
1321:
1322: ptoup((P)ARG0(arg),&f);
1323: ptoup((P)ARG1(arg),&xp);
1324: d = f->d;
1325:
1326: tab = (UP *)ALLOCA(d*sizeof(UP));
1327: powertabup(f,xp,tab);
1328: MKVECT(r,d); *rp = r;
1329: for ( i = 0; i < d; i++ )
1330: uptop(tab[i],(P *)&BDY(r)[i]);
1331: }
1332:
1.16 ! noro 1333: void Plazy_lm(NODE arg,Q *rp)
1.1 noro 1334: {
1335: lm_lazy = QTOS((Q)ARG0(arg));
1336: *rp = 0;
1337: }
1338:
1.16 ! noro 1339: void Pkmul(NODE arg,P *rp)
1.1 noro 1340: {
1341: P n1,n2;
1342:
1343: n1 = (P)ARG0(arg); n2 = (P)ARG1(arg);
1344: asir_assert(n1,O_P,"kmul");
1345: asir_assert(n2,O_P,"kmul");
1346: kmulp(CO,n1,n2,rp);
1347: }
1348:
1.16 ! noro 1349: void Pksquare(NODE arg,P *rp)
1.1 noro 1350: {
1351: P n1;
1352:
1353: n1 = (P)ARG0(arg);
1354: asir_assert(n1,O_P,"ksquare");
1355: ksquarep(CO,n1,rp);
1356: }
1357:
1.16 ! noro 1358: void Pktmul(NODE arg,P *rp)
1.1 noro 1359: {
1360: UP p1,p2,r;
1361:
1362: ptoup((P)ARG0(arg),&p1);
1363: ptoup((P)ARG1(arg),&p2);
1364: tkmulup(p1,p2,QTOS((Q)ARG2(arg))+1,&r);
1365: uptop(r,rp);
1366: }
1367:
1.16 ! noro 1368: void Pumul(NODE arg,P *rp)
1.1 noro 1369: {
1370: P a1,a2;
1371: UP p1,p2,r;
1372:
1373: a1 = (P)ARG0(arg); a2 = (P)ARG1(arg);
1374: if ( !a1 || !a2 || NUM(a1) || NUM(a2) )
1375: mulp(CO,a1,a2,rp);
1376: else {
1.16 ! noro 1377: if ( !uzpcheck((Obj)a1) || !uzpcheck((Obj)a2) || VR(a1) != VR(a2) )
1.1 noro 1378: error("umul : invalid argument");
1379: ptoup(a1,&p1);
1380: ptoup(a2,&p2);
1381: hybrid_mulup(0,p1,p2,&r);
1382: uptop(r,rp);
1383: }
1384: }
1385:
1.16 ! noro 1386: void Pusquare(NODE arg,P *rp)
1.1 noro 1387: {
1.16 ! noro 1388: UP p1,r;
1.1 noro 1389:
1390: ptoup((P)ARG0(arg),&p1);
1391: hybrid_squareup(0,p1,&r);
1392: uptop(r,rp);
1393: }
1394:
1.16 ! noro 1395: void Putmul(NODE arg,P *rp)
1.1 noro 1396: {
1397: UP p1,p2,r;
1398:
1399: ptoup((P)ARG0(arg),&p1);
1400: ptoup((P)ARG1(arg),&p2);
1401: hybrid_tmulup(0,p1,p2,QTOS((Q)ARG2(arg))+1,&r);
1402: uptop(r,rp);
1403: }
1404:
1.16 ! noro 1405: void Pumul_ff(NODE arg,Obj *rp)
1.1 noro 1406: {
1407: P a1,a2;
1408: UP p1,p2,r;
1409: P p;
1410:
1411: a1 = (P)ARG0(arg); a2 = (P)ARG1(arg);
1412: ptoup(a1,&p1);
1413: ptoup(a2,&p2);
1414: hybrid_mulup(current_ff,p1,p2,&r);
1415: uptop(r,&p);
1416: simp_ff((Obj)p,rp);
1417: }
1418:
1.16 ! noro 1419: void Pusquare_ff(NODE arg,Obj *rp)
1.1 noro 1420: {
1.16 ! noro 1421: UP p1,r;
1.1 noro 1422: P p;
1423:
1424: ptoup((P)ARG0(arg),&p1);
1425: hybrid_squareup(current_ff,p1,&r);
1426: uptop(r,&p);
1427: simp_ff((Obj)p,rp);
1428: }
1429:
1.16 ! noro 1430: void Putmul_ff(NODE arg,Obj *rp)
1.1 noro 1431: {
1432: UP p1,p2,r;
1433: P p;
1434:
1435: ptoup((P)ARG0(arg),&p1);
1436: ptoup((P)ARG1(arg),&p2);
1437: hybrid_tmulup(current_ff,p1,p2,QTOS((Q)ARG2(arg))+1,&r);
1438: uptop(r,&p);
1439: simp_ff((Obj)p,rp);
1440: }
1441:
1.16 ! noro 1442: void Phfmul_lm(NODE arg,P *rp)
1.1 noro 1443: {
1444: UP p1,p2,r;
1445: UP hi,lo,mid,t,s,p10,p11,p20,p21;
1446: unsigned int m,d;
1447: int i;
1448:
1449: ptoup((P)ARG0(arg),&p1);
1450: ptoup((P)ARG1(arg),&p2);
1451: d = MAX(p1->d,p2->d);
1452: for ( m = 1; m < d; m <<= 1 );
1453: if ( m > d )
1454: m >>= 1;
1455: if ( d-m < 10000 ) {
1456: decompup(p1,m,&p10,&p11); /* p1 = p11*x^m+p10 */
1457: decompup(p2,m,&p20,&p21); /* p2 = p21*x^m+p20 */
1458: fft_mulup_lm(p10,p20,&lo);
1459: kmulup(p11,p21,&hi);
1460: kmulup(p11,p20,&t); kmulup(p10,p21,&s); addup(t,s,&mid);
1461: r = UPALLOC(2*d);
1462: bzero((char *)COEF(r),(2*d+1)*sizeof(Num));
1463: if ( lo )
1464: bcopy((char *)COEF(lo),(char *)COEF(r),(DEG(lo)+1)*sizeof(Num));
1465: if ( hi )
1466: bcopy((char *)COEF(hi),(char *)(COEF(r)+2*m),(DEG(hi)+1)*sizeof(Num));
1467: for ( i = 2*d; i >= 0 && !COEF(r)[i]; i-- );
1468: if ( i < 0 )
1469: r = 0;
1470: else {
1471: DEG(r) = i;
1472: t = UPALLOC(DEG(mid)+m); DEG(t) = DEG(mid)+m;
1473: bzero((char *)COEF(t),(DEG(mid)+m+1)*sizeof(Num));
1474: bcopy((char *)COEF(mid),(char *)(COEF(t)+m),(DEG(mid)+1)*sizeof(Num));
1475: addup(t,r,&s);
1476: r = s;
1477: }
1478: } else
1479: fft_mulup_lm(p1,p2,&r);
1480: uptop(r,rp);
1481: }
1482:
1.16 ! noro 1483: void Pfmultest(NODE arg,LIST *rp)
1.1 noro 1484: {
1485: P p1,p2,r;
1486: int d1,d2;
1487: UM w1,w2,wr;
1488: unsigned int *f1,*f2,*fr,*w;
1489: int index,mod,root,d,maxint,i;
1490: int cond;
1491: Q prime;
1492: NODE n0,n1;
1493:
1494: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg); index = QTOS((Q)ARG2(arg));
1495: FFT_primes(index,&mod,&root,&d);
1496: maxint = 1<<d;
1497: d1 = UDEG(p1); d2 = UDEG(p2);
1498: if ( maxint < d1+d2+1 )
1499: *rp = 0;
1500: else {
1501: w1 = W_UMALLOC(d1); w2 = W_UMALLOC(d2);
1502: wr = W_UMALLOC(d1+d2);
1503: ptoum(mod,p1,w1); ptoum(mod,p2,w2);
1504: f1 = (unsigned int *)ALLOCA(maxint*sizeof(unsigned int));
1505: f2 = (unsigned int *)ALLOCA(maxint*sizeof(unsigned int));
1506: fr = (unsigned int *)ALLOCA(maxint*sizeof(unsigned int));
1507: w = (unsigned int *)ALLOCA(12*maxint*sizeof(unsigned int));
1508:
1509: for ( i = 0; i <= d1; i++ )
1510: f1[i] = (unsigned int)w1->c[i];
1511: for ( i = 0; i <= d2; i++ )
1512: f2[i] = (unsigned int)w2->c[i];
1513:
1514: cond = FFT_pol_product(d1,f1,d2,f2,fr,index,w);
1515: if ( cond )
1516: error("fmultest : ???");
1517: wr->d = d1+d2;
1518: for ( i = 0; i <= d1+d2; i++ )
1519: wr->c[i] = (unsigned int)fr[i];
1520: umtop(VR(p1),wr,&r);
1521: STOQ(mod,prime);
1522: MKNODE(n1,prime,0);
1523: MKNODE(n0,r,n1);
1524: MKLIST(*rp,n0);
1525: }
1526: }
1527:
1.16 ! noro 1528: void Pkmulum(NODE arg,P *rp)
1.1 noro 1529: {
1530: P p1,p2;
1531: int d1,d2,mod;
1532: UM w1,w2,wr;
1533:
1534: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg); mod = QTOS((Q)ARG2(arg));
1535: d1 = UDEG(p1); d2 = UDEG(p2);
1536: w1 = W_UMALLOC(d1); w2 = W_UMALLOC(d2);
1537: wr = W_UMALLOC(d1+d2);
1538: ptoum(mod,p1,w1); ptoum(mod,p2,w2);
1539: kmulum(mod,w1,w2,wr);
1540: umtop(VR(p1),wr,rp);
1541: }
1542:
1.16 ! noro 1543: void Pksquareum(NODE arg,P *rp)
1.1 noro 1544: {
1545: P p1;
1546: int d1,mod;
1547: UM w1,wr;
1548:
1549: p1 = (P)ARG0(arg); mod = QTOS((Q)ARG1(arg));
1550: d1 = UDEG(p1);
1551: w1 = W_UMALLOC(d1);
1552: wr = W_UMALLOC(2*d1);
1553: ptoum(mod,p1,w1);
1554: kmulum(mod,w1,w1,wr);
1555: umtop(VR(p1),wr,rp);
1556: }
1557:
1.16 ! noro 1558: void Ptracemod_gf2n(NODE arg,P *rp)
1.1 noro 1559: {
1560: UP g,f,r;
1561:
1562: ptoup((P)ARG0(arg),&g);
1563: ptoup((P)ARG1(arg),&f);
1564: tracemodup_gf2n(g,f,(Q)ARG2(arg),&r);
1565: uptop(r,rp);
1.2 noro 1566: }
1567:
1.16 ! noro 1568: void Pumul_specialmod(NODE arg,P *rp)
1.2 noro 1569: {
1570: P a1,a2;
1571: UP p1,p2,r;
1572: int i,nmod;
1573: int *modind;
1574: NODE t,n;
1575:
1576: a1 = (P)ARG0(arg); a2 = (P)ARG1(arg);
1577: if ( !a1 || !a2 || NUM(a1) || NUM(a2) )
1578: mulp(CO,a1,a2,rp);
1579: else {
1.16 ! noro 1580: if ( !uzpcheck((Obj)a1) || !uzpcheck((Obj)a2) || VR(a1) != VR(a2) )
1.2 noro 1581: error("umul_specialmod : invalid argument");
1582: ptoup(a1,&p1);
1583: ptoup(a2,&p2);
1584: n = BDY((LIST)ARG2(arg));
1585: nmod = length(n);
1586: modind = (int *)MALLOC_ATOMIC(nmod*sizeof(int));
1587: for ( i = 0, t = n; i < nmod; i++, t = NEXT(t) )
1588: modind[i] = QTOS((Q)BDY(t));
1589: fft_mulup_specialmod_main(p1,p2,0,modind,nmod,&r);
1590: uptop(r,rp);
1591: }
1592: }
1593:
1.16 ! noro 1594: void Pusquare_specialmod(NODE arg,P *rp)
1.2 noro 1595: {
1596: P a1;
1597: UP p1,r;
1598: int i,nmod;
1599: int *modind;
1600: NODE t,n;
1601:
1602: a1 = (P)ARG0(arg);
1603: if ( !a1 || NUM(a1) )
1604: mulp(CO,a1,a1,rp);
1605: else {
1.16 ! noro 1606: if ( !uzpcheck((Obj)a1) )
1.2 noro 1607: error("usquare_specialmod : invalid argument");
1608: ptoup(a1,&p1);
1609: n = BDY((LIST)ARG1(arg));
1610: nmod = length(n);
1611: modind = (int *)MALLOC_ATOMIC(nmod*sizeof(int));
1612: for ( i = 0, t = n; i < nmod; i++, t = NEXT(t) )
1613: modind[i] = QTOS((Q)BDY(t));
1614: fft_mulup_specialmod_main(p1,p1,0,modind,nmod,&r);
1615: uptop(r,rp);
1616: }
1617: }
1618:
1.16 ! noro 1619: void Putmul_specialmod(NODE arg,P *rp)
1.2 noro 1620: {
1621: P a1,a2;
1622: UP p1,p2,r;
1623: int i,nmod;
1624: int *modind;
1625: NODE t,n;
1626:
1627: a1 = (P)ARG0(arg); a2 = (P)ARG1(arg);
1628: if ( !a1 || !a2 || NUM(a1) || NUM(a2) )
1629: mulp(CO,a1,a2,rp);
1630: else {
1.16 ! noro 1631: if ( !uzpcheck((Obj)a1) || !uzpcheck((Obj)a2) || VR(a1) != VR(a2) )
1.2 noro 1632: error("utmul_specialmod : invalid argument");
1633: ptoup(a1,&p1);
1634: ptoup(a2,&p2);
1635: n = BDY((LIST)ARG3(arg));
1636: nmod = length(n);
1637: modind = (int *)MALLOC_ATOMIC(nmod*sizeof(int));
1638: for ( i = 0, t = n; i < nmod; i++, t = NEXT(t) )
1639: modind[i] = QTOS((Q)BDY(t));
1640: fft_mulup_specialmod_main(p1,p2,QTOS((Q)ARG2(arg))+1,modind,nmod,&r);
1641: uptop(r,rp);
1642: }
1.1 noro 1643: }