Line data    Source code

       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@gmail.com, http://libtom.org
      16             :  */
      17             : 
      18             : /* computes xR**-1 == x (mod N) via Montgomery Reduction */
      19             : int
      20           0 : 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           0 :   digs = n->used * 2 + 1;
      32           0 :   if ((digs < MP_WARRAY) &&
      33           0 :       n->used <
      34             :       (1 << ((CHAR_BIT * sizeof (mp_word)) - (2 * DIGIT_BIT)))) {
      35           0 :     return fast_mp_montgomery_reduce (x, n, rho);
      36             :   }
      37             : 
      38             :   /* grow the input as required */
      39           0 :   if (x->alloc < digs) {
      40           0 :     if ((res = mp_grow (x, digs)) != MP_OKAY) {
      41           0 :       return res;
      42             :     }
      43             :   }
      44           0 :   x->used = digs;
      45             : 
      46           0 :   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           0 :     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           0 :       tmpn = n->dp;
      65             : 
      66             :       /* alias for the digits of x [the input] */
      67           0 :       tmpx = x->dp + ix;
      68             : 
      69             :       /* set the carry to zero */
      70           0 :       u = 0;
      71             : 
      72             :       /* Multiply and add in place */
      73           0 :       for (iy = 0; iy < n->used; iy++) {
      74             :         /* compute product and sum */
      75           0 :         r       = ((mp_word)mu) * ((mp_word)*tmpn++) +
      76           0 :                   ((mp_word) u) + ((mp_word) * tmpx);
      77             : 
      78             :         /* get carry */
      79           0 :         u       = (mp_digit)(r >> ((mp_word) DIGIT_BIT));
      80             : 
      81             :         /* fix digit */
      82           0 :         *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           0 :       while (u) {
      89           0 :         *tmpx   += u;
      90           0 :         u        = *tmpx >> DIGIT_BIT;
      91           0 :         *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           0 :   mp_clamp(x);
     105           0 :   mp_rshd (x, n->used);
     106             : 
     107             :   /* if x >= n then x = x - n */
     108           0 :   if (mp_cmp_mag (x, n) != MP_LT) {
     109           0 :     return s_mp_sub (x, n, x);
     110             :   }
     111             : 
     112           0 :   return MP_OKAY;
     113             : }
     114             : #endif
     115             : 
     116             : /* $Source: /cvs/libtom/libtommath/bn_mp_montgomery_reduce.c,v $ */
     117             : /* $Revision: 1.4 $ */
     118             : /* $Date: 2006/12/28 01:25:13 $ */