s_mp_mul_high_comba.c

   1 #include "tommath_private.h"
   2 #ifdef S_MP_MUL_HIGH_COMBA_C
   3 /* LibTomMath, multiple-precision integer library -- Tom St Denis */
   4 /* SPDX-License-Identifier: Unlicense */
   5
   6 /* this is a modified version of s_mp_mul_comba that only produces
   7  * output digits *above* digs.  See the comments for s_mp_mul_comba
   8  * to see how it works.
   9  *
  10  * This is used in the Barrett reduction since for one of the multiplications
  11  * only the higher digits were needed.  This essentially halves the work.
  12  *
  13  * Based on Algorithm 14.12 on pp.595 of HAC.
  14  */
  15 mp_err s_mp_mul_high_comba(const mp_int *a, const mp_int *b, mp_int *c, int digs)
  16 {
  17    int     oldused, pa, ix;
  18    mp_err   err;
  19    mp_digit MP_ALLOC_WARRAY(W);
  20    mp_word  _W;
  21
  22    MP_CHECK_WARRAY(W);
  23
  24    if (digs < 0) {
  25       MP_FREE_WARRAY(W);
  26       return MP_VAL;
  27    }
  28
  29    /* grow the destination as required */
  30    pa = a->used + b->used;
  31    if ((err = mp_grow(c, pa)) != MP_OKAY) {
  32       MP_FREE_WARRAY(W);
  33       return err;
  34    }
  35
  36    /* number of output digits to produce */
  37    pa = a->used + b->used;
  38    _W = 0;
  39    for (ix = digs; ix < pa; ix++) {
  40       int      tx, ty, iy, iz;
  41
  42       /* get offsets into the two bignums */
  43       ty = MP_MIN(b->used-1, ix);
  44       tx = ix - ty;
  45
  46       /* this is the number of times the loop will iterate, essentially its
  47          while (tx++ < a->used && ty-- >= 0) { ... }
  48        */
  49       iy = MP_MIN(a->used-tx, ty+1);
  50
  51       /* execute loop */
  52       for (iz = 0; iz < iy; iz++) {
  53          _W += (mp_word)a->dp[tx + iz] * (mp_word)b->dp[ty - iz];
  54       }
  55
  56       /* store term */
  57       W[ix] = (mp_digit)_W & MP_MASK;
  58
  59       /* make next carry */
  60       _W = _W >> (mp_word)MP_DIGIT_BIT;
  61    }
  62
  63    /* setup dest */
  64    oldused  = c->used;
  65    c->used = pa;
  66
  67    for (ix = digs; ix < pa; ix++) {
  68       /* now extract the previous digit [below the carry] */
  69       c->dp[ix] = W[ix];
  70    }
  71
  72    /* clear unused digits [that existed in the old copy of c] */
  73    s_mp_zero_digs(c->dp + c->used, oldused - c->used);
  74
  75    mp_clamp(c);
  76    MP_FREE_WARRAY(W);
  77    return MP_OKAY;
  78 }
  79 #endif