如何产生 [0, 2147483647] 之间的随机数

一、简介

随机数是编程中经常要用到的东西，但很遗憾的是在windows下vc和MinGW中的 RAND_MAX 都是32767，

也就是说调用系统的rand()函数只能产生范围在[0, 32767]之间的随机数。

那如何得到范围达到[0, 2147483647]的随机数呢？

二、最直观的方法

最直接的想法显然是调用rand()生成更大区间的随机数，比如两个rand()相加或者相乘，但是，

两个rand()相加时，越大的数出现的概率越高，因为2=2+0=1+1=0+2，而5=0+5=1+4=2+3=3+2=4+1=5+0；

两个rand()相乘时，是永远不可能得到更大区间范围内的质数的，比如两个[0, 4]范围内的随机数相乘的所有组合都不可能得到7。

综上，似乎要寻找一种更靠谱的方法。

三、基于按位随机的随机数产生方法

首先我们知道，由一个较大的区间产生较小区间的随机数是十分方便的，只要把多余的数剔除掉，然后等分区间就可以了。

比如，由[0, 17]之间的随机数发生器A产生[0, 2]之间的随机数发生器B，只需要将 A 做一个映射：

[0, 4] --> 0

[5, 9] --> 1

[10, 14] --> 2

[15, 17] --> 丢弃

这样0、1、2出现的概率都是相等的，对于A而言都是5/18，对于B而言都是1/3。

特别的，给定任意一个随机数发生器A（产生[0, m]之间的随机数，m>0)，一定能用这个随机数发生器产生[0, 1]的随机数发生器B：

m=1时A即为[0, 1]随机数发生器，m>1时，

若m为奇数，则将A产生的数直接%2；若m为偶数，若A产生的数等于m则丢弃，其它情况直接%2。

这样就得到了[0, 1]随机数发生器B。

那么，无论我们需要的随机数范围有多大，只要按位随机出0或者1，那么这个数一定随机的。

比如需要[0, 2147483647]之间的随机数，那么只要按位随机出31个二进制位，那么就得到了[0, 2147483647]之间的随机数。

实现代码如下：

View Code

 1 View Code 
 2  /*生成随机数*/
 3  int randrange(int min, int max)
 4  {
 5      int range = max - min + 1;
 6  
 7      if (range < 1)
 8      {
 9          throw std::range_error("The upper limit is less than the lower limit!");
10          return 0;
11      }
12      //系统rand()函数产生的随机数范围是0到RAND_MAX(32767)
13      //这里产生0到2147483647之间的随机数
14      std::bitset<31> randomset;
15      for (int i = 0; i != 31; ++i)
16      {
17          randomset[i] = (rand() % 2 == 1);
18      }    
19      int randomnumber =  static_cast<int>(randomset.to_ulong());
20  
21      return (randomnumber % range + min);
22  }

看上去这个思路应该是没有问题的，因为每一位都是随机的，那整个数一定是随机的，但是测试结果不尽如人意。

测试代码：

View Code

 1 std::bitset<2147483648> randomSet;
 2 
 3 int main(void)
 4 {
 5     unsigned int seed = (unsigned int)time(0) + rand();
 6     long long maxCount = 21474836480L;
 7     srand(seed);
 8     randomSet.reset();
 9 
10     std::cout << randomSet.count() << std::endl;
11     std::cout << "seed = " << seed << std::endl;
12     
13     for (long long i = 0; i != maxCount; ++i)
14     {
15         int randNum = randRange(0, 2147483647);
16         randomSet[randNum] = true;
17     }
18         
19     std::cout << randomSet.count() << std::endl << std::endl;
20 
21     system("pause");
22     return 0;
23 }

测试结果是：

randomSet.count()等于131066

测试代码中for循环执行了21474836480次，如果随机算法 “足够好”，那么randomSet中每一位都会有10次机会被置成true。

然而有测试结果来看，randomSet中却只有131066位被置为ture，这与2147483648相差甚远，看来这个随机算法真是 “足够不好”。

因此我们可以“武断”的认为，基于rand()来按位随机是完全不可行的。

那为什么这种做法结果不对呢。。。。？？？？

四、系统的rand()函数是怎么实现的

其实这种按位随机的做法本没问题，但这基于一个假设：用于产生每位0或1的rand()必须是绝对随机的。

我对绝对随机（假设在[0, 99]范围内）的理解是：

1、要保证在随机次数足够大的情况下，随机出来的数在[0, 99]区间内服从均匀分布。这意味着如果[0, 99]的数不是等概率取到就不是绝对随机的。

2、每次产生的随机数必须是相互独立的。这意味着{0, 1, 2, ……, 99}这种有规律的数列虽然满足第一个条件，但它也不是绝对随机的。

