libc随机数生成器有缺陷？
Wol*_*ram 20 c random algorithm math glibc
考虑一种算法来测试在特定次数的尝试之后从一组N个唯一数字中挑选某个数字的概率(例如,N = 2,轮盘中的概率是什么(没有0),它需要X尝试黑赢？)
对此的正确分布是pow(1-1/N,X-1)*(1/N).
但是,当我使用以下代码对其进行测试时,在X = 31处始终存在深沟,独立于N,并且独立于种子.
这是一个内在的缺陷,由于PRNG的实施细节在使用中无法防止,这是一个真正的错误,还是我忽略了一些明显的东西？
// C
#include <sys/times.h>
#include <math.h>
#include <stdio.h>
int array[101];
void main(){
int nsamples=10000000;
double breakVal,diffVal;
int i,cnt;
// seed, but doesn't change anything
struct tms time;
srandom(times(&time));
// sample
for(i=0;i<nsamples;i++){
cnt=1;
do{
if((random()%36)==0) // break if 0 is chosen
break;
cnt++;
}while(cnt<100);
array[cnt]++;
}
// show distribution
for(i=1;i<100;i++){
breakVal=array[i]/(double)nsamples; // normalize
diffVal=breakVal-pow(1-1/36.,i-1)*1/36.; // difference to expected value
printf("%d %.12g %.12g\n",i,breakVal,diffVal);
}
}
使用libc6软件包2.15-0ubuntu20和Intel Core i5-2500 SandyBridge测试了最新的Xubuntu 12.10,但几年前我在一台较旧的Ubuntu机器上发现了这一点.
我也在Windows 7上使用Unity3D/Mono进行了测试(虽然不确定哪个Mono版本),这里使用System.Random时,X = 55时沟渠发生,而Unity内置的Unity.Random没有可见的沟渠(至少没有)对于X <100).
分布: 
差异: 
int*_*jay 11
这是因为glibc的random()功能不够随意.根据这个页面,对于返回的随机数random(),我们有:
oi = (oi-3 + oi-31) % 2^31
要么:
oi = (oi-3 + oi-31 + 1) % 2^31.
现在考虑,并假设上面的第一个等式是使用的等式(每个数字的概率为50%).现在如果和,那么机会小于1/36.这是因为50%的时间将少于2 ^ 31,当发生这种情况时,xi = oi % 36xi-31=0xi-3!=0xi=0oi-31 + oi-3
xi = oi % 36 = (oi-3 + oi-31) % 36 = oi-3 % 36 = xi-3,
这是非零的.这会导致您在0样本后看到31个样本的沟渠.
在这个实验中测量的是伯努利实验的成功试验之间的间隔,其中成功被定义random() mod k == 0为一些k(OP中的36个).不幸的是,实施random()意味着伯努利试验不具有统计独立性,这一事实更加严重.
我们将写出`random()' 的输出,我们注意到:rndiith
rndi = rndi-31 + rndi-3 概率为0.75
rndi = rndi-31 + rndi-3 + 1 概率为0.25
(见下面的证明大纲.)
让我们假设,我们目前正在关注.那么一定是这样的情况,因为否则我们会将周期计为长度.rndi-31 mod k == 0rndirndi-3 mod k ≠ 0k-3
但是(大部分时间).(mod k): rndi = rndi-31 + rndi-3 = rndi-3 ≠ 0
所以,目前的审判是没有统计学独立于以前的试验中,和31 日成功试用后是不太可能会比在公正的系列伯努利试验的成功.
使用线性同余生成器(通常不适用于random()算法)的通常建议是使用高阶位而不是低阶位,因为高阶位是"更随机"(即,与连续值的相关性较小).但是在这种情况下也不会起作用,因为上述身份同样适用于函数high log k bits的功能mod k == low log k bits.
事实上,我们可能期望线性同余发生器更好地工作,特别是如果我们使用输出的高阶位,因为尽管LCG在蒙特卡罗模拟中不是特别好,但它不会受到线性反馈的影响.random().
random 算法,默认情况下:
让我们state成为无符号长子的向量.使用种子,一些固定值和混合算法进行初始化.为简单起见,我们可以认为状态向量是无限的,尽管只使用了最后的31个值,所以它实际上是作为环形缓冲区实现的.state0...state30
生成 rndi: (Note: ⊕ is addition mod 232.)
statei = statei-31 ⊕ statei-3
rndi = (statei - (statei mod 2)) / 2
Now, note that:
(i + j) mod 2 = i mod 2 + j mod 2 if i mod 2 == 0 or j mod 2 == 0
(i + j) mod 2 = i mod 2 + j mod 2 - 2 if i mod 2 == 1 and j mod 2 == 1
If i and j are uniformly distributed, the first case will occur 75% of the time, and the second case 25%.
So, by substitution in the generation formula:
rndi = (statei-31 ⊕ statei-3 - ((statei-31 + statei-3) mod 2)) / 2
= ((statei-31 - (statei-31 mod 2)) ⊕ (statei-3 - (statei-3 mod 2))) / 2 or
= ((statei-31 - (statei-31 mod 2)) ⊕ (statei-3 - (statei-3 mod 2)) + 2) / 2
The two cases can be further reduced to:
rndi = rndi-31 ⊕ rndi-3
random() mod k == 0
如上所述,第一种情况发生在75%的时间,假设rnd i-31和rnd i-3独立地从均匀分布中得出(它们不是,但它是合理的第一近似).