From e3dc2d225ea40a7453ae19ef4d12cf4c001076bd Mon Sep 17 00:00:00 2001 From: mdw Date: Fri, 11 Aug 2000 21:34:59 +0000 Subject: [PATCH] New restartable interface to Maurer testing. --- maurer.c | 193 +++++++++++++++++++++++++++++++++++++++++++++------------------ maurer.h | 73 ++++++++++++++++++++---- rspit.c | 64 +++++++++++++-------- 3 files changed, 241 insertions(+), 89 deletions(-) diff --git a/maurer.c b/maurer.c index 1828015..9846f89 100644 --- a/maurer.c +++ b/maurer.c @@ -1,6 +1,6 @@ /* -*-c-*- * - * $Id: maurer.c,v 1.1 2000/06/17 11:29:49 mdw Exp $ + * $Id: maurer.c,v 1.2 2000/08/11 21:34:59 mdw Exp $ * * Maurer's universal statistical test * @@ -30,6 +30,9 @@ /*----- Revision history --------------------------------------------------* * * $Log: maurer.c,v $ + * Revision 1.2 2000/08/11 21:34:59 mdw + * New restartable interface to Maurer testing. + * * Revision 1.1 2000/06/17 11:29:49 mdw * Maurer's universal statistical test. * @@ -84,32 +87,115 @@ static unsigned more(bitscan *b, unsigned l) return (x); } -/* --- @maurer@ --- * +/* --- @maurer_init@ --- * * - * Arguments: @const octet *p@ = pointer to a buffer of data - * @size_t sz@ = size of the buffer + * Arguments: @maurer_ctx *m@ = pointer to a context to initialize * @unsigned l@ = number of bits to take at a time * + * Returns: Minimum recommended amount of data to test. + * + * Use: Initializes a Maurer testing context. Note that @l@ values + * larger than 16 are not supported, and values less than 6 can + * give bizarre results. + */ + +unsigned long maurer_init(maurer_ctx *m, unsigned l) +{ + size_t i; + + assert(((void)"(maurer_init): chunks larger than 16 bits not supported", + 0 < l && l <= 16)); + + m->l = l; + m->a = 0; + m->b = 0; + m->n = 0; + m->x = 0; + + m->t = xmalloc(sizeof(unsigned long) << l); + for (i = 0; i < (1 << l); i++) + m->t[i] = 0; + return (((1000ul << l) * l + 7)/8); +} + +/* --- @bits@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to testing context + * @const octet **p@ = address of a buffer pointer + * @const octet *q@ = limit of the buffer pointer + * @unsigned *x@ = address to store next chunk + * + * Returns: Nonzero if some more bits arrived. + * + * Use: Reads some bits from a buffer. + */ + +static int bits(maurer_ctx *m, const octet **p, const octet *q, unsigned *x) +{ + while (m->b < m->l) { + if (*p >= q) + return (0); + m->a |= U8(*(*p)++) << m->b; + m->b += 8; + } + *x = m->a & ((1 << m->l) - 1); + m->a >>= m->l; + m->b -= m->l; + m->n++; + return (1); +} + +/* --- @maurer_test@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context + * @const void *buf@ = pointer to a buffer of data + * @size_t sz@ = size of the buffer + * + * Returns: --- + * + * Use: Scans a buffer of data and updates the testing context. + */ + +void maurer_test(maurer_ctx *m, const void *buf, size_t sz) +{ + const octet *p = buf, *l = p + sz; + unsigned long q = 10 << m->l; + unsigned x; + + /* --- Initialize the table --- */ + + while (m->n < q) { + if (!bits(m, &p, l, &x)) + return; + m->t[x] = m->n; + } + + /* --- Update the statistic --- */ + + while (bits(m, &p, l, &x)) { + m->x += log(m->n - m->t[x]); + m->t[x] = m->n; + } +} + +/* --- @maurer_done@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context + * * Returns: The statistic %$Z_u = (X_u - \mu)/\sigma$%, which should be - * normally distributed with a mean of zero and standard - * deviation of 1. + * normally distributed with a mean of 0 and variance of 1. * - * Use: Performs Maurer's universal statistical test on a buffer of - * data. You must ensure that the buffer size @sz@ is larger - * than %$11 l \cdot 2^l$% bits long, and ideally much larger. + * Use: Returns the result of Maurer's universal statistical test. + * Also frees the memory allocated during initialization. * - * Note that, under Unix systems, this function requires the - * maths library. + * If insufficient data was received, the value @HUGE_VAL@ is + * returned. */ -double maurer(const octet *p, size_t sz, unsigned l) +double maurer_done(maurer_ctx *m) { - size_t *t; /* Table of occurrences */ - size_t q, k; /* Block counts */ - size_t i; /* Loop counter */ - bitscan b; /* Bitscanning context */ - double x; /* Statistic accumulator */ - double mu, c, sigma; /* Other useful things */ + double x, mu, c, sigma; + unsigned long q = 10 << m->l, k; /* --- Table for means and variances --- */ @@ -132,59 +218,56 @@ double maurer(const octet *p, size_t sz, unsigned l) { 15.167379 , 3.421 } }; - assert(((void)"bit chunk size out of range", 0 < l && l <= 16)); - - /* --- Initialize everything --- */ - - t = xmalloc(sizeof(size_t) << l); - for (i = 0; i < (1 << l); i++) - t[i] = 0; - - b.p = p; - b.a = 0; - b.b = 0; - - k = (sz * 8)/l; - q = 10 << l; - assert(((void)"buffer too small for test", k > q && k - q > (1 << l))); - - /* --- Set up the table --- */ - - for (i = 0; i < q; i++) { - unsigned a = more(&b, l); - t[a] = i; - } - - /* --- Now compute the statistic --- */ + /* --- Check for bogus requests --- */ - x = 0; - for (; i < k; i++) { - unsigned a = more(&b, l); - x += log(i - t[a]); - t[a] = i; + if (m->n < q) { + x = HUGE_VAL; + goto done; } - /* --- Twiddle at the end --- * + /* --- Do the maths --- * * * This produces %$X_u$%, which should be normally distributed with a mean * and variance we can compute. */ - k -= q; - x /= k * log(2); + k = m->n - q; + x = m->x/(k * log(2)); /* --- Now translate this into %$Z_u$% distributed as %$N(0, 1)$% --- */ - mu = vtab[l - 1].mu; - c = 0.7 - 0.8/l + (1.6 + 12.8/l) * pow(k, -4.0/l); - sigma = sqrt(c * c * vtab[l - 1].var_1 / k); - + mu = vtab[m->l - 1].mu; + c = 0.7 - 0.8/m->l + (1.6 + 12.8/m->l) * pow(k, -4.0/m->l); + sigma = sqrt(c * c * vtab[m->l - 1].var_1 / k); x = (x - mu)/sigma; /* --- Done (phew!) --- */ - xfree(t); +done: + xfree(m->t); return (x); } +/* --- @maurer@ --- * + * + * Arguments: @const octet *p@ = pointer to a buffer of data + * @size_t sz@ = size of the buffer + * @unsigned l@ = number of bits to take at a time + * + * Returns: The statistic %$Z_u$%. + * + * Use: Simple interface to Maurer's universal statistical test. For + * large %$l$% you should use the more complicated restartable + * interface. + */ + +double maurer(const octet *p, size_t sz, unsigned l) +{ + maurer_ctx m; + + maurer_init(&m, l); + maurer_test(&m, p, sz); + return (maurer_done(&m)); +} + /*----- That's all, folks -------------------------------------------------*/ diff --git a/maurer.h b/maurer.h index cd74e8c..904e4a2 100644 --- a/maurer.h +++ b/maurer.h @@ -1,6 +1,6 @@ /* -*-c-*- * - * $Id: maurer.h,v 1.1 2000/06/17 11:29:49 mdw Exp $ + * $Id: maurer.h,v 1.2 2000/08/11 21:34:59 mdw Exp $ * * Maurer's universal statistical test * @@ -30,6 +30,9 @@ /*----- Revision history --------------------------------------------------* * * $Log: maurer.h,v $ + * Revision 1.2 2000/08/11 21:34:59 mdw + * New restartable interface to Maurer testing. + * * Revision 1.1 2000/06/17 11:29:49 mdw * Maurer's universal statistical test. * @@ -46,24 +49,74 @@ #include +/*----- Data structures ---------------------------------------------------*/ + +typedef struct maurer_ctx { + unsigned l; /* Chunk size, in bits */ + uint32 a; /* Bitscanner accumulator */ + unsigned b; /* Number of bits in accumulator */ + unsigned long n; /* Number of chunks read so far */ + double x; /* Statistic value */ + unsigned long *t; /* Offset table */ +} maurer_ctx; + /*----- Functions provided ------------------------------------------------*/ +/* --- @maurer_init@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to a context to initialize + * @unsigned l@ = number of bits to take at a time + * + * Returns: Minimum recommended amount of data to test. + * + * Use: Initializes a Maurer testing context. Note that @l@ values + * larger than 16 are not supported, and values less than 6 can + * give bizarre results. + */ + +extern unsigned long maurer_init(maurer_ctx */*m*/, unsigned /*l*/); + +/* --- @maurer_test@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context + * @const void *buf@ = pointer to a buffer of data + * @size_t sz@ = size of the buffer + * + * Returns: --- + * + * Use: Scans a buffer of data and updates the testing context. + */ + +extern void maurer_test(maurer_ctx */*m*/, + const void */*buf*/, size_t /*sz*/); + +/* --- @maurer_done@ --- * + * + * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context + * + * Returns: The statistic %$Z_u = (X_u - \mu)/\sigma$%, which should be + * normally distributed with a mean of 0 and variance of 1. + * + * Use: Returns the result of Maurer's universal statistical test. + * Also frees the memory allocated during initialization. + * + * If insufficient data was received, the value @HUGE_VAL@ is + * returned. + */ + +extern double maurer_done(maurer_ctx */*m*/); + /* --- @maurer@ --- * * * Arguments: @const octet *p@ = pointer to a buffer of data * @size_t sz@ = size of the buffer * @unsigned l@ = number of bits to take at a time * - * Returns: The statistic %$Z_u = (X_u - \mu)/\sigma$%, which should be - * normally distributed with a mean of zero and standard - * deviation of 1. - * - * Use: Performs Maurer's universal statistical test on a buffer of - * data. You must ensure that the buffer size @sz@ is larger - * than %$11 l \cdot 2^l$% bits long, and ideally much larger. + * Returns: The statistic %$Z_u$%. * - * Note that, under Unix systems, this function requires the - * maths library. + * Use: Simple interface to Maurer's universal statistical test. For + * large %$l$% you should use the more complicated restartable + * interface. */ extern double maurer(const octet */*p*/, size_t /*sz*/, unsigned /*l*/); diff --git a/rspit.c b/rspit.c index 01ebf4d..4622d90 100644 --- a/rspit.c +++ b/rspit.c @@ -1,6 +1,6 @@ /* -*-c-*- * - * $Id: rspit.c,v 1.9 2000/08/04 23:24:15 mdw Exp $ + * $Id: rspit.c,v 1.10 2000/08/11 21:34:59 mdw Exp $ * * Spit out random numbers * @@ -30,6 +30,9 @@ /*----- Revision history --------------------------------------------------* * * $Log: rspit.c,v $ + * Revision 1.10 2000/08/11 21:34:59 mdw + * New restartable interface to Maurer testing. + * * Revision 1.9 2000/08/04 23:24:15 mdw * Add a timer and a discard option. * @@ -373,7 +376,7 @@ static int opt(void) char *p; unsigned long lo, hi; lo = strtoul(optarg, &p, 0); - if (*p == '-') + if (*p == '-' || *p == ',') hi = strtoul(p + 1, &p, 0); else hi = lo; @@ -1119,6 +1122,21 @@ static int genbuf(const void *buf, size_t sz, void *p) return (0); } +typedef struct genmaurer_ctx { + size_t n; + maurer_ctx *m; +} genmaurer_ctx; + +static int genmaurer(const void *buf, size_t sz, void *p) +{ + genmaurer_ctx *g = p; + size_t i; + + for (i = 0; i < g->n; i++) + maurer_test(&g->m[i], buf, sz); + return (0); +} + static int generate(grand *r, size_t outsz, int (*func)(const void *buf, size_t sz, void *p), void *p) @@ -1236,14 +1254,17 @@ static int generate(grand *r, size_t outsz, if (flags & f_progress) fputc('\n', stderr); if (flags & f_timer) { - double sec = (double)clk/CLOCKS_PER_SEC; - double bps = (outsz << 3)/sec; - char *kk; - - for (kk = kmg; bps > 1024 && kk[1]; kk++, bps /= 1024) - ; - fprintf(stderr, "generated %lu bytes in %g secs (%g %cb/s)\n", - (unsigned long)outsz, sec, bps, *kk); + fprintf(stderr, "generated %lu bytes ", (unsigned long)outsz); + if (!clk) + fputs("too quickly to measure\n", stderr); + else { + char *kk; + double sec = (double)clk/CLOCKS_PER_SEC; + double bps = (outsz << 3)/sec; + for (kk = kmg; bps > 1024 && kk[1]; kk++, bps /= 1024) + ; + fprintf(stderr, "in %g secs (%g %cb/s)\n", sec, bps, *kk); + } } return (0); } @@ -1330,13 +1351,10 @@ int main(int ac, char *av[]) /* --- Do Maurer's test --- */ if (flags & f_maurer) { - octet *buf; size_t bufsz; unsigned i; unsigned rc = 0; - unsigned f = 0, jj = 0; - double maxz = 0; - octet *p; + genmaurer_ctx g; static struct { double x; const char *sig; } sigtab[] = { { 3.2905, "1e-3" }, @@ -1346,24 +1364,21 @@ int main(int ac, char *av[]) { 0 , 0 } }; + g.n = maurer_hi - maurer_lo + 1; + g.m = xmalloc(g.n * sizeof(maurer_ctx)); + for (i = 0; i < g.n; i++) + maurer_init(&g.m[i], i + maurer_lo); bufsz = (100 * maurer_hi) << maurer_hi; - if ((buf = a_alloc(arena_global, bufsz)) == 0) - die(EXIT_FAILURE, "not enough memory for data buffer"); - p = buf; - generate(r, bufsz, genbuf, &p); + + generate(r, bufsz, genmaurer, &g); for (i = maurer_lo; i <= maurer_hi; i++) { - double z = maurer(buf, bufsz, i); + double z = maurer_done(&g.m[i - maurer_lo]); double zz = fabs(z); unsigned j; for (j = 0; sigtab[j].sig; j++) { if (zz > sigtab[j].x) { - if (zz > fabs(maxz)) { - maxz = z; - f = i; - jj = j; - } rc = EXIT_FAILURE; moan("failed, bits = %u, sig = %s, Z_u = %g", i, sigtab[j].sig, z); @@ -1374,6 +1389,7 @@ int main(int ac, char *av[]) fprintf(stderr, "bits = %u, Z_u = %g\n", i, z); } + xfree(g.m); return (rc); } -- 2.11.0