X-Git-Url: https://git.distorted.org.uk/u/mdw/catacomb/blobdiff_plain/ba6e6b64033b1f9de49feccb5c9cd438354481f7..0f00dc4c8eb47e67bc0f148c2dd109f73a451e0a:/lcrand.c diff --git a/lcrand.c b/lcrand.c deleted file mode 100644 index a108f13..0000000 --- a/lcrand.c +++ /dev/null @@ -1,306 +0,0 @@ -/* -*-c-*- - * - * $Id: lcrand.c,v 1.5 2004/04/08 01:36:15 mdw Exp $ - * - * Simple linear congruential generator - * - * (c) 1999 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. - */ - -/*----- Header files ------------------------------------------------------*/ - -#include -#include -#include -#include - -#include -#include - -#include "grand.h" -#include "lcrand.h" - -/*----- Magic numbers -----------------------------------------------------*/ - -/* --- The generator parameters --- */ - -#define P LCRAND_P /* Modulus */ -#define A LCRAND_A /* Multiplier (primitive mod @p@) */ -#define C LCRAND_C /* Additive constant */ - -/* --- Precomputed values for modular reduction --- */ - -#define D 5 /* %$p = 2^{32} - d$% */ - -/* --- Other useful bits --- */ - -#define P256 4294967040u /* Highest multiple of 256 < %$p$% */ - -/*----- Main code ---------------------------------------------------------*/ - -/* --- @lcrand@ --- * - * - * Arguments: @uint32 x@ = seed value - * - * Returns: New state of the generator. - * - * Use: Steps the generator. Returns %$ax + c \bmod p$%. - */ - -uint32 lcrand(uint32 x) -{ - uint32 a[2], xx[2]; - uint32 yy[2]; - - /* --- Unpack things into the arrays --- */ - - a[0] = U16(A); a[1] = U16(A >> 16); - xx[0] = U16(x); xx[1] = U16(x >> 16); - - /* --- Multiply everything together --- * - * - * This is plain old long multiplication, although it looks a bit strange. - * I set up the top and bottom partial products directly where they're - * supposed to be. The cross terms I add together, with the low 16 bits in - * @q@ and the high 32 bits in @p@. These I then add into the product. - */ - - { - uint32 p, q; - - yy[0] = a[0] * xx[0]; - yy[1] = a[1] * xx[1]; - - p = a[0] * xx[1]; - q = p + a[1] * xx[0]; - p = ((q < p) << 16) + (q >> 16); - q = U16(q) << 16; - - q += yy[0]; - if (q < yy[0]) - p++; - else - p += (q >> 16) >> 16; - yy[0] = q; - - yy[1] += p; - } - - /* --- Now reduce mod p --- * - * - * I'm using shifts and adds to do the multiply step here. This needs to - * be changed if @D@ ever becomes something other than 5. - */ - -#if D != 5 -# error "Change shift sequence!" -#endif - - { - uint32 q; - - q = yy[1]; - x = yy[0]; - - while (q) { - uint32 y, z; - y = q >> 30; - z = q << 2; - z += q; - if (z < q) - y++; - else - y += (q >> 16) >> 16; - q = y; - x += z; - if (x < z || x > P) - x -= P; - } - } - - /* --- Now add on the constant --- */ - - x += C; - if (x < C || x >= P) - x -= P; - - /* --- Done --- */ - - return (x); -} - -/* --- @lcrand_range@ --- * - * - * Arguments: @uint32 *x@ = pointer to seed value (updated) - * @uint32 m@ = limit allowable - * - * Returns: A uniformly distributed pseudorandom integer in the interval - * %$[0, m)$%. - */ - -uint32 lcrand_range(uint32 *x, uint32 m) -{ - uint32 xx = *x; - uint32 r = P - P % m; - do xx = lcrand(xx); while (xx >= r); - *x = xx; - return (xx % m); -} - -/*----- Generic interface -------------------------------------------------*/ - -typedef struct gctx { - grand r; - uint32 x; -} gctx; - -static void gdestroy(grand *r) -{ - gctx *g = (gctx *)r; - DESTROY(g); -} - -static int gmisc(grand *r, unsigned op, ...) -{ - gctx *g = (gctx *)r; - va_list ap; - int rc = 0; - va_start(ap, op); - - switch (op) { - case GRAND_CHECK: - switch (va_arg(ap, unsigned)) { - case GRAND_CHECK: - case GRAND_SEEDINT: - case GRAND_SEEDUINT32: - case GRAND_SEEDRAND: - rc = 1; - break; - default: - rc = 0; - break; - } - break; - case GRAND_SEEDINT: - g->x = va_arg(ap, unsigned); - break; - case GRAND_SEEDUINT32: - g->x = va_arg(ap, uint32); - break; - case GRAND_SEEDRAND: { - grand *rr = va_arg(ap, grand *); - uint32 x; - do x = rr->ops->word(rr); while (x >= P || x == LCRAND_FIXEDPT); - g->x = x; - } break; - default: - GRAND_BADOP; - break; - } - - va_end(ap); - return (rc); -} - -static uint32 graw(grand *r) -{ - gctx *g = (gctx *)r; - g->x = lcrand(g->x); - return (g->x); -} - -static octet gbyte(grand *r) -{ - gctx *g = (gctx *)r; - uint32 x = g->x; - do x = lcrand(x); while (x >= P256); - g->x = x; - return (x / (P256 / 256)); -} - -static uint32 grange(grand *r, uint32 l) -{ - gctx *g = (gctx *)r; - return (lcrand_range(&g->x, l)); -} - -static const grand_ops gops = { - "lcrand", - LCRAND_P, 0, - gmisc, gdestroy, - graw, gbyte, grand_word, grange, grand_fill -}; - -/* --- @lcrand_create@ --- * - * - * Arguments: @uint32 x@ = initial seed - * - * Returns: Pointer to a generic generator. - * - * Use: Constructs a generic generator interface over a linear - * congruential generator. - */ - -grand *lcrand_create(uint32 x) -{ - gctx *g = CREATE(gctx); - g->r.ops = &gops; - g->x = x; - return (&g->r); -} - -/*----- Test rig ----------------------------------------------------------*/ - -#ifdef TEST_RIG - -#include - -static int verify(dstr *v) -{ - uint32 x = *(uint32 *)v[0].buf; - uint32 y = *(uint32 *)v[1].buf; - uint32 z = lcrand(x); - int ok = 1; - if (y != z) { - fprintf(stderr, - "\n*** lcrand failed. lcrand(%lu) = %lu, expected %lu\n", - (unsigned long)x, (unsigned long)z, (unsigned long)y); - ok = 0; - } - return (ok); -} - -static test_chunk tests[] = { - { "lcrand", verify, { &type_uint32, &type_uint32, 0 } }, - { 0, 0, { 0 } } -}; - -int main(int argc, char *argv[]) -{ - test_run(argc, argv, tests, SRCDIR"/tests/lcrand"); - return (0); -} - -#endif - -/*----- That's all, folks -------------------------------------------------*/