← Back to branch summary

~ubuntu-branches/ubuntu/utopic/dropbear/utopic-proposed

« back to all changes in this revision

Viewing changes to libtommath/pre_gen/mpi.c

  • Committer: Bazaar Package Importer
  • Author(s): Matt Johnston
  • Date: 2005-12-08 19:20:21 UTC
  • mfrom: (1.2.2 upstream)
  • Revision ID: james.westby@ubuntu.com-20051208192021-nyp9rwnt77nsg6ty
Tags: 0.47-1
* New upstream release.
* SECURITY: Fix incorrect buffer sizing.

Show diffs side-by-side

added added

removed removed

Lines of Context:
libtommath/pre_gen/mpi.c
69
69
 * Based on slow invmod except this is optimized for the case where b is 
70
70
 * odd as per HAC Note 14.64 on pp. 610
71
71
 */
72
 
int
73
 
fast_mp_invmod (mp_int * a, mp_int * b, mp_int * c)
 
72
int fast_mp_invmod (mp_int * a, mp_int * b, mp_int * c)
74
73
{
75
74
  mp_int  x, y, u, v, B, D;
76
75
  int     res, neg;
87
86
 
88
87
  /* x == modulus, y == value to invert */
89
88
  if ((res = mp_copy (b, &x)) != MP_OKAY) {
90
 
    goto __ERR;
 
89
    goto LBL_ERR;
91
90
  }
92
91
 
93
92
  /* we need y = |a| */
94
 
  if ((res = mp_abs (a, &y)) != MP_OKAY) {
95
 
    goto __ERR;
 
93
  if ((res = mp_mod (a, b, &y)) != MP_OKAY) {
 
94
    goto LBL_ERR;
96
95
  }
97
96
 
98
97
  /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */
99
98
  if ((res = mp_copy (&x, &u)) != MP_OKAY) {
100
 
    goto __ERR;
 
99
    goto LBL_ERR;
101
100
  }
102
101
  if ((res = mp_copy (&y, &v)) != MP_OKAY) {
103
 
    goto __ERR;
 
102
    goto LBL_ERR;
104
103
  }
105
104
  mp_set (&D, 1);
106
105
 
109
108
  while (mp_iseven (&u) == 1) {
110
109
    /* 4.1 u = u/2 */
111
110
    if ((res = mp_div_2 (&u, &u)) != MP_OKAY) {
112
 
      goto __ERR;
 
111
      goto LBL_ERR;
113
112
    }
114
113
    /* 4.2 if B is odd then */
115
114
    if (mp_isodd (&B) == 1) {
116
115
      if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) {
117
 
        goto __ERR;
 
116
        goto LBL_ERR;
118
117
      }
119
118
    }
120
119
    /* B = B/2 */
121
120
    if ((res = mp_div_2 (&B, &B)) != MP_OKAY) {
122
 
      goto __ERR;
 
121
      goto LBL_ERR;
123
122
    }
124
123
  }
125
124
 
127
126
  while (mp_iseven (&v) == 1) {
128
127
    /* 5.1 v = v/2 */
129
128
    if ((res = mp_div_2 (&v, &v)) != MP_OKAY) {
130
 
      goto __ERR;
 
129
      goto LBL_ERR;
131
130
    }
132
131
    /* 5.2 if D is odd then */
133
132
    if (mp_isodd (&D) == 1) {
134
133
      /* D = (D-x)/2 */
135
134
      if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) {
136
 
        goto __ERR;
 
135
        goto LBL_ERR;
137
136
      }
138
137
    }
139
138
    /* D = D/2 */
140
139
    if ((res = mp_div_2 (&D, &D)) != MP_OKAY) {
141
 
      goto __ERR;
 
140
      goto LBL_ERR;
142
141
    }
143
142
  }
144
143
 
146
145
  if (mp_cmp (&u, &v) != MP_LT) {
147
146
    /* u = u - v, B = B - D */
148
147
    if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) {
149
 
      goto __ERR;
 
148
      goto LBL_ERR;
150
149
    }
151
150
 
152
151
    if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) {
153
 
      goto __ERR;
 
152
      goto LBL_ERR;
154
153
    }
155
154
  } else {
156
155
    /* v - v - u, D = D - B */
157
156
    if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) {
158
 
      goto __ERR;
 
157
      goto LBL_ERR;
159
158
    }
160
159
 
161
160
    if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) {
162
 
      goto __ERR;
 
161
      goto LBL_ERR;
163
162
    }
164
163
  }
165
164
 
173
172
  /* if v != 1 then there is no inverse */
174
173
  if (mp_cmp_d (&v, 1) != MP_EQ) {
175
174
    res = MP_VAL;
176
 
    goto __ERR;
 
175
    goto LBL_ERR;
177
176
  }
178
177
 
179
178
  /* b is now the inverse */
180
179
  neg = a->sign;
181
180
  while (D.sign == MP_NEG) {
182
181
    if ((res = mp_add (&D, b, &D)) != MP_OKAY) {
183
 
      goto __ERR;
 
182
      goto LBL_ERR;
184
183
    }
185
184
  }
186
185
  mp_exch (&D, c);
187
186
  c->sign = neg;
188
187
  res = MP_OKAY;
189
188
 
190
 
__ERR:mp_clear_multi (&x, &y, &u, &v, &B, &D, NULL);
 
189
LBL_ERR:mp_clear_multi (&x, &y, &u, &v, &B, &D, NULL);
191
190
  return res;
192
191
}
193
192
#endif
220
219
 *
221
220
 * Based on Algorithm 14.32 on pp.601 of HAC.
222
221
*/
223
 
int
224
 
fast_mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho)
 