那我们系统随机函数rand()是否满足绝对随机呢？

经过google，并查看vs2010和GLibC2.16.0中对rand()函数的定义可以知道，rand()能满足第一个条件却不满足第二个条件，下面看源码：

vs2010的 rand.c 默认在 “C:\Program Files\Microsoft Visual Studio 10.0\VC\crt\src”下，其源码如下：

rand.c

View Code

 1 void __cdecl srand (
 2         unsigned int seed
 3         )
 4 {
 5         _getptd()->_holdrand = (unsigned long)seed;
 6 }
 7 
 8 int __cdecl rand (
 9         void
10         )
11 {
12         _ptiddata ptd = _getptd();
13         /*!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*/
14         return( ((ptd->_holdrand = ptd->_holdrand * 214013L
15             + 2531011L) >> 16) & 0x7fff );
16 }

GLibC2.16.0在这里下载，其rand()函数的一套源码如下：

/stdlib/rand.c

View Code

1 /*   /stdlib/rand.c   */
2 int
3 rand ()
4 {
5   return (int) __random ();
6 }

/stdlib/random.c

View Code

 1 /*   /stdlib/random.c    */
 2 long int
 3 __random ()
 4 {
 5   int32_t retval;
 6 
 7   __libc_lock_lock (lock);
 8 
 9   (void) __random_r (&unsafe_state, &retval);
10 
11   __libc_lock_unlock (lock);
12 
13   return retval;
14 }
15 
16 weak_alias (__random, random)

/stdlib/random_r.c

View Code

 1 /*   /stdlib/random_r.c    */
 2 int
 3 __random_r (buf, result)
 4      struct random_data *buf;
 5      int32_t *result;
 6 {
 7   int32_t *state;
 8 
 9   if (buf == NULL || result == NULL)
10     goto fail;
11 
12   state = buf->state;
13 
14   if (buf->rand_type == TYPE_0)
15     {
16       int32_t val = state[0];
17       /*!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*/
18       val = ((state[0] * 1103515245) + 12345) & 0x7fffffff;
19       state[0] = val;
20       *result = val;
21     }
22   else
23     {
24       int32_t *fptr = buf->fptr;
25       int32_t *rptr = buf->rptr;
26       int32_t *end_ptr = buf->end_ptr;
27       int32_t val;
28 
29       val = *fptr += *rptr;
30       /* Chucking least random bit.  */
31       *result = (val >> 1) & 0x7fffffff;
32       ++fptr;
33       if (fptr >= end_ptr)
34     {
35       fptr = state;
36       ++rptr;
37     }
38       else
39     {
40       ++rptr;
41       if (rptr >= end_ptr)
42         rptr = state;
43     }
44       buf->fptr = fptr;
45       buf->rptr = rptr;
46     }
47   return 0;
48 
49  fail:
50   __set_errno (EINVAL);
51   return -1;
52 }
53 
54 weak_alias (__random_r, random_r)

注意代码中加了/*!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*/标记的语句，可以看到库函数产生随机数用的是线性同余法。

所谓线性同余法是基于 X(n+1) = (a * X(n) + c) % m 这样的公式递推产生随机数，其中a、c、m是常数，随机数种子也就是用以确定X(0)的数。

可以看到，在vs2010中，a=214013， c=2531011， m=32768(&0x7fff)，在GLibC中a=1103515245， c=12345， m=2147483648(&0x7fffffff)。

所以vs2010中产生的随机数范围是[0, 32767]，而GLibC中则可达到[0, 2147483647]。后面所说的rand()默认指vs2010中的rand()。

因为线性同余法的随机数是由递推公式产生的，虽然能够保证大样本下的近似均匀分布（后面会有测试），但显然不能满足连续随机数之间的相互独立。

一般而言，采用这种方法会用time(0)作为随机数种子，time(0)的返回值是从1970年1月1日00：00到当前时刻所经历的时间（单位秒）。

对于程序的运行而言，如果不要求程序开始运行的时刻（不受控制的），我们可以近似认为随机数种子是不可预知的，因此，即使我们知道a、c、m的值也

无法预测会产生什么样的随机数。实际中我们可以认为rand()产生的数是“随机数”正是基于以上的原因（不可预测性）。

然而如果我们知道a、c、m的值和X(n)，那么之后递推产生的所有的数都是可预知的。

所以，由于前面的方法需要连续调用rand()函数来随机出每一个二进制位，这本质上和调用一次rand()函数没有区别，除了第一次rand()之外后面都是有规律并直接取决于前一次结果。

那每次rand()都重置随机数种子呢？同样也是不行的，在编程实现上，除了程序开始运行的时刻是不可控（近似认为）的之外，后面的都是完全确定的，所以每次重置随机数

