'\" t .\" @(#)mwcrans.3m 1.2 98/02/04 SMI; .TH mwcrans 3M "28 Jan 1998" .SH NAME mwcrans \- multiply with carry pseudo-random number generators .SH SYNOPSIS .LP .B cc .RI "[ " flag " \|.\|.\|. ] " file " \|.\|.\|." .B \-lsunmath -lm .RI "[ " library " \|.\|.\|. ]" .LP .B #include .LP .B int i_mwcran_(void); .LP .B unsigned int u_mwcran_(void); .LP .B long i_lmwcran_(void); .LP .B unsigned long u_lmwcran_(void); .LP .B long long i_llmwcran_(void); .LP .B unsigned long long u_llmwcran_(void); .LP .B float r_mwcran_(void); .LP .B double d_mwcran_(void); .LP .BI "void i_mwcrans_(int *" x ", const int *" n , .BI "const int *" l ", const int *" u ");" .LP .BI "void u_mwcrans_(unsigned *" x ", const int *" n , .BI "const unsigned *" l ", const unsigned *" u ");" .LP .BI "void i_lmwcrans_(long *" x ", const int *" n , .BI "const long *" l ", const long *" u ");" .LP .BI "void u_lmwcrans_(unsigned long *" x ", const int *" n , .BI "const unsigned long *" l ", const unsigned long *" u ");" .LP .BI "void i_llmwcrans_(long long *" x ", const int *" n , .BI "const long long *" l ", const long long *" u ");" .LP .BI "void u_llmwcrans_(unsigned long long *" x ", const int *" n , .BI "const unsigned long long *" l ", const unsigned long long *" u ");" .LP .BI "void r_mwcrans_(float *" x ", const int *" n , .BI "const float *" l ", const float *" u ");" .LP .BI "void d_mwcrans_(double *" x ", const int *" n , .BI "const double *" l ", const double *" u ");" .LP .B void i_init_mwcrans_(void); .LP .BI "void smwcran_(const int *" seed ");" .LP .BI "void i_set_mwcrans_(const int *" p ");" .LP .BI "void i_get_mwcrans_(int *" p ");" .SH DESCRIPTION These functions generate sequences of pseudo-random numbers using the multiply-with-carry algorithm developed by Marsaglia. Given a multiplier \f2M\f1 and the initial value of a seed \f2X\f1 and carry \f2C\f1 (all 32 bit integers), the algorithm generates a new seed and carry as follows: .LP .nf 1. \f2Z\f1 = \f2X\f1*\f2M\f1 + \f2C\f1 (here \f2Z\f1 is a 64 bit integer) 2. new seed \f2X\f1 = lower 32 bits of \f2Z\f1 3. new carry \f2C\f1 = upper 32 bits of \f2Z\f1 .fi .LP The period of the sequence of seeds and carries is \f2M\f1*2**31 - 1, provided that both (\f2M\f1*2**32 - 1) and (\f2M\f1*2**31 - 1) are prime. .LP The functions listed above internally use two 32 bit multiply-with-carry generators denoted by \f3mwcran0\f1 and \f3mwcran1\f1. The multiplier used in \f3mwcran0\f1 is 526533 (0x808C5), and the multiplier used in \f3mwcran1\f1 is 557325 (0x8810D). The periods of both of these generators are approximately 2**50. The following descriptions refer to these generators. .LP \f3u_mwcran_(\|)\f1 invokes \f3mwcran0\f1 and returns a pseudo-random 32 bit unsigned integer between 0 and 2**32 - 1. .LP \f3i_mwcran_(\|)\f1 returns a pseudo-random 32 bit signed integer between 0 and 2**31 - 1 by calling \f3u_mwcran_(\|)\f1 and masking off the most significant bit of the result. .LP \f3u_llmwcran_(\|)\f1 invokes both \f3mwcran0\f1 and \f3mwcran1\f1 and concatenates 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. .LP \f3i_llmwcran_(\|)\f1 returns a pseudo-random 64 bit signed integer between 0 and 2**63 - 1 by calling \f3u_llmwcran_(\|)\f1 and masking off the most significant bit of the result. .LP The functions \f3i_lmwcran_(\|)\f1 and \f3u_lmwcran_(\|)\f1 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 \f3i_mwcran_(\|)\f1 and \f3u_mwcran_(\|)\f1, respectively; if a long int is 64 bits wide (the LP64 data model), these functions are the same as \f3i_llmwcran_(\|)\f1 and \f3u_llmwcran_(\|)\f1, respectively. .LP \f3r_mwcran_(\|)\f1 returns a pseudo-random single precision floating point number in the range [0,1). It produces this number by calling \f3mwcran0\f1 repeatedly to generate a sequence of random bits, which is then interpreted as a binary fraction 0.xxxxxxxxxxxxxxx... and truncated to the float format. The resulting sequence of floating point numbers is uniformly distributed in the range [0,1). .LP \f3d_mwcran_(\|)\f1 returns a pseudo-random double precision floating point number in the range [0,1). It produces this number by calling \f3mwcran0\f1 and \f3mwcran1\f1 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). .LP \f3i_mwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 repeatedly to generate *\f2n\f1 pseudo-random 32 bit integers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1. If [*\f2l\f1,*\f2u\f1] = [0,0x7fffffff], then \f3i_mwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3i_mwcran_(\|)\f1 *\f2n\f1 times. .LP \f3u_mwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 repeatedly to generate *\f2n\f1 pseudo-random 32 bit unsigned integers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1. If [*\f2l\f1,*\f2u\f1] = [0,0xffffffff], then \f3u_mwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3u_mwcran_(\|)\f1 *\f2n\f1 times. .LP \f3i_llmwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 and \f3mwcran1\f1 repeatedly to generate *\f2n\f1 pseudo-random 64 bit integers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1. If [*\f2l\f1,*\f2u\f1] = [0,0x7fffffffffffffff], then \f3i_llmwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3i_llmwcran_(\|)\f1 *\f2n\f1 times. .LP \f3u_llmwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 and \f3mwcran1\f1 repeatedly to generate *\f2n\f1 pseudo-random 64 bit unsigned integers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1. If [*\f2l\f1,*\f2u\f1] = [0,0xffffffffffffffff], then \f3u_llmwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3u_llmwcran_(\|)\f1 *\f2n\f1 times. .LP The functions \f3i_lmwcrans_(\|)\f1 and \f3u_lmwcrans_(\|)\f1 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 \f3i_mwcrans_(\|)\f1 and \f3u_mwcrans_(\|)\f1, respectively; if a long int is 64 bits wide (the LP64 data model), these functions are the same as \f3i_llmwcrans_(\|)\f1 and \f3u_llmwcrans_(\|)\f1, respectively. .LP \f3r_mwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 repeatedly to generate *\f2n\f1 pseudo-random single precision floating point numbers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1 (up to a truncation error). If *\f2l\f1 is zero and *\f2u\f1 is the largest single precision floating point number less than 1, then \f3r_mwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3r_mwcran_(\|)\f1 *\f2n\f1 times. .LP \f3d_mwcrans_(\f2x\f3, \f2n\f3, \f2l\f3, \f2u\f3)\f1 invokes \f3mwcran0\f1 and \f3mwcran1\f1 repeatedly to generate *\f2n\f1 pseudo-random double precision floating point numbers \f2x\f1[0], \f2x\f1[1], \|.\|.\|. \f2x\f1[*\f2n\f1 - 1] uniformly distributed between *\f2l\f1 and *\f2u\f1 (up to a truncation error). If *\f2l\f1 is zero and *\f2u\f1 is the largest double precision floating point number less than 1, then \f3d_mwcrans_(\|)\f1 yields the same *\f2n\f1 random numbers as would be generated by calling \f3d_mwcran_(\|)\f1 *\f2n\f1 times. .LP \f3i_init_mwcrans_(\|)\f1 initializes the seeds and carries for \f3mwcran0\f1 and \f3mwcran1\f1 to default values. .LP \f3smwcran_(\f2m\f3)\f1 initializes the seeds and carries for \f3mwcran0\f1 and \f3mwcran1\f1 with a scrambled value of *\f2m\f1 according to the following formulas. .LP .RS .nf For \f3mwcran0\f1: \f2X0\f1 = MWCRAN_SEED0 + (*\f2m\f1)*0x110005 \f2C0\f1 = MWCRAN_CARRY0 + (*\f2m\f1)*0x110005 .LP For \f3mwcran1\f1: \f2X1\f1 = MWCRAN_SEED1 + (*\f2m\f1)*0x100021 \f2C1\f1 = MWCRAN_CARRY1 + (*\f2m\f1)*0x100021 .fi .RE .LP Here MWCRAN_SEED0, MWCRAN_CARRY0, MWCRAN_SEED1, and MWCRAN_CARRY1 are the default initial values established by \f3i_init_mwcrans_(\|)\f1 for the seeds and carries of \f3mwcran0\f1 and \f3mwcran1\f1, respectively. In particular, calling \f3smwcran_(\f2m\f3)\f1 with *\f2m\f1 equal to zero is equivalent to calling \f3i_init_mwcrans_(\|)\f1. (This function provides a simple way to obtain different sequences of random numbers. For more precise control of the seeds and carries, use \f3i_set_mwcrans_(\|)\f1.) .LP \f3i_get_mwcrans_(\f2p\f3)\f1 and \f3i_set_mwcrans_(\f2p\f3)\f1 respectively get and set the state table of \f3mwcran0\f1 and \f3mwcran1\f1. The table is an array of four integers \f2p\f1[0]-\f2p\f1[3] containing the seeds and carries for both generators: \f2p\f1[0] and \f2p\f1[1] contain the seed and carry for \f3mwcran0\f1, and \f2p\f1[2] and \f2p\f1[3] contain the seed and carry for \f3mwcran1\f1. .SH EXAMPLE The following example verifies that \f3u_mwcran_(\|)\f1 and \f3u_mwcrans_(\|)\f1 are consistent and tests the uniformity of the distribution of the numbers they generate: .LP .RS .nf .vs 11p .ft B /* sample program to check u_mwcran*() */ #include 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>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>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>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; } .vs .fi .RE .LP The output from the preceding program reads: .LP .RS .nf .vs 11p .ft B 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. .vs .fi .RE .SH ATTRIBUTES See .BR attributes (5) for descriptions of the following attributes: .sp .ne 10 .TS box; cbp-1 | cbp-1 l | l . ATTRIBUTE TYPE ATTRIBUTE VALUE = Availability SPROsunms, SPROlang (32 bit libraries) SP64sunms, SP64lang (64 bit libraries) Interface Stability Evolving MT-Level MT-Safe .TE .SH "SEE ALSO" .BR addrans (3M), .BR lcrans (3M), .BR shufrans (3M), .BR attributes (5) .SH NOTES Because these functions share the internal generators \f3mwcran0\f1 and \f3mwcran1\f1, they are not independent of one another. For example, calling \f3i_init_mwcran_(\|)\f1 followed by \f3d_mwcran_(\|)\f1 will produce a different result than calling \f3i_init_mwcran_(\|)\f1 followed by \f3i_mwcran_(\|)\f1 and then \f3d_mwcran_(\|)\f1. (The generators \f2are\f1 independent of those in other threads, however.) .LP The random numbers generated by \f3mwcran0\f1 and \f3mwcran1\f1 pass Marsaglia's Diehard test.