222
int fast_mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho)
225
223
{
226
224
  int     ix, res, olduse;
227
225
  mp_word W[MP_WARRAY];
401
399
 * Based on Algorithm 14.12 on pp.595 of HAC.
402
400
 *
403
401
 */
404
 
int
405
 
fast_s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
 
402
int fast_s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
406
403
{
407
404
  int     olduse, res, pa, ix, iz;
408
405
  mp_digit W[MP_WARRAY];
420
417
 
421
418
  /* clear the carry */
422
419
  _W = 0;
423
 
  for (ix = 0; ix <= pa; ix++) { 
 
420
  for (ix = 0; ix < pa; ix++) { 
424
421
      int      tx, ty;
425
422
      int      iy;
426
423
      mp_digit *tmpx, *tmpy;
433
430
      tmpx = a->dp + tx;
434
431
      tmpy = b->dp + ty;
435
432
 
436
 
      /* this is the number of times the loop will iterrate, essentially its 
 
433
      /* this is the number of times the loop will iterrate, essentially 
437
434
         while (tx++ < a->used && ty-- >= 0) { ... }
438
435
       */
439
436
      iy = MIN(a->used-tx, ty+1);
450
447
      _W = _W >> ((mp_word)DIGIT_BIT);
451
448
  }
452
449
 
 
450
  /* store final carry */
 
451
  W[ix] = (mp_digit)(_W & MP_MASK);
 
452
 
453
453
  /* setup dest */
454
454
  olduse  = c->used;
455
 
  c->used = digs;
 
455
  c->used = pa;
456
456
 
457
457
  {
458
458
    register mp_digit *tmpc;
459
459
    tmpc = c->dp;
460
 
    for (ix = 0; ix < digs; ix++) {
 
460
    for (ix = 0; ix < pa+1; ix++) {
461
461
      /* now extract the previous digit [below the carry] */
462
462
      *tmpc++ = W[ix];
463
463
    }
501
501
 *
502
502
 * Based on Algorithm 14.12 on pp.595 of HAC.
503
503
 */
504
 
int
505
 
fast_s_mp_mul_high_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
 
504
int fast_s_mp_mul_high_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
506
505
{
507
506
  int     olduse, res, pa, ix, iz;
508
507
  mp_digit W[MP_WARRAY];
519
518
  /* number of output digits to produce */
520
519
  pa = a->used + b->used;
521
520
  _W = 0;
522
 
  for (ix = digs; ix <= pa; ix++) { 
 
521
  for (ix = digs; ix < pa; ix++) { 
523
522
      int      tx, ty, iy;
524
523
      mp_digit *tmpx, *tmpy;
525
524
 
547
546
      /* make next carry */
548
547
      _W = _W >> ((mp_word)DIGIT_BIT);
549
548
  }
 
549
  
 
550
  /* store final carry */
 
551
  W[ix] = (mp_digit)(_W & MP_MASK);
550
552
 
551
553
  /* setup dest */
552
554
  olduse  = c->used;
591
593
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
592
594
 */
593
595
 
594
 
/* fast squaring
595
 
 *
596
 
 * This is the comba method where the columns of the product
597
 
 * are computed first then the carries are computed.  This
598
 
 * has the effect of making a very simple inner loop that
599
 
 * is executed the most
600
 
 *
601
 
 * W2 represents the outer products and W the inner.
602
 
 *
603
 
 * A further optimizations is made because the inner
604
 
 * products are of the form "A * B * 2".  The *2 part does
605
 
 * not need to be computed until the end which is good
606
 
 * because 64-bit shifts are slow!
607
 
 *
608
 
 * Based on Algorithm 14.16 on pp.597 of HAC.
609
 
 *
610
 
 */
611
596
/* the jist of squaring...
612
 
 
613
 
you do like mult except the offset of the tmpx [one that starts closer to zero]
614
 
can't equal the offset of tmpy.  So basically you set up iy like before then you min it with
615
 
(ty-tx) so that it never happens.  You double all those you add in the inner loop
 
597
 * you do like mult except the offset of the tmpx [one that 
 
598
 * starts closer to zero] can't equal the offset of tmpy.  
 
599
 * So basically you set up iy like before then you min it with
 
600
 * (ty-tx) so that it never happens.  You double all those 
 
601
 * you add in the inner loop
616
602
 
617
603
After that loop you do the squares and add them in.
618
 
 
619
 
Remove W2 and don't memset W
620
 
 
621
604
*/
622
605
 
623
606
int fast_s_mp_sqr (mp_int * a, mp_int * b)
636
619
 
637
620
  /* number of output digits to produce */
638
621
  W1 = 0;
639
 
  for (ix = 0; ix <= pa; ix++) { 
 
622
  for (ix = 0; ix < pa; ix++) { 
640
623
      int      tx, ty, iy;
641
624
      mp_word  _W;
642
625
      mp_digit *tmpy;
652
635
      tmpx = a->dp + tx;
653
636
      tmpy = a->dp + ty;
654
637
 
655
 
      /* this is the number of times the loop will iterrate, essentially its 
 
638
      /* this is the number of times the loop will iterrate, essentially
656
639
         while (tx++ < a->used && ty-- >= 0) { ... }
657
640
       */
658
641
      iy = MIN(a->used-tx, ty+1);
677
660
      }
678
661
 
679
662
      /* store it */
680
 
      W[ix] = _W;
 
663
      W[ix] = (mp_digit)(_W & MP_MASK);
681
664
 
682
665
      /* make next carry */
683
666
      W1 = _W >> ((mp_word)DIGIT_BIT);
1539
1522
 
1540
1523
  mp_set(&tq, 1);
1541
1524
  n = mp_count_bits(a) - mp_count_bits(b);
1542
 
  if (((res = mp_copy(a, &ta)) != MP_OKAY) ||
1543
 
      ((res = mp_copy(b, &tb)) != MP_OKAY) || 
 
1525
  if (((res = mp_abs(a, &ta)) != MP_OKAY) ||
 
1526
      ((res = mp_abs(b, &tb)) != MP_OKAY) || 
1544
1527
      ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) ||
1545
1528
      ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) {
1546
 
      goto __ERR;
 
1529
      goto LBL_ERR;
1547
1530
  }
1548
1531
 
1549
1532
  while (n-- >= 0) {
1550
1533
     if (mp_cmp(&tb, &ta) != MP_GT) {
1551
1534
        if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) ||
1552
1535
            ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) {
1553
 
           goto __ERR;
 
1536
           goto LBL_ERR;
1554
1537
        }
1555
1538
     }
1556
1539
     if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) ||
1557
1540
         ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) {
1558
 
           goto __ERR;
 
1541
           goto LBL_ERR;
1559
1542
     }
1560
1543
  }
1561
1544
 
1564
1547
  n2 = (a->sign == b->sign ? MP_ZPOS : MP_NEG);
1565
1548
  if (c != NULL) {
1566
1549
     mp_exch(c, &q);
1567
 
     c->sign  = n2;
 
1550
     c->sign  = (mp_iszero(c) == MP_YES) ? MP_ZPOS : n2;
1568
1551
  }
1569
1552
  if (d != NULL) {
1570
1553
     mp_exch(d, &ta);
1571
 
     d->sign = n;
 
1554
     d->sign = (mp_iszero(d) == MP_YES) ? MP_ZPOS : n;
1572
1555
  }
1573
 
__ERR:
 
1556
LBL_ERR:
1574
1557
   mp_clear_multi(&ta, &tb, &tq, &q, NULL);
1575
1558
   return res;
1576
1559
}
1619
1602
  q.used = a->used + 2;
1620
1603
 
1621
1604
  if ((res = mp_init (&t1)) != MP_OKAY) {
1622
 
    goto __Q;
 
1605
    goto LBL_Q;
1623
1606
  }
1624
1607
 
1625
1608
  if ((res = mp_init (&t2)) != MP_OKAY) {
1626
 
    goto __T1;
 
1609
    goto LBL_T1;
1627
1610
  }
1628
1611
 
1629
1612
  if ((res = mp_init_copy (&x, a)) != MP_OKAY) {
1630
 
    goto __T2;
 
1613
    goto LBL_T2;
1631
1614
  }
1632
1615
 
1633
1616
  if ((res = mp_init_copy (&y, b)) != MP_OKAY) {
1634
 
    goto __X;
 
1617
    goto LBL_X;
1635
1618
  }
1636
1619
 
1637
1620
  /* fix the sign */
1643
1626
  if (norm < (int)(DIGIT_BIT-1)) {
1644
1627
     norm = (DIGIT_BIT-1) - norm;
1645
1628
     if ((res = mp_mul_2d (&x, norm, &x)) != MP_OKAY) {
1646
 
       goto __Y;
 
1629
       goto LBL_Y;
1647
1630
     }
1648
1631
     if ((res = mp_mul_2d (&y, norm, &y)) != MP_OKAY) {
1649
 
       goto __Y;
 
1632
       goto LBL_Y;
1650
1633
     }
1651
1634
  } else {
1652
1635
     norm = 0;
1658
1641
 
1659
1642
  /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */
1660
1643
  if ((res = mp_lshd (&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */
1661
 
    goto __Y;
 
1644
    goto LBL_Y;
1662
1645
  }
1663
1646
 
1664
1647
  while (mp_cmp (&x, &y) != MP_LT) {
1665
1648
    ++(q.dp[n - t]);
1666
1649
    if ((res = mp_sub (&x, &y, &x)) != MP_OKAY) {
1667
 
      goto __Y;
 
1650
      goto LBL_Y;
1668
1651
    }
1669
1652
  }
1670
1653
 
1706
1689
      t1.dp[1] = y.dp[t];
1707
1690
      t1.used = 2;
1708
1691
      if ((res = mp_mul_d (&t1, q.dp[i - t - 1], &t1)) != MP_OKAY) {
1709
 
        goto __Y;
 
1692
        goto LBL_Y;
1710
1693
      }
1711
1694
 
1712
1695
      /* find right hand */
1718
1701
 
1719
1702
    /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */
1720
1703
    if ((res = mp_mul_d (&y, q.dp[i - t - 1], &t1)) != MP_OKAY) {
1721
 
      goto __Y;
 
1704
      goto LBL_Y;
1722
1705
    }
1723
1706
 
1724
1707
    if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) {
1725
 
      goto __Y;
 
1708
      goto LBL_Y;
1726
1709
    }
1727
1710
 
1728
1711
    if ((res = mp_sub (&x, &t1, &x)) != MP_OKAY) {
1729
 
      goto __Y;
 
1712
      goto LBL_Y;
1730
1713
    }
1731
1714
 
1732
1715
    /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */
1733
1716
    if (x.sign == MP_NEG) {
1734
1717
      if ((res = mp_copy (&y, &t1)) != MP_OKAY) {
1735
 
        goto __Y;
 
1718
        goto LBL_Y;
1736
1719
      }
1737
1720
      if ((res = mp_lshd (&t1, i - t - 1)) != MP_OKAY) {
1738
 
        goto __Y;
 
1721
        goto LBL_Y;
1739
1722
      }
1740
1723
      if ((res = mp_add (&x, &t1, &x)) != MP_OKAY) {
1741
 
        goto __Y;
 
1724
        goto LBL_Y;
1742
1725
      }
1743
1726
 
1744
1727
      q.dp[i - t - 1] = (q.dp[i - t - 1] - 1UL) & MP_MASK;
1765
1748
 
1766
1749
  res = MP_OKAY;
1767
1750
 
1768
 
__Y:mp_clear (&y);
1769
 
__X:mp_clear (&x);
1770
 
__T2:mp_clear (&t2);
1771
 
__T1:mp_clear (&t1);
1772
 
__Q:mp_clear (&q);
 
1751
LBL_Y:mp_clear (&y);
 
1752
LBL_X:mp_clear (&x);
 
1753
LBL_T2:mp_clear (&t2);
 
1754
LBL_T1:mp_clear (&t1);
 
1755
LBL_Q:mp_clear (&q);
1773
1756
  return res;
1774
1757
}
1775
1758
 
2199
2182
 * Based on algorithm from the paper
2200
2183
 *
2201
2184
 * "Generating Efficient Primes for Discrete Log Cryptosystems"
2202
 
 *                 Chae Hoon Lim, Pil Loong Lee,
 
2185
 *                 Chae Hoon Lim, Pil Joong Lee,
2203
2186
 *          POSTECH Information Research Laboratories
2204
2187
 *
2205
2188
 * The modulus must be of a special format [see manual]
2457
2440
     return err;
2458
2441
#else 
2459
2442
     /* no invmod */
2460
 
     return MP_VAL
2461
 
#endif
2462
 
  }
 
2443
     return MP_VAL;
 
2444
#endif
 
2445
  }
 
2446
 
 
2447
/* modified diminished radix reduction */
 
2448
#if defined(BN_MP_REDUCE_IS_2K_L_C) && defined(BN_MP_REDUCE_2K_L_C)
 
2449
  if (mp_reduce_is_2k_l(P) == MP_YES) {
 
2450
     return s_mp_exptmod(G, X, P, Y, 1);
 
2451
  }
 
2452
#endif
2463
2453
 
2464
2454
#ifdef BN_MP_DR_IS_MODULUS_C
2465
2455
  /* is it a DR modulus? */
2466
2456
  dr = mp_dr_is_modulus(P);
2467
2457
#else
 
2458
  /* default to no */
2468
2459
  dr = 0;
2469
2460
#endif
2470
2461
 
2471
2462
#ifdef BN_MP_REDUCE_IS_2K_C
2472
 
  /* if not, is it a uDR modulus? */
 
2463
  /* if not, is it a unrestricted DR modulus? */
2473
2464
  if (dr == 0) {
2474
2465
     dr = mp_reduce_is_2k(P) << 1;
2475
2466
  }
2476
2467
#endif
2477
2468
    
2478
 
  /* if the modulus is odd or dr != 0 use the fast method */
 
2469
  /* if the modulus is odd or dr != 0 use the montgomery method */
2479
2470
#ifdef BN_MP_EXPTMOD_FAST_C
2480
2471
  if (mp_isodd (P) == 1 || dr !=  0) {
2481
2472
    return mp_exptmod_fast (G, X, P, Y, dr);
2483
2474
#endif
2484
2475
#ifdef BN_S_MP_EXPTMOD_C
2485
2476
    /* otherwise use the generic Barrett reduction technique */
2486
 
    return s_mp_exptmod (G, X, P, Y);
 
2477
    return s_mp_exptmod (G, X, P, Y, 0);
2487
2478
#else
2488
2479
    /* no exptmod for evens */
2489
2480
    return MP_VAL;
2529
2520
   #define TAB_SIZE 256
2530
2521
#endif
2531
2522
 
2532
 
int
2533
 
mp_exptmod_fast (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode)
 
2523
int mp_exptmod_fast (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode)
2534
2524
{
2535
2525
  mp_int  M[TAB_SIZE], res;
2536
2526
  mp_digit buf, mp;
2588
2578
#ifdef BN_MP_MONTGOMERY_SETUP_C     
2589
2579
     /* now setup montgomery  */
2590
2580
     if ((err = mp_montgomery_setup (P, &mp)) != MP_OKAY) {
2591
 
        goto __M;
 
2581
        goto LBL_M;
2592
2582
     }
2593
2583
#else
2594
2584
     err = MP_VAL;
2595
 
     goto __M;
 
2585
     goto LBL_M;
2596
2586
#endif
2597
2587
 
2598
2588
     /* automatically pick the comba one if available (saves quite a few calls/ifs) */
2608
2598
        redux = mp_montgomery_reduce;
2609
2599
#else
2610
2600
        err = MP_VAL;
2611
 
        goto __M;
 
2601
        goto LBL_M;
2612
2602
#endif
2613
2603
     }
2614
2604
  } else if (redmode == 1) {
2618
2608
     redux = mp_dr_reduce;
2619
2609
#else
2620
2610
     err = MP_VAL;
2621
 
     goto __M;
 
2611
     goto LBL_M;
2622
2612
#endif
2623
2613
  } else {
2624
2614
#if defined(BN_MP_REDUCE_2K_SETUP_C) && defined(BN_MP_REDUCE_2K_C)
2625
2615
     /* setup DR reduction for moduli of the form 2**k - b */
2626
2616
     if ((err = mp_reduce_2k_setup(P, &mp)) != MP_OKAY) {
2627
 
        goto __M;
 
2617
        goto LBL_M;
2628
2618
     }
2629
2619
     redux = mp_reduce_2k;
2630
2620
#else
2631
2621
     err = MP_VAL;
2632
 
     goto __M;
 
2622
     goto LBL_M;
2633
2623
#endif
2634
2624
  }
2635
2625
 
2636
2626
  /* setup result */
2637
2627
  if ((err = mp_init (&res)) != MP_OKAY) {
2638
 
    goto __M;
 
2628
    goto LBL_M;
2639
2629
  }
2640
2630
 
2641
2631
  /* create M table
2649
2639
#ifdef BN_MP_MONTGOMERY_CALC_NORMALIZATION_C
2650
2640
     /* now we need R mod m */
2651
2641
     if ((err = mp_montgomery_calc_normalization (&res, P)) != MP_OKAY) {
2652
 
       goto __RES;
 
2642
       goto LBL_RES;
2653
2643
     }
2654
2644
#else 
2655
2645
     err = MP_VAL;
2656
 
     goto __RES;
 
2646
     goto LBL_RES;
2657
2647
#endif
2658
2648
 
2659
2649
     /* now set M[1] to G * R mod m */
2660
2650
     if ((err = mp_mulmod (G, &res, P, &M[1])) != MP_OKAY) {
2661
 
       goto __RES;
 
2651
       goto LBL_RES;
2662
2652
     }
2663
2653
  } else {
2664
2654
     mp_set(&res, 1);
2665
2655
     if ((err = mp_mod(G, P, &M[1])) != MP_OKAY) {
2666
 
        goto __RES;
 
2656
        goto LBL_RES;
2667
2657
     }
2668
2658
  }
2669
2659
 
2670
2660
  /* compute the value at M[1<<(winsize-1)] by squaring M[1] (winsize-1) times */
2671
2661
  if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) {
2672
 
    goto __RES;
 
2662
    goto LBL_RES;
2673
2663
  }
2674
2664
 
2675
2665
  for (x = 0; x < (winsize - 1); x++) {
2676
2666
    if ((err = mp_sqr (&M[1 << (winsize - 1)], &M[1 << (winsize - 1)])) != MP_OKAY) {
2677
 
      goto __RES;
 
2667
      goto LBL_RES;
2678
2668
    }
2679
2669
    if ((err = redux (&M[1 << (winsize - 1)], P, mp)) != MP_OKAY) {
2680
 
      goto __RES;
 
2670
      goto LBL_RES;
2681
2671
    }
2682
2672
  }
2683
2673
 
2684
2674
  /* create upper table */
2685
2675
  for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) {
2686
2676
    if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) {
2687
 
      goto __RES;
 
2677
      goto LBL_RES;
2688
2678
    }
