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