Home | History | Annotate | Download | only in libtommath
      1 #include <tommath.h>
      2 #ifdef BN_MP_MONTGOMERY_REDUCE_C
      3 /* LibTomMath, multiple-precision integer library -- Tom St Denis
      4  *
      5  * LibTomMath is a library that provides multiple-precision
      6  * integer arithmetic as well as number theoretic functionality.
      7  *
      8  * The library was designed directly after the MPI library by
      9  * Michael Fromberger but has been written from scratch with
     10  * additional optimizations in place.
     11  *
     12  * The library is free for all purposes without any express
     13  * guarantee it works.
     14  *
     15  * Tom St Denis, tomstdenis (at) gmail.com, http://math.libtomcrypt.com
     16  */
     17 
     18 /* computes xR**-1 == x (mod N) via Montgomery Reduction */
     19 int
     20 mp_montgomery_reduce (mp_int * x, mp_int * n, mp_digit rho)
     21 {
     22   int     ix, res, digs;
     23   mp_digit mu;
     24 
     25   /* can the fast reduction [comba] method be used?
     26    *
     27    * Note that unlike in mul you're safely allowed *less*
     28    * than the available columns [255 per default] since carries
     29    * are fixed up in the inner loop.
     30    */
     31   digs = n->used * 2 + 1;
     32   if ((digs < MP_WARRAY) &&
     33       n->used <
     34       (1 << ((CHAR_BIT * sizeof (mp_word)) - (2 * DIGIT_BIT)))) {
     35     return fast_mp_montgomery_reduce (x, n, rho);
     36   }
     37 
     38   /* grow the input as required */
     39   if (x->alloc < digs) {
     40     if ((res = mp_grow (x, digs)) != MP_OKAY) {
     41       return res;
     42     }
     43   }
     44   x->used = digs;
     45 
     46   for (ix = 0; ix < n->used; ix++) {
     47     /* mu = ai * rho mod b
     48      *
     49      * The value of rho must be precalculated via
     50      * montgomery_setup() such that
     51      * it equals -1/n0 mod b this allows the
     52      * following inner loop to reduce the
     53      * input one digit at a time
     54      */
     55     mu = (mp_digit) (((mp_word)x->dp[ix]) * ((mp_word)rho) & MP_MASK);
     56 
     57     /* a = a + mu * m * b**i */
     58     {
     59       register int iy;
     60       register mp_digit *tmpn, *tmpx, u;
     61       register mp_word r;
     62 
     63       /* alias for digits of the modulus */
     64       tmpn = n->dp;
     65 
     66       /* alias for the digits of x [the input] */
     67       tmpx = x->dp + ix;
     68 
     69       /* set the carry to zero */
     70       u = 0;
     71 
     72       /* Multiply and add in place */
     73       for (iy = 0; iy < n->used; iy++) {
     74         /* compute product and sum */
     75         r       = ((mp_word)mu) * ((mp_word)*tmpn++) +
     76                   ((mp_word) u) + ((mp_word) * tmpx);
     77 
     78         /* get carry */
     79         u       = (mp_digit)(r >> ((mp_word) DIGIT_BIT));
     80 
     81         /* fix digit */
     82         *tmpx++ = (mp_digit)(r & ((mp_word) MP_MASK));
     83       }
     84       /* At this point the ix'th digit of x should be zero */
     85 
     86 
     87       /* propagate carries upwards as required*/
     88       while (u) {
     89         *tmpx   += u;
     90         u        = *tmpx >> DIGIT_BIT;
     91         *tmpx++ &= MP_MASK;
     92       }
     93     }
     94   }
     95 
     96   /* at this point the n.used'th least
     97    * significant digits of x are all zero
     98    * which means we can shift x to the
     99    * right by n.used digits and the
    100    * residue is unchanged.
    101    */
    102 
    103   /* x = x/b**n.used */
    104   mp_clamp(x);
    105   mp_rshd (x, n->used);
    106 
    107   /* if x >= n then x = x - n */
    108   if (mp_cmp_mag (x, n) != MP_LT) {
    109     return s_mp_sub (x, n, x);
    110   }
    111 
    112   return MP_OKAY;
    113 }
    114 #endif
    115 
    116 /* $Source: /cvs/libtom/libtommath/bn_mp_montgomery_reduce.c,v $ */
    117 /* $Revision: 1.3 $ */
    118 /* $Date: 2006/03/31 14:18:44 $ */
    119