chiark - git - mdw - catacomb/blob - mpbarrett.c

   1 /* -*-c-*-
   2  *
   3  * $Id: mpbarrett.c,v 1.4 2000/06/17 11:45:09 mdw Exp $
   4  *
   5  * Barrett modular reduction
   6  *
   7  * (c) 1999 Straylight/Edgeware
   8  */
   9
  10 /*----- Licensing notice --------------------------------------------------*
  11  *
  12  * This file is part of Catacomb.
  13  *
  14  * Catacomb is free software; you can redistribute it and/or modify
  15  * it under the terms of the GNU Library General Public License as
  16  * published by the Free Software Foundation; either version 2 of the
  17  * License, or (at your option) any later version.
  18  *
  19  * Catacomb is distributed in the hope that it will be useful,
  20  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  22  * GNU Library General Public License for more details.
  23  *
  24  * You should have received a copy of the GNU Library General Public
  25  * License along with Catacomb; if not, write to the Free
  26  * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
  27  * MA 02111-1307, USA.
  28  */
  29
  30 /*----- Revision history --------------------------------------------------*
  31  *
  32  * $Log: mpbarrett.c,v $
  33  * Revision 1.4  2000/06/17 11:45:09  mdw
  34  * Major memory management overhaul.  Added arena support.  Use the secure
  35  * arena for secret integers.  Replace and improve the MP management macros
  36  * (e.g., replace MP_MODIFY by MP_DEST).
  37  *
  38  * Revision 1.3  1999/12/12 15:08:52  mdw
  39  * Don't bother shifting %$q$% in @mpbarrett_reduce@, just skip the least
  40  * significant digits.
  41  *
  42  * Revision 1.2  1999/12/11 01:50:56  mdw
  43  * Improve initialization slightly.
  44  *
  45  * Revision 1.1  1999/12/10 23:21:59  mdw
  46  * Barrett reduction support: works with even moduli.
  47  *
  48  */
  49
  50 /*----- Header files ------------------------------------------------------*/
  51
  52 #include "mp.h"
  53 #include "mpbarrett.h"
  54
  55 /*----- Main code ---------------------------------------------------------*/
  56
  57 /* --- @mpbarrett_create@ --- *
  58  *
  59  * Arguments:   @mpbarrett *mb@ = pointer to Barrett reduction context
  60  *              @mp *m@ = modulus to work to
  61  *
  62  *
  63  * Returns:     ---
  64  *
  65  * Use:         Initializes a Barrett reduction context ready for use.
  66  */
  67
  68 void mpbarrett_create(mpbarrett *mb, mp *m)
  69 {
  70   mp *b;
  71
  72   /* --- Validate the arguments --- */
  73
  74   assert(((void)"Barrett modulus must be positive", (m->f & MP_NEG) == 0));
  75
  76   /* --- Compute %$\mu$% --- */
  77
  78   mp_shrink(m);
  79   mb->k = MP_LEN(m);
  80   mb->m = MP_COPY(m);
  81   b = mp_new(2 * mb->k + 1, 0);
  82   MPX_ZERO(b->v, b->vl - 1);
  83   b->vl[-1] = 1;
  84   mp_div(&b, 0, b, m);
  85   mb->mu = b;
  86 }
  87
  88 /* --- @mpbarrett_destroy@ --- *
  89  *
  90  * Arguments:   @mpbarrett *mb@ = pointer to Barrett reduction context
  91  *
  92  * Returns:     ---
  93  *
  94  * Use:         Destroys a Barrett reduction context releasing any resources
  95  *              claimed.
  96  */
  97
  98 void mpbarrett_destroy(mpbarrett *mb)
  99 {
 100   mp_drop(mb->m);
 101   mp_drop(mb->mu);
 102 }
 103
 104 /* --- @mpbarrett_reduce@ --- *
 105  *
 106  * Arguments:   @mpbarrett *mb@ = pointer to Barrett reduction context
 107  *              @mp *d@ = destination for result
 108  *              @mp *m@ = number to reduce
 109  *
 110  * Returns:     The residue of @m@ modulo the number in the reduction
 111  *              context.
 112  *
 113  * Use:         Performs an efficient modular reduction.  The argument is
 114  *              assumed to be positive.
 115  */
 116
 117 mp *mpbarrett_reduce(mpbarrett *mb, mp *d, mp *m)
 118 {
 119   mp *q;
 120   size_t k = mb->k;
 121
 122   /* --- Special case if @m@ is too small --- */
 123
 124   if (MP_LEN(m) < k) {
 125     m = MP_COPY(m);
 126     if (d)
 127       MP_DROP(d);
 128     return (m);
 129   }
 130
 131   /* --- First stage --- */
 132
 133   {
 134     mp qq;
 135     mp_build(&qq, m->v + (k - 1), m->vl);
 136     q = mp_mul(MP_NEW, &qq, mb->mu);
 137     if (MP_LEN(q) <= k) {
 138       m = MP_COPY(m);
 139       if (d)
 140         MP_DROP(d);
 141       return (m);
 142     }
 143   }
 144
 145   /* --- Second stage --- */
 146
 147   {
 148     mp *r;
 149     mpw *mvl;
 150
 151     MP_COPY(m);
 152     if (MP_LEN(m) <= k + 1)
 153       mvl = m->vl;
 154     else
 155       mvl = m->v + k + 1;
 156     r = mp_new(k + 1, (q->f | mb->m->f) & MP_BURN);
 157     mpx_umul(r->v, r->vl, q->v + k + 1, q->vl, mb->m->v, mb->m->vl);
 158     MP_DEST(d, k + 1, r->f);
 159     mpx_usub(d->v, d->vl, m->v, mvl, r->v, r->vl);
 160     d->f = (m->f | r->f) & MP_BURN;
 161     MP_DROP(r);
 162     MP_DROP(q);
 163     MP_DROP(m);
 164   }
 165
 166   /* --- Final stage --- */
 167
 168   MP_SHRINK(d);
 169   while (MPX_UCMP(d->v, d->vl, >=, mb->m->v, mb->m->vl))
 170     mpx_usub(d->v, d->vl, d->v, d->vl, mb->m->v, mb->m->vl);
 171
 172   MP_SHRINK(d);
 173   return (d);
 174 }
 175
 176 /* --- @mpbarrett_exp@ --- *
 177  *
 178  * Arguments:   @mpbarrett *mb@ = pointer to Barrett reduction context
 179  *              @mp *d@ = fake destination
 180  *              @mp *a@ = base
 181  *              @mp *e@ = exponent
 182  *
 183  * Returns:     Result, %$a^e \bmod m$%.
 184  */
 185
 186 mp *mpbarrett_exp(mpbarrett *mb, mp *d, mp *a, mp *e)
 187 {
 188   mpscan sc;
 189   mp *x = MP_ONE;
 190   mp *spare = MP_NEW;
 191
 192   a = MP_COPY(a);
 193   mp_scan(&sc, e);
 194   if (e->f & MP_BURN)
 195     spare = MP_NEWSEC;
 196   if (MP_STEP(&sc)) {
 197     size_t sq = 0;
 198     for (;;) {
 199       mp *dd;
 200       if (MP_BIT(&sc)) {
 201         while (sq) {
 202           dd = mp_sqr(spare, a);
 203           dd = mpbarrett_reduce(mb, dd, dd);
 204           spare = a; a = dd;
 205           sq--;
 206         }
 207         dd = mp_mul(spare, x, a);
 208         dd = mpbarrett_reduce(mb, dd, dd);
 209         spare = x; x = dd;
 210       }
 211       sq++;
 212       if (!MP_STEP(&sc))
 213         break;
 214     }
 215   }
 216
 217   MP_DROP(a);
 218   if (spare != MP_NEW)
 219     MP_DROP(spare);
 220   if (d != MP_NEW)
 221     MP_DROP(d);
 222   return (x);
 223 }
 224
 225 /*----- Test rig ----------------------------------------------------------*/
 226
 227 #ifdef TEST_RIG
 228
 229 static int vmod(dstr *v)
 230 {
 231   mp *x = *(mp **)v[0].buf;
 232   mp *n = *(mp **)v[1].buf;
 233   mp *r = *(mp **)v[2].buf;
 234   mp *s;
 235   mpbarrett mb;
 236   int ok = 1;
 237
 238   mpbarrett_create(&mb, n);
 239   s = mpbarrett_reduce(&mb, MP_NEW, x);
 240
 241   if (MP_CMP(s, !=, r)) {
 242     fputs("\n*** barrett reduction failure\n", stderr);
 243     fputs("x = ", stderr); mp_writefile(x, stderr, 10); fputc('\n', stderr);
 244     fputs("n = ", stderr); mp_writefile(n, stderr, 10); fputc('\n', stderr);
 245     fputs("r = ", stderr); mp_writefile(r, stderr, 10); fputc('\n', stderr);
 246     fputs("s = ", stderr); mp_writefile(s, stderr, 10); fputc('\n', stderr);
 247     ok = 0;
 248   }
 249
 250   mpbarrett_destroy(&mb);
 251   mp_drop(x);
 252   mp_drop(n);
 253   mp_drop(r);
 254   mp_drop(s);
 255   assert(mparena_count(MPARENA_GLOBAL) == 0);
 256   return (ok);
 257 }
 258
 259 static int vexp(dstr *v)
 260 {
 261   mp *m = *(mp **)v[0].buf;
 262   mp *a = *(mp **)v[1].buf;
 263   mp *b = *(mp **)v[2].buf;
 264   mp *r = *(mp **)v[3].buf;
 265   mp *mr;
 266   int ok = 1;
 267
 268   mpbarrett mb;
 269   mpbarrett_create(&mb, m);
 270
 271   mr = mpbarrett_exp(&mb, MP_NEW, a, b);
 272
 273   if (MP_CMP(mr, !=, r)) {
 274     fputs("\n*** barrett modexp failed", stderr);
 275     fputs("\n m = ", stderr); mp_writefile(m, stderr, 10);
 276     fputs("\n a = ", stderr); mp_writefile(a, stderr, 10);
 277     fputs("\n e = ", stderr); mp_writefile(b, stderr, 10);
 278     fputs("\n r = ", stderr); mp_writefile(r, stderr, 10);
 279     fputs("\nmr = ", stderr); mp_writefile(mr, stderr, 10);
 280     fputc('\n', stderr);
 281     ok = 0;
 282   }
 283
 284   mp_drop(m);
 285   mp_drop(a);
 286   mp_drop(b);
 287   mp_drop(r);
 288   mp_drop(mr);
 289   mpbarrett_destroy(&mb);
 290   assert(mparena_count(MPARENA_GLOBAL) == 0);
 291   return ok;
 292 }
 293
 294 static test_chunk tests[] = {
 295   { "mpbarrett-reduce", vmod, { &type_mp, &type_mp, &type_mp, 0 } },
 296   { "mpbarrett-exp", vexp, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } },
 297   { 0, 0, { 0 } }
 298 };
 299
 300 int main(int argc, char *argv[])
 301 {
 302   sub_init();
 303   test_run(argc, argv, tests, SRCDIR "/tests/mpbarrett");
 304   return (0);
 305 }
 306
 307 #endif
 308
 309 /*----- That's all, folks -------------------------------------------------*/