userapps/opensource/sshd/libtommath/bn_mp_reduce.c

   1 #include <tommath.h>
   2 #ifdef BN_MP_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@iahu.ca, http://math.libtomcrypt.org
  16  */
  17
  18 /* reduces x mod m, assumes 0 < x < m**2, mu is
  19  * precomputed via mp_reduce_setup.
  20  * From HAC pp.604 Algorithm 14.42
  21  */
  22 int mp_reduce (mp_int * x, mp_int * m, mp_int * mu)
  23 {
  24   mp_int  q;
  25   int     res, um = m->used;
  26
  27   /* q = x */
  28   if ((res = mp_init_copy (&q, x)) != MP_OKAY) {
  29     return res;
  30   }
  31
  32   /* q1 = x / b**(k-1)  */
  33   mp_rshd (&q, um - 1);
  34
  35   /* according to HAC this optimization is ok */
  36   if (((unsigned long) um) > (((mp_digit)1) << (DIGIT_BIT - 1))) {
  37     if ((res = mp_mul (&q, mu, &q)) != MP_OKAY) {
  38       goto CLEANUP;
  39     }
  40   } else {
  41 #ifdef BN_S_MP_MUL_HIGH_DIGS_C
  42     if ((res = s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) {
  43       goto CLEANUP;
  44     }
  45 #elif defined(BN_FAST_S_MP_MUL_HIGH_DIGS_C)
  46     if ((res = fast_s_mp_mul_high_digs (&q, mu, &q, um)) != MP_OKAY) {
  47       goto CLEANUP;
  48     }
  49 #else
  50     {
  51       res = MP_VAL;
  52       goto CLEANUP;
  53     }
  54 #endif
  55   }
  56
  57   /* q3 = q2 / b**(k+1) */
  58   mp_rshd (&q, um + 1);
  59
  60   /* x = x mod b**(k+1), quick (no division) */
  61   if ((res = mp_mod_2d (x, DIGIT_BIT * (um + 1), x)) != MP_OKAY) {
  62     goto CLEANUP;
  63   }
  64
  65   /* q = q * m mod b**(k+1), quick (no division) */
  66   if ((res = s_mp_mul_digs (&q, m, &q, um + 1)) != MP_OKAY) {
  67     goto CLEANUP;
  68   }
  69
  70   /* x = x - q */
  71   if ((res = mp_sub (x, &q, x)) != MP_OKAY) {
  72     goto CLEANUP;
  73   }
  74
  75   /* If x < 0, add b**(k+1) to it */
  76   if (mp_cmp_d (x, 0) == MP_LT) {
  77     mp_set (&q, 1);
  78     if ((res = mp_lshd (&q, um + 1)) != MP_OKAY)
  79       goto CLEANUP;
  80     if ((res = mp_add (x, &q, x)) != MP_OKAY)
  81       goto CLEANUP;
  82   }
  83
  84   /* Back off if it's too big */
  85   while (mp_cmp (x, m) != MP_LT) {
  86     if ((res = s_mp_sub (x, m, x)) != MP_OKAY) {
  87       goto CLEANUP;
  88     }
  89   }
  90
  91 CLEANUP:
  92   mp_clear (&q);
  93
  94   return res;
  95 }
  96 #endif