Annotation of OpenXM_contrib2/asir2000/engine/up_gf2n.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/engine/up_gf2n.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
! 2: #include "ca.h"
! 3: #include <math.h>
! 4:
! 5: extern int debug_up;
! 6: extern int up_lazy;
! 7: extern GEN_UP2 current_mod_gf2n;
! 8:
! 9: void squarep_gf2n(vl,n1,nr)
! 10: VL vl;
! 11: P n1;
! 12: P *nr;
! 13: {
! 14: UP b1,br;
! 15:
! 16: if ( !n1 )
! 17: *nr = 0;
! 18: else if ( OID(n1) == O_N )
! 19: mulp(vl,n1,n1,nr);
! 20: else {
! 21: ptoup(n1,&b1);
! 22: squareup_gf2n(b1,&br);
! 23: uptop(br,nr);
! 24: }
! 25: }
! 26:
! 27: void squareup_gf2n(n1,nr)
! 28: UP n1;
! 29: UP *nr;
! 30: {
! 31: UP r;
! 32: GF2N *c1,*c;
! 33: int i,d1,d;
! 34:
! 35: if ( !n1 )
! 36: *nr = 0;
! 37: else if ( !n1->d ) {
! 38: *nr = r = UPALLOC(0); r->d = 0;
! 39: squaregf2n((GF2N)n1->c[0],(GF2N *)(&r->c[0]));
! 40: } else {
! 41: d1 = n1->d;
! 42: d = 2*d1;
! 43: *nr = r = UPALLOC(d); r->d = d;
! 44: c1 = (GF2N *)n1->c; c = (GF2N *)r->c;
! 45: bzero((char *)c,(d+1)*sizeof(GF2N *));
! 46: for ( i = 0; i <= d1; i++ )
! 47: squaregf2n(c1[i],&c[2*i]);
! 48: }
! 49: }
! 50:
! 51: /* x^(2^n) mod f */
! 52:
! 53: void powermodup_gf2n(f,xp)
! 54: UP f;
! 55: UP *xp;
! 56: {
! 57: UP x,t,invf;
! 58: int k,n;
! 59: GF2N lm;
! 60: struct oEGT eg_sq,eg_rem,eg_mul,eg_inv,eg0,eg1,eg2;
! 61:
! 62: n = degup2(current_mod_gf2n->dense);
! 63: MKGF2N(ONEUP2,lm);
! 64: x = UPALLOC(1); x->d = 1; x->c[1] = (Num)lm;
! 65:
! 66: reverseup(f,f->d,&t);
! 67: invmodup(t,f->d,&invf);
! 68: for ( k = 0; k < n; k++ ) {
! 69: squareup_gf2n(x,&t);
! 70: rembymulup_special(t,f,invf,&x);
! 71: /* remup(t,f,&x); */
! 72: }
! 73: *xp = x;
! 74: }
! 75:
! 76: /* g^d mod f */
! 77:
! 78: void generic_powermodup_gf2n(g,f,d,xp)
! 79: UP g,f;
! 80: Q d;
! 81: UP *xp;
! 82: {
! 83: N e;
! 84: UP x,y,t,invf,s;
! 85: int k;
! 86: GF2N lm;
! 87:
! 88: e = NM(d);
! 89: MKGF2N(ONEUP2,lm);
! 90: y = UPALLOC(0); y->d = 0; y->c[0] = (Num)lm;
! 91: remup(g,f,&x);
! 92: if ( !x ) {
! 93: *xp = !d ? y : 0;
! 94: return;
! 95: } else if ( !x->d ) {
! 96: pwrup(x,d,xp);
! 97: return;
! 98: }
! 99: reverseup(f,f->d,&t);
! 100: invmodup(t,f->d,&invf);
! 101: for ( k = n_bits(e)-1; k >= 0; k-- ) {
! 102: squareup_gf2n(y,&t);
! 103: rembymulup_special(t,f,invf,&s);
! 104: y = s;
! 105: if ( e->b[k/32] & (1<<(k%32)) ) {
! 106: mulup(y,x,&t);
! 107: remup(t,f,&s);
! 108: y = s;
! 109: }
! 110: }
! 111: *xp = y;
! 112: }
! 113:
! 114: /* g+g^2+...+g^(2^(nd-1)) mod f; where e = deg(mod) */
! 115:
! 116: void tracemodup_gf2n(g,f,d,xp)
! 117: UP g,f;
! 118: Q d;
! 119: UP *xp;
! 120: {
! 121: UP x,t,s,u,invf;
! 122: int en,i;
! 123:
! 124: en = QTOS(d)*degup2(current_mod_gf2n->dense);
! 125: remup(g,f,&x);
! 126: if ( !x ) {
! 127: *xp = 0;
! 128: return;
! 129: }
! 130: reverseup(f,f->d,&t);
! 131: invmodup(t,f->d,&invf);
! 132: for ( i = 1, t = s = x; i < en; i++ ) {
! 133: squareup_gf2n(t,&u);
! 134: rembymulup_special(u,f,invf,&t);
! 135: addup(s,t,&u); s = u;
! 136: }
! 137: *xp = s;
! 138: }
! 139:
! 140: void tracemodup_gf2n_slow(g,f,d,xp)
! 141: UP g,f;
! 142: Q d;
! 143: UP *xp;
! 144: {
! 145: UP x,t,s,u;
! 146: int en,i;
! 147:
! 148: en = QTOS(d)*degup2(current_mod_gf2n->dense);
! 149: remup(g,f,&x);
! 150: if ( !x ) {
! 151: *xp = 0;
! 152: return;
! 153: }
! 154: for ( i = 1, t = s = x; i < en; i++ ) {
! 155: squareup_gf2n(t,&u);
! 156: remup(u,f,&t);
! 157: addup(s,t,&u); s = u;
! 158: }
! 159: *xp = s;
! 160: }
! 161:
! 162: static struct oEGT eg_trace_tab,eg_trace_mul;
! 163:
! 164: void tracemodup_gf2n_tab(g,f,d,xp)
! 165: UP g,f;
! 166: Q d;
! 167: UP *xp;
! 168: {
! 169: UP x0,x2,t,s,u;
! 170: int en,i;
! 171: UP *tab;
! 172: GF2N one;
! 173: struct oEGT eg1,eg2;
! 174:
! 175: en = QTOS(d)*degup2(current_mod_gf2n->dense);
! 176: remup(g,f,&t); g = t;
! 177: if ( !g ) {
! 178: *xp = 0;
! 179: return;
! 180: }
! 181:
! 182: MKGF2N(ONEUP2,one);
! 183: x0 = UPALLOC(0); x0->d = 0; x0->c[0] = (Num)one;
! 184: x2 = UPALLOC(2); x2->d = 2; x2->c[2] = (Num)one;
! 185:
! 186: tab = (UP *)ALLOCA(en*sizeof(UP));
! 187: tab[0] = x0;
! 188: remup(x2,f,&tab[1]);
! 189:
! 190: for ( i = 2; i < en; i++ ) {
! 191: mulup(tab[i-1],tab[1],&t); remup(t,f,&tab[i]);
! 192: }
! 193:
! 194: for ( i = 1, t = s = g; i < en; i++ ) {
! 195: square_rem_tab_up_gf2n(t,tab,&u); t = u;
! 196: addup(s,t,&u); s = u;
! 197: }
! 198: *xp = s;
! 199: }
! 200:
! 201: void square_rem_tab_up_gf2n(f,tab,rp)
! 202: UP f;
! 203: UP *tab;
! 204: UP *rp;
! 205: {
! 206: UP s,t,u,n;
! 207: Num *c;
! 208: int i,d;
! 209:
! 210: n = UPALLOC(0); n->d = 0;
! 211: if ( !f )
! 212: *rp = 0;
! 213: else {
! 214: d = f->d; c = f->c;
! 215: up_lazy = 1;
! 216: for ( i = 0, s = 0; i <= d; i++ ) {
! 217: squaregf2n((GF2N)c[i],(GF2N *)(&n->c[0]));
! 218: mulup(tab[i],n,&t); addup(s,t,&u); s = u;
! 219: }
! 220: up_lazy = 0;
! 221: simpup(s,rp);
! 222: }
! 223: }
! 224:
! 225: void powertabup_gf2n(f,xp,tab)
! 226: UP f;
! 227: UP xp;
! 228: UP *tab;
! 229: {
! 230: UP y,t,invf;
! 231: int i,d;
! 232: GF2N lm;
! 233:
! 234: d = f->d;
! 235: MKGF2N(ONEUP2,lm);
! 236: y = UPALLOC(0); y->d = 0; y->c[0] = (Num)lm;
! 237: tab[0] = y;
! 238: tab[1] = xp;
! 239:
! 240: reverseup(f,f->d,&t);
! 241: invmodup(t,f->d,&invf);
! 242:
! 243: for ( i = 2; i < d; i++ ) {
! 244: if ( debug_up )
! 245: fprintf(stderr,".");
! 246: if ( !(i%2) )
! 247: squareup_gf2n(tab[i/2],&t);
! 248: else
! 249: kmulup(tab[i-1],xp,&t);
! 250: rembymulup_special(t,f,invf,&tab[i]);
! 251: /* remup(t,f,&tab[i]); */
! 252: }
! 253: }
! 254:
! 255: void find_root_gf2n(f,r)
! 256: UP f;
! 257: GF2N *r;
! 258: {
! 259: UP g,ut,c,t,h,rem;
! 260: int n;
! 261: GF2N rn;
! 262: struct oEGT eg0,eg1,eg_trace;
! 263:
! 264: n = degup2(current_mod_gf2n->dense);
! 265: g = f;
! 266: while ( g->d > 1 ) {
! 267: ut = UPALLOC(1); ut->c[0] = 0;
! 268: randomgf2n(&rn);
! 269: if ( !rn )
! 270: continue;
! 271: ut->c[1] = (Num)rn; ut->d = 1;
! 272: tracemodup_gf2n_tab(ut,f,ONE,&c);
! 273: gcdup(c,g,&h);
! 274: if ( h->d && h->d < g->d ) {
! 275: if ( 2*h->d > g->d ) {
! 276: qrup(g,h,&t,&rem); g = t;
! 277: if ( rem )
! 278: error("find_root_gf2n : cannot happen");
! 279: } else
! 280: g = h;
! 281: }
! 282: monicup(g,&t); g = t;
! 283: printf("deg(g)=%d\n",g->d);
! 284: }
! 285: divgf2n((GF2N)g->c[0],(GF2N)g->c[1],r);
! 286: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>