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**32 - 1) と (M*2**31 - 1) が両方とも素数である場合、シードとキャリーのシーケンスの周期は M*2**31 - 1 です。
上記の関数では内部的に、mwcran0 と mwcran1 という 2 つの 32 ビットのキャリー付き乗算ジェネレータが使用されます。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 の両方を呼び出し、2 つの 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 を繰り返し呼び出して、*l から *u までで均一に分布された *n 個の 32 ビット整数の疑似乱数 x[0], x[1], ... x[*n - 1] を生成します。[*l,*u] = [0,0x7fffffff] の場合、i_mwcrans_() は、i_mwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数が生成されます。
u_mwcrans_(x, n, l, u) は、mwcran0 を繰り返し呼び出して、*l から *u までで均一に分布された *n 個の 32 ビット符号なし整数の疑似乱数 x[0], x[1], ... x[*n - 1] を生成します。[*l,*u] = [0,0xffffffff] の場合、u_mwcrans_() は、u_mwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数が生成されます。
i_llmwcrans_(x, n, l, u) は、mwcran0 と mwcran1 を繰り返し呼び出して、*l から *u までで均一に分布された *n 個の 64 ビット整数の疑似乱数 x[0], x[1], ... x[*n - 1] を生成します。[*l,*u] = [0,0x7fffffffffffffff] の場合、i_llmwcrans_() は、i_llmwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数が生成されます。
u_llmwcrans_(x, n, l, u) は、mwcran0 と mwcran1 を繰り返し呼び出して、*l から *u までで均一に分布された *n 個の 64 ビット符号なし整数の乱数 x[0], x[1], ... x[*n - 1] を生成します。[*l,*u] = [0,0xffffffffffffffff] の場合、u_llmwcrans_() は、u_llmwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数が生成されます。
関数 i_lmwcrans_() と u_lmwcrans_() は、ロング整数型の疑似乱数の配列を生成します。生成される数値の範囲は、ロング整数型の幅によって異なります。ロング整数が 32 ビット幅 (つまり、ILP32 データモデル) の場合、これらの関数はそれぞれ i_mwcrans_() と u_mwcrans_() と同じです。ロング整数が 64 ビット幅 (LP64 データモデル) の場合、これらの関数はそれぞれ i_llmwcrans_() と u_llmwcrans_() と同じです。
r_mwcrans_(x, n, l, u) は、mwcran0 を繰り返し呼び出して、*l から *u までで均一に分布された *n 個の単精度浮動小数点の疑似乱数 x[0], x[1], ... x[*n - 1] を生成します。*l がゼロで、*u が 1 よりも小さい最大単精度浮動小数点数である場合、r_mwcrans_() は、r_mwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数を生成します。
d_mwcrans_(x, n, l, u) は、mwcran0 と mwcran1 を繰り返し呼び出して、*l から *u まで (切り捨てエラーが発生するまで) で均一に分布された *n 個の倍精度浮動小数点の疑似乱数 x[0], x[1], ... x[*n - 1] を生成します。*l がゼロで、*u が 1 よりも小さい最大倍精度浮動小数点数である場合、d_mwcrans_() は、d_mwcran_() を *n 回呼び出したときに生成される場合と同じ *n 個の乱数を生成します。
i_init_mwcrans_() は、mwcran0 と mwcran1 のシードとキャリーをデフォルト値に初期化します。
smwcran_(m) は、次の公式に従って、*m のスクランブル値で mwcran0 and 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 は、それぞれ mwcran0 と mwcran1 のシードとキャリーに対して i_init_mwcrans_() を呼び出す際に確立されるデフォルトの初期値です。特に、ゼロと等しい *m で smwcran_(m) を呼び出すことは、i_init_mwcrans_() を呼び出すことと同等です。(この関数は、異なる乱数のシーケンスを取得するための簡単な方法です。より正確にシードとキャリーを制御するには、i_set_mwcrans_() を使用してください。)
i_get_mwcrans_(p) と i_set_mwcrans_(p) は、それぞれ mwcran0 と mwcran1 の状態テーブルを取得および設定します。このテーブルは、両方のジェネレータのシードとキャリーが含まれている 4 つの整数 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 試験に合格しています。