Man Page mwcrans.3m




NAME

     mwcrans - multiply with carry pseudo-random  number  genera-
     tors


SYNOPSIS

     cc [ flag ... ] file ...  -lsunmath -lm [ library ... ]

     #include <sunmath.h>

     int i_mwcran_(void);

     unsigned int u_mwcran_(void);

     long i_lmwcran_(void);

     unsigned long u_lmwcran_(void);

     long long i_llmwcran_(void);

     unsigned long long u_llmwcran_(void);

     float r_mwcran_(void);

     double d_mwcran_(void);

     void i_mwcrans_(int *x, const int *n, const  int  *l,  const
     int *u);

     void u_mwcrans_(unsigned *x, const int  *n,  const  unsigned
     *l, const unsigned *u);

     void i_lmwcrans_(long *x, const int *n, const long *l, const
     long *u);

     void u_lmwcrans_(unsigned  long  *x,  const  int  *n,  const
     unsigned long *l, const unsigned long *u);

     void i_llmwcrans_(long long *x, const  int  *n,  const  long
     long *l, const long long *u);

     void u_llmwcrans_(unsigned long long *x, const int *n, const
     unsigned long long *l, const unsigned long long *u);

     void r_mwcrans_(float *x, const  int  *n,  const  float  *l,
     const float *u);

     void d_mwcrans_(double *x, const int *n,  const  double  *l,
     const double *u);

     void i_init_mwcrans_(void);


     void smwcran_(const int *seed);

     void i_set_mwcrans_(const int *p);

     void i_get_mwcrans_(int *p);


