From gary Wed Jun  5 10:44:38 1996
From: fischer-michael@CS.YALE.EDU
Date: Tue, 4 Jun 96 18:44:20 EDT
Subject: Patches to mrandom-3.0
To: Clark Thomborson <cthombor@cs>

*** 3.1	1993/05/28 19:25:35
--- bentley.c	1996/06/04 14:51:51
***************
*** 60,66 ****
  long rngstate[RNGstatesize_2];
  {
    int i, ii;
!   long last, next, v[55];
    arr[0] = last = seed; 
    next = 1;
    for (i=1; i < 55; i++) {
--- 60,66 ----
  long rngstate[RNGstatesize_2];
  {
    int i, ii;
!   long last, next;
    arr[0] = last = seed; 
    next = 1;
    for (i=1; i < 55; i++) {
*** 3.1	1993/05/28 19:25:36
--- mrandom.c	1996/06/04 19:35:53
***************
*** 1,4 ****
! /* mrandom.c 3.1 5/28/93 */
  /*
   *		mrandom.c
   *
--- 1,4 ----
! /* mrandom.c 3.1x 5/31/96, derived from mrandom.c 3.1 5/28/93 */
  /*
   *		mrandom.c
   *
***************
*** 27,32 ****
--- 27,40 ----
   *         - Allow choice of algorithms in mrandomrv
   *         - Dynamic allocation of RNGs
   *
+  *      Version 3.1x:
+  *              Michael Fischer, May 1996
+  *         - Ported to NeXT; removed unused variables
+  *         - Converted function headers to ISO-C
+  *         - Fixed memory-overwrite bug in scanvector()
+  *         - Fixed state saving/restoring for 4.3bsd random()
+  *         - Cleaned up nextval() and restart_rng().
+  *
   *	This material is based upon work supported by the National
   *	Science Foundation under grant number MIP-9023238.  The
   *	Government has certain rights in this material.
***************
*** 47,57 ****
--- 55,71 ----
  /*****************/
  /* Include files */
  /*****************/
+ #ifdef NeXT
+ #include <stdlib.h>
+ #else
  #include <malloc.h> /* To allocate space for RNGs*/
+ #endif
  #include <math.h>   /* we need floor() */
  #include <stdio.h>  /* We need FILE */
  #include <strings.h>
+ #ifndef NeXT
  #include <values.h> /* We need MAXLONG */
+ #endif
  
  /*******************************************************/
  /* Header files for RNGs currently included in package */
***************
*** 190,197 ****
  /* fill_buffer                       */
  /* Fill an RNG buffer                */
  /*************************************/
! void fill_buffer(rng)
!      RNGdata *rng;
  {
    long i,j; /* Loop variables */
    
--- 204,210 ----
  /* fill_buffer                       */
  /* Fill an RNG buffer                */
  /*************************************/
! void fill_buffer(RNGdata *rng)
  {
    long i,j; /* Loop variables */
    
***************
*** 291,298 ****
    inc_counts(rng,n)
  
  /* Auxiliary routines for 32-bit generators */
! double longtodouble(l)
!      long l;
  {
    if (l > 0)
      return((double)l);
--- 304,310 ----
    inc_counts(rng,n)
  
  /* Auxiliary routines for 32-bit generators */
! double longtodouble(long l)
  {
    if (l > 0)
      return((double)l);
***************
*** 330,341 ****
  /* drandomrv       */
  /* Returns doubles */
  /*******************/
! double drandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      double v[];
  {
-   long i;
    double vdef[1],*vreturn;
    
    if (RNGrange > MAXLONG + 1.0) {
--- 342,349 ----
  /* drandomrv       */
  /* Returns doubles */
  /*******************/
! double drandomrv(RNGdata *rng, long n, double v[])
  {
    double vdef[1],*vreturn;
    
    if (RNGrange > MAXLONG + 1.0) {
***************
*** 352,363 ****
  /* frandomrv      */
  /* Returns floats */
  /******************/
! float frandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      float v[];
  {
-   long i;
    float vdef[1],*vreturn;
    
    if (RNGrange > MAXLONG + 1.0) {
--- 360,367 ----
  /* frandomrv      */
  /* Returns floats */
  /******************/
! float frandomrv(RNGdata *rng, long n, float v[])
  {
    float vdef[1],*vreturn;
    
    if (RNGrange > MAXLONG + 1.0) {
***************
*** 374,383 ****
  /* lrandomrv     */
  /* Returns longs */
  /*****************/
! long lrandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      long v[];
  {
    long vdef[1],*vreturn;
  
--- 378,384 ----
  /* lrandomrv     */
  /* Returns longs */
  /*****************/
! long lrandomrv(RNGdata *rng, long n, long v[])
  {
    long vdef[1],*vreturn;
  
***************
*** 397,406 ****
  /*  i.e. each generate from the underlying*/
  /*  RNG is used to generate one bit.      */
  /******************************************/
! int brandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      int v[];
  {
    long i,j;
    int vdef[1],*vreturn;
--- 398,404 ----
  /*  i.e. each generate from the underlying*/
  /*  RNG is used to generate one bit.      */
  /******************************************/
! int brandomrv(RNGdata *rng, long n, int v[])
  {
    long i,j;
    int vdef[1],*vreturn;
***************
*** 429,438 ****
  /*  is not 2^32.                          */
  /******************************************/
  /* Splitting of this fast bit stream is not currently supported. */
! int brandomrv_f(rng,n,v)
!      RNGdata *rng;
!      long n;
!      int v[];
  {
    long i,j,r;
    int vdef[1],*vreturn;
--- 427,433 ----
  /*  is not 2^32.                          */
  /******************************************/
  /* Splitting of this fast bit stream is not currently supported. */
! int brandomrv_f(RNGdata *rng, long n, int v[])
  {
    long i,j,r;
    int vdef[1],*vreturn;
***************
*** 462,475 ****
  /* dxrandomrv      */
  /* Returns doubles */
  /*******************/
! double dxrandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      double v[];
  {
!   long i,j;
    double vdef[1],*vreturn;
-   double d;
  
    xrandom_setup();
    if (RNGreturns == RET_DOUBLE) {
--- 457,466 ----
  /* dxrandomrv      */
  /* Returns doubles */
  /*******************/
! double dxrandomrv(RNGdata *rng, long n, double v[])
  {
!   long j;
    double vdef[1],*vreturn;
  
    xrandom_setup();
    if (RNGreturns == RET_DOUBLE) {
***************
*** 495,506 ****
  /* fxrandomrv     */
  /* Returns floats */
  /******************/
! float fxrandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      float v[];
  {
!   long i,j;
    float vdef[1],*vreturn;
    
    xrandom_setup();
--- 486,494 ----
  /* fxrandomrv     */
  /* Returns floats */
  /******************/
! float fxrandomrv(RNGdata *rng, long n, float v[])
  {
!   long j;
    float vdef[1],*vreturn;
    
    xrandom_setup();
***************
*** 527,538 ****
  /* lxrandomrv    */
  /* Returns longs */
  /*****************/
! long lxrandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      long v[];
  {
!   long i,j;
    long vdef[1],*vreturn;
    
    xrandom_setup();
--- 515,523 ----
  /* lxrandomrv    */
  /* Returns longs */
  /*****************/
! long lxrandomrv(RNGdata *rng, long n, long v[])
  {
!   long j;
    long vdef[1],*vreturn;
    
    xrandom_setup();
***************
*** 559,570 ****
  /* bxrandomrv                  */
  /* Return "high quality" bits. */
  /*******************************/
! int bxrandomrv(rng,n,v)
!      RNGdata *rng;
!      long n;
!      int v[];
  {
!   long i,j;
    int vdef[1],*vreturn;
    
    xrandom_setup();
--- 544,552 ----
  /* bxrandomrv                  */
  /* Return "high quality" bits. */
  /*******************************/
! int bxrandomrv(RNGdata *rng, long n, int v[])
  {
!   long j;
    int vdef[1],*vreturn;
    
    xrandom_setup();
***************
*** 597,610 ****
   *  of bits.
   * Splitting of this fast bit stream is not currently supported.
   */
! int bxrandomrv_f(rng,n,v)
!      RNGdata *rng;
!      long n;
!      int v[];
  {
!   long i,j;
    int vdef[1],*vreturn;
!   long value;
    long index;
    
    inc_counts(rng,n);
--- 579,589 ----
   *  of bits.
   * Splitting of this fast bit stream is not currently supported.
   */
! int bxrandomrv_f(RNGdata *rng, long n, int v[])
  {
!   long i;
    int vdef[1],*vreturn;
!   long value=0;
    long index;
    
    inc_counts(rng,n);
***************
*** 628,636 ****
  /* seed_rng                          */
  /* Seeds the underlying RNG          */
  /*************************************/
! void seed_rng(rng, seed)
!      RNGdata *rng;
!      long *seed;
  {
    long i;
    
--- 607,613 ----
  /* seed_rng                          */
  /* Seeds the underlying RNG          */
  /*************************************/
! void seed_rng(RNGdata *rng, long *seed)
  {
    long i;
    
***************
*** 654,662 ****
  /*  stderr if alg is out of range    */
  /* Return 1 otherwise                */
  /*************************************/
! int setRNGparams(rng,alg)
!      RNGdata *rng;
!      long alg;
  {
    if (alg<0 || alg>NUM_RNGS-1) /* alg is out of range */
      return(0);
--- 631,637 ----
  /*  stderr if alg is out of range    */
  /* Return 1 otherwise                */
  /*************************************/
! int setRNGparams(RNGdata *rng, long alg)
  {
    if (alg<0 || alg>NUM_RNGS-1) /* alg is out of range */
      return(0);
***************
*** 677,684 ****
  /* Check the integrity of the RNG.   */
  /* Return 1 if ok, 0 if not.         */
  /*************************************/
! int check_rng(rng)
!      RNGdata *rng;
  {
    /* Temporary variables */
    long statesize, seedsize;
--- 652,658 ----
  /* Check the integrity of the RNG.   */
  /* Return 1 if ok, 0 if not.         */
  /*************************************/
! int check_rng(RNGdata *rng)
  {
    /* Temporary variables */
    long statesize, seedsize;
***************
*** 724,732 ****
  /* Currently only supports RNGs with */
  /*  at most two seeds.               */
  /*************************************/
! char *describe_rng (rng,rngid)
!      RNGdata *rng;
!      char rngid[RNGIDSTRLEN];
  {
    if (!check_rng(rng)) {
      BADRNG_ERR;
--- 698,704 ----
  /* Currently only supports RNGs with */
  /*  at most two seeds.               */
  /*************************************/
! char *describe_rng (RNGdata *rng, char rngid[RNGIDSTRLEN])
  {
    if (!check_rng(rng)) {
      BADRNG_ERR;
***************
*** 750,758 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int split_rng(rng,new_value)
!      RNGdata *rng;
!      long new_value;
  {
    if (check_rng(rng)) {
      if (new_value < 0)
--- 722,728 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int split_rng(RNGdata *rng, long new_value)
  {
    if (check_rng(rng)) {
      if (new_value < 0)
***************
*** 776,784 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int mralg_rng(rng,new_value)
!      RNGdata *rng;
!      long new_value;
  {
    if (check_rng(rng)) {
      if (new_value == 0 || new_value == 1) {
--- 746,752 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int mralg_rng(RNGdata *rng, long new_value)
  {
    if (check_rng(rng)) {
      if (new_value == 0 || new_value == 1) {
***************
*** 802,809 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! double range_rng(rng)
!      RNGdata *rng;
  {
    if (check_rng(rng))
      return(RNGrange);
--- 770,776 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! double range_rng(RNGdata *rng)
  {
    if (check_rng(rng))
      return(RNGrange);
***************
*** 813,822 ****
  }
  
  /* Auxiliary routine for allocating memory for an RNG */
! RNGdata *allocate_rng(rng, alg, bufsize)
!      RNGdata *rng;
!      long alg, bufsize;
  {
    if (bufsize <= 0)
      bufsize=1;
    if ((rng = (RNGdata *)malloc(sizeof(RNGdata))) == 0)
--- 780,788 ----
  }
  
  /* Auxiliary routine for allocating memory for an RNG */
! RNGdata *allocate_rng(long alg, long bufsize)
  {
+   RNGdata *rng;
    if (bufsize <= 0)
      bufsize=1;
    if ((rng = (RNGdata *)malloc(sizeof(RNGdata))) == 0)
***************
*** 931,948 ****
   * replicating experiments or debugging, it is better to restart an
   * existing generator (stored in a statefile) than to initialize a new one.
   */
! RNGdata *init_rng(alg, mrandom_alg, seed, count1, count2, bufsize)
!      long alg;         /* The RNG algorithm to use */
!      long mrandom_alg; /* Algorithm to use for mrandom:
  			  0 = Thomborson's slow but unbasied conversion
  			  1 = (int) (m*drandom()) */
!      long *seed;          /* Seed */
!      long count1, count2; /* Initial counts */
!      long bufsize;        /* Size of buffer */
  {
    RNGdata *rng;
    double x[64]; /* Temporary vector */
-   long i;       /* Loop variable    */
    
    /* Are count1 and count2 in range? */
    if (count1 >= BILLION || count1 < 0 || count2 < 0) {
--- 897,914 ----
   * replicating experiments or debugging, it is better to restart an
   * existing generator (stored in a statefile) than to initialize a new one.
   */
! RNGdata *init_rng(
!      long alg,         /* The RNG algorithm to use */
!      long mrandom_alg, /* Algorithm to use for mrandom:
  			  0 = Thomborson's slow but unbasied conversion
  			  1 = (int) (m*drandom()) */
!      long *seed,          /* Seed */
!      long count1,	  /* Initial counts */ 
!      long count2,
!      long bufsize)        /* Size of buffer */
  {
    RNGdata *rng;
    double x[64]; /* Temporary vector */
    
    /* Are count1 and count2 in range? */
    if (count1 >= BILLION || count1 < 0 || count2 < 0) {
***************
*** 953,959 ****
      fprintf(stderr, "Warning: this initialization will take a LONG time!\n");
      fflush(stderr);
    }
!   if ((rng=allocate_rng(rng,alg,bufsize))==0)
      return(0); /* Allocate space for the RNG */
    seed_rng(rng, seed);       /* Seed the RNG */
    if (mralg_rng(rng, mrandom_alg) == 0)
--- 919,925 ----
      fprintf(stderr, "Warning: this initialization will take a LONG time!\n");
      fflush(stderr);
    }
!   if ((rng=allocate_rng(alg,bufsize))==0)
      return(0); /* Allocate space for the RNG */
    seed_rng(rng, seed);       /* Seed the RNG */
    if (mralg_rng(rng, mrandom_alg) == 0)
***************
*** 980,989 ****
  /* Returns 0 if kill failed    */
  /* Returns 1 otherwise         */
  /*******************************/
! int kill_rng(rng)
!      RNGdata *rng;
  {
    if (check_rng(rng)) {
      free(RNGstate);
      free(RNGseed);
      free(RNGbuffer.lbuf);
--- 946,959 ----
  /* Returns 0 if kill failed    */
  /* Returns 1 otherwise         */
  /*******************************/
! int kill_rng(RNGdata *rng)
  {
    if (check_rng(rng)) {
+ 
+     /* Special fixup for 4.3bsd random(): deactivate current 
+        state before freeing its memory */
+     if (RNGalg == 1) random_deactivate_state();
+ 
      free(RNGstate);
      free(RNGseed);
      free(RNGbuffer.lbuf);
***************
*** 1001,1018 ****
  /* Return the next lxrandomrv() output */
  /*  without disturbing the RNG state.  */
  /***************************************/
! long nextval(rng)
!      RNGdata *rng;
  {
!   static long *state=0;
!   long i, r, tcount1, tcount2, retval;
!   
!   if (state != 0)
!     free(state);
!   state = (long *)calloc(RNGstatesize,sizeof(long));
!   /* Check that this only allocates the first time */
    
    /* Preserve old RNG state */
    tcount1 = RNGcount1;
    tcount2 = RNGcount2;
    for (i=0; i<RNGstatesize; i++) {
--- 971,987 ----
  /* Return the next lxrandomrv() output */
  /*  without disturbing the RNG state.  */
  /***************************************/
! long nextval(RNGdata *rng)
  {
!   long *state=0;
!   long i, tcount1, tcount2, retval;
    
+   /* Special fixup for 4.3bsd random(): deactivate
+      state vector before saving it */
+   if (RNGalg == 1) random_deactivate_state();
+ 
    /* Preserve old RNG state */
+   state = (long *)calloc(RNGstatesize,sizeof(long));
    tcount1 = RNGcount1;
    tcount2 = RNGcount2;
    for (i=0; i<RNGstatesize; i++) {
***************
*** 1021,1048 ****
    
    /* Find the next value in this pseudorandom sequence */
    retval=lxrandomr(rng);
!   
    /* Restore the RNG state */
    RNGcount1 = tcount1;
    RNGcount2 = tcount2;
-   /* Special fixup for 4.3bsd random()'s hidden state variables */
-   if (RNGalg == 1) for (i=1; i<31; i++) {
-     r = random();
-   }
    for (i=0; i<RNGstatesize; i++) {
      RNGstate[i] = state[i];
    }
    
    return(retval);
  }
  
  /* printvector */
  /* auxiliary routine for save_rng */
! void printvector(fp,size,vector,type)
!      FILE *fp;
!      long size;
!      long *vector;
!      int type;
  {
    long i,j,printval;
    long windex,bindex;
--- 990,1018 ----
    
    /* Find the next value in this pseudorandom sequence */
    retval=lxrandomr(rng);
! 
!   /* Special fixup for 4.3bsd random(): deactivate
!      state vector before updating it */
!   if (RNGalg == 1) random_deactivate_state();
! 
    /* Restore the RNG state */
    RNGcount1 = tcount1;
    RNGcount2 = tcount2;
    for (i=0; i<RNGstatesize; i++) {
      RNGstate[i] = state[i];
    }
    
+   free(state);
    return(retval);
  }
  
  /* printvector */
  /* auxiliary routine for save_rng */
! void printvector(
!      FILE *fp,
!      long size,
!      long *vector,
!      int type)
  {
    long i,j,printval;
    long windex,bindex;
***************
*** 1077,1087 ****
  
  /* scanvector */
  /* auxiliary routine for restart_rng */
! void scanvector(fp,size,vector,type)
!      FILE *fp;
!      long size;
!      long *vector;
!      int type;
  {
    long i,j,scanval;
    long windex,bindex;
--- 1047,1057 ----
  
  /* scanvector */
  /* auxiliary routine for restart_rng */
! void scanvector(
!      FILE *fp,
!      long size,
!      long *vector,
!      int type)
  {
    long i,j,scanval;
    long windex,bindex;
***************
*** 1090,1096 ****
    bindex=0;
    i=0;
  
!   while (windex<size) {
      fscanf(fp,RNGfileLINEnn[type],&scanval);
      for (j=0; j < (type << 3); j++) {
        if (bindex-- == 0) {
--- 1060,1066 ----
    bindex=0;
    i=0;
  
!   while (bindex!=0 || windex+1<size) {
      fscanf(fp,RNGfileLINEnn[type],&scanval);
      for (j=0; j < (type << 3); j++) {
        if (bindex-- == 0) {
***************
*** 1114,1126 ****
  /* Return 0 if RNG couldn't be saved.        */
  /* Returns 1 otherwise.                      */
  /*********************************************/
! int save_rng(rng,filename)
!      RNGdata *rng;
!      char *filename;
  {
    FILE *fp;
-   long i;
-   unsigned short *sarray;
    RNGdata *temp_rng;
  
    if (!check_rng(rng))
--- 1084,1092 ----
  /* Return 0 if RNG couldn't be saved.        */
  /* Returns 1 otherwise.                      */
  /*********************************************/
! int save_rng(RNGdata *rng, char *filename)
  {
    FILE *fp;
    RNGdata *temp_rng;
  
    if (!check_rng(rng))
***************
*** 1166,1194 ****
  /*     - Declare a new RNG                                      */
  /*     - Kill a previously initialized RNG using kill_rng()     */
  /****************************************************************/
! RNGdata *restart_rng(filename)
!      char *filename;
  {
    FILE *fp;
!   unsigned short *sarray; /* for reading alg #4's statefile */
!   long i, m, r; /* temps */
    long alg, bufsize;
    int errflag; /* initially 0, becomes 1 if a problem is discovered */
    RNGdata *rng;
    
-   /* If mru_rng == 0, we assume random() hasn't been called yet,
-    * so its internal state is at its startup values.
-    * If mru_rng != 0 and alg=1, then we've used random() already, so we
-    * must reset its internal state to its startup values.
-    */
-   if (mru_rng != 0 && mru_rng->rngalg == 1) {
-     m = (mru_rng->rngcount1%31) + mru_rng->rngcount2*(BILLION%31);
-     /* note: the hidden state variables are counters mod 31 */
-     for ( ; (m%31) != 0; m++) {
-       r = random();
-     }
-   }
-   
    /* Restore counter values, retrieve original RNG seeds and current state */
    fp = fopen(filename, "r");
    if (!fp)
--- 1132,1145 ----
  /*     - Declare a new RNG                                      */
  /*     - Kill a previously initialized RNG using kill_rng()     */
  /****************************************************************/
! RNGdata *restart_rng(char *filename)
  {
    FILE *fp;
!   long i;			/* temps */
    long alg, bufsize;
    int errflag; /* initially 0, becomes 1 if a problem is discovered */
    RNGdata *rng;
    
    /* Restore counter values, retrieve original RNG seeds and current state */
    fp = fopen(filename, "r");
    if (!fp)
***************
*** 1197,1203 ****
    fscanf(fp, RNGfileLINE0, &alg);
    fscanf(fp,RNGname_a[alg]);
    fscanf(fp, RNGfileLINE1, &bufsize);
!   if ((rng=allocate_rng(rng, alg, bufsize))==0)
      return(0);
    fscanf(fp, RNGfileLINE2);
    scanvector(fp,RNGseedsize,RNGseed,RNGstatetype);
--- 1148,1154 ----
    fscanf(fp, RNGfileLINE0, &alg);
    fscanf(fp,RNGname_a[alg]);
    fscanf(fp, RNGfileLINE1, &bufsize);
!   if ((rng=allocate_rng(alg, bufsize))==0)
      return(0);
    fscanf(fp, RNGfileLINE2);
    scanvector(fp,RNGseedsize,RNGseed,RNGstatetype);
***************
*** 1217,1222 ****
--- 1168,1178 ----
      RNGdata *test_rng;
      test_rng = init_rng(RNGalg, RNGmrandom_alg, RNGseed,
  			RNGcount1, RNGcount2, RNGbuffer.size);
+ 
+     /* Special fixup for 4.3bsd random(): deactivate state vector
+        before comparing it */
+     if (test_rng->rngalg == 1) random_deactivate_state();
+ 
      /* See if we got to the same state */
      for (i=0; i<RNGstatesize; i++) {
        if ((test_rng->rngstate)[i] != RNGstate[i]) {
***************
*** 1224,1244 ****
        }
      }
      kill_rng(test_rng);
-   } else { /* Just update mru_rng */
-     /* first, some special hacks for alg 1 */
-     if (RNGalg == 1) {
-       /* tell random() we've got a 31-word statefile */
-       RNGstate[0] = MRNGSTATE0;
-       setstate(RNGstate);
-       /* and modify random()'s hidden state to correspond to the RNGcount */
-       m = RNGcount1%31 + RNGcount2*(BILLION%31);
-       for (i=0 ; i<(m%31); i++) {
-         r = random();
-       }
-     }
-     mru_rng = rng; /* remember this rng, for use by mrandom() and frandom() */
    }
!   
    /* see if RNG state looks ok */
    errflag = errflag || !check_rng(rng);
  
--- 1180,1188 ----
        }
      }
      kill_rng(test_rng);
    }
!   mru_rng = rng; /* remember this rng, for use by mrandom() and frandom() */
! 
    /* see if RNG state looks ok */
    errflag = errflag || !check_rng(rng);
  
***************
*** 1259,1266 ****
  /* Flush the contents of the RNG's buffers   */
  /* Returns 1 upon success, 0 upon failure    */
  /*********************************************/
! int flush_rng(rng)
!      RNGdata *rng;
  {
    if (check_rng(rng)) {
      if (RNGreturns == RET_LONG) {
--- 1203,1209 ----
  /* Flush the contents of the RNG's buffers   */
  /* Returns 1 upon success, 0 upon failure    */
  /*********************************************/
! int flush_rng(RNGdata *rng)
  {
    if (check_rng(rng)) {
      if (RNGreturns == RET_LONG) {
***************
*** 1309,1326 ****
   * bias-removal algorithm makes an excessive number of calls to the
   * underlying generator.
   */
! long mrandomrv(rng,m,n,v)
!      RNGdata *rng;
!      long m,n,v[];
  {
    long i,yield; /* temps */
    long ncalls; /* counts number of RNG calls during this mrandomrv call */
    long vdef[1],*vreturn; /* default value for v */
    char rngdesc[RNGIDSTRLEN]; /* for a describe_rng() call */
-   double d; /* temporary variable */
  
    /* we can avoid some calcs&tests if m and rng range are unchanged */ 
!   static long lastm=0,lastrangem1=0; /* m on last call */
    static double lastrange=0;         /* range on last call */
    static double maxv;  /* largest "useful" RNG output for this m */
    static double ifreq; /* multiplicity of RNG -> mrandomr() mapping */
--- 1252,1266 ----
   * bias-removal algorithm makes an excessive number of calls to the
   * underlying generator.
   */
! long mrandomrv(RNGdata *rng, long m, long n, long v[])
  {
    long i,yield; /* temps */
    long ncalls; /* counts number of RNG calls during this mrandomrv call */
    long vdef[1],*vreturn; /* default value for v */
    char rngdesc[RNGIDSTRLEN]; /* for a describe_rng() call */
  
    /* we can avoid some calcs&tests if m and rng range are unchanged */ 
!   static long lastm=0;		     /* m on last call */
    static double lastrange=0;         /* range on last call */
    static double maxv;  /* largest "useful" RNG output for this m */
    static double ifreq; /* multiplicity of RNG -> mrandomr() mapping */
***************
*** 1429,1434 ****
--- 1369,1375 ----
      return(*vreturn);
      break;
    default:
+     return 0;
      break;
    }
  }
*** 3.1	1993/05/28 19:25:37
--- mrandom.h	1996/05/31 21:42:43
***************
*** 1,4 ****
! /* mrandom.h 3.1 5/28/93 */
  /*
   *		mrandom.h
   *
--- 1,4 ----
! /* mrandom.h 3.1x 5/31/96, derived from mrandom.h 3.1 5/28/93 */
  /*
   *		mrandom.h
   *
***************
*** 27,32 ****
--- 27,40 ----
   *         - Allow choice of algorithms in mrandomrv
   *         - Dynamic allocation of RNGs
   *
+  *      Version 3.1x:
+  *              Michael Fischer, May 1996
+  *         - Ported to NeXT; removed unused variables
+  *         - Converted function headers to ISO-C
+  *         - Fixed memory-overwrite bug in scanvector()
+  *         - Fixed state saving/restoring for 4.3bsd random()
+  *         - Cleaned up nextval() and restart_rng().
+  *
   *	This material is based upon work supported by the National
   *	Science Foundation under grant number MIP-9023238.  The
   *	Government has certain rights in this material.
***************
*** 175,181 ****
  /* drandomrv       */
  /* Returns doubles */
  /*******************/
! double drandomrv(/*RNGdata *rng, long n, double v[]*/);
  #define drandomr(rng)  drandomrv(rng,0,0)
  #define drandomv(n,v)  drandomrv(0,n,v)
  #define drandom()      drandomrv(0,0,0)
--- 183,189 ----
  /* drandomrv       */
  /* Returns doubles */
  /*******************/
! double drandomrv(RNGdata *rng, long n, double v[]);
  #define drandomr(rng)  drandomrv(rng,0,0)
  #define drandomv(n,v)  drandomrv(0,n,v)
  #define drandom()      drandomrv(0,0,0)
***************
*** 184,190 ****
  /* frandomrv      */
  /* Returns floats */
  /******************/
! float frandomrv(/*RNGdata *rng, long n, float v[]*/);
  #define frandomr(rng)  frandomrv(rng,0,0)
  #define frandomv(n,v)  frandomrv(0,n,v)
  #define frandom()      frandomrv(0,0,0)
--- 192,198 ----
  /* frandomrv      */
  /* Returns floats */
  /******************/
! float frandomrv(RNGdata *rng, long n, float v[]);
  #define frandomr(rng)  frandomrv(rng,0,0)
  #define frandomv(n,v)  frandomrv(0,n,v)
  #define frandom()      frandomrv(0,0,0)
***************
*** 193,199 ****
  /* lrandomrv     */
  /* Returns longs */
  /*****************/
! long lrandomrv(/*RNGdata *rng, long n, long v[]*/);
  #define lrandomr(rng)  lrandomrv(rng,0,0)
  #define lrandomv(n,v)  lrandomrv(0,n,v)
  #define lrandom()      lrandomrv(0,0,0)
--- 201,207 ----
  /* lrandomrv     */
  /* Returns longs */
  /*****************/
! long lrandomrv(RNGdata *rng, long n, long v[]);
  #define lrandomr(rng)  lrandomrv(rng,0,0)
  #define lrandomv(n,v)  lrandomrv(0,n,v)
  #define lrandom()      lrandomrv(0,0,0)
***************
*** 204,210 ****
  /*  i.e. each generate from the underlying*/
  /*  RNG is used to generate one bit.      */
  /******************************************/
! int brandomrv(/*RNGdata *rng, long n, int v[]*/);
  #define brandomr(rng)  brandomrv(rng,0,0)
  #define brandomv(n,v)  brandomrv(0,n,v)
  #define brandom()      brandomrv(0,0,0)
--- 212,218 ----
  /*  i.e. each generate from the underlying*/
  /*  RNG is used to generate one bit.      */
  /******************************************/
! int brandomrv(RNGdata *rng, long n, int v[]);
  #define brandomr(rng)  brandomrv(rng,0,0)
  #define brandomv(n,v)  brandomrv(0,n,v)
  #define brandom()      brandomrv(0,0,0)
***************
*** 218,224 ****
  /*  is not 2^32.                          */
  /******************************************/
  /* Splitting of this fast bit stream is not currently supported. */
! int brandomrv_f(/*RNGdata *rng, long n, int v[]*/);
  #define brandomr_f(rng)  brandomrv_f(rng,0,0)
  #define brandomv_f(n,v)  brandomrv_f(0,n,v)
  #define brandom_f()      brandomrv_f(0,0,0)
--- 226,232 ----
  /*  is not 2^32.                          */
  /******************************************/
  /* Splitting of this fast bit stream is not currently supported. */
! int brandomrv_f(RNGdata *rng, long n, int v[]);
  #define brandomr_f(rng)  brandomrv_f(rng,0,0)
  #define brandomv_f(n,v)  brandomrv_f(0,n,v)
  #define brandom_f()      brandomrv_f(0,0,0)
***************
*** 230,236 ****
  /* dxrandomrv      */
  /* Returns doubles */
  /*******************/
! double dxrandomrv(/*RNGdata *rng, long n, double v[]*/);
  #define dxrandomr(rng) dxrandomrv(rng,0,0)
  #define dxrandomv(n,v) dxrandomrv(0,n,v)
  #define dxrandom()     dxrandomrv(0,0,0)
--- 238,244 ----
  /* dxrandomrv      */
  /* Returns doubles */
  /*******************/
! double dxrandomrv(RNGdata *rng, long n, double v[]);
  #define dxrandomr(rng) dxrandomrv(rng,0,0)
  #define dxrandomv(n,v) dxrandomrv(0,n,v)
  #define dxrandom()     dxrandomrv(0,0,0)
***************
*** 239,245 ****
  /* fxrandomrv     */
  /* Returns floats */
  /******************/
! float fxrandomrv(/*RNGdata *rng, long n, float v[]*/);
  #define fxrandomr(rng) fxrandomrv(rng,0,0)
  #define fxrandomv(n,v) fxrandomrv(0,n,v)
  #define fxrandom()     fxrandomrv(0,0,0)
--- 247,253 ----
  /* fxrandomrv     */
  /* Returns floats */
  /******************/
! float fxrandomrv(RNGdata *rng, long n, float v[]);
  #define fxrandomr(rng) fxrandomrv(rng,0,0)
  #define fxrandomv(n,v) fxrandomrv(0,n,v)
  #define fxrandom()     fxrandomrv(0,0,0)
***************
*** 248,254 ****
  /* lxrandomrv    */
  /* Returns longs */
  /*****************/
! long lxrandomrv(/*RNGdata *rng, long n, long v[]*/);
  #define lxrandomr(rng) lxrandomrv(rng,0,0)
  #define lxrandomv(n,v) lxrandomrv(0,n,v)
  #define lxrandom()     lxrandomrv(0,0,0)
--- 256,262 ----
  /* lxrandomrv    */
  /* Returns longs */
  /*****************/
! long lxrandomrv(RNGdata *rng, long n, long v[]);
  #define lxrandomr(rng) lxrandomrv(rng,0,0)
  #define lxrandomv(n,v) lxrandomrv(0,n,v)
  #define lxrandom()     lxrandomrv(0,0,0)
***************
*** 257,263 ****
  /* bxrandomrv                  */
  /* Return "high quality" bits. */
  /*******************************/
! int bxrandomrv(/*RNGdata *rng, long n, int v[]*/);
  #define bxrandomr(rng) bxrandomrv(rng,0,0)
  #define bxrandomv(n,v) bxrandomrv(0,n,v)
  #define bxrandom()     bxrandomrv(0,0,0)
--- 265,271 ----
  /* bxrandomrv                  */
  /* Return "high quality" bits. */
  /*******************************/
! int bxrandomrv(RNGdata *rng, long n, int v[]);
  #define bxrandomr(rng) bxrandomrv(rng,0,0)
  #define bxrandomv(n,v) bxrandomrv(0,n,v)
  #define bxrandom()     bxrandomrv(0,0,0)
***************
*** 275,281 ****
   *  of bits.
   * Splitting of this fast bit stream is not currently supported.
   */
! int bxrandomrv_f(/*RNGdata *rng, long n, int v[]*/);
  #define bxrandomr_f(rng) bxrandomrv_f(rng,0,0)
  #define bxrandomv_f(n,v) bxrandomrv_f(0,n,v)
  #define bxrandom_f()     bxrandomrv_f(0,0,0)
--- 283,289 ----
   *  of bits.
   * Splitting of this fast bit stream is not currently supported.
   */
! int bxrandomrv_f(RNGdata *rng, long n, int v[]);
  #define bxrandomr_f(rng) bxrandomrv_f(rng,0,0)
  #define bxrandomv_f(n,v) bxrandomrv_f(0,n,v)
  #define bxrandom_f()     bxrandomrv_f(0,0,0)
***************
*** 284,297 ****
  /* seed_rng                 */
  /* Seeds the underlying RNG */
  /****************************/
! void seed_rng(/*RNGdata *rng, long *seed*/);
  
  /*************************************/
  /* check_rng                         */
  /* Check the integrity of the RNG.   */
  /* Return 1 if ok, 0 if not.         */
  /*************************************/
! int check_rng(/*RNGdata *rng*/);
  
  /*************************************/
  /* describe_rng                      */
--- 292,305 ----
  /* seed_rng                 */
  /* Seeds the underlying RNG */
  /****************************/
! void seed_rng(RNGdata *rng, long *seed);
  
  /*************************************/
  /* check_rng                         */
  /* Check the integrity of the RNG.   */
  /* Return 1 if ok, 0 if not.         */
  /*************************************/
! int check_rng(RNGdata *rng);
  
  /*************************************/
  /* describe_rng                      */
***************
*** 307,313 ****
  /* Currently only supports RNGs with */
  /*  at most two seeds.               */
  /*************************************/
! char *describe_rng (/*RNGdata *rng, char rngid[RNGIDSTRLEN]*/);
  
  /*************************************/
  /* split_rng                         */
--- 315,321 ----
  /* Currently only supports RNGs with */
  /*  at most two seeds.               */
  /*************************************/
! char *describe_rng (RNGdata *rng, char rngid[RNGIDSTRLEN]);
  
  /*************************************/
  /* split_rng                         */
***************
*** 321,327 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int split_rng(/* RNGdata *rng, long new_value*/);
  
  /*************************************/
  /* mralg_rng                         */
--- 329,335 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int split_rng(RNGdata *rng, long new_value);
  
  /*************************************/
  /* mralg_rng                         */
***************
*** 332,338 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int mralg_rng(/* RNGdata *rng, long new_value*/);
  
  /*************************************/
  /* range_rng                         */
--- 340,346 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! int mralg_rng(RNGdata *rng, long new_value);
  
  /*************************************/
  /* range_rng                         */
***************
*** 343,349 ****
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! double range_rng(/*RNGdata *rng*/);
  
  /************/
  /* init_rng */
--- 351,357 ----
  /* Exits on error if rng has not     */
  /*  been initialized                 */
  /*************************************/
! double range_rng(RNGdata *rng);
  
  /************/
  /* init_rng */
***************
*** 426,433 ****
   * replicating experiments or debugging, it is better to restart an
   * existing generator (stored in a statefile) than to initialize a new one.
   */
! RNGdata *init_rng(/* long alg, long mrandom_alg, long *seed,
! 		     long count1, long count2, long bufsize*/);
  
  /*******************************/
  /* kill_rng                    */
--- 434,441 ----
   * replicating experiments or debugging, it is better to restart an
   * existing generator (stored in a statefile) than to initialize a new one.
   */
! RNGdata *init_rng( long alg, long mrandom_alg, long *seed,
! 		   long count1, long count2, long bufsize);
  
  /*******************************/
  /* kill_rng                    */
***************
*** 435,441 ****
  /* Returns 0 if kill failed    */
  /* Returns 1 otherwise         */
  /*******************************/
! int kill_rng(/*RNGdata *rng*/);
  
  /*********************************************/
  /* save_rng                                  */
--- 443,449 ----
  /* Returns 0 if kill failed    */
  /* Returns 1 otherwise         */
  /*******************************/
! int kill_rng(RNGdata *rng);
  
  /*********************************************/
  /* save_rng                                  */
***************
*** 443,449 ****
  /* Return 0 if RNG couldn't be saved.        */
  /* Returns 1 otherwise.                      */
  /*********************************************/
! int save_rng(/* RNGdata *rng, char *filename*/);
  
  /****************************************************************/
  /* restart_rng                                                  */
--- 451,457 ----
  /* Return 0 if RNG couldn't be saved.        */
  /* Returns 1 otherwise.                      */
  /*********************************************/
! int save_rng(RNGdata *rng, char *filenamed);
  
  /****************************************************************/
  /* restart_rng                                                  */
***************
*** 458,471 ****
  /*     - Declare a new RNG                                      */
  /*     - Kill a previously initialized RNG using kill_rng()     */
  /****************************************************************/
! RNGdata *restart_rng(/* char *filename*/);
  
  /*********************************************/
  /* flush_rng                                 */
  /* Flush the contents of the RNG's buffers   */
  /* Returns 1 upon success, 0 upon failure    */
  /*********************************************/
! int flush_rng(/*RNGdata *rng*/);
  
  /*************/
  /* mrandomrv */
--- 466,479 ----
  /*     - Declare a new RNG                                      */
  /*     - Kill a previously initialized RNG using kill_rng()     */
  /****************************************************************/
! RNGdata *restart_rng(char *filename);
  
  /*********************************************/
  /* flush_rng                                 */
  /* Flush the contents of the RNG's buffers   */
  /* Returns 1 upon success, 0 upon failure    */
  /*********************************************/
! int flush_rng(RNGdata *rng);
  
  /*************/
  /* mrandomrv */
***************
*** 496,502 ****
   * bias-removal algorithm makes an excessive number of calls to the
   * underlying generator.
   */
! long mrandomrv (/* RNGdata *rng, long m, long n, long v*/);
  #define mrandomr(rng,m)  mrandomrv(rng,m,0,0)
  #define mrandomv(m,n,v)  mrandomrv(0,m,n,v)
  #define mrandom(m)       mrandomrv(0,m,0,0)
--- 504,510 ----
   * bias-removal algorithm makes an excessive number of calls to the
   * underlying generator.
   */
! long mrandomrv (RNGdata *rng, long m, long n, long v[]);
  #define mrandomr(rng,m)  mrandomrv(rng,m,0,0)
  #define mrandomv(m,n,v)  mrandomrv(0,m,n,v)
  #define mrandom(m)       mrandomrv(0,m,0,0)
*** 3.1	1993/05/28 19:25:37
--- mrtest.c	1996/05/31 21:44:42
***************
*** 1,4 ****
! /* mrtest.c 3.1 5/28/93 */
  /*
   *		mrtest.c
   *
--- 1,4 ----
! /* mrtest.c 3.1x 5/31/96, derived from mrtest.c 3.1 5/28/93 */
  /*
   *		mrtest.c
   *
***************
*** 28,33 ****
--- 28,38 ----
   *      Robert Plotkin, May 1993
   *              -update for use with mrandom 3.0
   *
+  *      Michael Fischer, May 1996
+  *              -ported to NeXT: fixed includes, moved big arrays off stack
+  *              -minor bug fixes in error printouts
+  *              -changed function headers to ISO-C, eliminated unused vars
+  *
   *	Possible future additions:
   *		-more RNGs
   *		-faster (hash-based) large-m equidistribution test 
***************
*** 64,74 ****
   *	This code is neither copyrighted nor patented.
   */
  
  # include <sys/file.h> /* we use access() */
  # include "mrandom.h"
  # include "xsq.h" /* for interpretation of X-squared values */
  # include <math.h>
! # include <values.h> /* we need MAXLONG */
  
  # ifndef RNGFILENAME
  # define RNGFILENAME "RNGstatefile" /* where the RNG state is stored */
--- 69,87 ----
   *	This code is neither copyrighted nor patented.
   */
  
+ # include <stdio.h>
+ # include <stdlib.h>
+ # ifdef NeXT
+ # include <libc.h>     /* we use access() */
+ # else
  # include <sys/file.h> /* we use access() */
+ # endif
  # include "mrandom.h"
  # include "xsq.h" /* for interpretation of X-squared values */
  # include <math.h>
! # ifndef NeXT
! # include <values.h>   /* we need MAXLONG */
! # endif
  
  # ifndef RNGFILENAME
  # define RNGFILENAME "RNGstatefile" /* where the RNG state is stored */
***************
*** 137,147 ****
  
  # ifndef VECTOR
  /* slower code with runtime options */
! long random_value(rng,m,method)
! RNGdata *rng;
! long m,method;
  {
- 
    switch (method) {
      case NORMAL:
  	 return( mrandom(m) );
--- 150,157 ----
  
  # ifndef VECTOR
  /* slower code with runtime options */
! long random_value(RNGdata *rng, long m, long method)
  {
    switch (method) {
      case NORMAL:
  	 return( mrandom(m) );
***************
*** 220,228 ****
  		     */
    int test2tuples, test3tuples;
    /* various xsq stats */
!   double xsq,xsq2,xsq3, xsq2l, xsq3l;
    /* random generates: current, previous, previous^2, first, second */
!   long x, prevx, pprevx, firstx, secondx;
  
    /* var for max test */
    long maxx;
--- 230,238 ----
  		     */
    int test2tuples, test3tuples;
    /* various xsq stats */
!   double xsq,xsq2, xsq2l, xsq3l;
    /* random generates: current, previous, previous^2, first, second */
!   long x, prevx=0, pprevx=0, firstx=0, secondx=0;
  
    /* var for max test */
    long maxx;
***************
*** 235,246 ****
     * freq[i] = x[i], if testequidist == 3
     * freq[i] = undefined, otherwise
     */
!   long freq[EQUIDISTMAX];
    /* data for 2-tuple test */
!   long pairs[PAIRMAXM][PAIRMAXM];
    /* vars for 3-tuple test on low-order 3 bits */
    long lowpairs[8][8];
    long lowtrips[8][8][8];
  
    if(argc > 1) {
      for(i=1;i<argc;i++) {
--- 245,265 ----
     * freq[i] = x[i], if testequidist == 3
     * freq[i] = undefined, otherwise
     */
!   long* freq; 
    /* data for 2-tuple test */
!   long* pairs[PAIRMAXM];
    /* vars for 3-tuple test on low-order 3 bits */
    long lowpairs[8][8];
    long lowtrips[8][8][8];
+   int errflag = 0;
+ 
+   if ((freq = (long*)calloc(EQUIDISTMAX, sizeof(long))) == NULL) errflag++;
+   for (i=0; i<PAIRMAXM; i++)
+     if ((pairs[i] = (long*)calloc(PAIRMAXM, sizeof(long))) == NULL) errflag++;
+   if (errflag) {
+     fprintf(stderr,"Not enough memory\n");
+     exit(1);
+   }
  
    if(argc > 1) {
      for(i=1;i<argc;i++) {
***************
*** 259,265 ****
  	  count2 = 0;
  	  sscanf (&(argv[i][2]), "%ld,%ld,%ld,%ld,%ld,&ld,&ld",
  		 &alg, &seed1, &seed2, &count1, &count2);
- 	  seed[0]=seed1; seed[1]=seed2;
  	  reseeding = 1;
  	  if (seed1<0) {
  	    printf("Illegal value, %ld, for RNG seed.\n",seed1);
--- 278,283 ----
***************
*** 270,282 ****
  	    argerr((char*) argv[0]);
            }
  	  if (count1<0) {
! 	    printf("Illegal value, %ld, for number of times ");
! 	    printf(" to cycle rng before starting tests.\n",count1);
  	    argerr((char*) argv[0]);
            }
  	  if (count2<0) {
! 	    printf("Illegal value, %ld, for number of billions of times ");
! 	    printf(" to cycle rng before starting tests.\n",count1);
  	    argerr((char*) argv[0]);
            }
  	  break;
--- 288,300 ----
  	    argerr((char*) argv[0]);
            }
  	  if (count1<0) {
! 	    printf("Illegal value, %ld, for number of times "
! 	           " to cycle rng before starting tests.\n",count1);
  	    argerr((char*) argv[0]);
            }
  	  if (count2<0) {
! 	    printf("Illegal value, %ld, for number of billions of times "
! 	           " to cycle rng before starting tests.\n",count1);
  	    argerr((char*) argv[0]);
            }
  	  break;
***************
*** 341,346 ****
--- 359,366 ----
    }
    printf("\n");
  
+   seed[0]=seed1; seed[1]=seed2;
+ 
    dn = (double) n;
  
    if (m > 0) {
***************
*** 525,537 ****
        if (m < 0 && x < 0) { /* print an unsigned long */
          printf("  %.0f\n", (double)(x&MAXLONG)+(double)MAXLONG+1.0);
        } else if (m <= 10) {
!         printf("  %d\n", x);
        } else if (m <= 100) {
!         printf("  %2d\n", x);
        } else if (m <= 1000) {
!         printf("  %3d\n", x);
        } else {
!         printf("  %4d\n", x);
        }
      }
      /* now discard some generates: d=2^k, k>9, gives nrand48() trouble */ 
--- 545,557 ----
        if (m < 0 && x < 0) { /* print an unsigned long */
          printf("  %.0f\n", (double)(x&MAXLONG)+(double)MAXLONG+1.0);
        } else if (m <= 10) {
!         printf("  %ld\n", x);
        } else if (m <= 100) {
!         printf("  %2ld\n", x);
        } else if (m <= 1000) {
!         printf("  %3ld\n", x);
        } else {
!         printf("  %4ld\n", x);
        }
      }
      /* now discard some generates: d=2^k, k>9, gives nrand48() trouble */ 
*** 3.1	1993/05/28 19:25:41
--- rngs.h	1996/06/04 21:19:46
***************
*** 7,13 ****
  #endif
  
  /* RNGstate[0] for random(), table length 32 */
! #define MRNGSTATE0	3
  
  /* RNG state sizes */
  #define RNGstatesize_0 1
--- 7,26 ----
  #endif
  
  /* RNGstate[0] for random(), table length 32 */
! #define MRNGTYPE	3  /* bsd4.3 random() alg type (for 128-byte buffer) */
! #define MRNG_MAX_TYPES  5  /* number of types supported by bsd4.3 random()   */
!                            /* state[0]%MAX_TYPES should equal MRNGTYPE       */
! 
! /* Globals pointing to original and current state vectors
!  * for bsd4.3 random() */
! long *MRNGorigstate = NULL;
! long *MRNGstate = NULL;
! 
! /* Prototypes */
! long random(void);
! void *initstate(unsigned seed, void *state, int n);
! void *setstate(void* state);
! long nrand48(unsigned short xseed[3]);
  
  /* RNG state sizes */
  #define RNGstatesize_0 1
***************
*** 78,135 ****
  /* A trivial RNG, which can be set up to return a constant or an
   * arithmetic progression.  Can come in handy for debugging. 
   */
! long _triv(rng)
! RNGdata *rng;
  {
    RNGstate[0] += RNGseed[1];
    RNGstate[0] &= MAXLONG; /* clear sign bit */
    return(RNGstate[0]);
  }
  
! void _triv_seed(rng, seed)
! RNGdata *rng;
! long *seed;
  {
  RNGstate[0] = seed[0];
  }
  
! int _triv_check(rng)
! RNGdata *rng;
  {
  return(RNGstate[0] >= 0);
  }
  
  /* bsd4.3 random */
! long _random(rng)
! RNGdata *rng;
  {
! return(random());
  }
  
! void _random_seed(rng, seed)
! RNGdata *rng;
! long *seed;
  {
! /* Call the 4.3bsd random() initialization routine */
! initstate(seed[0], RNGstate, RNGstatesize*4);
  }
  
! int _random_check(rng)
! RNGdata *rng;
  {
! return(RNGstate[0] == MRNGSTATE0);
  }
  
  /* 4.3bsd nrand48 */
! long _nrand48(rng)
! RNGdata *rng;
  {
! return(nrand48(RNGstate));
  }
  
! void _nrand48_seed(rng, seed)
! RNGdata *rng;
! long *seed;
  {
  unsigned short *sarray;
  
--- 91,182 ----
  /* A trivial RNG, which can be set up to return a constant or an
   * arithmetic progression.  Can come in handy for debugging. 
   */
! long _triv(RNGdata *rng)
  {
    RNGstate[0] += RNGseed[1];
    RNGstate[0] &= MAXLONG; /* clear sign bit */
    return(RNGstate[0]);
  }
  
! void _triv_seed(RNGdata *rng, long *seed)
  {
  RNGstate[0] = seed[0];
  }
  
! int _triv_check(RNGdata *rng)
  {
  return(RNGstate[0] >= 0);
  }
  
+ /*****************/
  /* bsd4.3 random */
! /*****************/
! long _random(RNGdata *rng)
  {
!   if (MRNGstate != RNGstate) {
!     long *oldstate = setstate(RNGstate);
!     MRNGstate = RNGstate;
!     if (MRNGorigstate == NULL) MRNGorigstate = oldstate;
!   }
!   return(random());
  }
  
! void random_deactivate_state()
  {
!   if (MRNGorigstate != NULL) {
!     setstate(MRNGorigstate);
!     MRNGstate = MRNGorigstate;
!   }
  }
  
! void random_update_state() 
  {
!   if (MRNGorigstate != NULL) {
!     setstate(setstate(MRNGorigstate));
!   }
  }
  
+ void _random_seed(RNGdata *rng, long *seed)
+ {
+   /* Call the 4.3bsd random() initialization routine.  Assumes new
+    * state is different from old state.  Note that if the old and new
+    * state vectors are the same, the hidden pointers will be lost.
+    * This is because the new state[0] is first saved, then the hidden
+    * pointers are stored in the old state[0], then the new state
+    * becomes current, and finally the new pointers are set from the
+    * saved values.  */
+ 
+   if (MRNGorigstate == NULL) {	/* first time */
+     MRNGorigstate = initstate(seed[0], RNGstate, RNGstatesize*4);
+   }
+   else {
+     random_deactivate_state();
+     initstate(seed[0], RNGstate, RNGstatesize*4);
+   }
+   MRNGstate = RNGstate;
+   random_update_state();
+ }
+ 
+ int _random_check(RNGdata *rng)
+ {
+   int type;
+   int rear;
+ 
+   random_update_state();
+   type = RNGstate[0] % MRNG_MAX_TYPES;
+   rear = RNGstate[0] / MRNG_MAX_TYPES;
+   return type == MRNGTYPE && rear >= 0 && rear < RNGstatesize;
+ }
+ 
+ /******************/
  /* 4.3bsd nrand48 */
! /******************/
! long _nrand48(RNGdata *rng)
  {
! return(nrand48((unsigned short*)RNGstate));
  }
  
! void _nrand48_seed(RNGdata *rng, long *seed)
  {
  unsigned short *sarray;
  
***************
*** 141,148 ****
    sarray[3] = 0; /* stuff a 0 into the other half of RNGstate[1] */
  }
  
! int _nrand48_check(rng)
! RNGdata *rng;
  {
  unsigned short *sarray;
  
--- 188,194 ----
    sarray[3] = 0; /* stuff a 0 into the other half of RNGstate[1] */
  }
  
! int _nrand48_check(RNGdata *rng)
  {
  unsigned short *sarray;
  
***************
*** 151,174 ****
  }
  
  /* 4.3bsd rand */
! long _rand(rng)
! RNGdata *rng;
  {
  /* Very inefficient for many calls!!! */
  srand(RNGstate[0]); /* In case RNGs were switched */
  return(RNGstate[0]=rand());
  }
  
! void _rand_seed(rng, seed)
! RNGdata *rng;
! long *seed;
  {
  /* Use ran() */
  RNGstate[0] = seed[0]; /* The srand() call will be made later */
  }
  
! int _rand_check(rng)
! RNGdata *rng;
  {
  return(RNGstate[0] >= 0);
  }
--- 197,216 ----
  }
  
  /* 4.3bsd rand */
! long _rand(RNGdata *rng)
  {
  /* Very inefficient for many calls!!! */
  srand(RNGstate[0]); /* In case RNGs were switched */
  return(RNGstate[0]=rand());
  }
  
! void _rand_seed(RNGdata *rng, long *seed)
  {
  /* Use ran() */
  RNGstate[0] = seed[0]; /* The srand() call will be made later */
  }
  
! int _rand_check(RNGdata *rng)
  {
  return(RNGstate[0] >= 0);
  }
*** 3.1	1993/05/28 19:25:42
--- ultra.c	1996/06/04 14:53:12
***************
*** 95,101 ****
  
  long swb32fill(rng)
  RNGdata *rng;
! { long temp;
    swb32p = swb32;
    SWBfill(swb32,rng);
    swb32n = N-1;
--- 95,101 ----
  
  long swb32fill(rng)
  RNGdata *rng;
! {
    swb32p = swb32;
    SWBfill(swb32,rng);
    swb32n = N-1;
***************
*** 144,150 ****
  unsigned long *seed;
  {
    short i,j;
!   unsigned long tidbits;
    unsigned long congy, shrgx;
  
    congy=seed[0];
--- 144,150 ----
  unsigned long *seed;
  {
    short i,j;
!   unsigned long tidbits = 0;
    unsigned long congy, shrgx;
  
    congy=seed[0];
*** 3.1	1993/05/28 19:25:42
--- xsq.c	1996/05/31 15:41:28
***************
*** 8,13 ****
--- 8,14 ----
   *   correction.
   */
  
+ #include <stdio.h>
  #include <math.h>
  
  #define min(x,y)	( ((x)<(y)) ? (x) : (y) )
***************
*** 66,71 ****
--- 67,73 ----
   * otherwise uses "%.2g".  Always prints at least 2 significant digits,
   * never prints more than 3.
   */
+ void
  print2digits(f)
  double f;
  {
***************
*** 100,106 ****
  {
    double stddev;/* temp for calculating # of std devations */
    int i; /* loop counter */
-   int robust; /* true if we can do a 2-tailed test with p < 1e-6 */
    int accurate; /* true if the stddev estimate is accurate */ 
    double nrobust,ngiveup; /* min n to give (robust,adequate) test results */
    double smallxsq; /* temp for nrobust and ngiveup calculations */
--- 102,107 ----
***************
*** 195,205 ****
    print2digits(xsq);
    printf(", k = ");
    print2digits(k);
!   printf(", equivalent normal standard deviation ");
!   printf(accurate? "= " : "~ ");
    stddev = xsq_to_stddev(xsq,k,n,&accurate);
    print2digits(stddev);
!   printf("\n", stddev);
  
    if (fabs(stddev) > 6.1) {
      printf("\n");
--- 196,206 ----
    print2digits(xsq);
    printf(", k = ");
    print2digits(k);
!   printf(", equivalent normal standard deviation %s",
! 	 accurate? "= " : "~ ");
    stddev = xsq_to_stddev(xsq,k,n,&accurate);
    print2digits(stddev);
!   printf("\n");
  
    if (fabs(stddev) > 6.1) {
      printf("\n");