Annotation of OpenXM_contrib2/asir2018/builtin/ec.c, Revision 1.1
1.1 ! 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
! 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
! 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: *
! 48: * $OpenXM$
! 49: */
! 50: #include "ca.h"
! 51: #include "parse.h"
! 52: #include "inline.h"
! 53:
! 54: /*
! 55: * Elliptic curve related functions
! 56: *
! 57: * Argument specifications
! 58: * Point = vector(x,y,z) (projective representation)
! 59: * ECInfo = vector(a,b,p)
! 60: * a,b : integer for GF(p), GF2N for GF(2^n)
! 61: * p : integer for GF(p), polynomial for GF(2^n)
! 62: *
! 63: */
! 64:
! 65: struct oKeyIndexPair {
! 66: unsigned int key;
! 67: int index;
! 68: };
! 69:
! 70: void Psha1();
! 71: void Psha1_ec();
! 72: void Pecm_add_ff();
! 73: void Pecm_chsgn_ff();
! 74: void Pecm_sub_ff();
! 75: void Pecm_compute_all_key_homo_ff();
! 76: void Pnextvect1(),Psort_ktarray(),Pecm_find_match(),Pseparate_vect();
! 77: void Pecm_set_addcounter();
! 78: void Pecm_count_order();
! 79:
! 80: void ecm_add_ff(VECT,VECT,VECT,VECT *);
! 81: void ecm_add_gfp(VECT,VECT,VECT,VECT *);
! 82: void ecm_add_gf2n(VECT,VECT,VECT,VECT *);
! 83:
! 84: void ecm_chsgn_ff(VECT,VECT *);
! 85: void ecm_sub_ff(VECT,VECT,VECT,VECT *);
! 86:
! 87: void compute_all_key_homo_gfp(VECT *,int,unsigned int *);
! 88: void compute_all_key_homo_gf2n(VECT *,int,unsigned int *);
! 89:
! 90: unsigned int separate_vect(double *,int);
! 91: void ecm_find_match(unsigned int *,int,unsigned int *,int,LIST *);
! 92: int find_match(unsigned int,unsigned int *,int);
! 93:
! 94: void sort_ktarray(VECT,VECT,LIST *);
! 95: int comp_kip(struct oKeyIndexPair *,struct oKeyIndexPair *);
! 96:
! 97: int nextvect1(VECT,VECT);
! 98:
! 99: unsigned int ecm_count_order_gfp(unsigned int,unsigned int,unsigned int);
! 100: unsigned int ecm_count_order_gf2n(UP2,GF2N,GF2N);
! 101:
! 102: void sha_memory(unsigned char *,unsigned int,unsigned int *);
! 103:
! 104: unsigned int ecm_addcounter;
! 105: extern int current_ff,lm_lazy;
! 106:
! 107: struct ftab ec_tab[] = {
! 108: /* point addition */
! 109: {"ecm_add_ff",Pecm_add_ff,3},
! 110:
! 111: /* point change sign */
! 112: {"ecm_chsgn_ff",Pecm_chsgn_ff,1},
! 113:
! 114: /* point subtraction */
! 115: {"ecm_sub_ff",Pecm_sub_ff,3},
! 116:
! 117: /* key computation for sort and match */
! 118: {"ecm_compute_all_key_homo_ff",Pecm_compute_all_key_homo_ff,1},
! 119:
! 120: /* exhausitve search of rational points */
! 121: {"ecm_count_order",Pecm_count_order,1},
! 122:
! 123: /* common functions */
! 124: {"nextvect1",Pnextvect1,2},
! 125: {"sort_ktarray",Psort_ktarray,2},
! 126: {"separate_vect",Pseparate_vect,1},
! 127: {"ecm_find_match",Pecm_find_match,2},
! 128: {"ecm_set_addcounter",Pecm_set_addcounter,-1},
! 129: {"sha1",Psha1,-2},
! 130: #if 0
! 131: {"sha1_free",Psha1_free,1},
! 132: #endif
! 133: {0,0,0},
! 134: };
! 135:
! 136: void Psha1(NODE arg,Z *rp)
! 137: {
! 138: unsigned char *s;
! 139: unsigned int digest[5];
! 140: size_t bl;
! 141: Z z;
! 142: mpz_t r;
! 143:
! 144: asir_assert(ARG0(arg),O_N,"sha1_free");
! 145: z = (Z)ARG0(arg);
! 146: s = (unsigned char *)mpz_export(0,&bl,1,sizeof(unsigned char),1,0,BDY(z));
! 147: sha_memory(s,bl,digest);
! 148: mpz_init(r);
! 149: mpz_import(r,5,1,sizeof(int),0,0,digest);
! 150: MPZTOZ(r,*rp);
! 151: }
! 152:
! 153: #if 0
! 154: void Psha1_ec(arg,rp)
! 155: NODE arg;
! 156: Q *rp;
! 157: {
! 158: #include <fj_crypt.h>
! 159: SHS_CTX context;
! 160: unsigned char *s;
! 161: int i,j,l,bl,n;
! 162: unsigned int t;
! 163: N z,r;
! 164: extern int little_endian;
! 165:
! 166: asir_assert(ARG0(arg),O_N,"sha1");
! 167: z = NM((Q)ARG0(arg));
! 168: n = PL(z);
! 169: l = n_bits(z);
! 170: bl = (l+7)/8;
! 171: s = (unsigned char *)MALLOC(bl);
! 172: for ( i = 0, j = bl-1; i < n; i++ ) {
! 173: t = BD(z)[i];
! 174: if ( j >= 0 ) s[j--] = t&0xff; t>>=8;
! 175: if ( j >= 0 ) s[j--] = t&0xff; t>>=8;
! 176: if ( j >= 0 ) s[j--] = t&0xff; t>>=8;
! 177: if ( j >= 0 ) s[j--] = t;
! 178: }
! 179: SHSInit(&context);
! 180: SHSUpdate(&context,s,bl);
! 181: SHSFinal(&context);
! 182: r = NALLOC(5);
! 183: if ( little_endian )
! 184: for ( i = 0; i < 5; i++ ) {
! 185: t = context.digest[4-i];
! 186: BD(r)[i] = (t>>24)|((t&0xff0000)>>8)|((t&0xff00)<<8)|(t<<24);
! 187: } else
! 188: for ( i = 0; i < 5; i++ )
! 189: BD(r)[i] = context.digest[4-i];
! 190: for ( i = 4; i >= 0 && !BD(r)[i]; i-- );
! 191: if ( i < 0 )
! 192: *rp = 0;
! 193: else {
! 194: PL(r) = i+1;
! 195: NTOQ(r,1,*rp);
! 196: }
! 197: }
! 198: #endif
! 199:
! 200: void Pecm_count_order(NODE arg,Z *rp)
! 201: {
! 202: Z p;
! 203: UP2 d;
! 204: Obj a,b;
! 205: unsigned int p0,a0,b0,ord;
! 206: VECT ec;
! 207: Obj *vb;
! 208:
! 209: switch ( current_ff ) {
! 210: case FF_GFP:
! 211: getmod_lm(&p);
! 212: if ( z_bits((Q)p) > 32 )
! 213: error("ecm_count_order : ground field too large");
! 214: p0 = QTOS(p);
! 215: ec = (VECT)ARG0(arg);
! 216: vb = (Obj *)BDY(ec); simp_ff(vb[0],&a); simp_ff(vb[1],&b);
! 217: a0 = LMTOS((LM)a);
! 218: b0 = LMTOS((LM)b);
! 219: ord = ecm_count_order_gfp(p0,a0,b0);
! 220: UTOQ(ord,*rp);
! 221: break;
! 222: case FF_GF2N:
! 223: getmod_gf2n(&d);
! 224: if ( degup2(d) > 10 )
! 225: error("ecm_count_order : ground field too large");
! 226: ec = (VECT)ARG0(arg);
! 227: vb = (Obj *)BDY(ec); simp_ff(vb[0],&a); simp_ff(vb[1],&b);
! 228: ord = ecm_count_order_gf2n(d,(GF2N)a,(GF2N)b);
! 229: UTOQ(ord,*rp);
! 230: break;
! 231: default:
! 232: error("ecm_count_order : current_ff is not set");
! 233: }
! 234: }
! 235:
! 236: void Pecm_set_addcounter(NODE arg,Z *rp)
! 237: {
! 238: if ( arg )
! 239: ecm_addcounter = QTOS((Q)ARG0(arg));
! 240: UTOQ(ecm_addcounter,*rp);
! 241: }
! 242:
! 243: void Pecm_compute_all_key_homo_ff(NODE arg,VECT *rp)
! 244: {
! 245: unsigned int *ka;
! 246: int len,i;
! 247: VECT *pa;
! 248: VECT r,v;
! 249: LIST *vb;
! 250: USINT *b;
! 251:
! 252: v = (VECT)ARG0(arg);
! 253: len = v->len;
! 254: vb = (LIST *)v->body;
! 255: pa = (VECT *)ALLOCA(len*sizeof(VECT));
! 256: ka = (unsigned int *)ALLOCA(len*sizeof(unsigned int));
! 257: for ( i = 0; i < len; i++ )
! 258: pa[i] = (VECT)BDY(NEXT(BDY(vb[i])));
! 259: switch ( current_ff ) {
! 260: case FF_GFP:
! 261: compute_all_key_homo_gfp(pa,len,ka); break;
! 262: case FF_GF2N:
! 263: compute_all_key_homo_gf2n(pa,len,ka); break;
! 264: default:
! 265: error("ecm_compute_all_key_homo_ff : current_ff is not set");
! 266: }
! 267: MKVECT(r,len); *rp = r;
! 268: b = (USINT *)r->body;
! 269: for ( i = 0; i < len; i++ )
! 270: MKUSINT(b[i],ka[i]);
! 271: }
! 272:
! 273: void Psort_ktarray(NODE arg,LIST *rp)
! 274: {
! 275: sort_ktarray((VECT)ARG0(arg),(VECT)ARG1(arg),rp);
! 276: }
! 277:
! 278: void Pecm_add_ff(NODE arg,VECT *rp)
! 279: {
! 280: ecm_add_ff(ARG0(arg),ARG1(arg),ARG2(arg),rp);
! 281: }
! 282:
! 283: void Pecm_sub_ff(NODE arg,VECT *rp)
! 284: {
! 285: ecm_sub_ff(ARG0(arg),ARG1(arg),ARG2(arg),rp);
! 286: }
! 287:
! 288: void Pecm_chsgn_ff(NODE arg,VECT *rp)
! 289: {
! 290: ecm_chsgn_ff(ARG0(arg),rp);
! 291: }
! 292:
! 293: void Pnextvect1(NODE arg,Z *rp)
! 294: {
! 295: int index;
! 296:
! 297: index = nextvect1(ARG0(arg),ARG1(arg));
! 298: STOQ(index,*rp);
! 299: }
! 300:
! 301: /* XXX at least n < 32 must hold. What is the strict restriction for n ? */
! 302:
! 303: void Pseparate_vect(NODE arg,LIST *rp)
! 304: {
! 305: VECT v;
! 306: int n,i;
! 307: Z *b;
! 308: double *w;
! 309: unsigned int s;
! 310: NODE ns,nc,t,t1;
! 311: Z iq;
! 312: LIST ls,lc;
! 313:
! 314: v = (VECT)ARG0(arg);
! 315: n = v->len; b = (Z *)v->body;
! 316: w = (double *)ALLOCA(n*sizeof(double));
! 317: for ( i = 0; i < n; i++ )
! 318: w[i] = (double)QTOS(b[i]);
! 319: s = separate_vect(w,n);
! 320: ns = nc = 0;
! 321: for ( i = n-1; i >= 0; i-- )
! 322: if ( s & (1<<i) ) {
! 323: STOQ(i,iq); MKNODE(t,iq,ns); ns = t;
! 324: } else {
! 325: STOQ(i,iq); MKNODE(t,iq,nc); nc = t;
! 326: }
! 327: MKLIST(ls,ns); MKLIST(lc,nc);
! 328: MKNODE(t,lc,0); MKNODE(t1,ls,t);
! 329: MKLIST(*rp,t1);
! 330: }
! 331:
! 332: void Pecm_find_match(NODE arg,LIST *rp)
! 333: {
! 334: VECT g,b;
! 335: int ng,nb,i;
! 336: USINT *p;
! 337: unsigned int *kg,*kb;
! 338:
! 339: g = (VECT)ARG0(arg); ng = g->len;
! 340: kg = (unsigned int *)ALLOCA(ng*sizeof(unsigned int));
! 341: for ( i = 0, p = (USINT *)g->body; i < ng; i++ )
! 342: kg[i] = p[i]?BDY(p[i]):0;
! 343: b = (VECT)ARG1(arg); nb = b->len;
! 344: kb = (unsigned int *)ALLOCA(nb*sizeof(unsigned int));
! 345: for ( i = 0, p = (USINT *)b->body; i < nb; i++ )
! 346: kb[i] = p[i]?BDY(p[i]):0;
! 347: ecm_find_match(kg,ng,kb,nb,rp);
! 348: }
! 349:
! 350: void ecm_add_ff(VECT p1,VECT p2,VECT ec,VECT *pr)
! 351: {
! 352: if ( !p1 )
! 353: *pr = p2;
! 354: else if ( !p2 )
! 355: *pr = p1;
! 356: else {
! 357: switch ( current_ff ) {
! 358: case FF_GFP:
! 359: ecm_add_gfp(p1,p2,ec,pr); break;
! 360: case FF_GF2N:
! 361: ecm_add_gf2n(p1,p2,ec,pr); break;
! 362: default:
! 363: error("ecm_add_ff : current_ff is not set");
! 364: }
! 365: }
! 366: }
! 367:
! 368: /* ec = [AX,BC] */
! 369:
! 370: void ecm_add_gf2n(VECT p1,VECT p2,VECT ec,VECT *rp)
! 371: {
! 372: GF2N ax,bc,a0,a1,a2,b0,b1,b2;
! 373: GF2N a2b0,a0b2,a2b1,a1b2,a02,a04,a22,a24,a0a2,a0a22,a1a2;
! 374: GF2N t,s,u,r0,r1,r00,r01,r02,r002,r003,r02q;
! 375: VECT r;
! 376: GF2N *vb,*rb;
! 377:
! 378: ecm_addcounter++;
! 379: /* addition with O */
! 380: if ( !p1 ) {
! 381: *rp = p2;
! 382: return;
! 383: }
! 384: if ( !p2 ) {
! 385: *rp = p1;
! 386: return;
! 387: }
! 388: vb = (GF2N *)BDY(ec); ax = vb[0]; bc = vb[1];
! 389: vb = (GF2N *)BDY(p1); a0 = vb[0]; a1 = vb[1]; a2 = vb[2];
! 390: vb = (GF2N *)BDY(p2); b0 = vb[0]; b1 = vb[1]; b2 = vb[2];
! 391:
! 392: mulgf2n(a2,b0,&a2b0); mulgf2n(a0,b2,&a0b2);
! 393: if ( !cmpgf2n(a2b0,a0b2) ) {
! 394: mulgf2n(a2,b1,&a2b1);
! 395: mulgf2n(a1,b2,&a1b2);
! 396: if ( !cmpgf2n(a2b1,a1b2) ) {
! 397: if ( !a0 )
! 398: *rp = 0;
! 399: else {
! 400: squaregf2n(a0,&a02); squaregf2n(a02,&a04);
! 401: squaregf2n(a2,&a22); squaregf2n(a22,&a24);
! 402: mulgf2n(a0,a2,&a0a2); squaregf2n(a0a2,&a0a22);
! 403: mulgf2n(bc,a24,&t); addgf2n(a04,t,&r0);
! 404: mulgf2n(a04,a0a2,&t); mulgf2n(a1,a2,&a1a2);
! 405: addgf2n(a02,a1a2,&s); addgf2n(s,a0a2,&u);
! 406: mulgf2n(u,r0,&s); addgf2n(t,s,&r1);
! 407:
! 408: MKVECT(r,3); rb = (GF2N *)r->body;
! 409: mulgf2n(r0,a0a2,&rb[0]); rb[1] = r1; mulgf2n(a0a22,a0a2,&rb[2]);
! 410: *rp = r;
! 411: }
! 412: } else
! 413: *rp = 0;
! 414: } else {
! 415: mulgf2n(a1,b2,&a1b2); addgf2n(a0b2,a2b0,&r00);
! 416: mulgf2n(a2,b1,&t); addgf2n(a1b2,t,&r01); mulgf2n(a2,b2,&r02);
! 417: squaregf2n(r00,&r002); mulgf2n(r002,r00,&r003);
! 418:
! 419: addgf2n(r00,r01,&t); mulgf2n(t,r01,&s); mulgf2n(s,r02,&t);
! 420: if ( ax ) {
! 421: mulgf2n(r02,ax,&r02q);
! 422: addgf2n(t,r003,&s); mulgf2n(r02q,r002,&t); addgf2n(s,t,&r0);
! 423: } else
! 424: addgf2n(t,r003,&r0);
! 425:
! 426: mulgf2n(a0b2,r002,&t); addgf2n(t,r0,&s); mulgf2n(r01,s,&t);
! 427: mulgf2n(r002,a1b2,&s); addgf2n(r0,s,&u); mulgf2n(r00,u,&s);
! 428: addgf2n(t,s,&r1);
! 429:
! 430: MKVECT(r,3); rb = (GF2N *)r->body;
! 431: mulgf2n(r0,r00,&rb[0]); rb[1] = r1; mulgf2n(r003,r02,&rb[2]);
! 432: *rp = r;
! 433: }
! 434: }
! 435:
! 436: extern LM THREE_LM,FOUR_LM,EIGHT_LM;
! 437:
! 438: /* 0 < p < 2^16, 0 <= a,b < p */
! 439:
! 440: unsigned int ecm_count_order_gfp(unsigned int p,unsigned int a,unsigned int b)
! 441: {
! 442: unsigned int x,rhs,ord,t;
! 443:
! 444: for ( x = 0, ord = 1; x < p; x++ ) {
! 445: DMAR(x,x,a,p,t) /* t = x^2+a mod p */
! 446: DMAR(t,x,b,p,rhs) /* rhs = x*(x^2+a)+b mod p */
! 447: if ( !rhs )
! 448: ord++;
! 449: else if ( small_jacobi(rhs,p)==1 )
! 450: ord+=2;
! 451: }
! 452: return ord;
! 453: }
! 454:
! 455: unsigned int ecm_count_order_gf2n(UP2 d,GF2N a,GF2N b)
! 456: {
! 457: error("ecm_count_order_gf2n : not implemented yet");
! 458: /* NOTREACHED */
! 459: return 0;
! 460: }
! 461:
! 462: /* ec = [AX,BC] */
! 463:
! 464: void ecm_add_gfp(VECT p1,VECT p2,VECT ec,VECT *pr)
! 465: {
! 466: LM aa,bb,x1,y1,z1,x2,y2,z2,x1z2,v1,y1z2,u1,u2,v2,v3,z1z2;
! 467: LM v2x1z2,a1,x3,y3,z3,w1,s1,s2,s3,s1y1,b1,h1;
! 468: LM t,s,u;
! 469: LM *vb;
! 470: VECT r;
! 471:
! 472: ecm_addcounter++;
! 473: /* addition with O */
! 474: if( !p1 ) {
! 475: *pr = p2;
! 476: return;
! 477: }
! 478: if( !p2 ) {
! 479: *pr = p1;
! 480: return;
! 481: }
! 482:
! 483: /* set parameters */
! 484:
! 485: /* aa = ec[0]; bb = ec[1]; */
! 486: vb = (LM *)BDY(ec); aa = vb[0]; bb = vb[1];
! 487:
! 488: /* x1 = p1[0]; y1 = p1[1]; z1 = p1[2]; */
! 489: vb = (LM *)BDY(p1); x1 = vb[0]; y1 = vb[1]; z1 = vb[2];
! 490:
! 491: /* x2 = p2[0]; y2 = p2[1]; z2 = p2[2]; */
! 492: vb = (LM *)BDY(p2); x2 = vb[0]; y2 = vb[1]; z2 = vb[2];
! 493:
! 494: /* addition */
! 495:
! 496: /* x1z2 = (x1*z2) %p; */
! 497: mullm(x1,z2,&x1z2);
! 498:
! 499: /* v1 = (x2*z1-x1z2) %p; */
! 500: lm_lazy = 1;
! 501: mullm(x2,z1,&t); sublm(t,x1z2,&s);
! 502: lm_lazy = 0; simplm(s,&v1);
! 503:
! 504: /* y1z2 = (y1*z2) %p; */
! 505: mullm(y1,z2,&y1z2);
! 506:
! 507: /* u1 = (y2*z1-y1z2) %p; */
! 508: lm_lazy = 1;
! 509: mullm(y2,z1,&t); sublm(t,y1z2,&s);
! 510: lm_lazy = 0; simplm(s,&u1);
! 511:
! 512: if( v1 != 0 ) {
! 513: /* u2 = (u1*u1) %p; */
! 514: mullm(u1,u1,&u2);
! 515:
! 516: /* v2 = (v1*v1) %p; */
! 517: mullm(v1,v1,&v2);
! 518:
! 519: /* v3 = (v1*v2) %p; */
! 520: mullm(v1,v2,&v3);
! 521:
! 522: /* z1z2 = (z1*z2) %p; */
! 523: mullm(z1,z2,&z1z2);
! 524:
! 525: /* v2x1z2 = (v2*x1z2) %p; */
! 526: mullm(v2,x1z2,&v2x1z2);
! 527:
! 528: /* a1 = (u2*z1z2-v3-2*v2x1z2) %p; */
! 529: lm_lazy = 1;
! 530: mullm(u2,z1z2,&t); addlm(v2x1z2,v2x1z2,&s);
! 531: addlm(v3,s,&u); sublm(t,u,&s);
! 532: lm_lazy = 0; simplm(s,&a1);
! 533:
! 534: /* x3 = ( v1 * a1 ) %p; */
! 535: mullm(v1,a1,&x3);
! 536:
! 537: /* y3 = ( u1 * ( v2x1z2 - a1 ) - v3 * y1z2 ) %p; */
! 538: lm_lazy = 1;
! 539: sublm(v2x1z2,a1,&t); mullm(u1,t,&s); mullm(v3,y1z2,&u); sublm(s,u,&t);
! 540: lm_lazy = 0; simplm(t,&y3);
! 541:
! 542: /* z3 = ( v3 * z1z2 ) %p; */
! 543: mullm(v3,z1z2,&z3);
! 544: } else if( u1 == 0 ) {
! 545: /* w1 = (aa*z1*z1+3*x1*x1) %p; */
! 546: lm_lazy = 1;
! 547: mullm(z1,z1,&t); mullm(aa,t,&s);
! 548: mullm(x1,x1,&t); mullm(THREE_LM,t,&u); addlm(s,u,&t);
! 549: lm_lazy = 0; simplm(t,&w1);
! 550:
! 551: /* s1 = (y1*z1) %p; */
! 552: mullm(y1,z1,&s1);
! 553:
! 554: /* s2 = (s1*s1) %p; */
! 555: mullm(s1,s1,&s2);
! 556:
! 557: /* s3 = (s1*s2) %p; */
! 558: mullm(s1,s2,&s3);
! 559:
! 560: /* s1y1 = (s1*y1) %p; */
! 561: mullm(s1,y1,&s1y1);
! 562:
! 563: /* b1 = (s1y1*x1) %p; */
! 564: mullm(s1y1,x1,&b1);
! 565:
! 566: /* h1 = (w1*w1-8*b1) %p; */
! 567: lm_lazy = 1;
! 568: mullm(w1,w1,&t); mullm(EIGHT_LM,b1,&s); sublm(t,s,&u);
! 569: lm_lazy = 0; simplm(u,&h1);
! 570:
! 571: /* x3 = ( 2 * h1 * s1 ) %p; */
! 572: lm_lazy = 1;
! 573: mullm(h1,s1,&t); addlm(t,t,&s);
! 574: lm_lazy = 0; simplm(s,&x3);
! 575:
! 576: /* y3 = ( w1 * ( 4 * b1 - h1 ) - 8 * s1y1 * s1y1 ) %p; */
! 577: lm_lazy = 1;
! 578: mullm(FOUR_LM,b1,&t); sublm(t,h1,&s); mullm(w1,s,&u);
! 579: mullm(s1y1,s1y1,&t); mullm(EIGHT_LM,t,&s); sublm(u,s,&t);
! 580: lm_lazy = 0; simplm(t,&y3);
! 581:
! 582: /* z3 = ( 8 * s3 ) %p; */
! 583: mullm(EIGHT_LM,s3,&z3);
! 584: } else {
! 585: *pr = 0;
! 586: return;
! 587: }
! 588: if ( !z3 )
! 589: *pr = 0;
! 590: else {
! 591: MKVECT(r,3); *pr = r;
! 592: vb = (LM *)BDY(r); vb[0] = x3; vb[1] = y3; vb[2] = z3;
! 593: }
! 594: }
! 595:
! 596: void ecm_chsgn_ff(VECT p,VECT *pr)
! 597: {
! 598: Obj x,y,z;
! 599: LM tl;
! 600: GF2N tg;
! 601: Obj *vb;
! 602: VECT r;
! 603:
! 604: if( !p ) {
! 605: *pr = 0;
! 606: return;
! 607: }
! 608:
! 609: /* x = p[0]; y = p[1]; z = p[2]; */
! 610: vb = (Obj *)BDY(p); x = vb[0]; y = vb[1]; z = vb[2];
! 611: switch ( current_ff ) {
! 612: case FF_GFP:
! 613: if ( !y )
! 614: *pr = p;
! 615: else {
! 616: chsgnlm((LM)y,&tl);
! 617: MKVECT(r,3); *pr = r;
! 618: vb = (Obj *)BDY(r); vb[0] = x; vb[1] = (Obj)tl; vb[2] = z;
! 619: }
! 620: break;
! 621: case FF_GF2N:
! 622: addgf2n((GF2N)x,(GF2N)y,&tg);
! 623: MKVECT(r,3); *pr = r;
! 624: vb = (Obj *)BDY(r); vb[0] = x; vb[1] = (Obj)tg; vb[2] = z;
! 625: break;
! 626: default:
! 627: error("ecm_chsgn_ff : current_ff is not set");
! 628: }
! 629: }
! 630:
! 631: void ecm_sub_ff(VECT p1,VECT p2,VECT ec,VECT *pr)
! 632: {
! 633: VECT mp2;
! 634:
! 635: ecm_chsgn_ff(p2,&mp2);
! 636: ecm_add_ff(p1,mp2,ec,pr);
! 637: }
! 638:
! 639: /* tplist = [[t,p],...]; t:interger, p=[p0,p1]:point (vector) */
! 640:
! 641: int comp_kip(struct oKeyIndexPair *a,struct oKeyIndexPair *b)
! 642: {
! 643: unsigned int ka,kb;
! 644:
! 645: ka = a->key; kb = b->key;
! 646: if ( ka > kb )
! 647: return 1;
! 648: else if ( ka < kb )
! 649: return -1;
! 650: else
! 651: return 0;
! 652: }
! 653:
! 654: #define EC_GET_XZ(p,x,z) \
! 655: if ( !(p) ) {\
! 656: (x)=0; (z)=(LM)ONE;\
! 657: } else { \
! 658: LM *vb;\
! 659: vb = (LM *)BDY((VECT)(p));\
! 660: (x) = vb[0]; (z) = vb[2];\
! 661: }
! 662:
! 663: #define EC_GET_XZ_GF2N(p,x,z) \
! 664: if ( !(p) ) {\
! 665: (x)=0; (z)=(GF2N)ONE;\
! 666: } else { \
! 667: GF2N *vb;\
! 668: vb = (GF2N *)BDY((VECT)(p));\
! 669: (x) = vb[0]; (z) = vb[2];\
! 670: }
! 671:
! 672: void compute_all_key_homo_gfp(VECT *pa,int len,unsigned int *ka)
! 673: {
! 674: LM *b,*x,*z;
! 675: int i;
! 676: LM t,s,m;
! 677:
! 678: b = (LM *)ALLOCA((len+1)*sizeof(LM));
! 679: x = (LM *)ALLOCA(len*sizeof(LM));
! 680: z = (LM *)ALLOCA(len*sizeof(LM));
! 681: b[0] = ONELM;
! 682: for ( i = 1; i <= len; i++ ) {
! 683: EC_GET_XZ(pa[i-1],x[i-1],z[i-1]);
! 684: mullm(b[i-1],z[i-1],&b[i]);
! 685: }
! 686: /* b[0] = 1 */
! 687: divlm(b[0],b[len],&m);
! 688: for ( i = len-1; i >= 0; i-- ) {
! 689: mullm(m,b[i],&s); mullm(s,x[i],&t); s = t;
! 690: ka[i] = LMTOS(s); ka[i] |= 0x80000000;
! 691: mullm(m,z[i],&s); m = s;
! 692: }
! 693: }
! 694:
! 695: void compute_all_key_homo_gf2n(VECT *pa,int len,unsigned int *ka)
! 696: {
! 697: GF2N *b,*x,*z;
! 698: int i;
! 699: GF2N t,s,m;
! 700:
! 701: b = (GF2N *)ALLOCA((len+1)*sizeof(Q));
! 702: x = (GF2N *)ALLOCA(len*sizeof(Q));
! 703: z = (GF2N *)ALLOCA(len*sizeof(Q));
! 704: MKGF2N(ONEUP2,b[0]);
! 705: for ( i = 1; i <= len; i++ ) {
! 706: EC_GET_XZ_GF2N(pa[i-1],x[i-1],z[i-1]);
! 707: mulgf2n(b[i-1],z[i-1],&b[i]);
! 708: }
! 709: invgf2n(b[len],&m);
! 710: for ( i = len-1; i >= 0; i-- ) {
! 711: mulgf2n(m,b[i],&s); mulgf2n(s,x[i],&t); s = t;
! 712: ka[i] = s ? s->body->b[0] : 0; ka[i] |= 0x80000000;
! 713: mulgf2n(m,z[i],&s); m = s;
! 714: }
! 715: }
! 716:
! 717: unsigned int separate_vect(double *v,int n)
! 718: {
! 719: unsigned int max = 1<<n;
! 720: unsigned int i,j,i0;
! 721: double all,a,total,m;
! 722:
! 723: for ( i = 0, all = 1; i < (unsigned int)n; i++ )
! 724: all *= v[i];
! 725:
! 726: for ( i = 0, m = 0; i < max; i++ ) {
! 727: for ( a = 1, j = 0; j < (unsigned int)n; j++ )
! 728: if ( i & (1<<j) )
! 729: a *= v[j];
! 730: total = a+(all/a)*2;
! 731: if ( !m || total < m ) {
! 732: m = total;
! 733: i0 = i;
! 734: }
! 735: }
! 736: return i0;
! 737: }
! 738:
! 739: void ecm_find_match(unsigned int *g,int ng,unsigned int *b,int nb,LIST *r)
! 740: {
! 741: int i,j;
! 742: Z iq,jq;
! 743: NODE n0,n1,c0,c;
! 744: LIST l;
! 745:
! 746: for ( i = 0, c0 = 0; i < ng; i++ ) {
! 747: j = find_match(g[i],b,nb);
! 748: if ( j >= 0 ) {
! 749: STOQ(i,iq); STOQ(j,jq);
! 750: MKNODE(n1,jq,0); MKNODE(n0,iq,n1); MKLIST(l,n0);
! 751: NEXTNODE(c0,c);
! 752: BDY(c) = (pointer)l;
! 753: }
! 754: }
! 755: if ( c0 )
! 756: NEXT(c) = 0;
! 757: MKLIST(*r,c0);
! 758: }
! 759:
! 760: int find_match(unsigned int k,unsigned int *key,int n)
! 761: {
! 762: int s,e,m;
! 763:
! 764: for ( s = 0, e = n; (e-s) > 1; ) {
! 765: m = (s+e)/2;
! 766: if ( k==key[m] )
! 767: return m;
! 768: else if ( k > key[m] )
! 769: s = m;
! 770: else
! 771: e = m;
! 772: }
! 773: if ( k == key[s] )
! 774: return s;
! 775: else
! 776: return -1;
! 777: }
! 778:
! 779: int nextvect1(VECT vect,VECT bound)
! 780: {
! 781: int size,i,a;
! 782: Z *vb,*bb;
! 783:
! 784: size = vect->len;
! 785: vb = (Z *)vect->body;
! 786: bb = (Z *)bound->body;
! 787: for ( i = size-1; i >= 0; i-- )
! 788: if ( (a=QTOS(vb[i])) < QTOS(bb[i]) ) {
! 789: a++; STOQ(a,vb[i]);
! 790: break;
! 791: } else
! 792: vb[i] = 0;
! 793: return i;
! 794: }
! 795:
! 796: void sort_ktarray(VECT karray,VECT tarray,LIST *rp)
! 797: {
! 798: NODE r,r1;
! 799: int i,i0,k,len,same,tsame;
! 800: struct oKeyIndexPair *kip;
! 801: VECT key,value,v;
! 802: Q *tb,*samebuf;
! 803: USINT *kb;
! 804: Obj *svb;
! 805: USINT *skb;
! 806:
! 807: len = karray->len;
! 808: kb = (USINT *)karray->body;
! 809:
! 810: kip = (struct oKeyIndexPair *)ALLOCA(len*sizeof(struct oKeyIndexPair));
! 811: for ( i = 0; i < len; i++ ) {
! 812: kip[i].key = BDY(kb[i]); kip[i].index = i;
! 813: }
! 814: qsort((void *)kip,len,sizeof(struct oKeyIndexPair),
! 815: (int (*)(const void *,const void *))comp_kip);
! 816:
! 817: for ( same = tsame = i = i0 = 0, k = 1; i < len; i++, tsame++ )
! 818: if ( kip[i0].key != kip[i].key ) {
! 819: i0 = i; k++;
! 820: same = MAX(tsame,same);
! 821: tsame = 0;
! 822: }
! 823: same = MAX(tsame,same);
! 824: samebuf = (Q *)ALLOCA(same*sizeof(Q));
! 825:
! 826: MKVECT(key,k); skb = (USINT *)BDY(key);
! 827: MKVECT(value,k); svb = (Obj *)BDY(value);
! 828:
! 829: tb = (Q *)tarray->body;
! 830: for ( same = i = i0 = k = 0; i <= len; i++ ) {
! 831: if ( i == len || kip[i0].key != kip[i].key ) {
! 832: skb[k] = kb[kip[i0].index];
! 833: if ( same > 1 ) {
! 834: MKVECT(v,same);
! 835: bcopy((char *)samebuf,(char *)v->body,same*sizeof(Q));
! 836: svb[k] = (Obj)v;
! 837: } else
! 838: svb[k] = (Obj)samebuf[0];
! 839: i0 = i;
! 840: k++;
! 841: same = 0;
! 842: if ( i == len )
! 843: break;
! 844: }
! 845: samebuf[same++] = tb[kip[i].index];
! 846: }
! 847: MKNODE(r1,value,0); MKNODE(r,key,r1); MKLIST(*rp,r);
! 848: }
! 849:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to ec.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2018 / builtin