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, constg 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 *, 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);
这些函数使用 Marsaglia 开发的带进位相乘算法生成伪随机数序列。在提供了乘数 M 以及种子 X 和进位 C 的初始值(全部为 32 位整数)的情况下,该算法会按如下方式生成一个新种子和进位:
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
种子和进位的序列的周期为 M*2**31 - 1,前提是 (M*2**32 - 1) 和 (M*2**31 - 1) 都是素数。
上面列出的函数内部使用两个 32 位带进位相乘生成器,分别用 mwcran0 和 mwcran1 表示。mwcran0 中使用的乘数为 526533 (0x808C5),mwcran1 中使用的乘数为 557325 (0x8810D)。这两个生成器的周期都大约为 2**50。以下说明涉及这些生成器。
u_mwcran_() 调用 mwcran0 并返回介于 0 和 2**32 - 1 之间的伪随机 32 位无符号整数。
i_mwcran_() 通过调用 u_mwcran_() 并屏蔽结果的最高有效位,返回介于 0 和 2**31 - 1 之间的伪随机 32 位有符号整数。
u_llmwcran_() 同时调用 mwcran0 和 mwcran1 并将两个 32 位值连接在一起以返回介于 0 和 2**64 - 1 之间的伪随机 64 位无符号整数。此随机数生成器的周期大约为 2**100。
i_llmwcran_() 通过调用 u_llmwcran_() 并屏蔽结果的最高有效位,返回介于 0 和 2**63 - 1 之间的伪随机 64 位有符号整数。
函数 i_lmwcran_() 和 u_lmwcran_() 返回伪随机长整数。这些数的范围取决于长整数类型的宽度:如果长整数为 32 位宽(即 ILP32 数据模型),则这些函数分别与 i_mwcran_() 和 u_mwcran_() 相同;如果长整数为 64 位宽(LP64 数据模型),则这些函数分别与 i_llmwcran_() 和 u_llmwcran_() 相同。
r_mwcran_() 返回一个 [0,1) 范围中的伪随机单精度浮点数。此函数重复调用 mwcran0 生成一个随机位序列,然后将此序列解释为二进制小数 0.xxxxxxxxxxxxxxx... 并将其截断为浮点格式,从而生成此数。结果浮点数序列在 [0,1) 范围内均匀分布。
d_mwcran_() 返回一个 [0,1) 范围中的伪随机双精度浮点数。此函数重复调用 mwcran0 和 mwcran1 生成一个随机位序列,然后将此序列解释为二进制小数 0.xxxxxxxxxxxxxxx... 并将其截断为双精度格式,从而生成此数。结果浮点数序列在 [0,1) 范围内均匀分布。
i_mwcrans_(x, n, l, u) 重复调用 mwcran0 以生成 *n 个伪随机 32 位整数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u 之间均匀分布。如果 [*l,*u] = [0,0x7fffffff],则 i_mwcrans_() 获得的 *n 个随机数与通过调用 i_mwcran_() *n 次生成的随机数相同。
u_mwcrans_(x, n, l, u) 重复调用 mwcran0 以生成 *n 个伪随机 32 位无符号整数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u 之间均匀分布。如果 [*l,*u] = [0,0xffffffff],则 u_mwcrans_() 获得的 *n 个随机数与通过调用 u_mwcran_() *n 次生成的随机数相同。
i_llmwcrans_(x, n, l, u) 重复调用 mwcran0 和 mwcran1 以生成 *n 个伪随机 64 位整数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u 之间均匀分布。如果 [*l,*u] = [0,0x7fffffffffffffff],则 i_llmwcrans_() 获得的 *n 个随机数与通过调用 i_llmwcran_() *n 次生成的随机数相同。
u_llmwcrans_(x, n, l, u) 重复调用 mwcran0 和 mwcran1 以生成 *n 个伪随机 64 位无符号整数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u 之间均匀分布。如果 [*l,*u] = [0,0xffffffffffffffff],则 u_llmwcrans_() 获得的 *n 个随机数与通过调用 u_llmwcran_() *n 次生成的随机数相同。
函数 i_lmwcrans_() 和 u_lmwcrans_() 生成伪随机长整数的数组。它们生成的数的范围取决于长整数类型的宽度:如果长整数为 32 位宽(即 ILP32 数据模型),则这些函数分别与 i_mwcrans_() 和 u_mwcrans_() 相同;如果长整数为 64 位宽(LP64 数据模型),则这些函数分别与 i_llmwcrans_() 和 u_llmwcrans_() 相同。
r_mwcrans_(x, n, l, u) 重复调用 mwcran0 以生成 *n 个伪随机单精度浮点数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u(最高到发生截断错误)之间均匀分布。如果 *l 是零并且 *u 是小于 1 的最大单精度浮点数,则 r_mwcrans_() 获得的 *n 个随机数与通过调用 r_mwcran_() *n 次生成的随机数相同。
d_mwcrans_(x, n, l, u) 重复调用 mwcran0 和 mwcran1 以生成 *n 个伪随机双精度浮点数 x[0]、x[1]、...x[*n - 1] 在 *l 与 *u(最高到发生截断错误)之间均匀分布。如果 *l 是零并且 *u 是小于 1 的最大双精度浮点数,则 d_mwcrans_() 获得的 *n 个随机数与通过调用 d_mwcran_() *n 次生成的随机数相同。
i_init_mwcrans_() 将 mwcran0 和 mwcran1 的种子和进位初始化为缺省值。
smwcran_(m) 根据以下公式使用 *m 的杂乱值初始化 mwcran0 和 mwcran1 的种子和进位。
For mwcran0: X0 = MWCRAN_SEED0 + (*m)*0x110005 C0 = MWCRAN_CARRY0 + (*m)*0x110005 For mwcran1: X1 = MWCRAN_SEED1 + (*m)*0x100021 C1 = MWCRAN_CARRY1 + (*m)*0x100021
此处,MWCRAN_SEED0、MWCRAN_CARRY0、MWCRAN_SEED1 和 MWCRAN_CARRY1 是 i_init_mwcrans_() 分别为 mwcran0 和 mwcran1 的种子和进位建立的缺省初始值。需要特别指出的是,调用 smwcran_(m) 时 *m 等于零等效于调用 i_init_mwcrans_()。(此函数提供了一种获取不同随机数序列的简单方式。如需更加精确地控制种子和进位,可使用 i_set_mwcrans_()。)
i_get_mwcrans_(p) 和 i_set_mwcrans_(p) 分别获取和设置 mwcran0 和 mwcran1 的状态表。该表是四个整数 p[0]-p[3] 的数组,其中包含两个生成器的种子和进位:p[0] 和 p[1] 包含 mwcran0 的种子和进位,p[2] 和 p[3] 包含 mwcran1 的种子和进位。
以下示例验证 u_mwcran_() 和 u_mwcrans_() 是否一致,并测试它们生成的数的分布均匀性:
/* 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; }
上一个程序的输出如下:
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(5):
|
addrans(3M)、lcrans(3M)、shufrans(3M)、attributes(5)
因为这些函数共享内部生成器 mwcran0 和 mwcran1,因此它们不是互相独立的。例如,调用 i_init_mwcran_() 之后再调用 d_mwcran_() 产生的结果不同于依次调用 i_init_mwcran_()、i_mwcran_() 和 d_mwcran_()。(但是,生成器独立于其他线程中的那些生成器。)
mwcran0 和 mwcran1 生成的随机数通过 Marsaglia 的 Diehard 测试。