package/network/services/ead/src/tinysrp/t_conf.c

   1 /*
   2  * Copyright (c) 1997-1999  The Stanford SRP Authentication Project
   3  * All Rights Reserved.
   4  *
   5  * Permission is hereby granted, free of charge, to any person obtaining
   6  * a copy of this software and associated documentation files (the
   7  * "Software"), to deal in the Software without restriction, including
   8  * without limitation the rights to use, copy, modify, merge, publish,
   9  * distribute, sublicense, and/or sell copies of the Software, and to
  10  * permit persons to whom the Software is furnished to do so, subject to
  11  * the following conditions:
  12  *
  13  * The above copyright notice and this permission notice shall be
  14  * included in all copies or substantial portions of the Software.
  15  *
  16  * THE SOFTWARE IS PROVIDED "AS-IS" AND WITHOUT WARRANTY OF ANY KIND,
  17  * EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
  18  * WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
  19  *
  20  * IN NO EVENT SHALL STANFORD BE LIABLE FOR ANY SPECIAL, INCIDENTAL,
  21  * INDIRECT OR CONSEQUENTIAL DAMAGES OF ANY KIND, OR ANY DAMAGES WHATSOEVER
  22  * RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER OR NOT ADVISED OF
  23  * THE POSSIBILITY OF DAMAGE, AND ON ANY THEORY OF LIABILITY, ARISING OUT
  24  * OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
  25  *
  26  * In addition, the following conditions apply:
  27  *
  28  * 1. Any software that incorporates the SRP authentication technology
  29  *    must display the following acknowlegment:
  30  *    "This product uses the 'Secure Remote Password' cryptographic
  31  *     authentication system developed by Tom Wu (tjw@CS.Stanford.EDU)."
  32  *
  33  * 2. Any software that incorporates all or part of the SRP distribution
  34  *    itself must also display the following acknowledgment:
  35  *    "This product includes software developed by Tom Wu and Eugene
  36  *     Jhong for the SRP Distribution (http://srp.stanford.edu/srp/)."
  37  *
  38  * 3. Redistributions in source or binary form must retain an intact copy
  39  *    of this copyright notice and list of conditions.
  40  */
  41
  42 #include <stdio.h>
  43
  44 #include "t_defines.h"
  45 #include "t_pwd.h"
  46 #include "t_read.h"
  47 #include "bn.h"
  48 #include "bn_lcl.h"
  49 #include "bn_prime.h"
  50
  51 #define TABLE_SIZE      32
  52
  53 static int witness(BIGNUM *w, const BIGNUM *a, const BIGNUM *a1,
  54         const BIGNUM *a1_odd, int k, BN_CTX *ctx, BN_MONT_CTX *mont);
  55
  56 /*
  57  * This is the safe prime generation logic.
  58  * To generate a safe prime p (where p = 2q+1 and q is prime), we start
  59  * with a random odd q that is one bit shorter than the desired length
  60  * of p.  We use a simple 30-element sieve to filter the values of q
  61  * and consider only those that are 11, 23, or 29 (mod 30).  (If q were
  62  * anything else, either q or p would be divisible by 2, 3, or 5).
  63  * For the values of q that are left, we apply the following tests in
  64  * this order:
  65  *
  66  *   trial divide q
  67  *   let p = 2q + 1
  68  *   trial divide p
  69  *   apply Fermat test to q (2^q == 2 (mod q))
  70  *   apply Fermat test to p (2^p == 2 (mod p))
  71  *   apply real probablistic primality test to q
  72  *   apply real probablistic primality test to p
  73  *
  74  * A number that passes all these tests is considered a safe prime for
  75  * our purposes.  The tests are ordered this way for efficiency; the
  76  * slower tests are run rarely if ever at all.
  77  */
  78
  79 static int
  80 trialdiv(x)
  81      const BigInteger x;
  82 {
  83   static int primes[] = {               /* All odd primes < 256 */
  84       3,   5,   7,  11,  13,  17,  19,  23,  29,
  85      31,  37,  41,  43,  47,  53,  59,  61,  67,
  86      71,  73,  79,  83,  89,  97, 101, 103,
  87     107, 109, 113, 127, 131, 137, 139, 149, 151,
  88     157, 163, 167, 173, 179, 181, 191, 193, 197,
  89     199, 211, 223, 227, 229, 233, 239, 241, 251
  90   };
  91   static int nprimes = sizeof(primes) / sizeof(int);
  92   int i;
  93
  94   for(i = 0; i < nprimes; ++i) {
  95     if(BigIntegerModInt(x, primes[i]) == 0)
  96       return primes[i];
  97   }
  98   return 1;
  99 }
 100
 101 /* x + sieve30[x%30] == 11, 23, or 29 (mod 30) */
 102
 103 static int sieve30[] =
 104 {  11, 10,  9,  8,  7,  6,  5,  4,  3,  2,
 105     1, 12, 11, 10,  9,  8,  7,  6,  5,  4,
 106     3,  2,  1,  6,  5,  4,  3,  2,  1, 12
 107 };
 108
 109 /* Find a Sophie-Germain prime between "lo" and "hi".  NOTE: this is not
 110    a "safe prime", but the smaller prime.  Take 2q+1 to get the safe prime. */
 111
 112 static void
 113 sophie_germain(q, lo, hi)
 114      BigInteger q;              /* assumed initialized */
 115      const BigInteger lo;
 116      const BigInteger hi;
 117 {
 118   BigInteger m, p, r;
 119   char parambuf[MAXPARAMLEN];
 120   int foundprime = 0;
 121   int i, mod30;
 122
 123   m = BigIntegerFromInt(0);
 124   BigIntegerSub(m, hi, lo);
 125   i = (BigIntegerBitLen(m) + 7) / 8;
 126   t_random(parambuf, i);
 127   r = BigIntegerFromBytes(parambuf, i);
 128   BigIntegerMod(r, r, m);
 129
 130   BigIntegerAdd(q, r, lo);
 131   if(BigIntegerModInt(q, 2) == 0)
 132     BigIntegerAddInt(q, q, 1);          /* make q odd */
 133
 134   mod30 = BigIntegerModInt(q, 30);      /* mod30 = q % 30 */
 135
 136   BigIntegerFree(m);
 137   m = BigIntegerFromInt(2);                     /* m = 2 */
 138   p = BigIntegerFromInt(0);
 139
 140   while(BigIntegerCmp(q, hi) < 0) {
 141     if(trialdiv(q) < 2) {
 142       BigIntegerMulInt(p, q, 2);                        /* p = 2 * q */
 143       BigIntegerAddInt(p, p, 1);                /* p += 1 */
 144       if(trialdiv(p) < 2) {
 145         BigIntegerModExp(r, m, q, q);           /* r = 2^q % q */
 146         if(BigIntegerCmpInt(r, 2) == 0) {       /* if(r == 2) */
 147           BigIntegerModExp(r, m, p, p);         /* r = 2^p % p */
 148           if(BigIntegerCmpInt(r, 2) == 0) {     /* if(r == 2) */
 149             if(BigIntegerCheckPrime(q) && BigIntegerCheckPrime(p)) {
 150               ++foundprime;
 151               break;
 152             }
 153           }
 154         }
 155       }
 156     }
 157
 158     i = sieve30[mod30];
 159     BigIntegerAddInt(q, q, i);          /* q += i */
 160     mod30 = (mod30 + i) % 30;
 161   }
 162
 163   /* should wrap around on failure */
 164   if(!foundprime) {
 165     fprintf(stderr, "Prime generation failed!\n");
 166     exit(1);
 167   }
 168
 169   BigIntegerFree(r);
 170   BigIntegerFree(m);
 171   BigIntegerFree(p);
 172 }
 173
 174 _TYPE( struct t_confent * )
 175 t_makeconfent(tc, nsize)
 176      struct t_conf * tc;
 177      int nsize;
 178 {
 179   BigInteger n, g, q, t, u;
 180
 181   t = BigIntegerFromInt(0);
 182   u = BigIntegerFromInt(1);             /* u = 1 */
 183   BigIntegerLShift(t, u, nsize - 2);    /* t = 2^(nsize-2) */
 184   BigIntegerMulInt(u, t, 2);            /* u = 2^(nsize-1) */
 185
 186   q = BigIntegerFromInt(0);
 187   sophie_germain(q, t, u);
 188
 189   n = BigIntegerFromInt(0);
 190   BigIntegerMulInt(n, q, 2);
 191   BigIntegerAddInt(n, n, 1);
 192
 193   /* Look for a generator mod n */
 194   g = BigIntegerFromInt(2);
 195   while(1) {
 196     BigIntegerModExp(t, g, q, n);               /* t = g^q % n */
 197     if(BigIntegerCmpInt(t, 1) == 0)             /* if(t == 1) */
 198       BigIntegerAddInt(g, g, 1);                /* ++g */
 199     else
 200       break;
 201   }
 202   BigIntegerFree(t);
 203   BigIntegerFree(u);
 204   BigIntegerFree(q);
 205
 206   tc->tcbuf.modulus.data = tc->modbuf;
 207   tc->tcbuf.modulus.len = BigIntegerToBytes(n, tc->tcbuf.modulus.data);
 208   BigIntegerFree(n);
 209
 210   tc->tcbuf.generator.data = tc->genbuf;
 211   tc->tcbuf.generator.len = BigIntegerToBytes(g, tc->tcbuf.generator.data);
 212   BigIntegerFree(g);
 213
 214   tc->tcbuf.index = 1;
 215   return &tc->tcbuf;
 216 }
 217
 218 _TYPE( struct t_confent * )
 219 t_makeconfent_c(tc, nsize)
 220      struct t_conf * tc;
 221      int nsize;
 222 {
 223   BigInteger g, n, p, q, j, k, t, u;
 224   int psize, qsize;
 225
 226   psize = nsize / 2;
 227   qsize = nsize - psize;
 228
 229   t = BigIntegerFromInt(1);             /* t = 1 */
 230   u = BigIntegerFromInt(0);
 231   BigIntegerLShift(u, t, psize - 3);    /* u = t*2^(psize-3) = 2^(psize-3) */
 232   BigIntegerMulInt(t, u, 3);                    /* t = 3*u = 1.5*2^(psize-2) */
 233   BigIntegerAdd(u, u, t);                       /* u += t [u = 2^(psize-1)] */
 234   j = BigIntegerFromInt(0);
 235   sophie_germain(j, t, u);
 236
 237   k = BigIntegerFromInt(0);
 238   if(qsize != psize) {
 239     BigIntegerFree(t);
 240     t = BigIntegerFromInt(1);           /* t = 1 */
 241     BigIntegerLShift(u, t, qsize - 3);  /* u = t*2^(qsize-3) = 2^(qsize-3) */
 242     BigIntegerMulInt(t, u, 3);          /* t = 3*u = 1.5*2^(qsize-2) */
 243     BigIntegerAdd(u, u, t);             /* u += t [u = 2^(qsize-1)] */
 244   }
 245   sophie_germain(k, t, u);
 246
 247   p = BigIntegerFromInt(0);
 248   BigIntegerMulInt(p, j, 2);            /* p = 2 * j */
 249   BigIntegerAddInt(p, p, 1);            /* p += 1 */
 250
 251   q = BigIntegerFromInt(0);
 252   BigIntegerMulInt(q, k, 2);            /* q = 2 * k */
 253   BigIntegerAddInt(q, q, 1);            /* q += 1 */
 254
 255   n = BigIntegerFromInt(0);
 256   BigIntegerMul(n, p, q);               /* n = p * q */
 257   BigIntegerMul(u, j, k);               /* u = j * k */
 258
 259   BigIntegerFree(p);
 260   BigIntegerFree(q);
 261   BigIntegerFree(j);
 262   BigIntegerFree(k);
 263
 264   g = BigIntegerFromInt(2);             /* g = 2 */
 265
 266   /* Look for a generator mod n */
 267   while(1) {
 268     BigIntegerModExp(t, g, u, n);       /* t = g^u % n */
 269     if(BigIntegerCmpInt(t, 1) == 0)
 270       BigIntegerAddInt(g, g, 1);        /* ++g */
 271     else
 272       break;
 273   }
 274
 275   BigIntegerFree(u);
 276   BigIntegerFree(t);
 277
 278   tc->tcbuf.modulus.data = tc->modbuf;
 279   tc->tcbuf.modulus.len = BigIntegerToBytes(n, tc->tcbuf.modulus.data);
 280   BigIntegerFree(n);
 281
 282   tc->tcbuf.generator.data = tc->genbuf;
 283   tc->tcbuf.generator.len = BigIntegerToBytes(g, tc->tcbuf.generator.data);
 284   BigIntegerFree(g);
 285
 286   tc->tcbuf.index = 1;
 287   return &tc->tcbuf;
 288 }
 289
 290 _TYPE( struct t_confent * )
 291 t_newconfent(tc)
 292     struct t_conf * tc;
 293 {
 294   tc->tcbuf.index = 0;
 295   tc->tcbuf.modulus.data = tc->modbuf;
 296   tc->tcbuf.modulus.len = 0;
 297   tc->tcbuf.generator.data = tc->genbuf;
 298   tc->tcbuf.generator.len = 0;
 299   return &tc->tcbuf;
 300 }
 301
 302 _TYPE( void )
 303 t_putconfent(ent, fp)
 304      const struct t_confent * ent;
 305      FILE * fp;
 306 {
 307   char strbuf[MAXB64PARAMLEN];
 308
 309   fprintf(fp, "%d:%s:", ent->index,
 310           t_tob64(strbuf, ent->modulus.data, ent->modulus.len));
 311   fprintf(fp, "%s\n",
 312           t_tob64(strbuf, ent->generator.data, ent->generator.len));
 313 }
 314
 315 int
 316 BigIntegerBitLen(b)
 317      BigInteger b;
 318 {
 319   return BN_num_bits(b);
 320 }
 321
 322 int
 323 BigIntegerCheckPrime(n)
 324      BigInteger n;
 325 {
 326   BN_CTX * ctx = BN_CTX_new();
 327   int rv = BN_is_prime(n, 25, NULL, ctx, NULL);
 328   BN_CTX_free(ctx);
 329   return rv;
 330 }
 331
 332 unsigned int
 333 BigIntegerModInt(d, m)
 334      BigInteger d;
 335      unsigned int m;
 336 {
 337   return BN_mod_word(d, m);
 338 }
 339
 340 void
 341 BigIntegerMod(result, d, m)
 342      BigInteger result, d, m;
 343 {
 344   BN_CTX * ctx = BN_CTX_new();
 345   BN_mod(result, d, m, ctx);
 346   BN_CTX_free(ctx);
 347 }
 348
 349 void
 350 BigIntegerMul(result, m1, m2)
 351      BigInteger result, m1, m2;
 352 {
 353   BN_CTX * ctx = BN_CTX_new();
 354   BN_mul(result, m1, m2, ctx);
 355   BN_CTX_free(ctx);
 356 }
 357
 358 void
 359 BigIntegerLShift(result, x, bits)
 360      BigInteger result, x;
 361      unsigned int bits;
 362 {
 363   BN_lshift(result, x, bits);
 364 }
 365
 366 int BN_is_prime(const BIGNUM *a, int checks, void (*callback)(int,int,void *),
 367         BN_CTX *ctx_passed, void *cb_arg)
 368         {
 369         return BN_is_prime_fasttest(a, checks, callback, ctx_passed, cb_arg, 0);
 370         }
 371
 372 int BN_is_prime_fasttest(const BIGNUM *a, int checks,
 373                 void (*callback)(int,int,void *),
 374                 BN_CTX *ctx_passed, void *cb_arg,
 375                 int do_trial_division)
 376         {
 377         int i, j, ret = -1;
 378         int k;
 379         BN_CTX *ctx = NULL;
 380         BIGNUM *A1, *A1_odd, *check; /* taken from ctx */
 381         BN_MONT_CTX *mont = NULL;
 382         const BIGNUM *A = NULL;
 383
 384         if (checks == BN_prime_checks)
 385                 checks = BN_prime_checks_for_size(BN_num_bits(a));
 386
 387         /* first look for small factors */
 388         if (!BN_is_odd(a))
 389                 return(0);
 390         if (do_trial_division)
 391                 {
 392                 for (i = 1; i < NUMPRIMES; i++)
 393                         if (BN_mod_word(a, primes[i]) == 0)
 394                                 return 0;
 395                 if (callback != NULL) callback(1, -1, cb_arg);
 396                 }
 397
 398         if (ctx_passed != NULL)
 399                 ctx = ctx_passed;
 400         else
 401                 if ((ctx=BN_CTX_new()) == NULL)
 402                         goto err;
 403         BN_CTX_start(ctx);
 404
 405         /* A := abs(a) */
 406         if (a->neg)
 407                 {
 408                 BIGNUM *t;
 409                 if ((t = BN_CTX_get(ctx)) == NULL) goto err;
 410                 BN_copy(t, a);
 411                 t->neg = 0;
 412                 A = t;
 413                 }
 414         else
 415                 A = a;
 416         A1 = BN_CTX_get(ctx);
 417         A1_odd = BN_CTX_get(ctx);
 418         check = BN_CTX_get(ctx);
 419         if (check == NULL) goto err;
 420
 421         /* compute A1 := A - 1 */
 422         if (!BN_copy(A1, A))
 423                 goto err;
 424         if (!BN_sub_word(A1, 1))
 425                 goto err;
 426         if (BN_is_zero(A1))
 427                 {
 428                 ret = 0;
 429                 goto err;
 430                 }
 431
 432         /* write  A1  as  A1_odd * 2^k */
 433         k = 1;
 434         while (!BN_is_bit_set(A1, k))
 435                 k++;
 436         if (!BN_rshift(A1_odd, A1, k))
 437                 goto err;
 438
 439         /* Montgomery setup for computations mod A */
 440         mont = BN_MONT_CTX_new();
 441         if (mont == NULL)
 442                 goto err;
 443         if (!BN_MONT_CTX_set(mont, A, ctx))
 444                 goto err;
 445
 446         for (i = 0; i < checks; i++)
 447                 {
 448                 if (!BN_pseudo_rand(check, BN_num_bits(A1), 0, 0))
 449                         goto err;
 450                 if (BN_cmp(check, A1) >= 0)
 451                         if (!BN_sub(check, check, A1))
 452                                 goto err;
 453                 if (!BN_add_word(check, 1))
 454                         goto err;
 455                 /* now 1 <= check < A */
 456
 457                 j = witness(check, A, A1, A1_odd, k, ctx, mont);
 458                 if (j == -1) goto err;
 459                 if (j)
 460                         {
 461                         ret=0;
 462                         goto err;
 463                         }
 464                 if (callback != NULL) callback(1,i,cb_arg);
 465                 }
 466         ret=1;
 467 err:
 468         if (ctx != NULL)
 469                 {
 470                 BN_CTX_end(ctx);
 471                 if (ctx_passed == NULL)
 472                         BN_CTX_free(ctx);
 473                 }
 474         if (mont != NULL)
 475                 BN_MONT_CTX_free(mont);
 476
 477         return(ret);
 478         }
 479
 480 static int witness(BIGNUM *w, const BIGNUM *a, const BIGNUM *a1,
 481         const BIGNUM *a1_odd, int k, BN_CTX *ctx, BN_MONT_CTX *mont)
 482         {
 483         if (!BN_mod_exp_mont(w, w, a1_odd, a, ctx, mont)) /* w := w^a1_odd mod a */
 484                 return -1;
 485         if (BN_is_one(w))
 486                 return 0; /* probably prime */
 487         if (BN_cmp(w, a1) == 0)
 488                 return 0; /* w == -1 (mod a),  'a' is probably prime */
 489         while (--k)
 490                 {
 491                 if (!BN_mod_mul(w, w, w, a, ctx)) /* w := w^2 mod a */
 492                         return -1;
 493                 if (BN_is_one(w))
 494                         return 1; /* 'a' is composite, otherwise a previous 'w' would
 495                                    * have been == -1 (mod 'a') */
 496                 if (BN_cmp(w, a1) == 0)
 497                         return 0; /* w == -1 (mod a), 'a' is probably prime */
 498                 }
 499         /* If we get here, 'w' is the (a-1)/2-th power of the original 'w',
 500          * and it is neither -1 nor +1 -- so 'a' cannot be prime */
 501         return 1;
 502         }
 503
 504 int BN_mod_exp_mont(BIGNUM *rr, BIGNUM *a, const BIGNUM *p,
 505                     const BIGNUM *m, BN_CTX *ctx, BN_MONT_CTX *in_mont)
 506         {
 507         int i,j,bits,ret=0,wstart,wend,window,wvalue;
 508         int start=1,ts=0;
 509         BIGNUM *d,*r;
 510         BIGNUM *aa;
 511         BIGNUM val[TABLE_SIZE];
 512         BN_MONT_CTX *mont=NULL;
 513
 514         bn_check_top(a);
 515         bn_check_top(p);
 516         bn_check_top(m);
 517
 518         if (!(m->d[0] & 1))
 519                 {
 520                 return(0);
 521                 }
 522         bits=BN_num_bits(p);
 523         if (bits == 0)
 524                 {
 525                 BN_one(rr);
 526                 return(1);
 527                 }
 528         BN_CTX_start(ctx);
 529         d = BN_CTX_get(ctx);
 530         r = BN_CTX_get(ctx);
 531         if (d == NULL || r == NULL) goto err;
 532
 533         /* If this is not done, things will break in the montgomery
 534          * part */
 535
 536         if (in_mont != NULL)
 537                 mont=in_mont;
 538         else
 539                 {
 540                 if ((mont=BN_MONT_CTX_new()) == NULL) goto err;
 541                 if (!BN_MONT_CTX_set(mont,m,ctx)) goto err;
 542                 }
 543
 544         BN_init(&val[0]);
 545         ts=1;
 546         if (BN_ucmp(a,m) >= 0)
 547                 {
 548                 if (!BN_mod(&(val[0]),a,m,ctx))
 549                         goto err;
 550                 aa= &(val[0]);
 551                 }
 552         else
 553                 aa=a;
 554         if (!BN_to_montgomery(&(val[0]),aa,mont,ctx)) goto err; /* 1 */
 555
 556         window = BN_window_bits_for_exponent_size(bits);
 557         if (window > 1)
 558                 {
 559                 if (!BN_mod_mul_montgomery(d,&(val[0]),&(val[0]),mont,ctx)) goto err; /* 2 */
 560                 j=1<<(window-1);
 561                 for (i=1; i<j; i++)
 562                         {
 563                         BN_init(&(val[i]));
 564                         if (!BN_mod_mul_montgomery(&(val[i]),&(val[i-1]),d,mont,ctx))
 565                                 goto err;
 566                         }
 567                 ts=i;
 568                 }
 569
 570         start=1;        /* This is used to avoid multiplication etc
 571                          * when there is only the value '1' in the
 572                          * buffer. */
 573         wvalue=0;       /* The 'value' of the window */
 574         wstart=bits-1;  /* The top bit of the window */
 575         wend=0;         /* The bottom bit of the window */
 576
 577         if (!BN_to_montgomery(r,BN_value_one(),mont,ctx)) goto err;
 578         for (;;)
 579                 {
 580                 if (BN_is_bit_set(p,wstart) == 0)
 581                         {
 582                         if (!start)
 583                                 {
 584                                 if (!BN_mod_mul_montgomery(r,r,r,mont,ctx))
 585                                 goto err;
 586                                 }
 587                         if (wstart == 0) break;
 588                         wstart--;
 589                         continue;
 590                         }
 591                 /* We now have wstart on a 'set' bit, we now need to work out
 592                  * how bit a window to do.  To do this we need to scan
 593                  * forward until the last set bit before the end of the
 594                  * window */
 595                 j=wstart;
 596                 wvalue=1;
 597                 wend=0;
 598                 for (i=1; i<window; i++)
 599                         {
 600                         if (wstart-i < 0) break;
 601                         if (BN_is_bit_set(p,wstart-i))
 602                                 {
 603                                 wvalue<<=(i-wend);
 604                                 wvalue|=1;
 605                                 wend=i;
 606                                 }
 607                         }
 608
 609                 /* wend is the size of the current window */
 610                 j=wend+1;
 611                 /* add the 'bytes above' */
 612                 if (!start)
 613                         for (i=0; i<j; i++)
 614                                 {
 615                                 if (!BN_mod_mul_montgomery(r,r,r,mont,ctx))
 616                                         goto err;
 617                                 }
 618
 619                 /* wvalue will be an odd number < 2^window */
 620                 if (!BN_mod_mul_montgomery(r,r,&(val[wvalue>>1]),mont,ctx))
 621                         goto err;
 622
 623                 /* move the 'window' down further */
 624                 wstart-=wend+1;
 625                 wvalue=0;
 626                 start=0;
 627                 if (wstart < 0) break;
 628                 }
 629         if (!BN_from_montgomery(rr,r,mont,ctx)) goto err;
 630         ret=1;
 631 err:
 632         if ((in_mont == NULL) && (mont != NULL)) BN_MONT_CTX_free(mont);
 633         BN_CTX_end(ctx);
 634         for (i=0; i<ts; i++)
 635                 BN_clear_free(&(val[i]));
 636         return(ret);
 637         }
 638
 639 BN_ULONG BN_mod_word(const BIGNUM *a, BN_ULONG w)
 640         {
 641 #ifndef BN_LLONG
 642         BN_ULONG ret=0;
 643 #else
 644         BN_ULLONG ret=0;
 645 #endif
 646         int i;
 647
 648         w&=BN_MASK2;
 649         for (i=a->top-1; i>=0; i--)
 650                 {
 651 #ifndef BN_LLONG
 652                 ret=((ret<<BN_BITS4)|((a->d[i]>>BN_BITS4)&BN_MASK2l))%w;
 653                 ret=((ret<<BN_BITS4)|(a->d[i]&BN_MASK2l))%w;
 654 #else
 655                 ret=(BN_ULLONG)(((ret<<(BN_ULLONG)BN_BITS2)|a->d[i])%
 656                         (BN_ULLONG)w);
 657 #endif
 658                 }
 659         return((BN_ULONG)ret);
 660         }
 661
 662 static int bnrand(int pseudorand, BIGNUM *rnd, int bits, int top, int bottom)
 663         {
 664         unsigned char *buf=NULL;
 665         int ret=0,bit,bytes,mask;
 666
 667         if (bits == 0)
 668                 {
 669                 BN_zero(rnd);
 670                 return 1;
 671                 }
 672
 673         bytes=(bits+7)/8;
 674         bit=(bits-1)%8;
 675         mask=0xff<<bit;
 676
 677         buf=(unsigned char *)malloc(bytes);
 678         if (buf == NULL)
 679                 {
 680                 goto err;
 681                 }
 682
 683         /* make a random number and set the top and bottom bits */
 684         /* this ignores the pseudorand flag */
 685
 686         t_random(buf, bytes);
 687
 688         if (top)
 689                 {
 690                 if (bit == 0)
 691                         {
 692                         buf[0]=1;
 693                         buf[1]|=0x80;
 694                         }
 695                 else
 696                         {
 697                         buf[0]|=(3<<(bit-1));
 698                         buf[0]&= ~(mask<<1);
 699                         }
 700                 }
 701         else
 702                 {
 703                 buf[0]|=(1<<bit);
 704                 buf[0]&= ~(mask<<1);
 705                 }
 706         if (bottom) /* set bottom bits to whatever odd is */
 707                 buf[bytes-1]|=1;
 708         if (!BN_bin2bn(buf,bytes,rnd)) goto err;
 709         ret=1;
 710 err:
 711         if (buf != NULL)
 712                 {
 713                 memset(buf,0,bytes);
 714                 free(buf);
 715                 }
 716         return(ret);
 717         }
 718
 719 /* BN_pseudo_rand is the same as BN_rand, now. */
 720
 721 int     BN_pseudo_rand(BIGNUM *rnd, int bits, int top, int bottom)
 722         {
 723         return bnrand(1, rnd, bits, top, bottom);
 724         }
 725
 726 #define MONT_WORD /* use the faster word-based algorithm */
 727
 728 int BN_mod_mul_montgomery(BIGNUM *r, BIGNUM *a, BIGNUM *b,
 729                           BN_MONT_CTX *mont, BN_CTX *ctx)
 730         {
 731         BIGNUM *tmp,*tmp2;
 732         int ret=0;
 733
 734         BN_CTX_start(ctx);
 735         tmp = BN_CTX_get(ctx);
 736         tmp2 = BN_CTX_get(ctx);
 737         if (tmp == NULL || tmp2 == NULL) goto err;
 738
 739         bn_check_top(tmp);
 740         bn_check_top(tmp2);
 741
 742         if (a == b)
 743                 {
 744                 if (!BN_sqr(tmp,a,ctx)) goto err;
 745                 }
 746         else
 747                 {
 748                 if (!BN_mul(tmp,a,b,ctx)) goto err;
 749                 }
 750         /* reduce from aRR to aR */
 751         if (!BN_from_montgomery(r,tmp,mont,ctx)) goto err;
 752         ret=1;
 753 err:
 754         BN_CTX_end(ctx);
 755         return(ret);
 756         }
 757
 758 int BN_from_montgomery(BIGNUM *ret, BIGNUM *a, BN_MONT_CTX *mont,
 759              BN_CTX *ctx)
 760         {
 761         int retn=0;
 762
 763 #ifdef MONT_WORD
 764         BIGNUM *n,*r;
 765         BN_ULONG *ap,*np,*rp,n0,v,*nrp;
 766         int al,nl,max,i,x,ri;
 767
 768         BN_CTX_start(ctx);
 769         if ((r = BN_CTX_get(ctx)) == NULL) goto err;
 770
 771         if (!BN_copy(r,a)) goto err;
 772         n= &(mont->N);
 773
 774         ap=a->d;
 775         /* mont->ri is the size of mont->N in bits (rounded up
 776            to the word size) */
 777         al=ri=mont->ri/BN_BITS2;
 778
 779         nl=n->top;
 780         if ((al == 0) || (nl == 0)) { r->top=0; return(1); }
 781
 782         max=(nl+al+1); /* allow for overflow (no?) XXX */
 783         if (bn_wexpand(r,max) == NULL) goto err;
 784         if (bn_wexpand(ret,max) == NULL) goto err;
 785
 786         r->neg=a->neg^n->neg;
 787         np=n->d;
 788         rp=r->d;
 789         nrp= &(r->d[nl]);
 790
 791         /* clear the top words of T */
 792 #if 1
 793         for (i=r->top; i<max; i++) /* memset? XXX */
 794                 r->d[i]=0;
 795 #else
 796         memset(&(r->d[r->top]),0,(max-r->top)*sizeof(BN_ULONG));
 797 #endif
 798
 799         r->top=max;
 800         n0=mont->n0;
 801
 802 #ifdef BN_COUNT
 803         printf("word BN_from_montgomery %d * %d\n",nl,nl);
 804 #endif
 805         for (i=0; i<nl; i++)
 806                 {
 807 #ifdef __TANDEM
 808                 {
 809                    long long t1;
 810                    long long t2;
 811                    long long t3;
 812                    t1 = rp[0] * (n0 & 0177777);
 813                    t2 = 037777600000l;
 814                    t2 = n0 & t2;
 815                    t3 = rp[0] & 0177777;
 816                    t2 = (t3 * t2) & BN_MASK2;
 817                    t1 = t1 + t2;
 818                    v=bn_mul_add_words(rp,np,nl,(BN_ULONG) t1);
 819                 }
 820 #else
 821                 v=bn_mul_add_words(rp,np,nl,(rp[0]*n0)&BN_MASK2);
 822 #endif
 823                 nrp++;
 824                 rp++;
 825                 if (((nrp[-1]+=v)&BN_MASK2) >= v)
 826                         continue;
 827                 else
 828                         {
 829                         if (((++nrp[0])&BN_MASK2) != 0) continue;
 830                         if (((++nrp[1])&BN_MASK2) != 0) continue;
 831                         for (x=2; (((++nrp[x])&BN_MASK2) == 0); x++) ;
 832                         }
 833                 }
 834         bn_fix_top(r);
 835
 836         /* mont->ri will be a multiple of the word size */
 837 #if 0
 838         BN_rshift(ret,r,mont->ri);
 839 #else
 840         ret->neg = r->neg;
 841         x=ri;
 842         rp=ret->d;
 843         ap= &(r->d[x]);
 844         if (r->top < x)
 845                 al=0;
 846         else
 847                 al=r->top-x;
 848         ret->top=al;
 849         al-=4;
 850         for (i=0; i<al; i+=4)
 851                 {
 852                 BN_ULONG t1,t2,t3,t4;
 853
 854                 t1=ap[i+0];
 855                 t2=ap[i+1];
 856                 t3=ap[i+2];
 857                 t4=ap[i+3];
 858                 rp[i+0]=t1;
 859                 rp[i+1]=t2;
 860                 rp[i+2]=t3;
 861                 rp[i+3]=t4;
 862                 }
 863         al+=4;
 864         for (; i<al; i++)
 865                 rp[i]=ap[i];
 866 #endif
 867 #else /* !MONT_WORD */
 868         BIGNUM *t1,*t2;
 869
 870         BN_CTX_start(ctx);
 871         t1 = BN_CTX_get(ctx);
 872         t2 = BN_CTX_get(ctx);
 873         if (t1 == NULL || t2 == NULL) goto err;
 874
 875         if (!BN_copy(t1,a)) goto err;
 876         BN_mask_bits(t1,mont->ri);
 877
 878         if (!BN_mul(t2,t1,&mont->Ni,ctx)) goto err;
 879         BN_mask_bits(t2,mont->ri);
 880
 881         if (!BN_mul(t1,t2,&mont->N,ctx)) goto err;
 882         if (!BN_add(t2,a,t1)) goto err;
 883         BN_rshift(ret,t2,mont->ri);
 884 #endif /* MONT_WORD */
 885
 886         if (BN_ucmp(ret, &(mont->N)) >= 0)
 887                 {
 888                 BN_usub(ret,ret,&(mont->N));
 889                 }
 890         retn=1;
 891  err:
 892         BN_CTX_end(ctx);
 893         return(retn);
 894         }
 895
 896 void BN_MONT_CTX_init(BN_MONT_CTX *ctx)
 897         {
 898         ctx->ri=0;
 899         BN_init(&(ctx->RR));
 900         BN_init(&(ctx->N));
 901         BN_init(&(ctx->Ni));
 902         ctx->flags=0;
 903         }
 904
 905 BN_MONT_CTX *BN_MONT_CTX_new(void)
 906         {
 907         BN_MONT_CTX *ret;
 908
 909         if ((ret=(BN_MONT_CTX *)malloc(sizeof(BN_MONT_CTX))) == NULL)
 910                 return(NULL);
 911
 912         BN_MONT_CTX_init(ret);
 913         ret->flags=BN_FLG_MALLOCED;
 914         return(ret);
 915         }
 916
 917 void BN_MONT_CTX_free(BN_MONT_CTX *mont)
 918         {
 919         if(mont == NULL)
 920             return;
 921
 922         BN_free(&(mont->RR));
 923         BN_free(&(mont->N));
 924         BN_free(&(mont->Ni));
 925         if (mont->flags & BN_FLG_MALLOCED)
 926                 free(mont);
 927         }
 928
 929 int BN_MONT_CTX_set(BN_MONT_CTX *mont, const BIGNUM *mod, BN_CTX *ctx)
 930         {
 931         BIGNUM Ri,*R;
 932
 933         BN_init(&Ri);
 934         R= &(mont->RR);                                 /* grab RR as a temp */
 935         BN_copy(&(mont->N),mod);                        /* Set N */
 936
 937 #ifdef MONT_WORD
 938                 {
 939                 BIGNUM tmod;
 940                 BN_ULONG buf[2];
 941
 942                 mont->ri=(BN_num_bits(mod)+(BN_BITS2-1))/BN_BITS2*BN_BITS2;
 943                 BN_zero(R);
 944                 BN_set_bit(R,BN_BITS2);                 /* R */
 945
 946                 buf[0]=mod->d[0]; /* tmod = N mod word size */
 947                 buf[1]=0;
 948                 tmod.d=buf;
 949                 tmod.top=1;
 950                 tmod.dmax=2;
 951                 tmod.neg=mod->neg;
 952                                                         /* Ri = R^-1 mod N*/
 953                 if ((BN_mod_inverse(&Ri,R,&tmod,ctx)) == NULL)
 954                         goto err;
 955                 BN_lshift(&Ri,&Ri,BN_BITS2);            /* R*Ri */
 956                 if (!BN_is_zero(&Ri))
 957                         BN_sub_word(&Ri,1);
 958                 else /* if N mod word size == 1 */
 959                         BN_set_word(&Ri,BN_MASK2);  /* Ri-- (mod word size) */
 960                 BN_div(&Ri,NULL,&Ri,&tmod,ctx); /* Ni = (R*Ri-1)/N,
 961                                                  * keep only least significant word: */
 962                 mont->n0=Ri.d[0];
 963                 BN_free(&Ri);
 964                 }
 965 #else /* !MONT_WORD */
 966                 { /* bignum version */
 967                 mont->ri=BN_num_bits(mod);
 968                 BN_zero(R);
 969                 BN_set_bit(R,mont->ri);                 /* R = 2^ri */
 970                                                         /* Ri = R^-1 mod N*/
 971                 if ((BN_mod_inverse(&Ri,R,mod,ctx)) == NULL)
 972                         goto err;
 973                 BN_lshift(&Ri,&Ri,mont->ri);            /* R*Ri */
 974                 BN_sub_word(&Ri,1);
 975                                                         /* Ni = (R*Ri-1) / N */
 976                 BN_div(&(mont->Ni),NULL,&Ri,mod,ctx);
 977                 BN_free(&Ri);
 978                 }
 979 #endif
 980
 981         /* setup RR for conversions */
 982         BN_zero(&(mont->RR));
 983         BN_set_bit(&(mont->RR),mont->ri*2);
 984         BN_mod(&(mont->RR),&(mont->RR),&(mont->N),ctx);
 985
 986         return(1);
 987 err:
 988         return(0);
 989         }
 990
 991 BIGNUM *BN_value_one(void)
 992         {
 993         static BN_ULONG data_one=1L;
 994         static BIGNUM const_one={&data_one,1,1,0};
 995
 996         return(&const_one);
 997         }
 998
 999 /* solves ax == 1 (mod n) */