DESCRIPTION

     These functions generate sequences of pseudo-random  numbers
     using  the  multiply-with-carry  algorithm developed by Mar-
     saglia.  Given a multiplier M and the  initial  value  of  a
     seed X and carry C (all 32 bit integers), the algorithm gen-
     erates a new seed and carry as follows:

          1. Z = X*M + C  (here Z is a 64 bit integer)
          2. new seed  X = lower 32 bits of Z
          3. new carry C = upper 32 bits of Z

     The period of the sequence of seeds and carries is M*2**31 -
     1,  provided  that  both (M*2**32 - 1) and (M*2**31 - 1) are
     prime.

     The  functions  listed  above  internally  use  two  32  bit
     multiply-with-carry   generators   denoted  by  mwcran0  and
     mwcran1.   The  multiplier  used  in   mwcran0   is   526533
     (0x808C5),  and  the  multiplier  used  in mwcran1 is 557325
     (0x8810D).  The periods of  both  of  these  generators  are
     approximately  2**50.   The  following descriptions refer to
     these generators.

     u_mwcran_() invokes mwcran0 and returns a  pseudo-random  32
     bit unsigned integer between 0 and 2**32 - 1.

     i_mwcran_() returns a pseudo-random 32  bit  signed  integer
     between  0  and 2**31 - 1 by calling u_mwcran_() and masking
     off the most significant bit of the result.

     u_llmwcran_() invokes both mwcran0 and mwcran1 and concaten-
     ates  the  two  32  bit  values together to return a pseudo-
     random 64 bit unsigned integer between 0 and 2**64 - 1.  The
     period  of  this  random  number  generator is approximately
     2**100.

     i_llmwcran_() returns a pseudo-random 64 bit signed  integer
     between 0 and 2**63 - 1 by calling u_llmwcran_() and masking
     off the most significant bit of the result.

     The functions i_lmwcran_() and u_lmwcran_()  return  pseudo-
     random long integers.  The ranges of these numbers depend on
     the width of the long integer type: if a long int is 32 bits
     wide  (i.e.,  the ILP32 data model), these functions are the
     same as i_mwcran_() and u_mwcran_(), respectively; if a long
     int  is  64 bits wide (the LP64 data model), these functions
     are the same as  i_llmwcran_()  and  u_llmwcran_(),  respec-
     tively.

     r_mwcran_() returns a pseudo-random single precision  float-
     ing  point  number  in  the  range  [0,1).  It produces this
     number by calling mwcran0 repeatedly to generate a  sequence
     of  random bits, which is then interpreted as a binary frac-
     tion 0.xxxxxxxxxxxxxxx... and truncated to the float format.
     The  resulting  sequence  of  floating point numbers is uni-
     formly distributed in the range [0,1).

     d_mwcran_() returns a pseudo-random double precision  float-
     ing  point  number  in  the  range  [0,1).  It produces this
     number by calling mwcran0 and mwcran1 repeatedly to generate
     a  sequence  of  random bits, which is then interpreted as a
     binary fraction 0.xxxxxxxxxxxxxxx... and  truncated  to  the
     double  format.   The  resulting  sequence of floating point
     numbers is uniformly distributed in the range [0,1).

     i_mwcrans_(x, n, l, u) invokes mwcran0  repeatedly  to  gen-
     erate  *n pseudo-random 32 bit integers x[0], x[1], ... x[*n
     - 1] uniformly distributed between *l and *u. If  [*l,*u]  =
     [0,0x7fffffff],  then i_mwcrans_() yields the same *n random
     numbers as would be  generated  by  calling  i_mwcran_()  *n
     times.

     u_mwcrans_(x, n, l, u) invokes mwcran0  repeatedly  to  gen-
     erate  *n pseudo-random 32 bit unsigned integers x[0], x[1],
     ... x[*n - 1] uniformly distributed between *l  and  *u.  If
     [*l,*u]  = [0,0xffffffff], then u_mwcrans_() yields the same
     *n  random  numbers  as  would  be  generated   by   calling
     u_mwcran_() *n times.

     i_llmwcrans_(x, n, l, u) invokes mwcran0 and mwcran1 repeat-
     edly  to  generate  *n  pseudo-random  64 bit integers x[0],
     x[1], ... x[*n - 1] uniformly distributed between *l and *u.
     If  [*l,*u]  =  [0,0x7fffffffffffffff],  then i_llmwcrans_()
     yields the same *n random numbers as would be  generated  by
     calling i_llmwcran_() *n times.

     u_llmwcrans_(x, n, l, u) invokes mwcran0 and mwcran1 repeat-
     edly  to  generate *n pseudo-random 64 bit unsigned integers
     x[0], x[1], ... x[*n - 1] uniformly distributed  between  *l
     and   *u.   If   [*l,*u]   =   [0,0xffffffffffffffff],  then
     u_llmwcrans_() yields the same *n random numbers as would be
     generated by calling u_llmwcran_() *n times.

     The  functions  i_lmwcrans_()  and  u_lmwcrans_()   generate
     arrays  of  pseudo-random  long integers.  The ranges of the
     numbers they generate  depend  on  the  width  of  the  long
     integer type: if a long int is 32 bits wide (i.e., the ILP32
     data model), these functions are the  same  as  i_mwcrans_()
     and  u_mwcrans_(),  respectively;  if  a long int is 64 bits
     wide (the LP64 data model), these functions are the same  as
     i_llmwcrans_() and u_llmwcrans_(), respectively.

     r_mwcrans_(x, n, l, u) invokes mwcran0  repeatedly  to  gen-
     erate  *n  pseudo-random  single  precision  floating  point
     numbers x[0], x[1], ...  x[*n  -  1]  uniformly  distributed
     between  *l and *u (up to a truncation error). If *l is zero
     and *u is the largest single precision floating point number
     less  than  1,  then  r_mwcrans_() yields the same *n random
     numbers as would be  generated  by  calling  r_mwcran_()  *n
     times.

     d_mwcrans_(x, n, l, u) invokes mwcran0 and  mwcran1  repeat-
     edly  to generate *n pseudo-random double precision floating
     point numbers x[0], x[1], ...  x[*n - 1]  uniformly  distri-
     buted  between  *l and *u (up to a truncation error).  If *l
     is zero and *u is  the  largest  double  precision  floating
     point  number less than 1, then d_mwcrans_() yields the same
     *n  random  numbers  as  would  be  generated   by   calling
     d_mwcran_() *n times.

     i_init_mwcrans_() initializes  the  seeds  and  carries  for
     mwcran0 and mwcran1 to default values.

     smwcran_(m) initializes the seeds and  carries  for  mwcran0
     and  mwcran1  with  a scrambled value of *m according to the
     following formulas.

          For mwcran0:
               X0 = MWCRAN_SEED0 + (*m)*0x110005
               C0 = MWCRAN_CARRY0 + (*m)*0x110005

          For mwcran1:
               X1 = MWCRAN_SEED1 + (*m)*0x100021
               C1 = MWCRAN_CARRY1 + (*m)*0x100021

     Here   MWCRAN_SEED0,   MWCRAN_CARRY0,   MWCRAN_SEED1,    and
     MWCRAN_CARRY1  are the default initial values established by
     i_init_mwcrans_() for the seeds and carries of  mwcran0  and
     mwcran1,  respectively.   In particular, calling smwcran_(m)
     with  *m  equal   to   zero   is   equivalent   to   calling
     i_init_mwcrans_().   (This function provides a simple way to
     obtain different sequences of random numbers.  For more pre-
     cise    control    of    the    seeds   and   carries,   use
     i_set_mwcrans_().)

     i_get_mwcrans_(p) and i_set_mwcrans_(p) respectively get and
     set the state table of mwcran0 and mwcran1.  The table is an
     array of four integers p[0]-p[3] containing  the  seeds  and
     carries  for both generators: p[0] and p[1] contain the seed
     and carry for mwcran0, and p[2] and p[3]  contain  the  seed
     and carry for mwcran1.