2689
2679
    if ((err = redux (&M[x], P, mp)) != MP_OKAY) {
2690
 
      goto __RES;
 
2680
      goto LBL_RES;
2691
2681
    }
2692
2682
  }
2693
2683
 
2727
2717
    /* if the bit is zero and mode == 1 then we square */
2728
2718
    if (mode == 1 && y == 0) {
2729
2719
      if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
2730
 
        goto __RES;
 
2720
        goto LBL_RES;
2731
2721
      }
2732
2722
      if ((err = redux (&res, P, mp)) != MP_OKAY) {
2733
 
        goto __RES;
 
2723
        goto LBL_RES;
2734
2724
      }
2735
2725
      continue;
2736
2726
    }
2744
2734
      /* square first */
2745
2735
      for (x = 0; x < winsize; x++) {
2746
2736
        if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
2747
 
          goto __RES;
 
2737
          goto LBL_RES;
2748
2738
        }
2749
2739
        if ((err = redux (&res, P, mp)) != MP_OKAY) {
2750
 
          goto __RES;
 
2740
          goto LBL_RES;
2751
2741
        }
2752
2742
      }
2753
2743
 
2754
2744
      /* then multiply */
2755
2745
      if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) {
2756
 
        goto __RES;
 
2746
        goto LBL_RES;
2757
2747
      }
2758
2748
      if ((err = redux (&res, P, mp)) != MP_OKAY) {
2759
 
        goto __RES;
 
2749
        goto LBL_RES;
2760
2750
      }
2761
2751
 
2762
2752
      /* empty window and reset */
2771
2761
    /* square then multiply if the bit is set */
2772
2762
    for (x = 0; x < bitcpy; x++) {
2773
2763
      if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
2774
 
        goto __RES;
 
2764
        goto LBL_RES;
2775
2765
      }
2776
2766
      if ((err = redux (&res, P, mp)) != MP_OKAY) {
2777
 
        goto __RES;
 
2767
        goto LBL_RES;
2778
2768
      }
2779
2769
 
2780
2770
      /* get next bit of the window */
2782
2772
      if ((bitbuf & (1 << winsize)) != 0) {
2783
2773
        /* then multiply */
2784
2774
        if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) {
2785
 
          goto __RES;
 
2775
          goto LBL_RES;
2786
2776
        }
2787
2777
        if ((err = redux (&res, P, mp)) != MP_OKAY) {
2788
 
          goto __RES;
 
2778
          goto LBL_RES;
2789
2779
        }
2790
2780
      }
2791
2781
    }
2799
2789
      * of R.
2800
2790
      */
2801
2791
     if ((err = redux(&res, P, mp)) != MP_OKAY) {
2802
 
       goto __RES;
 
2792
       goto LBL_RES;
2803
2793
     }
2804
2794
  }
2805
2795
 
2806
2796
  /* swap res with Y */
2807
2797
  mp_exch (&res, Y);
2808
2798
  err = MP_OKAY;
2809
 
__RES:mp_clear (&res);
2810
 
__M:
 
2799
LBL_RES:mp_clear (&res);
 
2800
LBL_M:
2811
2801
  mp_clear(&M[1]);
2812
2802
  for (x = 1<<(winsize-1); x < (1 << winsize); x++) {
2813
2803
    mp_clear (&M[x]);
2881
2871
       if ((err = mp_copy(&t3, &v3)) != MP_OKAY)                                  { goto _ERR; }
2882
2872
   }
2883
2873
 
 
2874
   /* make sure U3 >= 0 */
 
2875
   if (u3.sign == MP_NEG) {
 
2876
      mp_neg(&u1, &u1);
 
2877
      mp_neg(&u2, &u2);
 
2878
      mp_neg(&u3, &u3);
 
2879
   }
 
2880
 
2884
2881
   /* copy result out */
2885
2882
   if (U1 != NULL) { mp_exch(U1, &u1); }
2886
2883
   if (U2 != NULL) { mp_exch(U2, &u2); }
3059
3056
  }
3060
3057
 
3061
3058
  if ((res = mp_init_copy (&v, b)) != MP_OKAY) {
3062
 
    goto __U;
 
3059
    goto LBL_U;
3063
3060
  }
3064
3061
 
3065
3062
  /* must be positive for the remainder of the algorithm */
3073
3070
  if (k > 0) {
3074
3071
     /* divide the power of two out */
3075
3072
     if ((res = mp_div_2d(&u, k, &u, NULL)) != MP_OKAY) {
3076
 
        goto __V;
 
3073
        goto LBL_V;
3077
3074
     }
3078
3075
 
3079
3076
     if ((res = mp_div_2d(&v, k, &v, NULL)) != MP_OKAY) {
3080
 
        goto __V;
 
3077
        goto LBL_V;
3081
3078
     }
3082
3079
  }
3083
3080
 
3084
3081
  /* divide any remaining factors of two out */
3085
3082
  if (u_lsb != k) {
3086
3083
     if ((res = mp_div_2d(&u, u_lsb - k, &u, NULL)) != MP_OKAY) {
3087
 
        goto __V;
 
3084
        goto LBL_V;
3088
3085
     }
3089
3086
  }
3090
3087
 
3091
3088
  if (v_lsb != k) {
3092
3089
     if ((res = mp_div_2d(&v, v_lsb - k, &v, NULL)) != MP_OKAY) {
3093
 
        goto __V;
 
3090
        goto LBL_V;
3094
3091
     }
3095
3092
  }
3096
3093
 
3103
3100
     
3104
3101
     /* subtract smallest from largest */
3105
3102
     if ((res = s_mp_sub(&v, &u, &v)) != MP_OKAY) {
3106
 
        goto __V;
 
3103
        goto LBL_V;
3107
3104
     }
3108
3105
     
3109
3106
     /* Divide out all factors of two */
3110
3107
     if ((res = mp_div_2d(&v, mp_cnt_lsb(&v), &v, NULL)) != MP_OKAY) {
3111
 
        goto __V;
 
3108
        goto LBL_V;
3112
3109
     } 
3113
3110
  } 
3114
3111
 
3115
3112
  /* multiply by 2**k which we divided out at the beginning */
3116
3113
  if ((res = mp_mul_2d (&u, k, c)) != MP_OKAY) {
3117
 
     goto __V;
 
3114
     goto LBL_V;
3118
3115
  }
3119
3116
  c->sign = MP_ZPOS;
3120
3117
  res = MP_OKAY;
3121
 
__V:mp_clear (&u);
3122
 
__U:mp_clear (&v);
 
3118
LBL_V:mp_clear (&u);
 
3119
LBL_U:mp_clear (&v);
3123
3120
  return res;
3124
3121
}
3125
3122
#endif
3555
3552
  }
3556
3553
 
3557
3554
  /* x = a, y = b */
