crypto/ec/ec_point.c

   1 /*
   2  *
   3  *      ec_point.c
   4  *
   5  *      Elliptic Curve Arithmetic Functions
   6  *
   7  *      Copyright (C) Lenka Fibikova 2000
   8  *
   9  *
  10  */
  11
  12 #include <stdio.h>
  13 #include <stdlib.h>
  14 #include <assert.h>
  15 #include <memory.h>
  16
  17 #include <openssl/bn.h>
  18
  19 #include "../bn/bn_mont2.h" /* XXX */
  20 #include "ec.h"
  21
  22
  23 EC_POINT *ECP_new()
  24         {
  25         EC_POINT *ret;
  26
  27         ret=(EC_POINT *)malloc(sizeof(EC_POINT));
  28         if (ret == NULL) return NULL;
  29         ret->X = BN_new();
  30         ret->Y = BN_new();
  31         ret->Z = BN_new();
  32         ret->is_in_mont = 0;
  33
  34         if (ret->X == NULL || ret->Y == NULL || ret->Z == NULL)
  35                 {
  36                 if (ret->X != NULL) BN_free(ret->X);
  37                 if (ret->Y != NULL) BN_free(ret->Y);
  38                 if (ret->Z != NULL) BN_free(ret->Z);
  39                 free(ret);
  40                 return(NULL);
  41                 }
  42         return(ret);
  43         }
  44
  45
  46 void ECP_clear_free(EC_POINT *P)
  47         {
  48         if (P == NULL) return;
  49
  50         P->is_in_mont = 0;
  51         if (P->X != NULL) BN_clear_free(P->X);
  52         if (P->Y != NULL) BN_clear_free(P->Y);
  53         if (P->Z != NULL) BN_clear_free(P->Z);
  54         free(P);
  55         }
  56
  57
  58 void ECP_clear_free_precompute(ECP_PRECOMPUTE *prec)
  59         {
  60         int i;
  61         int max;
  62
  63         if (prec == NULL) return;
  64         if (prec->Pi != NULL)
  65                 {
  66                 max = 1;
  67                 max <<= (prec->r - 1);
  68
  69                 for (i = 0; i < max; i++)
  70                         {
  71                         if (prec->Pi[i] != NULL) ECP_clear_free(prec->Pi[i]);
  72                         }
  73                 }
  74         free(prec);
  75         }
  76
  77
  78 int ECP_is_on_ec(EC_POINT *P, EC *E, BN_CTX *ctx)
  79         {
  80         BIGNUM *n0, *n1, *n2, *p;
  81         int Pnorm;
  82         int ret = -1;
  83
  84         assert(P != NULL);
  85         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
  86
  87         assert(E != NULL);
  88         assert(E->A != NULL && E->B != NULL && E->p != NULL);
  89
  90         assert(ctx != NULL);
  91
  92         assert(!P->is_in_mont);
  93
  94         if (ECP_is_infty(P)) return 1;
  95
  96         BN_CTX_start(ctx);
  97         n0 = BN_CTX_get(ctx);
  98         n1 = BN_CTX_get(ctx);
  99         n2 = BN_CTX_get(ctx);
 100         if (n2 == NULL)
 101                 goto err;
 102
 103         p = E->p;
 104
 105         Pnorm = (ECP_is_norm(P));
 106
 107         if (!Pnorm)
 108                 {
 109                 if (!BN_mod_mul(n0, P->Z, P->Z, p, ctx)) goto err;
 110                 if (!BN_mod_mul(n1, n0, n0, p, ctx)) goto err;
 111                 if (!BN_mod_mul(n2, n0, n1, p, ctx)) goto err;
 112                 }
 113
 114         if (!BN_mod_mul(n0, P->X, P->X, p, ctx)) goto err;
 115         if (!BN_mod_mul(n0, n0, P->X, p, ctx)) goto err;
 116
 117         if (Pnorm)
 118                 {
 119                 if (!BN_mod_mul(n1, P->X, E->A, p, ctx)) goto err;
 120                 }
 121         else
 122                 {
 123                 if (!BN_mod_mul(n1, n1, P->X, p, ctx)) goto err;
 124                 if (!BN_mod_mul(n1, n1, E->A, p, ctx)) goto err;
 125                 }
 126         if (!BN_mod_add(n0, n0, n1, p, ctx)) goto err;
 127
 128         if (Pnorm)
 129                 {
 130                 if (!BN_mod_add(n0, n0, E->B, p, ctx)) goto err;
 131                 }
 132         else
 133                 {
 134                 if (!BN_mod_mul(n2, n2, E->B,  p, ctx)) goto err;
 135                 if (!BN_mod_add(n0, n0, n2, p, ctx)) goto err;
 136                 }
 137
 138         if (!BN_mod_mul(n1, P->Y, P->Y, p, ctx)) goto err;
 139
 140         if (BN_cmp(n0, n1))
 141                 ret = 0;
 142         else
 143                 ret = 1;
 144
 145 err:
 146         BN_CTX_end(ctx);
 147         return ret;
 148         }
 149
 150
 151 EC_POINT *ECP_generate(BIGNUM *x, BIGNUM *z,EC *E, BN_CTX *ctx)
 152 /* x == NULL || z = 0  -> point of infinity     */
 153 /* z == NULL || z = 1  -> normalized            */
 154         {
 155         BIGNUM *n0, *n1;
 156         EC_POINT *ret = NULL;
 157         int Pnorm, Pinfty, X0, A0;
 158
 159         assert(E != NULL);
 160         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
 161
 162         assert(ctx != NULL);
 163
 164         Pinfty = (x == NULL);
 165         Pnorm = (z == NULL);
 166         if (!Pnorm)
 167                 {
 168                 Pnorm = BN_is_one(z);
 169                 Pinfty = (Pinfty || BN_is_zero(z));
 170                 }
 171
 172         if (Pinfty)
 173                 {
 174                 if ((ret = ECP_new()) == NULL) return NULL;
 175                 if (!BN_zero(ret->Z))
 176                         {
 177                         ECP_clear_free(ret);
 178                         return NULL;
 179                         }
 180                 return ret;
 181                 }
 182
 183         X0 = BN_is_zero(x);
 184         A0 = BN_is_zero(E->A);
 185
 186         if ((ret = ECP_new()) == NULL) return NULL;
 187
 188         ret->is_in_mont = 0;
 189
 190         BN_CTX_start(ctx);
 191         n0 = BN_CTX_get(ctx);
 192         n1 = BN_CTX_get(ctx);
 193         if (n1 == NULL) goto err;
 194
 195         if (!BN_zero(n0)) goto err;
 196         if (!BN_zero(n1)) goto err;
 197
 198         if (!X0)
 199                 {
 200                 if (!BN_mod_sqr(n0, x, E->p, ctx)) goto err;
 201                 if (!BN_mod_mul(n0, n0, x, E->p, ctx)) goto err;        /* x^3 */
 202                 }
 203
 204         if (!X0 && !A0)
 205                 {
 206                 if (!BN_mod_mul(n1, E->A, x, E->p, ctx)) goto err;      /* Ax */
 207                 if (!BN_mod_add(n0, n0, n1, E->p, ctx)) goto err;       /* x^3 + Ax */
 208                 }
 209
 210         if (!BN_is_zero(E->B))
 211                 if (!BN_mod_add(n0, n0, E->B, E->p, ctx)) goto err;     /* x^3 + Ax +B */
 212
 213         if (!BN_mod_sqrt(ret->Y, n0, E->p, ctx)) goto err;
 214         if (BN_copy(ret->X, x) == NULL) goto err;
 215
 216         if (Pnorm)
 217                 {
 218                 if (!BN_one(ret->Z)) goto err;
 219                 }
 220         else
 221                 {
 222                 if (BN_copy(ret->Z, z) == NULL) goto err;
 223                 if (!BN_mod_sqr(n0, z, E->p, ctx)) goto err;
 224                 if (!BN_mod_mul(ret->X, ret->X, n0, E->p, ctx)) goto err;
 225                 if (!BN_mod_mul(n0, n0, z, E->p, ctx)) goto err;
 226                 if (!BN_mod_mul(ret->Y, ret->Y, n0, E->p, ctx)) goto err;
 227                 }
 228
 229 #ifdef TEST
 230         if (!ECP_is_on_ec(ret, E, ctx)) goto err;
 231 #endif
 232
 233         BN_CTX_end(ctx);
 234         return ret;
 235
 236 err:
 237         if (ret != NULL) ECP_clear_free(ret);
 238         BN_CTX_end(ctx);
 239         return NULL;
 240         }
 241
 242
 243 int ECP_ecp2bin(EC_POINT *P, unsigned char *to, int form)
 244 /* form =       1 ... compressed
 245                 2 ... uncompressed
 246                 3 ... hybrid */
 247         {
 248         int bytes, bx, by;
 249
 250         assert (P != NULL);
 251         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
 252         assert (!P->is_in_mont);
 253         assert (ECP_is_norm(P) || ECP_is_infty(P));
 254         assert (to != NULL);
 255         assert (form > 0 && form < 4);
 256
 257         if (BN_is_zero(P->Z))
 258                 {
 259                 to[0] = 0;
 260                 return 1;
 261                 }
 262
 263         bx = BN_num_bytes(P->X);
 264         if (form == 1 ) bytes = bx + 1;
 265         else
 266                 {
 267                 by = BN_num_bytes(P->Y);
 268                 bytes = (bx > by ? bx : by);
 269                 bytes = bytes * 2 + 1;
 270                 }
 271         memset(to, 0, bytes);
 272
 273         switch (form)
 274                 {
 275         case 1: to[0] = 2;      break;
 276         case 2: to[0] = 4;      break;
 277         case 3: to[0] = 6;      break;
 278                 }
 279         if (form != 2) to[0] += BN_is_bit_set(P->Y, 0);
 280
 281
 282         if ((BN_bn2bin(P->X, to + 1)) != bx) return 0;
 283         if (form != 1)
 284                 {
 285                 if ((BN_bn2bin(P->Y, to + bx + 1)) != by) return 0;
 286                 }
 287
 288         return bytes;
 289         }
 290
 291
 292 int ECP_bin2ecp(unsigned char *from, int len, EC_POINT *P, EC *E, BN_CTX *ctx)
 293         {
 294         int y;
 295         BIGNUM *x;
 296         EC_POINT *pp;
 297
 298         assert (E != NULL);
 299         assert (E->A != NULL && E->B != NULL && E->p != NULL);
 300         assert (!E->is_in_mont);
 301
 302         assert (ctx != NULL);
 303         assert (from != NULL);
 304         assert (P != NULL);
 305         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
 306
 307         if (len == 1 && from[0] != 0) return 0;
 308
 309         if (len == 0 || len == 1)
 310                 {
 311                 if (!BN_zero(P->Z)) return 0;
 312                 return 1;
 313                 }
 314
 315         switch (from[0])
 316                 {
 317         case 2:
 318         case 3:
 319                 y = from[0] - 2;
 320                 if ((x = BN_new()) == NULL) return 0;
 321                 if (BN_bin2bn(from + 1, len - 1, x) == NULL) return 0;
 322
 323                 pp = ECP_generate(x, NULL, E, ctx);
 324                 BN_clear_free(x);
 325                 if (pp == NULL) return 0;
 326
 327                 ECP_copy(P, pp);
 328                 ECP_clear_free(pp);
 329
 330                 if (BN_is_bit_set(P->Y, 0) != y)
 331                         if (!BN_sub(P->Y, E->p, P->Y)) return 0;
 332                 break;
 333
 334         case 4:
 335         case 6:
 336         case 7:
 337                 y = (len - 1)/2;
 338                 if (BN_bin2bn(from + 1, y, P->X) == NULL) return 0;
 339                 if (BN_bin2bn(from + y + 1, y, P->Y) == NULL) return 0;
 340                 if (!BN_set_word(P->Z, 1)) return 0;
 341                 break;
 342
 343         default:
 344                 assert(0);
 345
 346                 }
 347
 348         if (!ECP_is_on_ec(P, E, ctx)) return 0;
 349         return 1;
 350         }
 351
 352
 353 int ECP_normalize(EC_POINT *P, EC *E, BN_CTX *ctx)
 354         {
 355         BIGNUM *z, *zm;
 356
 357         assert (P != NULL);
 358         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
 359
 360         assert (E != NULL);
 361         assert (E->A != NULL && E->B != NULL && E->p != NULL);
 362
 363         assert (ctx != NULL);
 364
 365         if (ECP_is_norm(P)) return 1;
 366         if (ECP_is_infty(P)) return 0;
 367
 368         if ((zm = BN_mod_inverse(P->Z, P->Z, E->p, ctx)) == NULL) return 0;
 369
 370         assert(!P->is_in_mont);
 371
 372
 373         BN_CTX_start(ctx);
 374         z = BN_CTX_get(ctx);
 375         if (z == NULL) goto err;
 376
 377         if (!BN_mod_mul(z, zm, zm, E->p, ctx)) goto err;
 378         if (!BN_mod_mul(P->X, P->X, z, E->p, ctx)) goto err;
 379
 380         if (!BN_mod_mul(z, z, zm, E->p, ctx)) goto err;
 381         if (!BN_mod_mul(P->Y, P->Y, z, E->p, ctx)) goto err;
 382
 383         if (!BN_one(P->Z)) goto err;
 384
 385         if (zm != NULL) BN_clear_free(zm);
 386
 387         BN_CTX_end(ctx);
 388         return 1;
 389
 390 err:
 391         if (zm != NULL) BN_clear_free(zm);
 392         BN_CTX_end(ctx);
 393         return 0;
 394         }
 395
 396
 397 int ECP_copy(EC_POINT *R, EC_POINT *P)
 398         {
 399         assert(P != NULL);
 400         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 401
 402         assert(R != NULL);
 403         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
 404
 405         if (BN_copy(R->X, P->X) == NULL) return 0;
 406         if (BN_copy(R->Y, P->Y) == NULL) return 0;
 407         if (BN_copy(R->Z, P->Z) == NULL) return 0;
 408         R->is_in_mont = P->is_in_mont;
 409
 410         return 1;
 411         }
 412
 413
 414 EC_POINT *ECP_dup(EC_POINT *P)
 415         {
 416         EC_POINT *ret;
 417
 418         ret = ECP_new();
 419         if (ret == NULL) return NULL;
 420
 421         if (!ECP_copy(ret, P))
 422                 {
 423                 ECP_clear_free(ret);
 424                 return(NULL);
 425                 }
 426
 427         return(ret);
 428         }
 429
 430
 431 EC_POINT *ECP_minus(EC_POINT *P, BIGNUM *p) /* mont || non-mont */
 432         {
 433         EC_POINT *ret;
 434
 435         assert(P != NULL);
 436         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 437
 438         assert(p != NULL);
 439
 440         assert(BN_cmp(P->Y, p) < 0);
 441
 442         ret = ECP_dup(P);
 443         if (ret == NULL) return NULL;
 444
 445         if (BN_is_zero(ret->Y)) return ret;
 446
 447         if (!BN_sub(ret->Y, p, ret->Y))
 448                 {
 449                 ECP_clear_free(ret);
 450                 return NULL;
 451                 }
 452
 453         return ret;
 454         }
 455
 456
 457 #ifdef SIMPLE
 458 int ECP_cmp(EC_POINT *P, EC_POINT *Q, BIGNUM *p, BN_CTX *ctx)
 459 /* return values:
 460         -2 ... error
 461          0 ... P = Q
 462         -1 ... P = -Q
 463          1 ... else
 464 */
 465         {
 466         BIGNUM *n0, *n1, *n2, *n3, *n4;
 467         int Pnorm, Qnorm;
 468
 469         assert(P != NULL);
 470         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 471
 472         assert(Q != NULL);
 473         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
 474
 475         assert(p != NULL);
 476         assert(ctx != NULL);
 477
 478         assert(!P->is_in_mont);
 479         assert(!Q->is_in_mont);
 480
 481         if (ECP_is_infty(P) && ECP_is_infty(Q)) return 0;
 482         if (ECP_is_infty(P) || ECP_is_infty(Q)) return 1;
 483
 484
 485         Pnorm = (ECP_is_norm(P));
 486         Qnorm = (ECP_is_norm(Q));
 487
 488         BN_CTX_start(ctx);
 489         n0 = BN_CTX_get(ctx);
 490         n1 = BN_CTX_get(ctx);
 491         n2 = BN_CTX_get(ctx);
 492         n3 = BN_CTX_get(ctx);
 493         n4 = BN_CTX_get(ctx);
 494         if (n4 == NULL) goto err;
 495
 496         if (Qnorm)
 497                 {
 498                 if (BN_copy(n1, P->X) == NULL) goto err;                        /* L1 = x_p */
 499                 if (BN_copy(n2, P->Y) == NULL) goto err;                        /* L2 = y_p */
 500                 }
 501         else
 502                 {
 503                 if (!BN_sqr(n0, Q->Z, ctx)) goto err;
 504                 if (!BN_mod_mul(n1, P->X, n0, p, ctx)) goto err;        /* L1 = x_p * z_q^2 */
 505
 506                 if (!BN_mod_mul(n0, n0, Q->Z, p, ctx)) goto err;
 507                 if (!BN_mod_mul(n2, P->Y, n0, p, ctx)) goto err;        /* L2 = y_p * z_q^3 */
 508                 }
 509
 510         if (Pnorm)
 511                 {
 512                 if (BN_copy(n3, Q->X) == NULL) goto err;                        /* L3 = x_q */
 513                 if (BN_copy(n4, Q->Y) == NULL) goto err;                        /* L4 = y_q */
 514                 }
 515         else
 516                 {
 517                 if (!BN_sqr(n0, P->Z, ctx)) goto err;
 518                 if (!BN_mod_mul(n3, Q->X, n0, p, ctx)) goto err;        /* L3 = x_q * z_p^2 */
 519
 520                 if (!BN_mod_mul(n0, n0, P->Z, p, ctx)) goto err;
 521                 if (!BN_mod_mul(n4, Q->Y, n0, p, ctx)) goto err;        /* L4 = y_q * z_p^3 */
 522                 }
 523
 524         if (!BN_mod_sub(n0, n1, n3, p, ctx)) goto err;                  /* L5 = L1 - L3 */
 525
 526         if (!BN_is_zero(n0))
 527                 {
 528                 BN_CTX_end(ctx);
 529                 return 1;
 530                 }
 531
 532         if (!BN_mod_sub(n0, n2, n4, p, ctx)) goto err;                  /* L6 = L2 - L4 */
 533
 534         if (!BN_is_zero(n0))
 535                 {
 536                 BN_CTX_end(ctx);
 537                 return -1;
 538                 }
 539
 540         BN_CTX_end(ctx);
 541         return 0;
 542
 543 err:
 544         BN_CTX_end(ctx);
 545         return -2;
 546         }
 547
 548
 549 int ECP_double(EC_POINT *R, EC_POINT *P, EC *E, BN_CTX *ctx)
 550 /* R <- 2P (on E) */
 551         {
 552         BIGNUM *n0, *n1, *n2, *n3, *p;
 553         int Pnorm, A0;
 554
 555         assert(P != NULL);
 556         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 557
 558         assert(R != NULL);
 559         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
 560
 561         assert(E != NULL);
 562         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
 563
 564         assert(ctx != NULL);
 565
 566         assert(!P->is_in_mont);
 567
 568         if (ECP_is_infty(P))
 569                 {
 570                 if (!BN_zero(R->Z)) return 0;
 571                 return 1;
 572                 }
 573
 574         Pnorm = (ECP_is_norm(P));
 575         A0 = (BN_is_zero(E->A));
 576
 577         BN_CTX_start(ctx);
 578         n0 = BN_CTX_get(ctx);
 579         n1 = BN_CTX_get(ctx);
 580         n2 = BN_CTX_get(ctx);
 581         n3 = BN_CTX_get(ctx);
 582         if (n3 == NULL) goto err;
 583
 584         p = E->p;
 585
 586         /* L1 */
 587         if (Pnorm || A0)
 588                 {
 589                 if (!BN_mod_sqr(n1, P->X, p, ctx)) goto err;
 590                 if (!BN_mul_word(n1, 3)) goto err;
 591                 if (!A0)                                   /* if A = 0: L1 = 3 * x^2 + a * z^4 = 3 * x ^2    */
 592                         if (!BN_mod_add(n1, n1, E->A, p, ctx)) goto err; /* L1 = 3 * x^2 + a * z^4 = 3 * x^2 + a */
 593                 }
 594         else
 595                 {
 596                 if (!BN_mod_sqr(n0, P->Z, p, ctx)) goto err;
 597                 if (!BN_mod_mul(n0, n0, n0, p, ctx)) goto err;
 598                 if (!BN_mod_mul(n0, n0, E->A, p, ctx)) goto err;
 599                 if (!BN_mod_sqr(n1, P->X, p, ctx)) goto err;
 600                 if (!BN_mul_word(n1, 3)) goto err;
 601                 if (!BN_mod_add(n1, n1, n0, p, ctx)) goto err;          /* L1 = 3 * x^2 + a * z^4 */
 602                 }
 603
 604         /* Z */
 605         if (Pnorm)
 606                 {
 607                 if (BN_copy(n0, P->Y) == NULL) goto err;
 608                 }
 609         else
 610                 {
 611                 if (!BN_mod_mul(n0, P->Y, P->Z, p, ctx)) goto err;
 612                 }
 613         if (!BN_lshift1(n0, n0)) goto err;
 614         if (!BN_smod(R->Z, n0, p, ctx)) goto err;                               /* Z = 2 * y * z */
 615
 616         /* L2 */
 617         if (!BN_mod_sqr(n3, P->Y, p, ctx)) goto err;
 618         if (!BN_mod_mul(n2, P->X, n3, p, ctx)) goto err;
 619         if (!BN_lshift(n2, n2, 2)) goto err;
 620         if (!BN_smod(n2, n2, p, ctx)) goto err;                                 /* L2 = 4 * x * y^2 */
 621
 622         /* X */
 623         if (!BN_lshift1(n0, n2)) goto err;
 624         if (!BN_mod_sqr(R->X, n1, p, ctx)) goto err;
 625         if (!BN_mod_sub(R->X, R->X, n0, p, ctx)) goto err;              /* X = L1^2 - 2 * L2 */
 626
 627         /* L3 */
 628         if (!BN_mod_sqr(n0, n3, p, ctx)) goto err;
 629         if (!BN_lshift(n3, n0, 3)) goto err;
 630         if (!BN_smod(n3, n3, p, ctx)) goto err;                                 /* L3 = 8 * y^4 */
 631
 632         /* Y */
 633         if (!BN_mod_sub(n0, n2, R->X, p, ctx)) goto err;
 634         if (!BN_mod_mul(n0, n1, n0, p, ctx)) goto err;
 635         if (!BN_mod_sub(R->Y, n0, n3, p, ctx)) goto err;                /* Y = L1 * (L2 - X) - L3 */
 636
 637
 638 #ifdef TEST
 639         if (!ECP_is_on_ec(R, E, ctx)) return 0;
 640 #endif
 641
 642         BN_CTX_end(ctx);
 643         return 1;
 644
 645 err:
 646         BN_CTX_end(ctx);
 647         return 0;
 648         }
 649
 650
 651 int ECP_add(EC_POINT *R, EC_POINT *P, EC_POINT *Q, EC *E, BN_CTX *ctx)
 652 /* R <- P + Q (on E) */
 653         {
 654         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *n6, *p;
 655         int Pnorm, Qnorm;
 656
 657         assert(P != NULL);
 658         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 659
 660         assert(Q != NULL);
 661         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
 662
 663         assert(R != NULL);
 664         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
 665
 666         assert(E != NULL);
 667         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
 668         assert(!BN_is_zero(E->h));;
 669
 670         assert(ctx != NULL);
 671
 672         assert(!P->is_in_mont);
 673         assert(!Q->is_in_mont);
 674
 675         if (P == Q) return ECP_double(R, P, E, ctx);
 676
 677         if (ECP_is_infty(P)) return ECP_copy(R, Q);
 678         if (ECP_is_infty(Q)) return ECP_copy(R, P);
 679
 680         Pnorm = (ECP_is_norm(P));
 681         Qnorm = (ECP_is_norm(Q));
 682
 683         BN_CTX_start(ctx);
 684         n0 = BN_CTX_get(ctx);
 685         n1 = BN_CTX_get(ctx);
 686         n2 = BN_CTX_get(ctx);
 687         n3 = BN_CTX_get(ctx);
 688         n4 = BN_CTX_get(ctx);
 689         n5 = BN_CTX_get(ctx);
 690         n6 = BN_CTX_get(ctx);
 691         if (n6 == NULL) goto err;
 692
 693         p = E->p;
 694
 695         /* L1; L2 */
 696         if (Qnorm)
 697                 {
 698                 if (BN_copy(n1, P->X) == NULL) goto err;         /* L1 = x_p */
 699                 if (BN_copy(n2, P->Y) == NULL) goto err;         /* L2 = y_p */
 700                 }
 701         else
 702                 {
 703                 if (!BN_sqr(n0, Q->Z, ctx)) goto err;
 704                 if (!BN_mod_mul(n1, P->X, n0, p, ctx)) goto err; /* L1 = x_p * z_q^2 */
 705
 706                 if (!BN_mod_mul(n0, n0, Q->Z, p, ctx)) goto err;
 707                 if (!BN_mod_mul(n2, P->Y, n0, p, ctx)) goto err; /* L2 = y_p * z_q^3 */
 708                 }
 709
 710         /* L3; L4 */
 711         if (Pnorm)
 712                 {
 713                 if (BN_copy(n3, Q->X) == NULL) goto err;         /* L3 = x_q */
 714                 if (BN_copy(n4, Q->Y) == NULL) goto err;         /* L4 = y_q */
 715                 }
 716         else
 717                 {
 718                 if (!BN_sqr(n0, P->Z, ctx)) goto err;
 719                 if (!BN_mod_mul(n3, Q->X, n0, p, ctx)) goto err; /* L3 = x_q * z_p^2 */
 720
 721                 if (!BN_mod_mul(n0, n0, P->Z, p, ctx)) goto err;
 722                 if (!BN_mod_mul(n4, Q->Y, n0, p, ctx)) goto err; /* L4 = y_q * z_p^3 */
 723                 }
 724
 725         /* L5; L6 */
 726         if (!BN_mod_sub(n5, n1, n3, p, ctx)) goto err;          /* L5 = L1 - L3 */
 727         if (!BN_mod_sub(n6, n2, n4, p, ctx)) goto err;          /* L6 = L2 - L4 */
 728
 729         /* pata */
 730         if (BN_is_zero(n5))
 731                 {
 732                 if (BN_is_zero(n6))     /* P = Q => P + Q = 2P */
 733                         {
 734                         BN_CTX_end(ctx);
 735                         return ECP_double(R, P, E, ctx);
 736                         }
 737                 else                             /* P = -Q => P + Q = \infty */
 738                         {
 739                         BN_CTX_end(ctx);
 740                         if (!BN_zero(R->Z)) return 0;
 741                         return 1;
 742                         }
 743                 }
 744
 745         /* L7; L8 */
 746         if (!BN_mod_add(n1, n1, n3, p, ctx)) goto err;          /* L7 = L1 + L3 */
 747         if (!BN_mod_add(n2, n2, n4, p, ctx)) goto err;          /* L8 = L2 + L4 */
 748
 749         /* Z */
 750         if (Pnorm)
 751                 {
 752                 if (BN_copy(n0, Q->Z) == NULL) goto err;
 753                 }
 754         else
 755                 {
 756                 if (!BN_mod_mul(n0, P->Z, Q->Z, p, ctx)) goto err;
 757                 }
 758         if (!BN_mod_mul(R->Z, n0, n5, p, ctx)) goto err;        /* Z = z_p * z_q * L_5 */
 759
 760         /* X */
 761         if (!BN_mod_sqr(n0, n6, p, ctx)) goto err;
 762         if (!BN_mod_sqr(n4, n5, p, ctx)) goto err;
 763         if (!BN_mod_mul(n3, n1, n4, p, ctx)) goto err;
 764         if (!BN_mod_sub(R->X, n0, n3, p, ctx)) goto err;        /* X = L6^2 - L5^2 * L7 */
 765
 766         /* L9 */
 767         if (!BN_lshift1(n0, R->X)) goto err;
 768         if (!BN_mod_sub(n0, n3, n0, p, ctx)) goto err;          /* L9 = L5^2 * L7 - 2X */
 769
 770         /* Y */
 771         if (!BN_mod_mul(n0, n0, n6, p, ctx)) goto err;
 772         if (!BN_mod_mul(n5, n4, n5, p, ctx)) goto err;
 773         if (!BN_mod_mul(n1, n2, n5, p, ctx)) goto err;
 774         if (!BN_mod_sub(n0, n0, n1, p, ctx)) goto err;
 775         if (!BN_mod_mul(R->Y, n0, E->h, p, ctx)) goto err;      /* Y = (L6 * L9 - L8 * L5^3) / 2 */
 776
 777
 778
 779 #ifdef TEST
 780         if (!ECP_is_on_ec(R, E, ctx)) return 0;
 781 #endif
 782
 783         BN_CTX_end(ctx);
 784         return 1;
 785
 786 err:
 787         BN_CTX_end(cxt);
 788         return 0;
 789         }
 790
 791
 792 ECP_PRECOMPUTE *ECP_precompute(int r, EC_POINT *P, EC *E, BN_CTX *ctx)
 793         {
 794         ECP_PRECOMPUTE *ret;
 795         EC_POINT *P2;
 796         int i, max;
 797
 798         assert(r > 2);
 799         assert(!P->is_in_mont);
 800         assert(!E->is_in_mont);
 801
 802         ret=(ECP_PRECOMPUTE *)malloc(sizeof(ECP_PRECOMPUTE));
 803         if (ret == NULL) return NULL;
 804
 805         max = 1;
 806         max <<= (r - 1);
 807
 808         ret->r = 0;
 809
 810         ret->Pi=(EC_POINT **)malloc(sizeof(EC_POINT *) * max);
 811         if (ret->Pi == NULL) goto err;
 812
 813
 814         /* P2 = [2]P */
 815         if ((P2 = ECP_new()) == NULL) goto err;
 816         if (!ECP_double(P2, P, E, ctx)) goto err;
 817
 818         /* P_0 = P */
 819         if((ret->Pi[0] = ECP_dup(P)) == NULL) goto err;
 820
 821
 822         /* P_i = P_(i-1) + P2 */
 823         for (i = 1; i < max; i++)
 824                 {
 825                 if ((ret->Pi[i] = ECP_new()) == NULL) goto err;
 826
 827                 if (!ECP_add(ret->Pi[i], P2, ret->Pi[i - 1], E, ctx)) goto err;
 828                 }
 829
 830         ret->r = r;
 831         ECP_clear_free(P2);
 832
 833         return ret;
 834
 835 err:
 836         ECP_clear_free(P2);
 837         ECP_clear_free_precompute(ret);
 838         return NULL;
 839         }
 840
 841
 842 int ECP_multiply(EC_POINT *R, BIGNUM *k, ECP_PRECOMPUTE *prec, EC *E, BN_CTX *ctx)
 843 /* R = [k]P */
 844         {
 845         int j;
 846         int t, nextw, h, r;
 847
 848         assert(R != NULL);
 849         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
 850
 851         assert(E != NULL);
 852         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
 853
 854         assert(k != NULL);
 855         assert(!k->neg);
 856
 857         assert(ctx != NULL);
 858         assert(prec != NULL);
 859
 860         assert(!E->is_in_mont);
 861
 862         if (BN_is_zero(k))
 863                 {
 864                 if (!BN_zero(R->Z)) return 0;
 865                 R->is_in_mont = 0;
 866                 return 1;
 867                 }
 868
 869
 870         j = BN_num_bits(k);
 871         j--;
 872
 873         r = prec->r;
 874
 875         if (!BN_zero(R->Z)) return 0;
 876         R->is_in_mont = 0;
 877
 878         while(j >= 0)
 879                 {
 880                 if (!BN_is_bit_set(k, j))
 881                         {
 882                         if (!ECP_double(R, R, E, ctx)) return 0;
 883                         j--;
 884                         }
 885                 else
 886                         {
 887                         nextw = j - r;
 888                         if (nextw < -1) nextw = -1;
 889                         t = nextw + 1;
 890                         while(!BN_is_bit_set(k, t))
 891                                 t++;
 892
 893                         if (!ECP_double(R, R, E, ctx)) return 0;
 894
 895                         j--;
 896                         if (j < t) h = 0;
 897                         else
 898                                 {
 899                                 h = 1;
 900                                 for(; j > t; j--)
 901                                         {
 902                                         h <<= 1;
 903                                         if (BN_is_bit_set(k, j)) h++;
 904                                         if (!ECP_double(R, R, E, ctx)) return 0;
 905                                         }
 906                                 if (!ECP_double(R, R, E, ctx)) return 0;
 907                                 j--;
 908                                 }
 909
 910                         if (!ECP_add(R, R, prec->Pi[h], E, ctx)) return 0;
 911
 912                         for (; j > nextw; j--)
 913                                 {
 914                                 if (!ECP_double(R, R, E, ctx)) return 0;
 915                                 }
 916
 917                         }
 918                 }
 919
 920         return 1;
 921         }
 922
 923 #endif /* SIMPLE */
 924
 925
 926 #ifdef MONTGOMERY
 927
 928 int ECP_to_montgomery(EC_POINT *P, BN_MONTGOMERY *mont, BN_CTX *ctx)
 929         {
 930         assert(P != NULL);
 931         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 932
 933         assert(mont != NULL);
 934         assert(mont->p != NULL);
 935
 936         assert(ctx != NULL);
 937
 938         if (P->is_in_mont) return 1;
 939
 940         if (!BN_lshift(P->X, P->X, mont->R_num_bits)) return 0;
 941         if (!BN_mod(P->X, P->X, mont->p, ctx)) return 0;
 942
 943         if (!BN_lshift(P->Y, P->Y, mont->R_num_bits)) return 0;
 944         if (!BN_mod(P->Y, P->Y, mont->p, ctx)) return 0;
 945
 946         if (!BN_lshift(P->Z, P->Z, mont->R_num_bits)) return 0;
 947         if (!BN_mod(P->Z, P->Z, mont->p, ctx)) return 0;
 948
 949         P->is_in_mont = 1;
 950         return 1;
 951         }
 952
 953
 954 int ECP_from_montgomery(EC_POINT *P, BN_MONTGOMERY *mont, BN_CTX *ctx)
 955         {
 956
 957         assert(P != NULL);
 958         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 959
 960         assert(mont != NULL);
 961         assert(mont->p != NULL);
 962
 963         assert(ctx != NULL);
 964
 965         if (!P->is_in_mont) return 1;
 966
 967         if (!BN_mont_red(P->X, mont)) return 0;
 968         if (!BN_mont_red(P->Y, mont)) return 0;
 969         if (!BN_mont_red(P->Z, mont)) return 0;
 970
 971         P->is_in_mont = 0;
 972         return 1;
 973         }
 974
 975
 976 int ECP_mont_cmp(EC_POINT *P, EC_POINT *Q, BN_MONTGOMERY *mont, BN_CTX *ctx)
 977 /* return values:
 978         -2 ... error
 979          0 ... P = Q
 980         -1 ... P = -Q
 981          1 ... else
 982 */
 983         {
 984         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *p;
 985
 986         assert(P != NULL);
 987         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
 988
 989         assert(Q != NULL);
 990         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
 991
 992         assert(mont != NULL);
 993         assert(mont->p != NULL);
 994
 995         assert(ctx != NULL);
 996
 997         if (!P->is_in_mont)
 998                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
 999
1000         if (!Q->is_in_mont)
1001                 if (!ECP_to_montgomery(Q, mont, ctx)) return 0;
1002
1003
1004         if (ECP_is_infty(P) && ECP_is_infty(Q)) return 0;
1005         if (ECP_is_infty(P) || ECP_is_infty(Q)) return 1;
1006
1007
1008         BN_CTX_start(ctx);
1009         n0 = BN_CTX_get(ctx);
1010         n1 = BN_CTX_get(ctx);
1011         n2 = BN_CTX_get(ctx);
1012         n3 = BN_CTX_get(ctx);
1013         n4 = BN_CTX_get(ctx);
1014         n5 = BN_CTX_get(ctx);
1015         if (n5 == 0) goto err;
1016
1017
1018         p = mont->p;
1019
1020
1021         if (!BN_mont_mod_mul(n5, Q->Z, Q->Z, mont)) goto err;
1022         if (!BN_mont_mod_mul(n1, P->X, n5, mont)) goto err;     /* L1 = x_p * z_q^2 */
1023
1024         if (!BN_mont_mod_mul(n0, n5, Q->Z, mont)) goto err;
1025         if (!BN_mont_mod_mul(n2, P->Y, n0, mont)) goto err;     /* L2 = y_p * z_q^3 */
1026
1027         if (!BN_mont_mod_mul(n5, P->Z, P->Z, mont)) goto err;
1028         if (!BN_mont_mod_mul(n3, Q->X, n5, mont)) goto err;     /* L3 = x_q * z_p^2 */
1029
1030         if (!BN_mont_mod_mul(n0, n5, P->Z, mont)) goto err;
1031         if (!BN_mont_mod_mul(n4, Q->Y, n0, mont)) goto err;     /* L4 = y_q * z_p^3 */
1032
1033
1034         if (!BN_mod_sub_quick(n0, n1, n3, p)) goto err;                 /* L5 = L1 - L3 */
1035
1036         if (!BN_is_zero(n0))
1037                 {
1038                 BN_CTX_end(ctx);
1039                 return 1;
1040                 }
1041
1042         if (!BN_mod_sub_quick(n0, n2, n4, p)) goto err;                 /* L6 = L2 - L4 */
1043
1044         if (!BN_is_zero(n0))
1045                 {
1046                 BN_CTX_end(ctx);
1047                 return -1;
1048                 }
1049
1050         BN_CTX_end(ctx);
1051         return 0;
1052
1053 err:
1054         BN_CTX_end(ctx);
1055         return -2;
1056         }
1057
1058
1059 int ECP_mont_double(EC_POINT *R, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1060 /* R <- 2P (on E) */
1061         {
1062         BIGNUM *n0, *n1, *n2, *n3, *p;
1063
1064         assert(P != NULL);
1065         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1066
1067         assert(R != NULL);
1068         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1069
1070         assert(E != NULL);
1071         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1072
1073         assert(ctx != NULL);
1074
1075         if (!P->is_in_mont)
1076                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1077
1078         if (!E->is_in_mont)
1079                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1080
1081         R->is_in_mont = 1;
1082
1083         if (ECP_is_infty(P))
1084                 {
1085                 if (!BN_zero(R->Z)) return 0;
1086                 return 1;
1087                 }
1088
1089
1090         BN_CTX_start(ctx);
1091         n0 = BN_CTX_get(ctx);
1092         n1 = BN_CTX_get(ctx);
1093         n2 = BN_CTX_get(ctx);
1094         n3 = BN_CTX_get(ctx);
1095         if (n3 == 0) goto err;
1096
1097         p = E->p;
1098
1099         /* L1 */
1100         if (!BN_mont_mod_mul(n0, P->Z, P->Z, mont)) goto err;
1101         if (!BN_mont_mod_mul(n2, n0, n0, mont)) goto err;
1102         if (!BN_mont_mod_mul(n0, n2, E->A, mont)) goto err;
1103         if (!BN_mont_mod_mul(n1, P->X, P->X, mont)) goto err;
1104         if (!BN_mod_lshift1_quick(n2, n1, p)) goto err;
1105         if (!BN_mod_add_quick(n1, n1, n2, p)) goto err;
1106         if (!BN_mod_add_quick(n1, n1, n0, p)) goto err;         /* L1 = 3 * x^2 + a * z^4 */
1107
1108         /* Z */
1109         if (!BN_mont_mod_mul(n0, P->Y, P->Z, mont)) goto err;
1110         if (!BN_mod_lshift1_quick(R->Z, n0, p)) goto err;               /* Z = 2 * y * z */
1111
1112         /* L2 */
1113         if (!BN_mont_mod_mul(n3, P->Y, P->Y, mont)) goto err;
1114         if (!BN_mont_mod_mul(n2, P->X, n3, mont)) goto err;
1115         if (!BN_mod_lshift_quick(n2, n2, 2, p)) goto err;               /* L2 = 4 * x * y^2 */
1116
1117         /* X */
1118         if (!BN_mod_lshift1_quick(n0, n2, p)) goto err;
1119         if (!BN_mont_mod_mul(R->X, n1, n1, mont)) goto err;
1120         if (!BN_mod_sub_quick(R->X, R->X, n0, p)) goto err;     /* X = L1^2 - 2 * L2 */
1121
1122         /* L3 */
1123         if (!BN_mont_mod_mul(n0, n3, n3, mont)) goto err;
1124         if (!BN_mod_lshift_quick(n3, n0, 3, p)) goto err;               /* L3 = 8 * y^4 */
1125
1126
1127         /* Y */
1128         if (!BN_mod_sub_quick(n2, n2, R->X, p)) goto err;
1129         if (!BN_mont_mod_mul(n0, n1, n2, mont)) goto err;
1130         if (!BN_mod_sub_quick(R->Y, n0, n3, p)) goto err;               /* Y = L1 * (L2 - X) - L3 */
1131
1132         BN_CTX_end(ctx);
1133         return 1;
1134
1135 err:
1136         BN_CTX_end(ctx);
1137         return 0;
1138         }
1139
1140
1141 int ECP_mont_add(EC_POINT *R, EC_POINT *P, EC_POINT *Q, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1142 /* R <- P + Q (on E) */
1143         {
1144         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *n6, *p;
1145
1146         assert(P != NULL);
1147         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1148
1149         assert(Q != NULL);
1150         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
1151
1152         assert(R != NULL);
1153         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1154
1155         assert(E != NULL);
1156         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1157         assert(!BN_is_zero(E->h));;
1158
1159         assert(ctx != NULL);
1160
1161         if (!Q->is_in_mont)
1162                 if (!ECP_to_montgomery(Q, mont, ctx)) return 0;
1163
1164         if (!P->is_in_mont)
1165                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1166
1167         if (!E->is_in_mont)
1168                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1169
1170         if (P == Q) return ECP_mont_double(R, P, E, mont, ctx);
1171
1172         if (ECP_is_infty(P)) return ECP_copy(R, Q);
1173         if (ECP_is_infty(Q)) return ECP_copy(R, P);
1174
1175
1176         BN_CTX_start(ctx);
1177         n0 = BN_CTX_get(ctx);
1178         n1 = BN_CTX_get(ctx);
1179         n2 = BN_CTX_get(ctx);
1180         n3 = BN_CTX_get(ctx);
1181         n4 = BN_CTX_get(ctx);
1182         n5 = BN_CTX_get(ctx);
1183         n6 = BN_CTX_get(ctx);
1184         if (n6 == NULL) goto err;
1185
1186
1187         p = E->p;
1188
1189         R->is_in_mont = 1;
1190
1191         /* L1; L2 */
1192         if (!BN_mont_mod_mul(n6, Q->Z, Q->Z, mont)) goto err;
1193         if (!BN_mont_mod_mul(n1, P->X, n6, mont)) goto err;     /* L1 = x_p * z_q^2 */
1194
1195         if (!BN_mont_mod_mul(n0, n6, Q->Z, mont)) goto err;
1196         if (!BN_mont_mod_mul(n2, P->Y, n0, mont)) goto err;     /* L2 = y_p * z_q^3 */
1197
1198
1199         /* L3; L4 */
1200         if (!BN_mont_mod_mul(n6, P->Z, P->Z, mont)) goto err;
1201         if (!BN_mont_mod_mul(n3, Q->X, n6, mont)) goto err;     /* L3 = x_q * z_p^2 */
1202
1203         if (!BN_mont_mod_mul(n0, n6, P->Z, mont)) goto err;
1204         if (!BN_mont_mod_mul(n4, Q->Y, n0, mont)) goto err;     /* L4 = y_q * z_p^3 */
1205
1206
1207         /* L5; L6 */
1208         if (!BN_mod_sub_quick(n5, n1, n3, p)) goto err;                 /* L5 = L1 - L3 */
1209         if (!BN_mod_sub_quick(n6, n2, n4, p)) goto err;                 /*L6 = L2 - L4 */
1210
1211
1212         /* pata */
1213         if (BN_is_zero(n5))
1214                 {
1215                 if (BN_is_zero(n6))  /* P = Q => P + Q = 2P */
1216                         {
1217                         BN_CTX_end(ctx);
1218                         return ECP_mont_double(R, P, E, mont, ctx);
1219                         }
1220                 else                             /* P = -Q => P + Q = \infty */
1221                         {
1222                         BN_CTX_end(ctx);
1223                         if (!BN_zero(R->Z)) return 0;
1224                         return 1;
1225                         }
1226                 }
1227
1228         /* L7; L8 */
1229         if (!BN_mod_add_quick(n1, n1, n3, p)) goto err;                 /* L7 = L1 + L3 */
1230         if (!BN_mod_add_quick(n2, n2, n4, p)) goto err;                 /* L8 = L2 + L4 */
1231
1232
1233         /* Z */
1234         if (!BN_mont_mod_mul(n0, P->Z, Q->Z, mont)) goto err;
1235         if (!BN_mont_mod_mul(R->Z, n0, n5, mont)) goto err;     /* Z = z_p * z_q * L_5 */
1236
1237
1238         /* X */
1239         if (!BN_mont_mod_mul(n0, n6, n6, mont)) goto err;
1240         if (!BN_mont_mod_mul(n4, n5, n5, mont)) goto err;
1241         if (!BN_mont_mod_mul(n3, n1, n4, mont)) goto err;
1242         if (!BN_mod_sub_quick(R->X, n0, n3, p)) goto err;                       /* X = L6^2 - L5^2 * L7 */
1243
1244
1245         /* L9 */
1246         if (!BN_mod_lshift1_quick(n0, R->X, p)) goto err;
1247         if (!BN_mod_sub_quick(n3, n3, n0, p)) goto err;                 /* L9 = L5^2 * L7 - 2X */
1248
1249
1250         /* Y */
1251         if (!BN_mont_mod_mul(n0, n3, n6, mont)) goto err;
1252         if (!BN_mont_mod_mul(n6, n4, n5, mont)) goto err;
1253         if (!BN_mont_mod_mul(n1, n2, n6, mont)) goto err;
1254         if (!BN_mod_sub_quick(n0, n0, n1, p)) goto err;
1255         if (!BN_mont_mod_mul(R->Y, n0, E->h, mont)) goto err;   /* Y = (L6 * L9 - L8 * L5^3) / 2 */
1256
1257
1258         BN_CTX_end(ctx);
1259         return 1;
1260
1261 err:
1262         BN_CTX_end(ctx);
1263         return 0;
1264         }
1265
1266
1267 ECP_PRECOMPUTE *ECP_mont_precompute(int r, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1268         {
1269         ECP_PRECOMPUTE *ret;
1270         EC_POINT *P2;
1271         int i, max;
1272
1273         assert(r > 2);
1274         assert(r < sizeof(unsigned int) * 8 - 1);
1275
1276         assert(mont != NULL);
1277         assert(mont->p != NULL);
1278
1279         if (!P->is_in_mont)
1280                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1281
1282         if (!E->is_in_mont)
1283                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1284
1285         ret=(ECP_PRECOMPUTE *)malloc(sizeof(ECP_PRECOMPUTE));
1286         if (ret == NULL) return NULL;
1287
1288         max = 1;
1289         max <<= (r - 1);
1290
1291         ret->r = 0;
1292
1293         ret->Pi=(EC_POINT **)malloc(sizeof(EC_POINT *) * max);
1294         if (ret->Pi == NULL) goto err;
1295
1296
1297         /* P2 = [2]P */
1298         if ((P2 = ECP_new()) == NULL) goto err;
1299         if (!ECP_mont_double(P2, P, E, mont, ctx)) goto err;
1300
1301         /* P_0 = P */
1302         if((ret->Pi[0] = ECP_dup(P)) == NULL) goto err;
1303
1304
1305         /* P_i = P_(i-1) + P2 */
1306         for (i = 1; i < max; i++)
1307                 {
1308                 if ((ret->Pi[i] = ECP_new()) == NULL) goto err;
1309                 if (!ECP_mont_add(ret->Pi[i], P2, ret->Pi[i - 1], E, mont, ctx)) goto err;
1310                 }
1311
1312         ret->r = r;
1313         ECP_clear_free(P2);
1314
1315         return ret;
1316
1317 err:
1318         ECP_clear_free(P2);
1319         ECP_clear_free_precompute(ret);
1320         return NULL;
1321         }
1322
1323
1324 int ECP_mont_multiply(EC_POINT *R, BIGNUM *k, ECP_PRECOMPUTE *prec, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1325 /* R = [k]P   P = prec->Pi[0]*/
1326         {
1327         int j;
1328         int t, nextw, h, r;
1329
1330         assert(R != NULL);
1331         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1332
1333         assert(E != NULL);
1334         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1335
1336         assert(k != NULL);
1337         assert(!k->neg);
1338
1339         assert(ctx != NULL);
1340         assert(prec != NULL);
1341
1342         assert(mont != NULL);
1343         assert(mont->p != NULL);
1344
1345         if (!E->is_in_mont)
1346                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1347
1348
1349         if (BN_is_zero(k))
1350                 {
1351                 if (!BN_zero(R->Z)) return 0;
1352                 R->is_in_mont = 1;
1353                 return 1;
1354                 }
1355
1356         j = BN_num_bits(k);
1357         j--;
1358
1359         r = prec->r;
1360
1361         if (!BN_zero(R->Z)) return 0;
1362         R->is_in_mont = 1;
1363
1364         while(j >= 0)
1365                 {
1366                 if (!BN_is_bit_set(k, j))
1367                         {
1368                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1369                         j--;
1370                         }
1371                 else
1372                         {
1373                         nextw = j - r;
1374                         if (nextw < -1) nextw = -1;
1375                         t = nextw + 1;
1376                         while(!BN_is_bit_set(k, t))
1377                                 t++;
1378
1379                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1380
1381                         j--;
1382                         if (j < t) h = 0;
1383                         else
1384                                 {
1385                                 h = 1;
1386                                 for(; j > t; j--)
1387                                         {
1388                                         h <<= 1;
1389                                         if (BN_is_bit_set(k, j)) h++;
1390                                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1391                                         }
1392                                 if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1393                                 j--;
1394                                 }
1395
1396                         if (!ECP_mont_add(R, R, prec->Pi[h], E, mont, ctx)) return 0;
1397
1398                         for (; j > nextw; j--)
1399                                 {
1400                                 if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1401                                 }
1402
1403                         }
1404                 }
1405
1406         return 1;
1407         }
1408
1409
1410 int ECP_mont_multiply2(EC_POINT *R, BIGNUM *k, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1411 /* R = [k]P */
1412         {
1413         int j, hj, kj;
1414         BIGNUM *h;
1415         EC_POINT *mP;
1416
1417         assert(R != NULL);
1418         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1419
1420         assert(P != NULL);
1421         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1422
1423         assert(E != NULL);
1424         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1425
1426         assert(k != NULL);
1427         assert(!k->neg);
1428
1429         assert(ctx != NULL);
1430
1431         assert(mont != NULL);
1432         assert(mont->p != NULL);
1433
1434         if (!E->is_in_mont)
1435                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1436
1437         if (!P->is_in_mont)
1438                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1439
1440
1441         if (BN_is_zero(k))
1442                 {
1443                 if (!BN_zero(R->Z)) return 0;
1444                 R->is_in_mont = 1;
1445                 return 1;
1446                 }
1447
1448         if ((h = BN_dup(k)) == NULL) return 0;
1449
1450         if (!BN_lshift1(h, h)) goto err;
1451         if (!BN_add(h, h, k)) goto err;
1452
1453         if (!ECP_copy(R, P)) goto err;
1454         if ((mP = ECP_mont_minus(P, mont)) == NULL) goto err;
1455
1456         for(j = BN_num_bits(h) - 2; j > 0; j--)
1457                 {
1458                 if (!ECP_mont_double(R, R, E, mont, ctx)) goto err;
1459                 kj = BN_is_bit_set(k, j);
1460                 hj = BN_is_bit_set(h, j);
1461                 if (hj == 1 && kj == 0)
1462                         if (!ECP_mont_add(R, R, P, E, mont, ctx)) goto err;
1463                 if (hj == 0 && kj == 1)
1464                         if (!ECP_mont_add(R, R, mP, E, mont, ctx)) goto err;
1465                 }
1466
1467         if (h != NULL) BN_free(h);
1468         if (mP != NULL) ECP_clear_free(mP);
1469         return 1;
1470
1471 err:
1472         if (h != NULL) BN_free(h);
1473         if (mP != NULL) ECP_clear_free(mP);
1474         return 0;
1475         }
1476
1477 #endif /* MONTGOMERY */