EXAMPLE

     The  following  example  verifies   that   u_mwcran_()   and
     u_mwcrans_()  are consistent and tests the uniformity of the
     distribution of the numbers they generate:

          /* sample program to check u_mwcran*() */

          #include <sunmath.h>

          int main() {
               unsigned x[1000],y[1000],i,lb,ub;
               int n,hex1[16],hex2[16],hex3[16],seed,j;
               double t1,t2,t3,v;

               seed = 40;
               /* generate 1000 random unsigned ints with seed initialized to 40 */
               smwcran_(&seed);
               for (i=0;i<1000;i++) x[i] = u_mwcran_();

               /* generate the same 1000 random unsigned ints by u_mwcrans_() */
               n = 1000; lb = 0; ub = 0xffffffff;
               smwcran_(&seed);    /* reset seed */
               u_mwcrans_(y,&n,&lb,&ub);

               /* verify that x and y are indeed equal */
               for (n=0,i=0;i<1000;i++) n |= (x[i]-y[i]);
               if(n==0) printf("Array x is equal to array y.\n\n");
               else     printf("***Error: array x is not equal to array y.\n\n");

               /* Check uniformity by counting the appearance of hexadecimal digits
                  of 1000 random numbers. The expected value is 500. The counting is
                  repeated three times with different seed values. */

               /* initialize the counting array hex1[],hex2[],hex3[] */
               n = 1000; for (i=0;i<16;i++) hex1[i]=hex2[i]=hex3[i]=0;

               /* count the hex digits of 1000 random numbers with seed=1 */
               seed = 1; smwcran_(&seed);
               u_mwcrans_(x,&n,&lb,&ub);
               for (i=0;i<n;i++) {
                   j = x[i]; hex1[j&0xf] += 1; hex1[(j>>4)&0xf] += 1;
                   hex1[(j>>8)&0xf] += 1; hex1[(j>>12)&0xf] += 1;
                   hex1[(j>>16)&0xf] += 1; hex1[(j>>20)&0xf] += 1;
                   hex1[(j>>24)&0xf] += 1; hex1[(j>>28)&0xf] += 1;
               }

               /* count the hex digits of 1000 random number with seed=2 */
               seed = 2; smwcran_(&seed);
               u_mwcrans_(x,&n,&lb,&ub);
               for (i=0;i<n;i++) {
                   j = x[i]; hex2[j&0xf] += 1; hex2[(j>>4)&0xf] += 1;
                   hex2[(j>>8)&0xf] += 1; hex2[(j>>12)&0xf] += 1;
                   hex2[(j>>16)&0xf] += 1; hex2[(j>>20)&0xf] += 1;
                   hex2[(j>>24)&0xf] += 1; hex2[(j>>28)&0xf] += 1;
               }

               /* count the hex digits of 1000 random number with seed=3 */
               seed = 3; smwcran_(&seed);
               u_mwcrans_(x,&n,&lb,&ub);
               for (i=0;i<n;i++) {
                   j = x[i]; hex3[j&0xf] += 1; hex3[(j>>4)&0xf] += 1;
                   hex3[(j>>8)&0xf] += 1; hex3[(j>>12)&0xf] += 1;
                   hex3[(j>>16)&0xf] += 1; hex3[(j>>20)&0xf] += 1;
                   hex3[(j>>24)&0xf] += 1; hex3[(j>>28)&0xf] += 1;
               }

               /* Compute the Chi-square of each test */
               t1 = t2 = t3 = 0.0;
               for (i=0;i<16;i++) {
                   v  = hex1[i]-500; t1 += v*v/500.0;
                   v  = hex2[i]-500; t2 += v*v/500.0;
                   v  = hex3[i]-500; t3 += v*v/500.0;
               }

               /* print out result */
               printf("Expected value of each hex digit's appearance is 500.\n");
               printf("Observed result with seed=1,2, and 3:\n");
               printf("      Hex digit    Number of appearances with\n");
               printf("                   seed=1     seed=2    seed=3\n");
               for (i=0;i<16;i++) {
                   printf("           %01X:       %4d       %4d      %4d\n",
                   i,hex1[i],hex2[i],hex3[i]);
               }
               printf("                   ----------------------------\n");
               printf("Chi-square value%7.2g   %8.2g   %8.2g\n",t1,t2,t3);
               printf("Note:  A reasonable range of the Chi-square value is\n");
               printf("       within [7.26, 25.00] which corresponds to the 5\n");
               printf("       percent and 95 percent points of the Chi-square\n");
               printf("       distribution with degree of freedom equal to 15.\n");

               return 0;
          }

     The output from the preceding program reads:

          Array x is equal to array y.

          Expected value of each hex digit's appearance is 500.
          Observed result with seed=1,2, and 3:
                Hex digit    Number of appearances with
                             seed=1     seed=2    seed=3
                     0:        514        493       521
                     1:        529        507       480
                     2:        477        495       493
                     3:        495        541       517
                     4:        518        504       486
                     5:        496        464       484
                     6:        467        488       484
                     7:        511        487       540
                     8:        517        499       525
                     9:        500        489       490
                     A:        492        506       511
                     B:        475        504       482
                     C:        499        504       504
                     D:        485        514       493
                     E:        520        531       515
                     F:        505        474       475
                             ----------------------------
          Chi-square value    9.5         11         11
          Note:  A reasonable range of the Chi-square value is
                 within [7.26, 25.00] which corresponds to the 5
                 percent and 95 percent points of the Chi-square
                 distribution with degree of freedom equal to 15.


ATTRIBUTES

     See attributes(5) for descriptions of the  following  attri-
     butes:

     ______________________________________________________________
    |   ATTRIBUTE TYPE   |             ATTRIBUTE VALUE            |
    |____________________|_________________________________________|_
    | Availability       |  SPROsunms, SPROlang (32 bit libraries)|
    |                    |  SP64sunms, SP64lang (64 bit libraries)|
    | Interface Stability|  Evolving                              |
    | MT-Level           |  MT-Safe                               |
    |____________________|________________________________________|


SEE ALSO

     addrans(3M), lcrans(3M), shufrans(3M), attributes(5)


NOTES

     Because  these  functions  share  the  internal   generators
     mwcran0  and  mwcran1,  they  are  not  independent  of  one
     another.  For example, calling i_init_mwcran_() followed  by
     d_mwcran_()  will  produce  a  different result than calling
     i_init_mwcran_()   followed   by   i_mwcran_()   and    then
     d_mwcran_().   (The  generators  are independent of those in
     other threads, however.)

     The random numbers generated by  mwcran0  and  mwcran1  pass
     Marsaglia's Diehard test.