3558
 
  if ((res = mp_copy (a, &x)) != MP_OKAY) {
3559
 
    goto __ERR;
 
3555
  if ((res = mp_mod(a, b, &x)) != MP_OKAY) {
 
3556
      goto LBL_ERR;
3560
3557
  }
3561
3558
  if ((res = mp_copy (b, &y)) != MP_OKAY) {
3562
 
    goto __ERR;
 
3559
    goto LBL_ERR;
3563
3560
  }
3564
3561
 
3565
3562
  /* 2. [modified] if x,y are both even then return an error! */
3566
3563
  if (mp_iseven (&x) == 1 && mp_iseven (&y) == 1) {
3567
3564
    res = MP_VAL;
3568
 
    goto __ERR;
 
3565
    goto LBL_ERR;
3569
3566
  }
3570
3567
 
3571
3568
  /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */
3572
3569
  if ((res = mp_copy (&x, &u)) != MP_OKAY) {
3573
 
    goto __ERR;
 
3570
    goto LBL_ERR;
3574
3571
  }
3575
3572
  if ((res = mp_copy (&y, &v)) != MP_OKAY) {
3576
 
    goto __ERR;
 
3573
    goto LBL_ERR;
3577
3574
  }
3578
3575
  mp_set (&A, 1);
3579
3576
  mp_set (&D, 1);
3583
3580
  while (mp_iseven (&u) == 1) {
3584
3581
    /* 4.1 u = u/2 */
3585
3582
    if ((res = mp_div_2 (&u, &u)) != MP_OKAY) {
3586
 
      goto __ERR;
 
3583
      goto LBL_ERR;
3587
3584
    }
3588
3585
    /* 4.2 if A or B is odd then */
3589
3586
    if (mp_isodd (&A) == 1 || mp_isodd (&B) == 1) {
3590
3587
      /* A = (A+y)/2, B = (B-x)/2 */
3591
3588
      if ((res = mp_add (&A, &y, &A)) != MP_OKAY) {
3592
 
         goto __ERR;
 
3589
         goto LBL_ERR;
3593
3590
      }
3594
3591
      if ((res = mp_sub (&B, &x, &B)) != MP_OKAY) {
3595
 
         goto __ERR;
 
3592
         goto LBL_ERR;
3596
3593
      }
3597
3594
    }
3598
3595
    /* A = A/2, B = B/2 */
3599
3596
    if ((res = mp_div_2 (&A, &A)) != MP_OKAY) {
3600
 
      goto __ERR;
 
3597
      goto LBL_ERR;
3601
3598
    }
3602
3599
    if ((res = mp_div_2 (&B, &B)) != MP_OKAY) {
3603
 
      goto __ERR;
 
3600
      goto LBL_ERR;
3604
3601
    }
3605
3602
  }
3606
3603
 
3608
3605
  while (mp_iseven (&v) == 1) {
3609
3606
    /* 5.1 v = v/2 */
3610
3607
    if ((res = mp_div_2 (&v, &v)) != MP_OKAY) {
3611
 
      goto __ERR;
 
3608
      goto LBL_ERR;
3612
3609
    }
3613
3610
    /* 5.2 if C or D is odd then */
3614
3611
    if (mp_isodd (&C) == 1 || mp_isodd (&D) == 1) {
3615
3612
      /* C = (C+y)/2, D = (D-x)/2 */
3616
3613
      if ((res = mp_add (&C, &y, &C)) != MP_OKAY) {
3617
 
         goto __ERR;
 
3614
         goto LBL_ERR;
3618
3615
      }
3619
3616
      if ((res = mp_sub (&D, &x, &D)) != MP_OKAY) {
3620
 
         goto __ERR;
 
3617
         goto LBL_ERR;
3621
3618
      }
3622
3619
    }
3623
3620
    /* C = C/2, D = D/2 */
3624
3621
    if ((res = mp_div_2 (&C, &C)) != MP_OKAY) {
3625
 
      goto __ERR;
 
3622
      goto LBL_ERR;
3626
3623
    }
3627
3624
    if ((res = mp_div_2 (&D, &D)) != MP_OKAY) {
3628
 
      goto __ERR;
 
3625
      goto LBL_ERR;
3629
3626
    }
3630
3627
  }
3631
3628
 
3633
3630
  if (mp_cmp (&u, &v) != MP_LT) {
3634
3631
    /* u = u - v, A = A - C, B = B - D */
3635
3632
    if ((res = mp_sub (&u, &v, &u)) != MP_OKAY) {
3636
 
      goto __ERR;
 
3633
      goto LBL_ERR;
3637
3634
    }
3638
3635
 
3639
3636
    if ((res = mp_sub (&A, &C, &A)) != MP_OKAY) {
3640
 
      goto __ERR;
 
3637
      goto LBL_ERR;
3641
3638
    }
3642
3639
 
3643
3640
    if ((res = mp_sub (&B, &D, &B)) != MP_OKAY) {
3644
 
      goto __ERR;
 
3641
      goto LBL_ERR;
3645
3642
    }
3646
3643
  } else {
3647
3644
    /* v - v - u, C = C - A, D = D - B */
3648
3645
    if ((res = mp_sub (&v, &u, &v)) != MP_OKAY) {
3649
 
      goto __ERR;
 
3646
      goto LBL_ERR;
3650
3647
    }
3651
3648
 
3652
3649
    if ((res = mp_sub (&C, &A, &C)) != MP_OKAY) {
3653
 
      goto __ERR;
 
3650
      goto LBL_ERR;
3654
3651
    }
3655
3652
 
3656
3653
    if ((res = mp_sub (&D, &B, &D)) != MP_OKAY) {
3657
 
      goto __ERR;
 
3654
      goto LBL_ERR;
3658
3655
    }
3659
3656
  }
3660
3657
 
3667
3664
  /* if v != 1 then there is no inverse */
3668
3665
  if (mp_cmp_d (&v, 1) != MP_EQ) {
3669
3666
    res = MP_VAL;
3670
 
    goto __ERR;
 
3667
    goto LBL_ERR;
3671
3668
  }
3672
3669
 
3673
3670
  /* if its too low */
3674
3671
  while (mp_cmp_d(&C, 0) == MP_LT) {
3675
3672
      if ((res = mp_add(&C, b, &C)) != MP_OKAY) {
3676
 
         goto __ERR;
 
3673
         goto LBL_ERR;
3677
3674
      }
3678
3675
  }
3679
3676
  
3680
3677
  /* too big */
3681
3678
  while (mp_cmp_mag(&C, b) != MP_LT) {
3682
3679
      if ((res = mp_sub(&C, b, &C)) != MP_OKAY) {
3683
 
         goto __ERR;
 
3680
         goto LBL_ERR;
3684
3681
      }
3685
3682
  }
3686
3683
  
3687
3684
  /* C is now the inverse */
3688
3685
  mp_exch (&C, c);
3689
3686
  res = MP_OKAY;
3690
 
__ERR:mp_clear_multi (&x, &y, &u, &v, &A, &B, &C, &D, NULL);
 
3687
LBL_ERR:mp_clear_multi (&x, &y, &u, &v, &A, &B, &C, &D, NULL);
3691
3688
  return res;
3692
3689
}
3693
3690
#endif
3856
3853
  }
3857
3854
 
3858
3855
  if ((res = mp_init (&p1)) != MP_OKAY) {
3859
 
    goto __A1;
 
3856
    goto LBL_A1;
3860
3857
  }
3861
3858
 
3862
3859
  /* divide out larger power of two */
3863
3860
  k = mp_cnt_lsb(&a1);
3864
3861
  if ((res = mp_div_2d(&a1, k, &a1, NULL)) != MP_OKAY) {
3865
 
     goto __P1;
 
3862
     goto LBL_P1;
3866
3863
  }
3867
3864
 
3868
3865
  /* step 4.  if e is even set s=1 */
3890
3887
  } else {
3891
3888
    /* n1 = n mod a1 */
3892
3889
    if ((res = mp_mod (p, &a1, &p1)) != MP_OKAY) {
3893
 
      goto __P1;
 
3890
      goto LBL_P1;
3894
3891
    }
3895
3892
    if ((res = mp_jacobi (&p1, &a1, &r)) != MP_OKAY) {
3896
 
      goto __P1;
 
3893
      goto LBL_P1;
3897
3894
    }
3898
3895
    *c = s * r;
3899
3896
  }
3900
3897
 
3901
3898
  /* done */
3902
3899
  res = MP_OKAY;
3903
 
__P1:mp_clear (&p1);
3904
 
__A1:mp_clear (&a1);
 
3900
LBL_P1:mp_clear (&p1);
 
3901
LBL_A1:mp_clear (&a1);
3905
3902
  return res;
3906
3903
}
3907
3904
#endif
4227
4224
 
4228
4225
  /* t1 = get the GCD of the two inputs */
4229
4226
  if ((res = mp_gcd (a, b, &t1)) != MP_OKAY) {
4230
 
    goto __T;
 
4227
    goto LBL_T;
4231
4228
  }
4232
4229
 
4233
4230
  /* divide the smallest by the GCD */
4234
4231
  if (mp_cmp_mag(a, b) == MP_LT) {
4235
4232
     /* store quotient in t2 such that t2 * b is the LCM */
4236
4233
     if ((res = mp_div(a, &t1, &t2, NULL)) != MP_OKAY) {
4237
 
        goto __T;
 
4234
        goto LBL_T;
4238
4235
     }
4239
4236
     res = mp_mul(b, &t2, c);
4240
4237
  } else {
4241
4238
     /* store quotient in t2 such that t2 * a is the LCM */
4242
4239
     if ((res = mp_div(b, &t1, &t2, NULL)) != MP_OKAY) {
4243
 
        goto __T;
 
4240
        goto LBL_T;
4244
4241
     }
4245
4242
     res = mp_mul(a, &t2, c);
4246
4243
  }
4248
4245
  /* fix the sign to positive */
4249
4246
  c->sign = MP_ZPOS;
4250
4247
 
4251
 
__T:
 
4248
LBL_T:
4252
4249
  mp_clear_multi (&t1, &t2, NULL);
4253
4250
  return res;
4254
4251
}
4402
4399
  }
4403
4400
 
4404
4401
  /* if the modulus is larger than the value than return */
