New generic exponentation code. Includes sliding-window simultaneous
authormdw <mdw>
Sat, 16 Jun 2001 13:00:59 +0000 (13:00 +0000)
committermdw <mdw>
Sat, 16 Jun 2001 13:00:59 +0000 (13:00 +0000)
exponentiation.

exp.c [new file with mode: 0644]
exp.h [new file with mode: 0644]

diff --git a/exp.c b/exp.c
new file mode 100644 (file)
index 0000000..ed11485
--- /dev/null
+++ b/exp.c
@@ -0,0 +1,80 @@
+/* -*-c-*-
+ *
+ * $Id: exp.c,v 1.1 2001/06/16 13:00:59 mdw Exp $
+ *
+ * Generalized exponentiation
+ *
+ * (c) 2001 Straylight/Edgeware
+ */
+
+/*----- Licensing notice --------------------------------------------------* 
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ * 
+ * Catacomb is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Library General Public License for more details.
+ * 
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+/*----- Revision history --------------------------------------------------* 
+ *
+ * $Log: exp.c,v $
+ * Revision 1.1  2001/06/16 13:00:59  mdw
+ * New generic exponentation code.  Includes sliding-window simultaneous
+ * exponentiation.
+ *
+ */
+
+/*----- Header files ------------------------------------------------------*/
+
+#define EXP_TYPE /* Hack */
+#include "exp.h"
+
+/*----- Main code ---------------------------------------------------------*/
+
+/* --- @exp_simulnext@ --- *
+ *
+ * Arguments:  @exp_simul *e@ = pointer to state structure
+ *             @size_t x@ = a current accumulator
+ *
+ * Returns:    The next column of bits.
+ *
+ * Use:                Scans the next column of bits for a simultaneous
+ *             exponentiation.
+ */
+
+size_t exp_simulnext(exp_simul *e, size_t x)
+{
+  size_t i;
+
+  /* --- Move to the next word along --- */
+
+  if (!e->b) {
+    e->o--;
+    for (i = 0; i < e->n; i++)
+      e->s[i].w = e->o < e->s[i].len ? e->s[i].v[e->o] : 0;
+    e->b = MPW_BITS;
+  }
+
+  /* --- Scan out a column of bits --- */
+
+  for (i = 0; i < e->n; i++) {
+    x = (x << 1) | ((e->s[i].w >> (MPW_BITS - 1)) & 1u);
+    e->s[i].w <<= 1;
+  }      
+  e->b--;
+  return (x);
+}
+
+/*----- That's all, folks -------------------------------------------------*/
diff --git a/exp.h b/exp.h
new file mode 100644 (file)
index 0000000..6cfdfd8
--- /dev/null
+++ b/exp.h
@@ -0,0 +1,408 @@
+/* -*-c-*-
+ *
+ * $Id: exp.h,v 1.1 2001/06/16 13:00:59 mdw Exp $
+ *
+ * Generalized exponentiation
+ *
+ * (c) 2001 Straylight/Edgeware
+ */
+
+/*----- Licensing notice --------------------------------------------------* 
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ * 
+ * Catacomb is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Library General Public License for more details.
+ * 
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+/*----- Revision history --------------------------------------------------* 
+ *
+ * $Log: exp.h,v $
+ * Revision 1.1  2001/06/16 13:00:59  mdw
+ * New generic exponentation code.  Includes sliding-window simultaneous
+ * exponentiation.
+ *
+ */
+
+#ifdef CATACOMB_EXP_H
+#  error "Multiple inclusion of <catacomb/exp.h>"
+#endif
+
+#define CATACOMB_EXP_H
+
+#ifdef __cplusplus
+  extern "C" {
+#endif
+
+/*----- Header files ------------------------------------------------------*/
+
+#include <stddef.h>
+
+#include <mLib/alloc.h>
+
+#ifndef CATACOMB_MP_H
+#  include "mp.h"
+#endif
+
+/*----- Data structures ---------------------------------------------------*/
+
+typedef struct exp_simulscan {
+  mpw w;
+  size_t len;
+  const mpw *v;
+} exp_simulscan;  
+
+typedef struct exp_simul {
+  unsigned b;
+  size_t o, n;
+  exp_simulscan *s;
+} exp_simul;
+
+/*----- Macros provided ---------------------------------------------------*/
+
+/* --- Parameters --- */
+
+#ifndef EXP_WINSZ                      /* Sliding window size */
+#  define EXP_WINSZ 4                  /* Predefine if you need to */
+#endif
+
+/* --- These are determined from the window size --- */
+
+#define EXP_TABSZ (1 << EXP_WINSZ)
+#define EXP_THRESH (((MPW_BITS / EXP_WINSZ) << 2) + 1)
+
+/* --- Required operations --- *
+ *
+ * The macros here are independent of the underlying group elements.  You
+ * must provide the necessary group operations and other definitions.  The
+ * group operation is assumed to be written multiplicatively.
+ *
+ * @EXP_TYPE@                  The type of a group element, e.g., @mp *@.
+ *
+ * @EXP_COPY(d, x)@            Makes @d@ be a copy of @x@.
+ *
+ * @EXP_DROP(x)@               Discards the element @x@, reclaiming any
+ *                             memory it used.
+ *
+ * @EXP_MUL(a, x)@             Multiplies @a@ by @x@ (writing the result
+ *                             back to @a@).
+ *
+ * @EXP_SQR(a)@                        Multiplies @a@ by itself.
+ *
+ * @EXP_SETMUL(d, x, y)@       Sets @d@ to be the product of @x@ and @y@.
+ *                             The value @d@ has not been initialized.
+ *
+ * @EXP_SETSQR(d, x)@          Sets @d@ to be the square of @x@.
+ *
+ * Only @EXP_TYPE@, @EXP_MUL@ and @EXP_SQR@ are required for simple
+ * exponentation.  Sliding window and simultaneous exponentation require all
+ * of the operations.
+ */
+
+#ifndef EXP_TYPE
+#  error "EXP_TYPE not defined for <catacomb/exp.h>"
+#endif
+
+/* --- @EXP_SIMPLE@ --- *
+ *
+ * Arguments:  @a@ = the result object, initially a multiplicative identity
+ *             @g@ = the object to exponentiate
+ *             @x@ = the exponent, as a multiprecision integer
+ *
+ * Use:                Performs a simple left-to-right exponentiation.  At the end
+ *             of the code, the answer is left in @a@; @g@ and @x@ are
+ *             unchanged.
+ */
+
+#define EXP_SIMPLE(a, g, x) do {                                       \
+  mpscan sc;                                                           \
+  unsigned sq = 0;                                                     \
+                                                                       \
+  /* --- Begin scanning --- */                                         \
+                                                                       \
+  mp_rscan(&sc, x);                                                    \
+  if (!MP_RSTEP(&sc))                                                  \
+    goto exp_simple_exit;                                              \
+  while (!MP_RBIT(&sc))                                                        \
+    MP_RSTEP(&sc);                                                     \
+                                                                       \
+  /* --- Do the main body of the work --- */                           \
+                                                                       \
+  for (;;) {                                                           \
+    EXP_MUL(a, g);                                                     \
+    sq = 0;                                                            \
+    for (;;) {                                                         \
+      if (!MP_RSTEP(&sc))                                              \
+       goto exp_simple_done;                                           \
+      sq++;                                                            \
+      if (MP_RBIT(&sc))                                                        \
+       break;                                                          \
+    }                                                                  \
+    while (sq--) EXP_SQR(a);                                           \
+  }                                                                    \
+                                                                       \
+  /* --- Do a final round of squaring --- */                           \
+                                                                       \
+exp_simple_done:                                                       \
+  while (sq--) EXP_SQR(a);                                             \
+exp_simple_exit:;                                                      \
+} while (0)
+
+/* --- @EXP_WINDOW@ --- *
+ *
+ * Arguments:  @a@ = the result object, initially a multiplicative identity
+ *             @g@ = the object to exponentiate
+ *             @x@ = the exponent, as a multiprecision integer
+ *
+ * Use:                Performs a sliding-window exponentiation.  At the end of the
+ *             code, the answer is left in @a@; @g@ and @x@ are unchanged.
+ */
+
+#define EXP_WINDOW(a, g, x) do {                                       \
+  EXP_TYPE *v;                                                         \
+  EXP_TYPE g2;                                                         \
+  unsigned i, sq = 0;                                                  \
+  mpscan sc;                                                           \
+                                                                       \
+  /* --- Get going --- */                                              \
+                                                                       \
+  mp_rscan(&sc, x);                                                    \
+  if (!MP_RSTEP(&sc))                                                  \
+    goto exp_window_exit;                                              \
+                                                                       \
+  /* --- Do the precomputation --- */                                  \
+                                                                       \
+  EXP_SETSQR(g2, g);                                                   \
+  v = xmalloc(EXP_TABSZ * sizeof(EXP_TYPE));                           \
+  EXP_COPY(v[0], g);                                                   \
+  for (i = 1; i < EXP_TABSZ; i++)                                      \
+    EXP_SETMUL(v[i], v[i - 1], g2);                                    \
+  EXP_DROP(g2);                                                                \
+                                                                       \
+  /* --- Skip top-end zero bits --- *                                  \
+   *                                                                   \
+   * If the initial step worked, there must be a set bit somewhere, so \
+   * keep stepping until I find it.                                    \
+   */                                                                  \
+                                                                       \
+  while (!MP_RBIT(&sc))                                                        \
+    MP_RSTEP(&sc);                                                     \
+                                                                       \
+  /* --- Now for the main work --- */                                  \
+                                                                       \
+  for (;;) {                                                           \
+    unsigned l = 1;                                                    \
+    unsigned z = 0;                                                    \
+                                                                       \
+    /* --- The next bit is set, so read a window index --- *           \
+     *                                                                 \
+     * Reset @i@ to zero and increment @sq@.  Then, until either I read        \
+     * @WINSZ@ bits or I run out of bits, scan in a bit: if it's clear,        \
+     * bump the @z@ counter; if it's set, push a set bit into @i@,     \
+     * shift it over by @z@ bits, bump @sq@ by @z + 1@ and clear @z@.  \
+     * By the end of this palaver, @i@ is an index to the precomputed  \
+     * value in @v@.                                                   \
+     */                                                                        \
+                                                                       \
+    i = 0;                                                             \
+    sq++;                                                              \
+    while (l < EXP_WINSZ && MP_RSTEP(&sc)) {                           \
+      l++;                                                             \
+      if (!MP_RBIT(&sc))                                               \
+       z++;                                                            \
+      else {                                                           \
+       i = ((i << 1) | 1) << z;                                        \
+       sq += z + 1;                                                    \
+       z = 0;                                                          \
+      }                                                                        \
+    }                                                                  \
+                                                                       \
+    /* --- Do the squaring --- *                                       \
+     *                                                                 \
+     * Remember that @sq@ carries over from the zero-skipping stuff    \
+     * below.                                                          \
+     */                                                                        \
+                                                                       \
+    while (sq--) EXP_SQR(a);                                           \
+                                                                       \
+    /* --- Do the multiply --- */                                      \
+                                                                       \
+    EXP_MUL(a, v[i]);                                                  \
+                                                                       \
+    /* --- Now grind along through the rest of the bits --- */         \
+                                                                       \
+    sq = z;                                                            \
+    for (;;) {                                                         \
+      if (!MP_RSTEP(&sc))                                              \
+       goto exp_window_done;                                           \
+      if (MP_RBIT(&sc))                                                        \
+       break;                                                          \
+      sq++;                                                            \
+    }                                                                  \
+  }                                                                    \
+                                                                       \
+  /* --- Do a final round of squaring --- */                           \
+                                                                       \
+exp_window_done:                                                       \
+  while (sq--) EXP_SQR(a);                                             \
+  for (i = 0; i < EXP_TABSZ; i++)                                      \
+    EXP_DROP(v[i]);                                                    \
+  xfree(v);                                                            \
+exp_window_exit:;                                                      \
+} while (0)
+
+/* --- @EXP_SIMUL@ --- *
+ *
+ * Arguments:  @a@ = the result object, initially a multiplicative identity
+ *             @f@ = pointer to a vector of base/exp pairs
+ *             @n@ = the number of base/exp pairs
+ *
+ * Use:                Performs a simultaneous sliding-window exponentiation.  The
+ *             @f@ table is an array of structures containing members @base@
+ *             of type @EXP_TYPE@, and @exp@ of type @mp *@.
+ */
+
+#define EXP_SIMUL(a, f, n) do {                                                \
+  size_t i, j, jj, k;                                                  \
+  size_t vn = 1 << (EXP_WINSZ * n), m = (1 << n) - 1;                  \
+  EXP_TYPE *v = xmalloc(vn * sizeof(EXP_TYPE));                                \
+  exp_simul e;                                                         \
+  unsigned sq = 0;                                                     \
+                                                                       \
+  /* --- Fill in the precomputed table --- */                          \
+                                                                       \
+  j = 1;                                                               \
+  for (i = 0; i < n; i++) {                                            \
+    EXP_COPY(v[j], f[n - 1 - i].base);                                 \
+    j <<= 1;                                                           \
+  }                                                                    \
+  k = n * EXP_WINSZ;                                                   \
+  jj = 1;                                                              \
+  for (; i < k; i++) {                                                 \
+    EXP_SETSQR(v[j], v[jj]);                                           \
+    j <<= 1; jj <<= 1;                                                 \
+  }                                                                    \
+  for (i = 1; i < vn; i <<= 1) {                                       \
+    for (j = 1; j < i; j++)                                            \
+      EXP_SETMUL(v[j + i], v[j], v[i]);                                        \
+  }                                                                    \
+                                                                       \
+  /* --- Set up the bitscanners --- *                                  \
+   *                                                                   \
+   * Got to use custom scanners, to keep them all in sync.             \
+   */                                                                  \
+                                                                       \
+  e.n = n;                                                             \
+  e.b = 0;                                                             \
+  e.s = xmalloc(n * sizeof(*e.s));                                     \
+  e.o = 0;                                                             \
+  for (i = 0; i < n; i++) {                                            \
+    MP_SHRINK(f[i].exp);                                               \
+    e.s[i].len = MP_LEN(f[i].exp);                                     \
+    e.s[i].v = f[i].exp->v;                                            \
+    if (e.s[i].len > e.o)                                              \
+      e.o = e.s[i].len;                                                        \
+  }                                                                    \
+                                                                       \
+  /* --- Skip as far as a nonzero column in the exponent matrix --- */ \
+                                                                       \
+  do {                                                                 \
+    if (!e.o && !e.b)                                                  \
+      goto exp_simul_done;                                             \
+    i = exp_simulnext(&e, 0);                                          \
+  } while (!(i & m));                                                  \
+                                                                       \
+  /* --- Now for the main work --- */                                  \
+                                                                       \
+  for (;;) {                                                           \
+    unsigned l = 1;                                                    \
+    unsigned z = 0;                                                    \
+                                                                       \
+    /* --- Just read a nonzero column, so read a window index --- *    \
+     *                                                                 \
+     * Clear high bits of @i@ and increment @sq@.  Then, until either I        \
+     * read @WINSZ@ columns or I run out, scan in a column and append  \
+     * it to @i@.  If it's zero, bump the @z@ counter; if it's nonzero,        \
+     * bump @sq@ by @z + 1@ and clear @z@.  By the end of this palaver,        \
+     * @i@ is an index to the precomputed value in @v@, followed by    \
+     * @n * z@ zero bits.                                              \
+     */                                                                        \
+                                                                       \
+    sq++;                                                              \
+    while (l < EXP_WINSZ && (e.o || e.b)) {                            \
+      l++;                                                             \
+      i = exp_simulnext(&e, i);                                                \
+      if (!(i & m))                                                    \
+       z++;                                                            \
+      else {                                                           \
+       sq += z + 1;                                                    \
+       z = 0;                                                          \
+      }                                                                        \
+    }                                                                  \
+                                                                       \
+    /* --- Do the squaring --- *                                       \
+     *                                                                 \
+     * Remember that @sq@ carries over from the zero-skipping stuff    \
+     * below.                                                          \
+     */                                                                        \
+                                                                       \
+    while (sq--) EXP_SQR(a);                                           \
+                                                                       \
+    /* --- Do the multiply --- */                                      \
+                                                                       \
+    i >>= (z * n);                                                     \
+    EXP_MUL(a, v[i]);                                                  \
+                                                                       \
+    /* --- Now grind along through the rest of the bits --- */         \
+                                                                       \
+    sq = z;                                                            \
+    for (;;) {                                                         \
+      if (!e.o && !e.b)                                                        \
+       goto exp_simul_done;                                            \
+      if ((i = exp_simulnext(&e, 0)) != 0)                             \
+       break;                                                          \
+      sq++;                                                            \
+    }                                                                  \
+  }                                                                    \
+                                                                       \
+  /* --- Do a final round of squaring --- */                           \
+                                                                       \
+exp_simul_done:                                                                \
+  while (sq--) EXP_SQR(a);                                             \
+  for (i = 1; i < vn; i++)                                     \
+    EXP_DROP(v[i]);                                                    \
+  xfree(v);                                                            \
+} while (0)
+
+/*----- Functions provided ------------------------------------------------*/
+
+/* --- @exp_simulnext@ --- *
+ *
+ * Arguments:  @exp_simul *e@ = pointer to state structure
+ *             @size_t x@ = a current accumulator
+ *
+ * Returns:    The next column of bits.
+ *
+ * Use:                Scans the next column of bits for a simultaneous
+ *             exponentiation.
+ */
+
+extern size_t exp_simulnext(exp_simul */*e*/, size_t /*x*/);
+
+/*----- That's all, folks -------------------------------------------------*/
+
+#ifdef __cplusplus
+  }
+#endif