1000 BIGNUM *BN_mod_inverse(BIGNUM *in, BIGNUM *a, const BIGNUM *n, BN_CTX *ctx)
1001         {
1002         BIGNUM *A,*B,*X,*Y,*M,*D,*R=NULL;
1003         BIGNUM *T,*ret=NULL;
1004         int sign;
1005
1006         bn_check_top(a);
1007         bn_check_top(n);
1008
1009         BN_CTX_start(ctx);
1010         A = BN_CTX_get(ctx);
1011         B = BN_CTX_get(ctx);
1012         X = BN_CTX_get(ctx);
1013         D = BN_CTX_get(ctx);
1014         M = BN_CTX_get(ctx);
1015         Y = BN_CTX_get(ctx);
1016         if (Y == NULL) goto err;
1017
1018         if (in == NULL)
1019                 R=BN_new();
1020         else
1021                 R=in;
1022         if (R == NULL) goto err;
1023
1024         BN_zero(X);
1025         BN_one(Y);
1026         if (BN_copy(A,a) == NULL) goto err;
1027         if (BN_copy(B,n) == NULL) goto err;
1028         sign=1;
1029
1030         while (!BN_is_zero(B))
1031                 {
1032                 if (!BN_div(D,M,A,B,ctx)) goto err;
1033                 T=A;
1034                 A=B;
1035                 B=M;
1036                 /* T has a struct, M does not */
1037
1038                 if (!BN_mul(T,D,X,ctx)) goto err;
1039                 if (!BN_add(T,T,Y)) goto err;
1040                 M=Y;
1041                 Y=X;
1042                 X=T;
1043                 sign= -sign;
1044                 }
1045         if (sign < 0)
1046                 {
1047                 if (!BN_sub(Y,n,Y)) goto err;
1048                 }
1049
1050         if (BN_is_one(A))
1051                 { if (!BN_mod(R,Y,n,ctx)) goto err; }
1052         else
1053                 {
1054                 goto err;
1055                 }
1056         ret=R;
1057 err:
1058         if ((ret == NULL) && (in == NULL)) BN_free(R);
1059         BN_CTX_end(ctx);
1060         return(ret);
1061         }
1062
1063 int BN_set_bit(BIGNUM *a, int n)
1064         {
1065         int i,j,k;
1066
1067         i=n/BN_BITS2;
1068         j=n%BN_BITS2;
1069         if (a->top <= i)
1070                 {
1071                 if (bn_wexpand(a,i+1) == NULL) return(0);
1072                 for(k=a->top; k<i+1; k++)
1073                         a->d[k]=0;
1074                 a->top=i+1;
1075                 }
1076
1077         a->d[i]|=(((BN_ULONG)1)<<j);
1078         return(1);
1079         }
1080