4405
 
  if (b > (int) (a->used * DIGIT_BIT)) {
 
4402
  if (b >= (int) (a->used * DIGIT_BIT)) {
4406
4403
    res = mp_copy (a, c);
4407
4404
    return res;
4408
4405
  }
4484
4481
  /* how many bits of last digit does b use */
4485
4482
  bits = mp_count_bits (b) % DIGIT_BIT;
4486
4483
 
4487
 
 
4488
4484
  if (b->used > 1) {
4489
4485
     if ((res = mp_2expt (a, (b->used - 1) * DIGIT_BIT + bits - 1)) != MP_OKAY) {
4490
4486
        return res;
4983
4979
    u       = (mp_digit) (r >> ((mp_word) DIGIT_BIT));
4984
4980
  }
4985
4981
 
4986
 
  /* store final carry [if any] */
 
4982
  /* store final carry [if any] and increment ix offset  */
4987
4983
  *tmpc++ = u;
 
4984
  ++ix;
4988
4985
 
4989
4986
  /* now zero digits above the top */
4990
4987
  while (ix++ < olduse) {
5085
5082
  }
5086
5083
 
5087
5084
  if ((res = mp_init (&t2)) != MP_OKAY) {
5088
 
    goto __T1;
 
5085
    goto LBL_T1;
5089
5086
  }
5090
5087
 
5091
5088
  if ((res = mp_init (&t3)) != MP_OKAY) {
5092
 
    goto __T2;
 
5089
    goto LBL_T2;
5093
5090
  }
5094
5091
 
5095
5092
  /* if a is negative fudge the sign but keep track */
5102
5099
  do {
5103
5100
    /* t1 = t2 */
5104
5101
    if ((res = mp_copy (&t2, &t1)) != MP_OKAY) {
5105
 
      goto __T3;
 
5102
      goto LBL_T3;
5106
5103
    }
5107
5104
 
5108
5105
    /* t2 = t1 - ((t1**b - a) / (b * t1**(b-1))) */
5109
5106
    
5110
5107
    /* t3 = t1**(b-1) */
5111
5108
    if ((res = mp_expt_d (&t1, b - 1, &t3)) != MP_OKAY) {   
5112
 
      goto __T3;
 
5109
      goto LBL_T3;
5113
5110
    }
5114
5111
 
5115
5112
    /* numerator */
5116
5113
    /* t2 = t1**b */
5117
5114
    if ((res = mp_mul (&t3, &t1, &t2)) != MP_OKAY) {    
5118
 
      goto __T3;
 
5115
      goto LBL_T3;
5119
5116
    }
5120
5117
 
5121
5118
    /* t2 = t1**b - a */
5122
5119
    if ((res = mp_sub (&t2, a, &t2)) != MP_OKAY) {  
5123
 
      goto __T3;
 
5120
      goto LBL_T3;
5124
5121
    }
5125
5122
 
5126
5123
    /* denominator */
5127
5124
    /* t3 = t1**(b-1) * b  */
5128
5125
    if ((res = mp_mul_d (&t3, b, &t3)) != MP_OKAY) {    
5129
 
      goto __T3;
 
5126
      goto LBL_T3;
5130
5127
    }
5131
5128
 
5132
5129
    /* t3 = (t1**b - a)/(b * t1**(b-1)) */
5133
5130
    if ((res = mp_div (&t2, &t3, &t3, NULL)) != MP_OKAY) {  
5134
 
      goto __T3;
 
5131
      goto LBL_T3;
5135
5132
    }
5136
5133
 
5137
5134
    if ((res = mp_sub (&t1, &t3, &t2)) != MP_OKAY) {
5138
 
      goto __T3;
 
5135
      goto LBL_T3;
5139
5136
    }
5140
5137
  }  while (mp_cmp (&t1, &t2) != MP_EQ);
5141
5138
 
5142
5139
  /* result can be off by a few so check */
5143
5140
  for (;;) {
5144
5141
    if ((res = mp_expt_d (&t1, b, &t2)) != MP_OKAY) {
5145
 
      goto __T3;
 
5142
      goto LBL_T3;
5146
5143
    }
5147
5144
 
5148
5145
    if (mp_cmp (&t2, a) == MP_GT) {
5149
5146
      if ((res = mp_sub_d (&t1, 1, &t1)) != MP_OKAY) {
5150
 
         goto __T3;
 
5147
         goto LBL_T3;
5151
5148
      }
5152
5149
    } else {
5153
5150
      break;
5165
5162
 
5166
5163
  res = MP_OKAY;
5167
5164
 
5168
 
__T3:mp_clear (&t3);
5169
 
__T2:mp_clear (&t2);
5170
 
__T1:mp_clear (&t1);
 
5165
LBL_T3:mp_clear (&t3);
 
5166
LBL_T2:mp_clear (&t2);
 
5167
LBL_T1:mp_clear (&t1);
5171
5168
  return res;
5172
5169
}
5173
5170
#endif
5196
5193
int mp_neg (mp_int * a, mp_int * b)
5197
5194
{
5198
5195
  int     res;
5199
 
  if ((res = mp_copy (a, b)) != MP_OKAY) {
5200
 
    return res;
 
5196
  if (a != b) {
 
5197
     if ((res = mp_copy (a, b)) != MP_OKAY) {
 
5198
        return res;
 
5199
     }
5201
5200
  }
 
5201
 
5202
5202
  if (mp_iszero(b) != MP_YES) {
5203
5203
     b->sign = (a->sign == MP_ZPOS) ? MP_NEG : MP_ZPOS;
 
5204
  } else {
 
5205
     b->sign = MP_ZPOS;
5204
5206
  }
 
5207
 
5205
5208
  return MP_OKAY;
5206
5209
}
5207
5210
#endif
5304
5307
 
5305
5308
  /* compute t = b**a mod a */
5306
5309
  if ((err = mp_exptmod (b, a, a, &t)) != MP_OKAY) {
5307
 
    goto __T;
 
5310
    goto LBL_T;
5308
5311
  }
5309
5312
 
5310
5313
  /* is it equal to b? */
5313
5316
  }
5314
5317
 
5315
5318
  err = MP_OKAY;
5316
 
__T:mp_clear (&t);
 
5319
LBL_T:mp_clear (&t);
5317
5320
  return err;
5318
5321
}
5319
5322
#endif
5352
5355
  *result = MP_NO;
5353
5356
 
5354
5357
  for (ix = 0; ix < PRIME_SIZE; ix++) {
5355
 
    /* what is a mod __prime_tab[ix] */
5356
 
    if ((err = mp_mod_d (a, __prime_tab[ix], &res)) != MP_OKAY) {
 
5358
    /* what is a mod LBL_prime_tab[ix] */
 
5359
    if ((err = mp_mod_d (a, ltm_prime_tab[ix], &res)) != MP_OKAY) {
5357
5360
      return err;
5358
5361
    }
5359
5362
 
5410
5413
 
5411
5414
  /* is the input equal to one of the primes in the table? */
5412
5415
  for (ix = 0; ix < PRIME_SIZE; ix++) {
5413
 
      if (mp_cmp_d(a, __prime_tab[ix]) == MP_EQ) {
 
5416
      if (mp_cmp_d(a, ltm_prime_tab[ix]) == MP_EQ) {
5414
5417
         *result = 1;
5415
5418
         return MP_OKAY;
5416
5419
      }
5433
5436
 
5434
5437
  for (ix = 0; ix < t; ix++) {
5435
5438
    /* set the prime */
5436
 
    mp_set (&b, __prime_tab[ix]);
 
5439
    mp_set (&b, ltm_prime_tab[ix]);
5437
5440
 
5438
5441
    if ((err = mp_prime_miller_rabin (a, &b, &res)) != MP_OKAY) {
5439
 
      goto __B;
 
5442
      goto LBL_B;
5440
5443
    }
5441
5444
 
5442
5445
    if (res == MP_NO) {
5443
 
      goto __B;
 
5446
      goto LBL_B;
5444
5447
    }
5445
5448
  }
5446
5449
 
5447
5450
  /* passed the test */
5448
5451
  *result = MP_YES;
5449
 
__B:mp_clear (&b);
 
5452
LBL_B:mp_clear (&b);
5450
5453
  return err;
5451
5454
}
5452
5455
#endif
5496
5499
    return err;
5497
5500
  }
5498
5501
  if ((err = mp_sub_d (&n1, 1, &n1)) != MP_OKAY) {
5499
 
    goto __N1;
 
5502
    goto LBL_N1;
5500
5503
  }
5501
5504
 
5502
5505
  /* set 2**s * r = n1 */
5503
5506
  if ((err = mp_init_copy (&r, &n1)) != MP_OKAY) {
5504
 
    goto __N1;
 
5507
    goto LBL_N1;
5505
5508
  }
5506
5509
 
5507
5510
  /* count the number of least significant bits
5511
5514
 
5512
5515
  /* now divide n - 1 by 2**s */
5513
5516
  if ((err = mp_div_2d (&r, s, &r, NULL)) != MP_OKAY) {
5514
 
    goto __R;
 
5517
    goto LBL_R;
5515
5518
  }
5516
5519
 
5517
5520
  /* compute y = b**r mod a */
5518
5521
  if ((err = mp_init (&y)) != MP_OKAY) {
5519
 
    goto __R;
 
5522
    goto LBL_R;
5520
5523
  }
5521
5524
  if ((err = mp_exptmod (b, &r, a, &y)) != MP_OKAY) {
5522
 
    goto __Y;
 
5525
    goto LBL_Y;
5523
5526
  }
5524
5527
 
5525
5528
  /* if y != 1 and y != n1 do */
5528
5531
    /* while j <= s-1 and y != n1 */
5529
5532
    while ((j <= (s - 1)) && mp_cmp (&y, &n1) != MP_EQ) {
5530
5533
      if ((err = mp_sqrmod (&y, a, &y)) != MP_OKAY) {
5531
 
         goto __Y;
 
5534
         goto LBL_Y;
5532
5535
      }
5533
5536
 
5534
5537
      /* if y == 1 then composite */
5535
5538
      if (mp_cmp_d (&y, 1) == MP_EQ) {
5536
 
         goto __Y;
 
5539
         goto LBL_Y;
5537
5540
      }
5538
5541
 
5539
5542
      ++j;
5541
5544
 
5542
5545
    /* if y != n1 then composite */
5543
5546
    if (mp_cmp (&y, &n1) != MP_EQ) {
5544
 
      goto __Y;
 
5547
      goto LBL_Y;
5545
5548
    }
5546
5549
  }
5547
5550
 
5548
5551
  /* probably prime now */
5549
5552
  *result = MP_YES;
5550
 
__Y:mp_clear (&y);
5551
 
__R:mp_clear (&r);
5552
 
__N1:mp_clear (&n1);
 
5553
LBL_Y:mp_clear (&y);
 
5554
LBL_R:mp_clear (&r);
 
5555
LBL_N1:mp_clear (&n1);
5553
5556
  return err;
5554
5557
}
5555
5558
#endif
5594
5597
   a->sign = MP_ZPOS;
5595
5598
 
5596
5599
   /* simple algo if a is less than the largest prime in the table */
5597
 
   if (mp_cmp_d(a, __prime_tab[PRIME_SIZE-1]) == MP_LT) {
 
5600
   if (mp_cmp_d(a, ltm_prime_tab[PRIME_SIZE-1]) == MP_LT) {
5598
5601
      /* find which prime it is bigger than */
5599
5602
      for (x = PRIME_SIZE - 2; x >= 0; x--) {
5600
 
          if (mp_cmp_d(a, __prime_tab[x]) != MP_LT) {
 
5603
          if (mp_cmp_d(a, ltm_prime_tab[x]) != MP_LT) {
5601
5604
             if (bbs_style == 1) {
5602
5605
                /* ok we found a prime smaller or
5603
5606
                 * equal [so the next is larger]
5605
5608
                 * however, the prime must be
5606
5609
                 * congruent to 3 mod 4
5607
5610
                 */
5608
 
                if ((__prime_tab[x + 1] & 3) != 3) {
 
5611
                if ((ltm_prime_tab[x + 1] & 3) != 3) {
5609
5612
                   /* scan upwards for a prime congruent to 3 mod 4 */
5610
5613
                   for (y = x + 1; y < PRIME_SIZE; y++) {
5611
 
                       if ((__prime_tab[y] & 3) == 3) {
5612
 
                          mp_set(a, __prime_tab[y]);
 
5614
                       if ((ltm_prime_tab[y] & 3) == 3) {
 
5615
                          mp_set(a, ltm_prime_tab[y]);
5613
5616
                          return MP_OKAY;
5614
5617
                       }
5615
5618
                   }
5616
5619
                }
5617
5620
             } else {
5618
 
                mp_set(a, __prime_tab[x + 1]);
 
5621
                mp_set(a, ltm_prime_tab[x + 1]);
5619
5622
                return MP_OKAY;
5620
5623
             }
5621
5624
          }
5653
5656
 
5654
5657
   /* generate the restable */
5655
5658
   for (x = 1; x < PRIME_SIZE; x++) {
5656
 
      if ((err = mp_mod_d(a, __prime_tab[x], res_tab + x)) != MP_OKAY) {
 
5659
      if ((err = mp_mod_d(a, ltm_prime_tab[x], res_tab + x)) != MP_OKAY) {
5657
5660
         return err;
5658
5661
      }
5659
5662
   }
5679
5682
             res_tab[x] += kstep;
5680
5683
 
5681
5684
             /* subtract the modulus [instead of using division] */
5682
 
             if (res_tab[x] >= __prime_tab[x]) {
5683
 
                res_tab[x]  -= __prime_tab[x];
 
5685
             if (res_tab[x] >= ltm_prime_tab[x]) {
 
5686
                res_tab[x]  -= ltm_prime_tab[x];
5684
5687
             }
5685
5688
 
5686
5689
             /* set flag if zero */
5692
5695
 
5693
5696
      /* add the step */
5694
5697
      if ((err = mp_add_d(a, step, a)) != MP_OKAY) {
5695
 
         goto __ERR;
 
5698
         goto LBL_ERR;
5696
5699
      }
5697
5700
 
5698
5701
      /* if didn't pass sieve and step == MAX then skip test */
5702
5705
 
5703
5706
      /* is this prime? */
5704
5707
      for (x = 0; x < t; x++) {
5705
 
          mp_set(&b, __prime_tab[t]);
 
5708
          mp_set(&b, ltm_prime_tab[t]);
5706
5709
          if ((err = mp_prime_miller_rabin(a, &b, &res)) != MP_OKAY) {
5707
 
             goto __ERR;
 
5710
             goto LBL_ERR;
5708
5711
          }
5709
5712
          if (res == MP_NO) {
5710
5713
             break;
5717
5720
   }
5718
5721
 
5719
5722
   err = MP_OKAY;
5720
 
__ERR:
 
5723
LBL_ERR:
5721
5724
   mp_clear(&b);
5722
5725
   return err;
5723
5726
}
5828
5831
   }
5829
5832
 
5830
5833
   /* calc the byte size */
5831
 
   bsize = (size>>3)+(size&7?1:0);
 
5834
   bsize = (size>>3) + ((size&7)?1:0);
5832
5835
 
5833
5836
   /* we need a buffer of bsize bytes */
5834
5837
   tmp = OPT_CAST(unsigned char) XMALLOC(bsize);
5837
5840
   }
5838
5841
 
5839
5842
   /* calc the maskAND value for the MSbyte*/
5840
 
   maskAND = 0xFF >> (8 - (size & 7));
 
5843
   maskAND = ((size&7) == 0) ? 0xFF : (0xFF >> (8 - (size & 7)));
5841
5844
 
5842
5845
   /* calc the maskOR_msb */
5843
5846
   maskOR_msb        = 0;
5844
 
   maskOR_msb_offset = (size - 2) >> 3;
 
5847
   maskOR_msb_offset = ((size & 7) == 1) ? 1 : 0;
5845
5848
   if (flags & LTM_PRIME_2MSB_ON) {
5846
5849
      maskOR_msb     |= 1 << ((size - 2) & 7);
5847
5850
   } else if (flags & LTM_PRIME_2MSB_OFF) {
5848
5851
      maskAND        &= ~(1 << ((size - 2) & 7));
5849
 
   }
 
5852
   } 
5850
5853
 
5851
5854
   /* get the maskOR_lsb */
5852
 
   maskOR_lsb         = 0;
 
5855
   maskOR_lsb         = 1;
5853
5856
   if (flags & LTM_PRIME_BBS) {
5854
5857
      maskOR_lsb     |= 3;
5855
5858
   }
5943
5946
    return MP_VAL;
5944
5947
  }
5945
5948
 
5946
 
  /* init a copy of the input */
5947
 
  if ((res = mp_init_copy (&t, a)) != MP_OKAY) {
5948
 
    return res;
 
5949
  if (mp_iszero(a) == MP_YES) {
 
5950
     *size = 2;
 
5951
    return MP_OKAY;
5949
5952
  }
5950
5953
 
5951
5954
  /* digs is the digit count */
5952
5955
  digs = 0;
5953
5956
 
5954
5957
  /* if it's negative add one for the sign */
5955
 
  if (t.sign == MP_NEG) {
 
5958
  if (a->sign == MP_NEG) {
5956
5959
    ++digs;
5957
 
    t.sign = MP_ZPOS;
5958
 
  }
 
5960
  }
 
5961
 
 
5962
  /* init a copy of the input */
 
5963
  if ((res = mp_init_copy (&t, a)) != MP_OKAY) {
 
5964
    return res;
 
5965
  }
 
5966
 
 
5967
  /* force temp to positive */
 
5968
  t.sign = MP_ZPOS; 
5959
5969
 
5960
5970
  /* fetch out all of the digits */
5961
 
  while (mp_iszero (&t) == 0) {
 
5971
  while (mp_iszero (&t) == MP_NO) {
5962
5972
    if ((res = mp_div_d (&t, (mp_digit) radix, &t, &d)) != MP_OKAY) {
5963
5973
      mp_clear (&t);
5964
5974
      return res;
6032
6042
 
6033
6043
  /* first place a random non-zero digit */
6034
6044
  do {
6035
 
    d = ((mp_digit) abs (rand ()));
 
6045
    d = ((mp_digit) abs (rand ())) & MP_MASK;
6036
6046
  } while (d == 0);
6037
6047
 
6038
6048
  if ((res = mp_add_d (a, d, a)) != MP_OKAY) {
6039
6049
    return res;
6040
6050
  }
6041
6051
 
6042
 
  while (digits-- > 0) {
 
6052
  while (--digits > 0) {
6043
6053
    if ((res = mp_lshd (a, 1)) != MP_OKAY) {
6044
6054
      return res;
6045
6055
    }
6074
6084
 */
6075
6085
 
6076
6086
/* read a string [ASCII] in a given radix */
6077
 
int mp_read_radix (mp_int * a, char *str, int radix)
 
6087
int mp_read_radix (mp_int * a, const char *str, int radix)
6078
6088
{
6079
6089
  int     y, res, neg;
6080
6090
  char    ch;
6257
6267
 * precomputed via mp_reduce_setup.
6258
6268
 * From HAC pp.604 Algorithm 14.42
6259
6269
 */
6260
 
int
6261
 
mp_reduce (mp_int * x, mp_int * m, mp_int * mu)
 
6270
int mp_reduce (mp_int * x, mp_int * m, mp_int * mu)
6262
6271
{
6263
6272
  mp_int  q;
6264
6273
  int     res, um = m->used;
6278
6287
    }
6279
6288
  } else {
6280
6289
#ifdef BN_S_MP_MUL_HIGH_DIGS_C
6281
 
    if ((res = s_mp_mul_high_digs (&q, mu, &q, um - 1)) != MP_OKAY) {
 
6290
    if ((res = s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) {
6282
6291
      goto CLEANUP;
6283
6292
    }
6284
6293
#elif defined(BN_FAST_S_MP_MUL_HIGH_DIGS_C)
6285
 
    if ((res = fast_s_mp_mul_high_digs (&q, mu, &q, um - 1)) != MP_OKAY) {
 
6294
    if ((res = fast_s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) {
6286
6295
      goto CLEANUP;
6287
6296
    }
6288
6297
#else 
6355
6364
 */
6356
6365
 
6357
6366
/* reduces a modulo n where n is of the form 2**p - d */
6358
 
int
6359
 
mp_reduce_2k(mp_int *a, mp_int *n, mp_digit d)
 
6367
int mp_reduce_2k(mp_int *a, mp_int *n, mp_digit d)
6360
6368
{
6361
6369
   mp_int q;
6362
6370
   int    p, res;
6398
6406
 
6399
6407
/* End: bn_mp_reduce_2k.c */
6400
6408
 
 
6409
/* Start: bn_mp_reduce_2k_l.c */
 
6410
#include <tommath.h>
 
6411
#ifdef BN_MP_REDUCE_2K_L_C
 
6412
/* LibTomMath, multiple-precision integer library -- Tom St Denis
 
6413
 *
 
6414
 * LibTomMath is a library that provides multiple-precision
 
6415
 * integer arithmetic as well as number theoretic functionality.
 
6416
 *
 
6417
 * The library was designed directly after the MPI library by
 
6418
 * Michael Fromberger but has been written from scratch with
 
6419
 * additional optimizations in place.
 
6420
 *
 
6421
 * The library is free for all purposes without any express
 
6422
 * guarantee it works.
 
6423
 *
 
6424
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
 
6425
 */
 
6426
 
 
6427
/* reduces a modulo n where n is of the form 2**p - d 
 
6428
   This differs from reduce_2k since "d" can be larger
 
6429
   than a single digit.
 
6430
*/
 
6431
int mp_reduce_2k_l(mp_int *a, mp_int *n, mp_int *d)
 
6432
{
 
6433
   mp_int q;
 
6434
   int    p, res;
 
6435
   
 
6436
   if ((res = mp_init(&q)) != MP_OKAY) {
 
6437
      return res;
 
6438
   }
 
6439
   
 
6440
   p = mp_count_bits(n);    
 
6441
top:
 
6442
   /* q = a/2**p, a = a mod 2**p */
 
6443
   if ((res = mp_div_2d(a, p, &q, a)) != MP_OKAY) {
 
6444
      goto ERR;
 
6445
   }
 
6446
   
 
6447
   /* q = q * d */
 
6448
   if ((res = mp_mul(&q, d, &q)) != MP_OKAY) { 
 
6449
      goto ERR;
 
6450
   }
 
6451
   
 
6452
   /* a = a + q */
 
6453
   if ((res = s_mp_add(a, &q, a)) != MP_OKAY) {
 
6454
      goto ERR;
 
6455
   }
 
6456
   
 
6457
   if (mp_cmp_mag(a, n) != MP_LT) {
 
6458
      s_mp_sub(a, n, a);
 
6459
      goto top;
 
6460
   }
 
6461
   
 
6462
ERR:
 
6463
   mp_clear(&q);
 
6464
   return res;
 
6465
}
 
6466
 
 
6467
#endif
 
6468
 
 
6469
/* End: bn_mp_reduce_2k_l.c */
 
6470
 
6401
6471
/* Start: bn_mp_reduce_2k_setup.c */
6402
6472
#include <tommath.h>
6403
6473
#ifdef BN_MP_REDUCE_2K_SETUP_C
6417
6487
 */
6418
6488
 
6419
6489
/* determines the setup value */
6420
 
int 
6421
 
mp_reduce_2k_setup(mp_int *a, mp_digit *d)
 
6490
int mp_reduce_2k_setup(mp_int *a, mp_digit *d)
6422
6491
{
6423
6492
   int res, p;
6424
6493
   mp_int tmp;
6446
6515
 
6447
6516
/* End: bn_mp_reduce_2k_setup.c */
6448
6517
 
 
6518
/* Start: bn_mp_reduce_2k_setup_l.c */
 
6519
#include <tommath.h>
 
6520
#ifdef BN_MP_REDUCE_2K_SETUP_L_C
 
6521
/* LibTomMath, multiple-precision integer library -- Tom St Denis
 
6522
 *
 
6523
 * LibTomMath is a library that provides multiple-precision
 
6524
 * integer arithmetic as well as number theoretic functionality.
 
6525
 *
 
6526
 * The library was designed directly after the MPI library by
 
6527
 * Michael Fromberger but has been written from scratch with
 
6528
 * additional optimizations in place.
 
6529
 *
 
6530
 * The library is free for all purposes without any express
 
6531
 * guarantee it works.
 
6532
 *
 
6533
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
 
6534
 */
 
6535
 
 
6536
/* determines the setup value */
 
6537
int mp_reduce_2k_setup_l(mp_int *a, mp_int *d)
 
6538
{
 
6539
   int    res;
 
6540
   mp_int tmp;
 
6541
   
 
6542
   if ((res = mp_init(&tmp)) != MP_OKAY) {
 
6543
      return res;
 
6544
   }
 
6545
   
 
6546
   if ((res = mp_2expt(&tmp, mp_count_bits(a))) != MP_OKAY) {
 
6547
      goto ERR;
 
6548
   }
 
6549
   
 
6550
   if ((res = s_mp_sub(&tmp, a, d)) != MP_OKAY) {
 
6551
      goto ERR;
 
6552
   }
 
6553
   
 
6554
ERR:
 
6555
   mp_clear(&tmp);
 
6556
   return res;
 
6557
}
 
6558
#endif
 
6559
 
 
6560
/* End: bn_mp_reduce_2k_setup_l.c */
 
6561
 
6449
6562
/* Start: bn_mp_reduce_is_2k.c */
6450
6563
#include <tommath.h>
6451
6564
#ifdef BN_MP_REDUCE_IS_2K_C
6471
6584
   mp_digit iz;
6472
6585
   
6473
6586
   if (a->used == 0) {
6474
 
      return 0;
 
6587
      return MP_NO;
6475
6588
   } else if (a->used == 1) {
6476
 
      return 1;
 
6589
      return MP_YES;
6477
6590
   } else if (a->used > 1) {
6478
6591
      iy = mp_count_bits(a);
6479
6592
      iz = 1;
6482
6595
      /* Test every bit from the second digit up, must be 1 */
6483
6596
      for (ix = DIGIT_BIT; ix < iy; ix++) {
6484
6597
          if ((a->dp[iw] & iz) == 0) {
6485
 
             return 0;
 
6598
             return MP_NO;
6486
6599
          }
6487
6600
          iz <<= 1;
6488
6601
          if (iz > (mp_digit)MP_MASK) {
6491
6604
          }
6492
6605
      }
6493
6606
   }
6494
 
   return 1;
 
6607
   return MP_YES;
6495
6608
}
6496
6609
 
6497
6610
#endif
6498
6611
 
6499
6612
/* End: bn_mp_reduce_is_2k.c */
6500
6613
 
 
6614
/* Start: bn_mp_reduce_is_2k_l.c */
 
6615
#include <tommath.h>
 
6616
#ifdef BN_MP_REDUCE_IS_2K_L_C
 
6617
/* LibTomMath, multiple-precision integer library -- Tom St Denis
 
6618
 *
 
6619
 * LibTomMath is a library that provides multiple-precision
 
6620
 * integer arithmetic as well as number theoretic functionality.
 
6621
 *
 
6622
 * The library was designed directly after the MPI library by
 
6623
 * Michael Fromberger but has been written from scratch with
 
6624
 * additional optimizations in place.
 
6625
 *
 
6626
 * The library is free for all purposes without any express
 
6627
 * guarantee it works.
 
6628
 *
 
6629
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
 
6630
 */
 
6631
 
 
6632
/* determines if reduce_2k_l can be used */
 
6633
int mp_reduce_is_2k_l(mp_int *a)
 
6634
{
 
6635
   int ix, iy;
 
6636
   
 
6637
   if (a->used == 0) {
 
6638
      return MP_NO;
 
6639
   } else if (a->used == 1) {
 
6640
      return MP_YES;
 
6641
   } else if (a->used > 1) {
 
6642
      /* if more than half of the digits are -1 we're sold */
 
6643
      for (iy = ix = 0; ix < a->used; ix++) {
 
6644
          if (a->dp[ix] == MP_MASK) {
 
6645
              ++iy;
 
6646
          }
 
6647
      }
 
6648
      return (iy >= (a->used/2)) ? MP_YES : MP_NO;
 
6649
      
 
6650
   }
 
6651
   return MP_NO;
 
6652
}
 
6653
 
 
6654
#endif
 
6655
 
 
6656
/* End: bn_mp_reduce_is_2k_l.c */
 
6657
 
6501
6658
/* Start: bn_mp_reduce_setup.c */
6502
6659
#include <tommath.h>
6503
6660
#ifdef BN_MP_REDUCE_SETUP_C
7132
7289
 */
7133
7290
 
7134
7291
/* store in signed [big endian] format */
7135
 
int
7136
 
mp_to_signed_bin (mp_int * a, unsigned char *b)
 
7292
int mp_to_signed_bin (mp_int * a, unsigned char *b)
7137
7293
{
7138
7294
  int     res;
7139
7295
 
7147
7303
 
7148
7304
/* End: bn_mp_to_signed_bin.c */
7149
7305
 
 
7306
/* Start: bn_mp_to_signed_bin_n.c */
 
7307
#include <tommath.h>
 
7308
#ifdef BN_MP_TO_SIGNED_BIN_N_C
 
7309
/* LibTomMath, multiple-precision integer library -- Tom St Denis
 
7310
 *
 
7311
 * LibTomMath is a library that provides multiple-precision
 
7312
 * integer arithmetic as well as number theoretic functionality.
 
7313
 *
 
7314
 * The library was designed directly after the MPI library by
 
7315
 * Michael Fromberger but has been written from scratch with
 
7316
 * additional optimizations in place.
 
7317
 *
 
7318
 * The library is free for all purposes without any express
 
7319
 * guarantee it works.
 
7320
 *
 
7321
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
 
7322
 */
 
7323
 
 
7324
/* store in signed [big endian] format */
 
7325
int mp_to_signed_bin_n (mp_int * a, unsigned char *b, unsigned long *outlen)
 
7326
{
 
7327
   if (*outlen < (unsigned long)mp_signed_bin_size(a)) {
 
7328
      return MP_VAL;
 
7329
   }
 
7330
   *outlen = mp_signed_bin_size(a);
 
7331
   return mp_to_signed_bin(a, b);
 
7332
}
 
7333
#endif
 
7334
 
 
7335
/* End: bn_mp_to_signed_bin_n.c */
 
7336
 
7150
7337
/* Start: bn_mp_to_unsigned_bin.c */
7151
7338
#include <tommath.h>
7152
7339
#ifdef BN_MP_TO_UNSIGNED_BIN_C
7166
7353
 */
7167
7354
 
7168
7355
/* store in unsigned [big endian] format */
7169
 
int
7170
 
mp_to_unsigned_bin (mp_int * a, unsigned char *b)
 
7356
int mp_to_unsigned_bin (mp_int * a, unsigned char *b)
7171
7357
{
7172
7358
  int     x, res;
7173
7359
  mp_int  t;
7196
7382
 
7197
7383
/* End: bn_mp_to_unsigned_bin.c */
7198
7384
 
 
7385
/* Start: bn_mp_to_unsigned_bin_n.c */
 
7386
#include <tommath.h>
 
7387
#ifdef BN_MP_TO_UNSIGNED_BIN_N_C
 
7388
/* LibTomMath, multiple-precision integer library -- Tom St Denis
 
7389
 *
 
7390
 * LibTomMath is a library that provides multiple-precision
 
7391
 * integer arithmetic as well as number theoretic functionality.
 
7392
 *
 
7393
 * The library was designed directly after the MPI library by
 
7394
 * Michael Fromberger but has been written from scratch with
 
7395
 * additional optimizations in place.
 
7396
 *
 
7397
 * The library is free for all purposes without any express
 
7398
 * guarantee it works.
 
7399
 *
 
7400
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
 
7401
 */
 
7402
 
 
7403
/* store in unsigned [big endian] format */
 
7404
int mp_to_unsigned_bin_n (mp_int * a, unsigned char *b, unsigned long *outlen)
 
7405
{
 
7406
   if (*outlen < (unsigned long)mp_unsigned_bin_size(a)) {
 
7407
      return MP_VAL;
 
7408
   }
 
7409
   *outlen = mp_unsigned_bin_size(a);
 
7410
   return mp_to_unsigned_bin(a, b);
 
7411
}
 
7412
#endif
 
7413
 
 
7414
/* End: bn_mp_to_unsigned_bin_n.c */
 
7415
 
7199
7416
/* Start: bn_mp_toom_mul.c */
7200
7417
#include <tommath.h>
7201
7418
#ifdef BN_MP_TOOM_MUL_C
7216
7433
 
7217
7434
/* multiplication using the Toom-Cook 3-way algorithm 
7218
7435
 *
7219
 
 * Much more complicated than Karatsuba but has a lower asymptotic running time of 
7220
 
 * O(N**1.464).  This algorithm is only particularly useful on VERY large
7221
 
 * inputs (we're talking 1000s of digits here...).
 
7436
 * Much more complicated than Karatsuba but has a lower 
 
7437
 * asymptotic running time of O(N**1.464).  This algorithm is 
 
7438
 * only particularly useful on VERY large inputs 
 
7439
 * (we're talking 1000s of digits here...).
7222
7440
*/
7223
7441
int mp_toom_mul(mp_int *a, mp_int *b, mp_int *c)
7224
7442
{
7888
8106
 */
7889
8107
 
7890
8108
/* get the size for an unsigned equivalent */
7891
 
int
7892
 
mp_unsigned_bin_size (mp_int * a)
 
8109
int mp_unsigned_bin_size (mp_int * a)
7893
8110
{
7894
8111
  int     size = mp_count_bits (a);
7895
8112
  return (size / 8 + ((size & 7) != 0 ? 1 : 0));
7938
8155
  }
7939
8156
 
7940
8157
  for (ix = 0; ix < px; ix++) {
7941
 
 
 
8158
     t.dp[ix] ^= x->dp[ix];
7942
8159
  }
7943
8160
  mp_clamp (&t);
7944
8161
  mp_exch (c, &t);
7968
8185
 */
7969
8186
 
7970
8187
/* set to zero */
7971
 
void
7972
 
mp_zero (mp_int * a)
 
8188
void mp_zero (mp_int * a)
7973
8189
{
 
8190
  int       n;
 
8191
  mp_digit *tmp;
 
8192
 
7974
8193
  a->sign = MP_ZPOS;
7975
8194
  a->used = 0;
7976
 
  memset (a->dp, 0, sizeof (mp_digit) * a->alloc);
 
8195
 
 
8196
  tmp = a->dp;
 
8197
  for (n = 0; n < a->alloc; n++) {
 
8198
     *tmp++ = 0;
 
8199
  }
7977
8200
}
7978
8201
#endif
7979
8202
 
7996
8219
 *
7997
8220
 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
7998
8221
 */
7999
 
const mp_digit __prime_tab[] = {
 
8222
const mp_digit ltm_prime_tab[] = {
8000
8223
  0x0002, 0x0003, 0x0005, 0x0007, 0x000B, 0x000D, 0x0011, 0x0013,
8001
8224
  0x0017, 0x001D, 0x001F, 0x0025, 0x0029, 0x002B, 0x002F, 0x0035,
8002
8225
  0x003B, 0x003D, 0x0043, 0x0047, 0x0049, 0x004F, 0x0053, 0x0059,
8212
8435
   #define TAB_SIZE 256
8213
8436
#endif
8214
8437
 
8215
 
int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y)
 
8438
int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode)
8216
8439
{
8217
8440
  mp_int  M[TAB_SIZE], res, mu;
8218
8441
  mp_digit buf;
8219
8442
  int     err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize;
 
8443
  int (*redux)(mp_int*,mp_int*,mp_int*);
8220
8444
 
8221
8445
  /* find window size */
8222
8446
  x = mp_count_bits (X);
8261
8485
 
8262
8486
  /* create mu, used for Barrett reduction */
8263
8487
  if ((err = mp_init (&mu)) != MP_OKAY) {
8264
 
    goto __M;
8265
 
  }
8266
 
  if ((err = mp_reduce_setup (&mu, P)) != MP_OKAY) {
8267
 
    goto __MU;
8268
 
  }
 
8488
    goto LBL_M;
 
8489
  }
 
8490
  
 
8491
  if (redmode == 0) {
 
8492
     if ((err = mp_reduce_setup (&mu, P)) != MP_OKAY) {
 
8493
        goto LBL_MU;
 
8494
     }
 
8495
     redux = mp_reduce;
 
8496
  } else {
 
8497
     if ((err = mp_reduce_2k_setup_l (P, &mu)) != MP_OKAY) {
 
8498
        goto LBL_MU;
 
8499
     }
 
8500
     redux = mp_reduce_2k_l;
 
8501
  }    
8269
8502
 
8270
8503
  /* create M table
8271
8504
   *
8276
8509
   * computed though accept for M[0] and M[1]
8277
8510
   */
8278
8511
  if ((err = mp_mod (G, P, &M[1])) != MP_OKAY) {
8279
 
    goto __MU;
 
8512
    goto LBL_MU;
8280
8513
  }
8281
8514
 
8282
8515
  /* compute the value at M[1<<(winsize-1)] by squaring 
8283
8516
   * M[1] (winsize-1) times 
8284
8517
   */
8285
8518
  if ((err = mp_copy (&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) {
8286
 
    goto __MU;
 
8519
    goto LBL_MU;
8287
8520
  }
8288
8521
 
8289
8522
  for (x = 0; x < (winsize - 1); x++) {
 
8523
    /* square it */
8290
8524
    if ((err = mp_sqr (&M[1 << (winsize - 1)], 
8291
8525
                       &M[1 << (winsize - 1)])) != MP_OKAY) {
8292
 
      goto __MU;
 
8526
      goto LBL_MU;
8293
8527
    }
8294
 
    if ((err = mp_reduce (&M[1 << (winsize - 1)], P, &mu)) != MP_OKAY) {
8295
 
      goto __MU;
 
8528
 
 
8529
    /* reduce modulo P */
 
8530
    if ((err = redux (&M[1 << (winsize - 1)], P, &mu)) != MP_OKAY) {
 
8531
      goto LBL_MU;
8296
8532
    }
8297
8533
  }
8298
8534
 
8301
8537
   */
8302
8538
  for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) {
8303
8539
    if ((err = mp_mul (&M[x - 1], &M[1], &M[x])) != MP_OKAY) {
8304
 
      goto __MU;
 
8540
      goto LBL_MU;
8305
8541
    }
8306
 
    if ((err = mp_reduce (&M[x], P, &mu)) != MP_OKAY) {
8307
 
      goto __MU;
 
8542
    if ((err = redux (&M[x], P, &mu)) != MP_OKAY) {
 
8543
      goto LBL_MU;
8308
8544
    }
8309
8545
  }
8310
8546
 
8311
8547
  /* setup result */
8312
8548
  if ((err = mp_init (&res)) != MP_OKAY) {
8313
 
    goto __MU;
 
8549
    goto LBL_MU;
8314
8550
  }
8315
8551
  mp_set (&res, 1);
8316
8552
 
8350
8586
    /* if the bit is zero and mode == 1 then we square */
8351
8587
    if (mode == 1 && y == 0) {
8352
8588
      if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
8353
 
        goto __RES;
 
8589
        goto LBL_RES;
8354
8590
      }
8355
 
      if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) {
8356
 
        goto __RES;
 
8591
      if ((err = redux (&res, P, &mu)) != MP_OKAY) {
 
8592
        goto LBL_RES;
8357
8593
      }
8358
8594
      continue;
8359
8595
    }
8367
8603
      /* square first */
8368
8604
      for (x = 0; x < winsize; x++) {
8369
8605
        if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
8370
 
          goto __RES;
 
8606
          goto LBL_RES;
8371
8607
        }
8372
 
        if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) {
8373
 
          goto __RES;
 
8608
        if ((err = redux (&res, P, &mu)) != MP_OKAY) {
 
8609
          goto LBL_RES;
8374
8610
        }
8375
8611
      }
8376
8612
 
8377
8613
      /* then multiply */
8378
8614
      if ((err = mp_mul (&res, &M[bitbuf], &res)) != MP_OKAY) {
8379
 
        goto __RES;
 
8615
        goto LBL_RES;
8380
8616
      }
8381
 
      if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) {
8382
 
        goto __RES;
 
8617
      if ((err = redux (&res, P, &mu)) != MP_OKAY) {
 
8618
        goto LBL_RES;
8383
8619
      }
8384
8620
 
8385
8621
      /* empty window and reset */
8394
8630
    /* square then multiply if the bit is set */
8395
8631
    for (x = 0; x < bitcpy; x++) {
8396
8632
      if ((err = mp_sqr (&res, &res)) != MP_OKAY) {
8397
 
        goto __RES;
 
8633
        goto LBL_RES;
8398
8634
      }
8399
 
      if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) {
8400
 
        goto __RES;
 
8635
      if ((err = redux (&res, P, &mu)) != MP_OKAY) {
 
8636
        goto LBL_RES;
8401
8637
      }
8402
8638
 
8403
8639
      bitbuf <<= 1;
8404
8640
      if ((bitbuf & (1 << winsize)) != 0) {
8405
8641
        /* then multiply */
8406
8642
        if ((err = mp_mul (&res, &M[1], &res)) != MP_OKAY) {
8407
 
          goto __RES;
 
8643
          goto LBL_RES;
8408
8644
        }
8409
 
        if ((err = mp_reduce (&res, P, &mu)) != MP_OKAY) {
8410
 
          goto __RES;
 
8645
        if ((err = redux (&res, P, &mu)) != MP_OKAY) {
 
8646
          goto LBL_RES;
8411
8647
        }
8412
8648
      }
8413
8649
    }
8415
8651
 
8416
8652
  mp_exch (&res, Y);
8417
8653
  err = MP_OKAY;
8418
 
__RES:mp_clear (&res);
8419
 
__MU:mp_clear (&mu);
8420
 
__M:
 
8654
LBL_RES:mp_clear (&res);
 
8655
LBL_MU:mp_clear (&mu);
 
8656
LBL_M:
8421
8657
  mp_clear(&M[1]);
8422
8658
  for (x = 1<<(winsize-1); x < (1 << winsize); x++) {
8423
8659
    mp_clear (&M[x]);
8450
8686
 * HAC pp. 595, Algorithm 14.12  Modified so you can control how 
8451
8687
 * many digits of output are created.
8452
8688
 */
8453
 
int
8454
 
s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
 
8689
int s_mp_mul_digs (mp_int * a, mp_int * b, mp_int * c, int digs)
8455
8690
{
8456
8691
  mp_int  t;
8457
8692
  int     res, pa, pb, ix, iy;
8619
8854
 */
8620
8855
 
8621
8856
/* low level squaring, b = a*a, HAC pp.596-597, Algorithm 14.16 */
8622
 
int
8623
 
s_mp_sqr (mp_int * a, mp_int * b)
 
8857
int s_mp_sqr (mp_int * a, mp_int * b)
8624
8858
{
8625
8859
  mp_int  t;
8626
8860
  int     res, ix, iy, pa;
8797
9031
 CPU                    /Compiler     /MUL CUTOFF/SQR CUTOFF
8798
9032
-------------------------------------------------------------
8799
9033
 Intel P4 Northwood     /GCC v3.4.1   /        88/       128/LTM 0.32 ;-)
 
9034
 AMD Athlon64           /GCC v3.4.4   /        74/       124/LTM 0.34
8800
9035
 
8801
9036
*/
8802
9037
 
8803
 
int     KARATSUBA_MUL_CUTOFF = 88,      /* Min. number of digits before Karatsuba multiplication is used. */
8804
 
        KARATSUBA_SQR_CUTOFF = 128,     /* Min. number of digits before Karatsuba squaring is used. */
 
9038
int     KARATSUBA_MUL_CUTOFF = 74,      /* Min. number of digits before Karatsuba multiplication is used. */
 
9039
        KARATSUBA_SQR_CUTOFF = 124,     /* Min. number of digits before Karatsuba squaring is used. */
8805
9040
        
8806
9041
        TOOM_MUL_CUTOFF      = 350,      /* no optimal values of these are known yet so set em high */
8807
9042
        TOOM_SQR_CUTOFF      = 400;