博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
C++随机数--——生成任意范围内等概率随机数“足够好”的做法
阅读量:5825 次
发布时间:2019-06-18

本文共 2999 字,大约阅读时间需要 9 分钟。

很容易想到一个简单的做法

x = rand () % RANGE; /* Poor! return a random number in [0,RANGE) */

 

这种做法有三个问题:

1 rand()等概率返回[0, RAND_MAX]间这RAND_MAX+1个数中的任意一个,如果RAND_MAX+1不能被RANGE整除,那么这个[0, RANGE)分布就blatantly不均匀了。举个例子RAND_MAX=32767,目标区间为[0,10),32767%10=7,则(7,10)间数的取值概率小于[0,7]间数。

2 rand()的伪随机性。大多数随机数发生器产生随机数的低位bits都distressingly不随机(有一个short cyclical pattern或者bits间有依赖关系)。

3 RAND_MAX在<stdlib.h>中定义,其值最小为32767(2^15-1),最大为2147483647(2^31-1),假如RAND_MAX+1 < RANGE呢?

 

先讨论RANGE<=RAND_MAX+1的情况

方法1

这种方法相当于取高位bits

(int)(rand() / (RAND_MAX + 1.0) * RANGE)

但是这种方法解决了问题2,却并没有解决问题1,因为如果RAND_MAX+1不能被RANGE整除,不管按照什么规则,RAND_MAX+1个数始终是无法均分放到RANGE个格子里的,输入数据量较大时这个bias就会比较明显。

 

方法2

去除x=rand()所产生最大的(RAND_MAX+1)%RANGE个数,再(int)x/(int)((RANG_MAX+1.0)/RANGE)

int myrandom_range(int start,int end)//[start,end)
{
int N = end - start;
int bound = (RAND_MAX + 1) / N * N;
int r;
do
{
r = rand();
} while(r >= bound);
return start + int(r / (RAND_MAX + 1.0) * N);
}

这种方法解决了问题1&2,但是这个做法在RAND_MAX%RANGE (< RANGE)比较大的时候,可能会很浪费时间(极端情况是接近1/2概率需要重试),一般情况下还是够用的。

Let p = (RAND_MAX % RANGE) / (RAND_MAX + 1.0). This is the probability that any given call to rand() will require a retry. Note that this value is maximized when RANGE*2 = RAND_MAX+3, and which will yield a value of p roughly equal to 1/2.

  1. The average number of times that rand() will be called is 1/(1-p) (with a worst case of about 2).
  2. The probability that at least n calls to rand() will be required beyond the first one is pn.

So for most values of p which will be much less than 1/32 say, performance should not be a concern at all.

 

再讨论RANGE>RAND_MAX+1的情况

现在我们来考虑一个相对简单的问题,RANGE能被RAND_MAX+1整除的情况,例如数据范围[0,2^32)。

方法1加强版:

前面提到,RAND_MAX其值最小为(2^15-1),最大为(2^31-1),按最坏的情况算至少需要调用三次rand(),才能获得32位随机数,我们采取方法1中取高位的做法,三次rand()分别取低15位中的前10,11,11位。

inline unsigned __int32 rand32()
{
return ((rand()&0x00007FE0)>>5) + ((rand()&0x00007FF0)<<6) + ((rand()&0x00007FF0)<<17);
}

 

方法3:

可以考虑分段抽样,分成[n/(RNAD_MAX+1)]段,先等概率得到段再得到每段内的某个元素,这样分段也类似地有一个尾数问题,不是每次都刚好分到整数段,一定或多或少有一个余数段,这部分的值如何选取?

选到余数段的数据拿出来选取,先进行一次选到余数段概率的事件发生,然后进行单独选取:
r = N % (RAND_MAX+1); //余数
if ( happened( (double)r/N ) )//选到余数段的概率
result = N-r+myrandom(r); // myrandom可以用情况1中的代码实现
else
result = rand()+myrandom(N/(RAND_MAX+1))*(RAND_MAX+1); // 如果选不到余数段再进行分段选取

const double MinProb=1.0/(RAND_MAX+1);
bool happened(double probability)//probability 0~1
{
if(probability<=0)
{
return false;
}
if(probability
{
return rand()==0&&happened(probability*(RAND_MAX+1));
}
if(rand()<=probability*(RAND_MAX+1))
{
return true;
}
return false;
}
 
long myrandom(long n)//产生0~n-1之间的等概率随机数
{
t=0;
if(n<=RAND_MAX)
{
long R=RAND_MAX-(RAND_MAX+1)%n;//尾数
t = rand();
while ( t > r )
{
t = rand();
}
return t % n;
}
else
{
long r = n%(RAND_MAX+1);//余数
if( happened( (double)r/n ) )//取到余数的概率
{
return n-r+myrandom(r);
}
else
{
return rand()+myrandom(n/(RAND_MAX+1))*(RAND_MAX+1);
}
}
}

 

最后,使用之前,别忘了设随机种子srand((unsigned int)time(0))哦。

 

如果你懒得自己写,干脆用, produces integers in the range [0, 232-1].

 

参考:

转载于:https://www.cnblogs.com/wei-li/archive/2012/10/04/2711629.html

你可能感兴趣的文章
asp.net怎样在URL中使用中文、空格、特殊字符
查看>>
路由器发布服务器
查看>>
实现跨交换机VLAN间的通信
查看>>
jquery中的data-icon和data-role
查看>>
python例子
查看>>
环境变量(总结)
查看>>
ios之UILabel
查看>>
Java基础之String,StringBuilder,StringBuffer
查看>>
1月9日学习内容整理:爬虫基本原理
查看>>
安卓中数据库的搭建与使用
查看>>
AT3908 Two Integers
查看>>
渐变色文字
查看>>
C++ 0X 新特性实例(比较常用的) (转)
查看>>
node生成自定义命令(yargs/commander)
查看>>
各种非算法模板
查看>>
node-express项目的搭建并通过mongoose操作MongoDB实现增删改查分页排序(四)
查看>>
如何创建Servlet
查看>>
.NET 设计规范--.NET约定、惯用法与模式-2.框架设计基础
查看>>
win7 64位+Oracle 11g 64位下使用 PL/SQL Developer 的解决办法
查看>>
BZOJ1997:[HNOI2010]PLANAR——题解
查看>>