chiark - git - mdw - catacomb/blob - pub/bbs-jump.c

   1 /* -*-c-*-
   2  *
   3  * Jumping around a BBS sequence
   4  *
   5  * (c) 1999 Straylight/Edgeware
   6  */
   7
   8 /*----- Licensing notice --------------------------------------------------*
   9  *
  10  * This file is part of Catacomb.
  11  *
  12  * Catacomb is free software; you can redistribute it and/or modify
  13  * it under the terms of the GNU Library General Public License as
  14  * published by the Free Software Foundation; either version 2 of the
  15  * License, or (at your option) any later version.
  16  *
  17  * Catacomb is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20  * GNU Library General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Library General Public
  23  * License along with Catacomb; if not, write to the Free
  24  * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
  25  * MA 02111-1307, USA.
  26  */
  27
  28 /*----- Header files ------------------------------------------------------*/
  29
  30 #include "bbs.h"
  31 #include "mp.h"
  32 #include "mpbarrett.h"
  33 #include "mpcrt.h"
  34 #include "mpint.h"
  35
  36 /*----- Main code ---------------------------------------------------------*/
  37
  38 /* --- @jump@ --- *
  39  *
  40  * Arguments:   @bbs *b@ = pointer to BBS generator context
  41  *              @const bbs_priv *bp@ = pointer to BBS modulus factors
  42  *              @mp *n@ = number of steps to move
  43  *              @mp *px@ = exponent mod @p@ for a one-step jump
  44  *              @mp *qx@ = exponent mod @q@ for a one-step jump
  45  *
  46  * Returns:     ---
  47  *
  48  * Use:         Jumps a BBS context a certain number of places (assuming the
  49  *              arguments are right).
  50  *
  51  *              Let the BBS modulus be %$n = pq$% and the current residue be
  52  *              %$x$%.  Then the computations performed are:
  53  *
  54  *                * Calculate %$x_p = x \bmod p$% and %$x_q = x \bmod q$%.
  55  *
  56  *                * Determine %$e_p = px^n \bmod (p - 1)$% and similarly
  57  *                  %$e_q = qx^n \bmod (p - 1)$%.
  58  *
  59  *                * Calculate %$x_p' = x_p^{e_p} \bmod p$% and
  60  *                  %$x_q' = x_q^{e_q} \bmod q$%.
  61  *
  62  *                * Combine %$x_p'$% and %$x_q'$% using the Chinese Remainder
  63  *                  Theorem.
  64  *
  65  *              If you want to step the generator forwards, simply set
  66  *              %$px = qx = 2$%.  If you want to step backwards, make
  67  *              %$px = (p + 1)/4$% and %$qx = (q + 1)/4$%.  Note that, if
  68  *              %$x$% is a quadratic residue mod $%p$%, then
  69  *
  70  *              %$(x^2) ^ {(p + 1)/4}$%
  71  *                %${} = x^{(p + 1)/2}$%
  72  *                %${} = x \cdot x^{(p - 1)/2}$%
  73  *                %${} = x$%
  74  *
  75  *              Simple, no?  (Note that the division works because
  76  *              %$p \equiv 3 \pmod 4$%.)
  77  */
  78
  79 static void jump(bbs *b, const bbs_priv *bp, mp *n,
  80                  mp *px, mp *qx)
  81 {
  82   mp *ep, *eq;
  83   mp *v[2] = { MP_NEW, MP_NEW };
  84
  85   /* --- First work out the exponents --- */
  86
  87   {
  88     mpbarrett mb;
  89     mp *m;
  90
  91     m = mp_sub(MP_NEW, bp->p, MP_ONE);
  92     mpbarrett_create(&mb, m);
  93     ep = mpbarrett_exp(&mb, MP_NEW, px, n);
  94     mpbarrett_destroy(&mb);
  95     if (qx == px)
  96       eq = MP_COPY(ep);
  97     else {
  98       m = mp_sub(m, bp->q, MP_ONE);
  99       mpbarrett_create(&mb, m);
 100       eq = mpbarrett_exp(&mb, MP_NEW, qx, n);
 101       mpbarrett_destroy(&mb);
 102     }
 103
 104     mp_drop(m);
 105   }
 106
 107   /* --- Now calculate the residues of @x@ --- */
 108
 109   mp_div(0, &v[0], b->x, bp->p);
 110   mp_div(0, &v[1], b->x, bp->q);
 111
 112   /* --- Exponentiate --- */
 113
 114   {
 115     mpbarrett mb;
 116
 117     mpbarrett_create(&mb, bp->p);
 118     v[0] = mpbarrett_exp(&mb, v[0], v[0], ep);
 119     mpbarrett_destroy(&mb);
 120
 121     mpbarrett_create(&mb, bp->q);
 122     v[1] = mpbarrett_exp(&mb, v[1], v[1], eq);
 123     mpbarrett_destroy(&mb);
 124
 125     mp_drop(ep);
 126     mp_drop(eq);
 127   }
 128
 129   /* --- Sort out the result using the Chinese Remainder Theorem --- */
 130
 131   {
 132     mpcrt_mod mv[2];
 133     mpcrt c;
 134     int i;
 135
 136     mv[0].m = MP_COPY(bp->p);
 137     mv[1].m = MP_COPY(bp->q);
 138     for (i = 0; i < 2; i++)
 139       mv[i].n = mv[i].ni = mv[i].nni = MP_NEW;
 140     mpcrt_create(&c, mv, 2, b->mb.m);
 141     b->x = mpcrt_solve(&c, b->x, v);
 142     mpcrt_destroy(&c);
 143   }
 144
 145   /* --- Tidy away --- */
 146
 147   mp_drop(v[0]);
 148   mp_drop(v[1]);
 149   b->r = b->x->v[0];
 150   b->b = b->k;
 151 }
 152
 153 /* --- @bbs_{ff,rew}{,n}@ --- *
 154  *
 155  * Arguments:   @bbs *b@ = pointer to a BBS generator state
 156  *              @const bbs_priv *bp@ = pointer to BBS modulus factors
 157  *              @mp *n@, @unsigned long n@ = number of steps to make
 158  *
 159  * Returns:     ---
 160  *
 161  * Use:         `Fast-forwards' or rewinds a Blum-Blum-Shub generator by @n@
 162  *              steps.  The @...n@ versions take an @unsigned long@ argument;
 163  *              the non-@...n@ versions a multiprecision integer.  If @n@ is
 164  *              negative then the generator is stepped in the reverse
 165  *              direction.
 166  */
 167
 168 static void ff(bbs *b, const bbs_priv *bp, mp *n)
 169   { jump(b, bp, n, MP_TWO, MP_TWO); }
 170
 171 static void rew(bbs *b, const bbs_priv *bp, mp *n)
 172 {
 173   mp *px = mp_lsr(MP_NEW, bp->p, 2);
 174   mp *qx = mp_lsr(MP_NEW, bp->q, 2);
 175   px = mp_add(px, px, MP_ONE);
 176   qx = mp_add(qx, qx, MP_ONE);
 177   jump(b, bp, n, px, qx);
 178   mp_drop(px);
 179   mp_drop(qx);
 180 }
 181
 182 void bbs_ff(bbs *b, const bbs_priv *bp, mp *n)
 183 {
 184   if (!MP_NEGP(n))
 185     ff(b, bp, n);
 186   else {
 187     n = mp_neg(MP_NEW, n);
 188     rew(b, bp, n);
 189     mp_drop(n);
 190   }
 191 }
 192
 193 void bbs_ffn(bbs *b, const bbs_priv *bp, unsigned long n)
 194   { mp *nn = mp_fromulong(MP_NEW, n); ff(b, bp, nn); mp_drop(nn); }
 195
 196 void bbs_rew(bbs *b, const bbs_priv *bp, mp *n)
 197 {
 198   if (!MP_NEGP(n))
 199     rew(b, bp, n);
 200   else {
 201     n = mp_neg(MP_NEW, n);
 202     ff(b, bp, n);
 203     mp_drop(n);
 204   }
 205 }
 206
 207 void bbs_rewn(bbs *b, const bbs_priv *bp, unsigned long n)
 208   { mp *nn = mp_fromulong(MP_NEW, n); bbs_rew(b, bp, nn); mp_drop(nn); }
 209
 210 /*----- Test rig ----------------------------------------------------------*/
 211
 212 #ifdef TEST_RIG
 213
 214 static int verify(dstr *v)
 215 {
 216   bbs_priv bp;
 217   bbs b;
 218   mp *x;
 219   unsigned long n;
 220   int ok = 1;
 221   uint32 p, q, r;
 222   int i;
 223
 224   bp.p = *(mp **)v[0].buf;
 225   bp.q = *(mp **)v[1].buf;
 226   bp.n = mp_mul(MP_NEW, bp.p, bp.q);
 227   x = *(mp **)v[2].buf;
 228   n = *(unsigned long *)v[3].buf;
 229
 230   bbs_create(&b, bp.n, x);
 231   p = bbs_bits(&b, 32);
 232
 233   bbs_seed(&b, x);
 234   for (i = 0; i < n; i++)
 235     bbs_step(&b);
 236   q = bbs_bits(&b, 32);
 237   bbs_wrap(&b);
 238
 239   bbs_rewn(&b, &bp, n + (32 + b.k - 1) / b.k);
 240   r = bbs_bits(&b, 32);
 241
 242   if (r != p) {
 243     fputs("\n*** bbs rewind failure\n", stderr);
 244     fputs("p = ", stderr); mp_writefile(bp.p, stderr, 10); fputc('\n', stderr);
 245     fputs("q = ", stderr); mp_writefile(bp.q, stderr, 10); fputc('\n', stderr);
 246     fputs("n = ", stderr); mp_writefile(bp.n, stderr, 10); fputc('\n', stderr);
 247     fputs("x = ", stderr); mp_writefile(x, stderr, 10); fputc('\n', stderr);
 248     fprintf(stderr, "stepped %lu back\n", n + (32 + b.k - 1) / b.k);
 249     fprintf(stderr, "expected output = %08lx, found %08lx\n",
 250             (unsigned long)p, (unsigned long)r);
 251     ok = 0;
 252   }
 253
 254   bbs_seed(&b, x);
 255   bbs_ffn(&b, &bp, n);
 256   r = bbs_bits(&b, 32);
 257
 258   if (q != r) {
 259     fputs("\n*** bbs fastforward failure\n", stderr);
 260     fputs("p = ", stderr); mp_writefile(bp.p, stderr, 10); fputc('\n', stderr);
 261     fputs("q = ", stderr); mp_writefile(bp.q, stderr, 10); fputc('\n', stderr);
 262     fputs("n = ", stderr); mp_writefile(bp.n, stderr, 10); fputc('\n', stderr);
 263     fputs("x = ", stderr); mp_writefile(x, stderr, 10); fputc('\n', stderr);
 264     fprintf(stderr, "stepped %lu back\n", n + (32 + b.k - 1) / b.k);
 265     fprintf(stderr, "expected output = %08lx, found %08lx\n",
 266             (unsigned long)q, (unsigned long)r);
 267     ok = 0;
 268   }
 269
 270   bbs_destroy(&b);
 271   mp_drop(bp.p);
 272   mp_drop(bp.q);
 273   mp_drop(bp.n);
 274   mp_drop(x);
 275
 276   assert(mparena_count(MPARENA_GLOBAL) == 0);
 277   return (ok);
 278 }
 279
 280 static test_chunk tests[] = {
 281   { "bbs-jump", verify, { &type_mp, &type_mp, &type_mp, &type_ulong, 0 } },
 282   { 0, 0, { 0 } }
 283 };
 284
 285 int main(int argc, char *argv[])
 286 {
 287   sub_init();
 288   test_run(argc, argv, tests, SRCDIR "/t/bbs");
 289   return (0);
 290 }
 291
 292 #endif
 293
 294 /*----- That's all, folks -------------------------------------------------*/