chiark - git - mdw - catacomb/blob - mpmont-mexp.c

   1 /* -*-c-*-
   2  *
   3  * $Id$
   4  *
   5  * Multiple simultaneous exponentiations
   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 /*----- Header files ------------------------------------------------------*/
  31
  32 #include "mp.h"
  33 #include "mpmont.h"
  34
  35 #define EXP_WINSZ 3
  36 #include "mpmont-exp.h"
  37
  38 /*----- Main code ---------------------------------------------------------*/
  39
  40 /* --- @mpmont_mexpr@ --- *
  41  *
  42  * Arguments:   @mpmont *mm@ = pointer to Montgomery reduction context
  43  *              @mp *d@ = fake destination
  44  *              @const mp_expfactor *f@ = pointer to array of factors
  45  *              @size_t n@ = number of factors supplied
  46  *
  47  * Returns:     If the bases are %$g_0, g_1, \ldots, g_{n-1}$% and the
  48  *              exponents are %$e_0, e_1, \ldots, e_{n-1}$% then the result
  49  *              is:
  50  *
  51  *              %$g_0^{e_0} g_1^{e_1} \ldots g_{n-1}^{e_{n-1}} \bmod m$%
  52  *
  53  *              except that the %$g_i$% and result are in Montgomery form.
  54  */
  55
  56 static mp *mexpr(mpmont *mm, mp *d, mp_expfactor *f, size_t n)
  57 {
  58   mp *a = MP_COPY(mm->r);
  59   mp *spare = MP_NEW;
  60   size_t i;
  61
  62   for (i = 0; i < n; i++) {
  63     mp *t;
  64     if (f[i].exp->f & MP_BURN)
  65       spare = MP_NEWSEC;
  66     if (MP_NEGP(f[i].exp)) {
  67       t = mpmont_reduce(mm, f[i].base, f[i].base);
  68       t = mp_modinv(t, t, mm->m);
  69       f[i].base = mpmont_mul(mm, t, t, mm->r2);
  70     }
  71   }
  72   EXP_SIMUL(a, f, n);
  73   mp_drop(d);
  74   mp_drop(spare);
  75   for (i = 0; i < n; i++)
  76     MP_DROP(f[i].base);
  77   xfree(f);
  78   return (a);
  79 }
  80
  81 mp *mpmont_mexpr(mpmont *mm, mp *d, const mp_expfactor *f, size_t n)
  82 {
  83   mp_expfactor *ff = xmalloc(n * sizeof(mp_expfactor));
  84   size_t i;
  85
  86   for (i = 0; i < n; i++) {
  87     ff[i].base = MP_COPY(f[i].base);
  88     ff[i].exp = f[i].exp;
  89   }
  90   return (mexpr(mm, d, ff, n));
  91 }
  92
  93 /* --- @mpmont_mexp@ --- *
  94  *
  95  * Arguments:   @mpmont *mm@ = pointer to Montgomery reduction context
  96  *              @mp *d@ = fake destination
  97  *              @const mp_expfactor *f@ = pointer to array of factors
  98  *              @size_t n@ = number of factors supplied
  99  *
 100  * Returns:     Product of bases raised to exponents, all mod @m@.
 101  *
 102  * Use:         Convenient interface over @mpmont_mexpr@.
 103  */
 104
 105 mp *mpmont_mexp(mpmont *mm, mp *d, const mp_expfactor *f, size_t n)
 106 {
 107   mp_expfactor *ff = xmalloc(n * sizeof(mp_expfactor));
 108   size_t i;
 109
 110   for (i = 0; i < n; i++) {
 111     ff[i].base = mpmont_mul(mm, MP_NEW, f[i].base, mm->r2);
 112     ff[i].exp = f[i].exp;
 113   }
 114   d = mexpr(mm, d, ff, n);
 115   return (mpmont_reduce(mm, d, d));
 116 }
 117
 118 /*----- Test rig ----------------------------------------------------------*/
 119
 120 #ifdef TEST_RIG
 121
 122 #include <mLib/testrig.h>
 123
 124 static int verify(size_t n, dstr *v)
 125 {
 126   mp *m = *(mp **)v[0].buf;
 127   mp_expfactor *f = xmalloc(n * sizeof(*f));
 128   mp *r, *rr;
 129   size_t i, j;
 130   mpmont mm;
 131   int ok = 1;
 132
 133   j = 1;
 134   for (i = 0; i < n; i++) {
 135     f[i].base = *(mp **)v[j++].buf;
 136     f[i].exp = *(mp **)v[j++].buf;
 137   }
 138
 139   rr = *(mp **)v[j].buf;
 140   mpmont_create(&mm, m);
 141   r = mpmont_mexp(&mm, MP_NEW, f, n);
 142   if (!MP_EQ(r, rr)) {
 143     fputs("\n*** mexp failed\n", stderr);
 144     fputs("m = ", stderr); mp_writefile(m, stderr, 10);
 145     for (i = 0; i < n; i++) {
 146       fprintf(stderr, "\ng_%u = ", i);
 147       mp_writefile(f[i].base, stderr, 10);
 148       fprintf(stderr, "\ne_%u = ", i);
 149       mp_writefile(f[i].exp, stderr, 10);
 150     }
 151     fputs("\nr = ", stderr); mp_writefile(r, stderr, 10);
 152     fputs("\nR = ", stderr); mp_writefile(rr, stderr, 10);
 153     fputc('\n', stderr);
 154     ok = 0;
 155   }
 156
 157   for (i = 0; i < n; i++) {
 158     MP_DROP(f[i].base);
 159     MP_DROP(f[i].exp);
 160   }
 161   MP_DROP(m);
 162   MP_DROP(r);
 163   MP_DROP(rr);
 164   mpmont_destroy(&mm);
 165   assert(mparena_count(MPARENA_GLOBAL) == 0);
 166   return (ok);
 167 }
 168
 169 static int t1(dstr *v) { return verify(1, v); }
 170 static int t2(dstr *v) { return verify(2, v); }
 171 static int t3(dstr *v) { return verify(3, v); }
 172 static int t4(dstr *v) { return verify(4, v); }
 173 static int t5(dstr *v) { return verify(5, v); }
 174
 175 static test_chunk tests[] = {
 176   { "mexp-1", t1, { &type_mp,
 177                     &type_mp, &type_mp,
 178                     &type_mp, 0 } },
 179   { "mexp-2", t2, { &type_mp,
 180                     &type_mp, &type_mp,
 181                     &type_mp, &type_mp,
 182                     &type_mp, 0 } },
 183   { "mexp-3", t3, { &type_mp,
 184                     &type_mp, &type_mp,
 185                     &type_mp, &type_mp,
 186                     &type_mp, &type_mp,
 187                     &type_mp, 0 } },
 188   { "mexp-4", t4, { &type_mp,
 189                     &type_mp, &type_mp,
 190                     &type_mp, &type_mp,
 191                     &type_mp, &type_mp,
 192                     &type_mp, &type_mp,
 193                     &type_mp, 0 } },
 194   { "mexp-5", t5, { &type_mp,
 195                     &type_mp, &type_mp,
 196                     &type_mp, &type_mp,
 197                     &type_mp, &type_mp,
 198                     &type_mp, &type_mp,
 199                     &type_mp, &type_mp,
 200                     &type_mp, 0 } },
 201   { 0, 0, { 0 } }
 202 };
 203
 204 int main(int argc, char *argv[])
 205 {
 206   sub_init();
 207   test_run(argc, argv, tests, SRCDIR "/tests/mpmont");
 208   return (0);
 209 }
 210
 211 #endif
 212
 213 /*----- That's all, folks -------------------------------------------------*/