种子是在编程实现上要么一直srand()，而time(0)的值在1s内是不变的，要么每次Sleep(rand())一段时间，而这个rand()的结果已经由上次的随机数种子确定了。

综上所述，若要利用rand()产生我们自己的随机数，那么产生每个数的时候rand()最多只能用一次。

遗憾的是，从小空间（比如[0, 32767]）到大空间（比如[0, 2147483647]）的映射不可能是一一映射，所以用rand()产生[0, 2147483647]范围内的随机数至少

要连续调用两次甚至更多的rand()。基于前面的结论，调用系统函数rand()产生[0, 2147483647]范围内均匀分布的随机数是不可能的。

再次遗憾的是，在MinGW（可以近似认为是windows下的GCC）下RAND_MAX仍然是32767，好吧，只要把GLibC中的rand()代码手动拿来用了。

五、GLibC中rand()函数的简单实现

下面是一个简单的版本，忽略了临界区访问控制等代码（如果不需要并行执行的话）。

实现代码：

View Code

 1 unsigned int latestRandNum = 1;
 2 
 3 /*设置随机数种子*/
 4 void srandGLibC(unsigned int seed)
 5 {
 6     latestRandNum = (seed == 0 ? 1 : seed);
 7 }
 8 
 9 /*生成[0, 2147483647]之间的随机数*/
10 int randGLibC(void)
11 {
12     int randNum = ((latestRandNum * 1103515245) + 12345) & 0x7fffffff;
13     latestRandNum = randNum;
14     
15     return (latestRandNum);    
16 }
17 
18 /*在指定范围内生成随机数*/
19 int randRange(int min, int max)
20 {
21     unsigned int range = max - min + 1;
22 
23     if (range == 1 || range > 2147483648)
24     {
25         return min;
26     }
27     else if (range < 1)
28     {
29         throw std::range_error("The upper limit is less than the lower limit!");
30         return -1;
31     }
32 
33     return (randGLibC() % range + min);
34 }

测试代码：

View Code

 1 std::bitset<2147483648> randomSet;
 2 
 3 int main(void)
 4 {
 5     unsigned int maxCount = 2147483648L;
 6 
 7     for (int count = 0; count != 10; ++count)
 8     {
 9         unsigned int seed = (unsigned int)time(0) + randGLibC();
10         srandGLibC(seed);
11         randomSet.reset();
12         std::cout << randomSet.count() << std::endl;
13         std::cout << "seed = " << seed << std::endl;
14         
15         for (long long i = 0; i != maxCount; ++i)
16         {
17             int randNum = randRange(0, 2147483647);
18             randomSet[randNum] = true;
19         }
20 
21         std::cout << randomSet.count() << std::endl << std::endl;
22     }
23 
24     system("pause");
25     return 0;
26 }

测试结果：

View Code

 1 0
 2 seed = 2456847069
 3 2147483648
 4 
 5 0
 6 seed = 3328132247
 7 2147483648
 8 
 9 0
10 seed = 2984931043
11 2147483648
12 
13 0
14 seed = 1708563292
15 2147483648
16 
17 0
18 seed = 2508745454
19 2147483648
20 
21 0
22 seed = 2926070301
23 2147483648
24 
25 0
26 seed = 2299831591
27 2147483648
28 
29 0
30 seed = 2616227439
31 2147483648
32 
33 0
34 seed = 3271328796
35 2147483648
36 
37 0
38 seed = 2417430222
39 2147483648

一共进行了10次测试，每次用不同的随机数种子连续在[0, 2147483647]区间内rand 2147483648次。

而结果是，每次测试中这2147483648个数都正好各被随机到一次。

而把实现代码中

int randNum = ((latestRandNum * 1103515245) + 12345) & 0x7fffffff;

改成

int randNum = ((latestRandNum * 1103515245) + 12344) & 0x7fffffff;

同样用不同的随机数种子进行10次测试，测试结果是：

View Code

 1 0
 2 seed = 2456847067
 3 536870912
 4 
 5 0
 6 seed = 3268585404
 7 268435456
 8 
 9 0
10 seed = 1792994210
11 134217728
12 
13 0
14 seed = 2662066861
15 536870912
16 
17 0
18 seed = 2955665000
19 268435456
20 
21 0
22 seed = 2643173644
23 268435456
24 
25 0
26 seed = 2339863270
27 33554432
28 
29 0
30 seed = 2380245690
31 134217728
32 
33 0
34 seed = 2764481028
35 268435456
36 
37 0
38 seed = 1370651278
39 67108864

由上面两次测试结果可知，对于线性同余法的递推公式 X(n+1) = (a * X(n) + c) % m

当m=2147483648时，a=1103515245, c=12345是一组能在[0, 2147483647]保证均匀分布的参数，至于为什么，我不知道。

------------------------------- End -----------------------------