Improve analysis section.

[storin] / sac.c

/* -*-c-*-
 *
 * $Id: sac.c,v 1.1 2000/05/21 11:28:30 mdw Exp $
 *
 * Testing for strict avalanche
 *
 * (c) 2000 Mark Wooding
 */

/*----- Licensing notice --------------------------------------------------* 
 *
 * Copyright (c) 2000 Mark Wooding
 * All rights reserved.
 * 
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are
 * met:
 * 
 * 1. Redistributions of source code must retain the above copyright
 *    notice, this list of conditions and the following disclaimer.
 * 
 * 2, Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in the
 *    documentation and/or other materials provided with the distribution.
 * 
 * 3. The name of the authors may not be used to endorse or promote
 *    products derived from this software without specific prior written
 *    permission.
 * 
 * THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESS OR IMPLIED
 * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
 * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.  IN
 * NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
 * INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
 * ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 * POSSIBILITY OF SUCH DAMAGE.
 * 
 * Instead of accepting the above terms, you may redistribute and/or modify
 * this software under the terms of either the GNU General Public License,
 * or the GNU Library General Public License, published by the Free
 * Software Foundation; either version 2 of the License, or (at your
 * option) any later version.
 */

/*----- Revision history --------------------------------------------------* 
 *
 * $Log: sac.c,v $
 * Revision 1.1  2000/05/21 11:28:30  mdw
 * Initial check-in.
 *
 */

/*----- Header files ------------------------------------------------------*/

#include <math.h>
#include <stdio.h>
#include <stdlib.h>

#include "bits.h"

#include "fibrand.h"
#include "matrix.h"
#include "storin-tab.h"

/*----- Static variables --------------------------------------------------*/

static fibrand r;

/*----- The constant matrix -----------------------------------------------*/

static const uint24 m[] = STORIN_M;

/*----- Magic numbers -----------------------------------------------------*/

#define PROBES 16384
#define ROUNDS 3
#define ROT(x) ((x) ^ ((x) >> 12))
/* #define ROT(x) ROL24(x, 7) */

/*----- Main code ---------------------------------------------------------*/

static octet w[256];

#define HAMWEIGHT(x) (w[U8((x) >> 16)] + w[U8((x) >> 8)] + w[U8(x)])

static void haminit(void) {
  unsigned i;
  for (i = 0; i < 256; i++) {
    unsigned ww = 0;
    unsigned j;
    for (j = 0; j < 8; j++) {
      if (i & (1 << j))
	ww++;
    }
    w[i] = ww;
  }
}

static void eblk(uint24 *x)
{
  uint24 t[4];
  unsigned i, j;

  for (i = 0; i < ROUNDS; i++) {
    /* No key mixing */
    matmul(t, m, x, 4, 4, 1);
    for (j = 0; j < 4; j++)
      x[j] = ROT(t[j]);
  }
}

typedef struct stats {
  double sx;
  double sx2;
  double n;
} stats;

static void probe(unsigned bit, uint24 *delta, struct stats *ps)
{
  uint24 x[4], y[4];
  unsigned i, j;
  struct stats s = { 0, 0, 0 };

  for (i = 0; i < PROBES; i++) {
    unsigned h;
    for (j = 0; j < 4; j++) {
      x[j] = U24(fibrand_step(&r));
      y[j] = x[j] ^ delta[j];
    }
    eblk(x);
    eblk(y);
    h = 0;
    for (j = 0; j < 4; j++) {
      uint24 ww = x[j] ^ y[j];
      h += HAMWEIGHT(ww);
    }
/*     printf("%06x %06x %06x %06x -> %06x %06x %06x %06x\n", */
/* 	   delta[0], delta[1], delta[2], delta[3], */
/* 	   x[0] ^ y[0], x[1] ^ y[1], x[2] ^ y[2], x[3] ^ y[3]); */
    s.sx += h; ps->sx += h;
    s.sx2 += h * h; ps->sx2 += h * h;
    s.n++; ps->n++;
  }

  {
    double mean = s.sx / s.n;
    double var = s.sx2 / s.n - mean * mean;
    printf("bit %u: mean = %g, sd = %g\n", bit, mean, sqrt(var));
  }
}

int main(void)
{
  uint24 delta[4] = { 0 };
  unsigned i;
  struct stats s = { 0, 0, 0 };
  double mean, var;

  haminit();
  fibrand_lcseed(&r, 0);

  for (i = 0; i < 96; i++) {
    uint24 *dd = delta + i / 24;
    uint24 m = 1 << (i % 24);
    *dd ^= m;
    probe(i, delta, &s);
    *dd ^= m;
  }

  mean = s.sx / s.n;
  var = s.sx2 / s.n - mean * mean;
  printf("\nsummary: mean = %g, sd = %g\n", mean, sqrt(var));

  return (0);
}

/*----- That's all, folks -------------------------------------------------*/