跳转至

随机函数

概述

要想使用随机化技巧,前提条件是能够快速生成随机数.本文将介绍生成随机数的常见方法.

随机数与伪随机数

说一个单独的数是「随机数」是无意义的,所以以下我们都默认讨论「随机数列」,即使提到「随机数」,指的也是「随机数列中的一个元素」.

现有的计算机的运算过程都是确定性的,因此,仅凭借算法来生成真正 不可预测不可重复 的随机数列是不可能的.

然而在绝大部分情况下,我们都不需要如此强的随机性,而只需要所生成的数列在统计学上具有随机数列的种种特征(比如均匀分布、互相独立等等).这样的数列即称为 伪随机数 序列.

随机数与伪随机数在实际生活和算法中的应用举例:

  • 抽样调查时往往只需使用伪随机数.这是因为我们本就只关心统计特征.
  • 网络安全中往往要用到(比刚刚提到的伪随机数)更强的随机数.这是因为攻击者可能会利用可预测性做文章.
  • OI/ICPC 中用到的随机算法,基本都只需要伪随机数.这是因为,这些算法往往是通过引入随机数来把概率引入复杂度分析,从而降低复杂度.这本质上依然只利用了随机数的统计特征.
  • 某些随机算法(例如 Moser 算法)用到了随机数的熵相关的性质,因此必须使用真正的随机数.

随机数生成方法

rand

用于生成伪随机数,缺点是比较慢,使用时需要 #include<cstdlib>

调用 rand() 函数会返回一个 [0,RAND_MAX] 中的随机非负整数,其中 RAND_MAX 是标准库中的一个宏,在 Linux 系统下 RAND_MAX 等于 2311.可以用取模来限制所生成的数的大小.

使用 rand() 需要一个随机数种子,可以使用 srand(seed) 函数来将随机种子更改为 seed,当然不初始化也是可以的.

同一程序使用相同的 seed 两次运行,在同一机器、同一编译器下,随机出的结果将会是相同的.

有一个选择是使用当前系统时间来作为随机种子:srand(time(nullptr))

Warning

Windows 系统下 rand() 返回值的取值范围为 [0,215)(即 RAND_MAX 等于 2151),当需要生成的数不小于 215 时建议使用 (rand() << 15 | rand()) 来生成更大的随机数.

关于 rand()rand()%n 的随机性:

  • C/C++ 标准并未关于 rand() 所生成随机数的任何方面的质量做任何规定.
  • GCC 对 rand() 所采用的实现方式,保证了分布的均匀性等基本性质,但具有低位周期长度短等明显缺陷.(例如在笔者的机器上,rand()%2 所生成的序列的周期长约 2106
  • 即使假设 rand() 是均匀随机的,rand()%n 也不能保证均匀性,因为 [0,n) 中的每个数在 0%n,1%n,...,RAND_MAX%n 中的出现次数可能不相同.严格保证均匀性的做法可参考 Daniel Lemire, Fast Random Integer Generation in an Interval

预定义随机数生成器

定义了数个特别的流行算法.如没有特别说明,均定义于头文件 <random>

Warning

预定义随机数生成器于 C++11 标准3开始使用.

梅森缠绕器

梅森缠绕器(Mersenne Twister)由松本与西村于 1998 年提出1,因其中一种实现 MT19937 具有梅森素数 M19937=2199371 这么长的周期而得名.

从 C++11 开始,模板类 std::mersenne_twister_engine 实现基于如上方法的随机数生成器.因其使用起来十分复杂,通常实际使用的是该模板类的特化:std::mt19937std::mt19937_64.其优点是随机数质量高,且速度比 rand() 快很多.

mt19937 基于 32 位梅森缠绕器,由松本与西村于 1998 年提出1,是一个随机数生成器类,效用同 rand(),随机数的范围同 unsigned int 类型的取值范围.使用时用其定义一个随机数生成器即可:mt19937 myrand(seed)seed 可不填,不填 seed 则会使用默认随机种子.其重载了 operator(),需要生成随机数时调用 myrand() 即可返回一个随机数.

另一个类似的生成器是 mt19937_64,基于 64 位梅森缠绕器,由松本与西村于 2000 年提出,使用方式同 mt19937,但随机数范围扩大到了 unsigned long long 类型的取值范围.

代码示例
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
#include <ctime>
#include <iostream>
#include <random>

using namespace std;

int main() {
  mt19937 myrand(SEED);  // 将 SEED 替换为任意整数或换为 time(nullptr)
  cout << myrand() << endl;
  return 0;
}

线性同余随机数生成器

线性同余随机数生成器(简称 LCG)由 Thomson 和 Rotenberg 于 1958 年提出.其计算公式如下,其中 A,C,M 为预定义常数:

si{seedi=0si1×A+Ci1(modM)

从 C++11 开始,模板类 std::linear_congruential_engine 实现基于如上方法的随机数生成器,模板参数如下:

1
2
template <class UIntType, UIntType A, UIntType C, UIntType M>
class linear_congruential_engine;

这里 UIntType 表示 s 的类型,随后是三个常数.

随后 Lewis、Goodman 及 Miller 于 1969 提出,取 A=75=16807C=0M=2311=2147483647 时,此算法可以达到该模数下近乎最大的周期 2312(因为 16807 是该模数下的 原根)且具有相对不错的统计特征.

之后于 1988 年,该参数的 LCG 由 Park 与 Miller 采纳为「最小标准」,只因其实现极其简单、易懂、高效、且质量相对不错.最后于 1993 年,Park、Miller 和 Stockmeyer 将参数 A 改为 48271,成为较新的「最小标准」.两个版本的「最小标准」作为 linear_congruential_engine 的特化,也都于 C++11 中预定义,分别为 minstd_rand0minstd_rand.具体而言:

  • 对于 minstd_rand0si 的类型为 std::uint_fast32_tA16807C0M2147483647

  • 对于 minstd_randsi 的类型为 std::uint_fast32_tA48271C0M2147483647

使用方法也很简单.首先定义 std::minstd_rand myrand(seed);,其中 seed 为种子,不填时默认为 1.然后调用 myrand() 即可获得一个随机数.

如果需要自定义该形式的随机数生成器,且 A,C,M 不是常数,则还可以考虑手写实现.该方法实现难度低,但生成的随机序列周期长度较短(周期最大为 M,但大多数情况下都会比 M 短).

参考实现
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
#include <iostream>
using namespace std;

struct myrand {
  int A, C, M, x;

  myrand(int A, int C, int M) {
    this->A = A;
    this->C = C;
    this->M = M;
    this->x = 0;
  }

  // 生成随机序列的下一个随机数
  int next() { return x = ((long long)A * x + C) % M; }
};

myrand rnd(3, 5, 97);  // 初始化一个随机数生成器

int main() {
  int x = rnd.next();
  cout << x << endl;
  return 0;
}

随机数引擎适配器

在预定义随机数生成器的基础上,可以对其进行二次封装得到新的随机数生成器.具体类名和使用方法请参见 伪随机数生成——随机数引擎适配器 的列表.

下面的代码使用 std::independent_bits_engine 包装 std::minstd_rand 的输出.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
#include <cstdint>
#include <iostream>
#include <random>

using std::cout;

int main() {
  std::independent_bits_engine<std::minstd_rand, 32, std::uint_fast32_t> rng;
  for (int i = 0; i < 10; ++i) {
    cout << rng() << " ";
  }
  return 0;
}

非确定随机数的均匀分布整数随机数生成器

random_device 是一个基于硬件的均匀分布随机数生成器,在熵池耗尽 前可以高速生成随机数.该类在 C++11 定义,需要 random 头文件.由于熵池耗尽后性能急剧下降,所以建议用此方法生成 mt19937 等伪随机数的种子,而不是直接生成.

random_device 是非确定的均匀随机位生成器,尽管若不支持非确定随机数生成,则允许实现用伪随机数引擎实现.目前笔者尚未接到报告称 NOIP 评测机不支持基于硬件的均匀分布随机数生成.但出于保守考虑,建议使用该算法生成随机数种子.

参考代码如下.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
#include <iostream>
#include <map>
#include <random>
#include <string>

int main() {
  std::random_device rd;
  std::map<int, int> hist;
  std::uniform_int_distribution<int> dist(0, 9);
  for (int n = 0; n < 20000; ++n) {
    ++hist[dist(rd)];
  }
  for (auto p : hist) {
    std::cout << p.first << " : " << std::string(p.second / 100, '*') << '\n';
  }
  return 0;
}

可能的输出如下.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
0 : ********************
1 : *******************
2 : ********************
3 : ********************
4 : ********************
5 : *******************
6 : ********************
7 : ********************
8 : *******************
9 : ********************

随机数分布

这里介绍的是要求生成的随机数遵从某一分布的随机数生成器,如 离散均匀分布伯努利分布二项分布几何分布标准正态(高斯)分布

具体类名请参见 伪随机数生成——随机数分布 的列表.

下面的程序模拟了一个六面体骰子.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
#include <iostream>
#include <random>
using std::cout;
using std::mt19937;

using std::uniform_int_distribution;
int main() {
  mt19937 gen(SEED);  // 播种标准 mersenne_twister_engine
  uniform_int_distribution<> dis(1, 6);

  for (int n = 0; n < 10; ++n)
    // 用 dis 变换 gen 所生成的随机 unsigned int 到 [1, 6] 中的 int
    cout << dis(gen) << ' ';
  cout << '\n';
  return 0;
}

其他实现方法

有的时候我们需要实现自己的随机数生成器.下面是一些常用的随机数生成方法.

Xorshift

Xorshift 系列随机数生成器由 George Marsaglia 于 2003 年提出,其主要基于异或一个数的移位这一种操作.在算法竞赛中常用的版本有 xorshift32xorshift64 两种:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
#include <cstdint>

using u32 = std::uint32_t;
using u64 = std::uint64_t;

class xorshift32 {
 private:
  u32 x;

 public:
  xorshift32() : x(2463534242) {}

  explicit xorshift32(u32 s) : x(s) {}

  u32 operator()() {
    x ^= x << 13;
    x ^= x >> 17;
    x ^= x << 5;
    return x;
  }
};

class xorshift64 {
 private:
  u64 x;

 public:
  xorshift64() : x(88172645463325252) {}

  explicit xorshift64(u64 s) : x(s) {}

  u64 operator()() {
    x ^= x << 13;
    x ^= x >> 7;
    x ^= x << 17;
    return x;
  }
};

使用方法与前面介绍的大多数随机数生成器类似:首先定义随机数生成器 xorshift32 rng32(seed)xorshift64 rng64(seed),其中 seed 为种子,不填时使用默认种子;之后调用 rng32()rng64() 即可获得随机数.这两个随机数生成器在种子不为 0 时分别具有 23212641 的周期,因其实现极其简单,常将其和种子下发给选手生成大范围的数据以减少 I/O 开销.

SplitMix

SplitMix 系列随机数生成器由 Stelle,Lee 和 Flood 于 2014 年提出,并应用于 Java 8 的 java.util.SplittableRandom 类中.在算法竞赛界最常用的版本为 splitmix64

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
#include <cstdint>

using u64 = std::uint64_t;

class splitmix64 {
 private:
  u64 x;

 public:
  splitmix64() : x(0) {}

  explicit splitmix64(u64 s) : x(s) {}

  u64 operator()() {
    u64 z = (x += 0x9e3779b97f4a7c15);
    z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
    z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
    return z ^ (z >> 31);
  }
};

其最主要的用途是对一个已经求出的哈希值进行二次哈希,以防止被人通过精心构造的数据卡哈希表卡到超时:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
#include <cstdint>

using u64 = std::uint64_t;

struct MyHash {
  u64 mix64(u64 x) const {
    x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9;
    x = (x ^ (x >> 27)) * 0x94d049bb133111eb;
    return x ^ (x >> 31);
  }

  size_t operator()(size_t x) const {
    // 推荐使用 <chrono> 库的
    // std::chrono::steady_clock::now().time_since_epoch().count() 作为
    // FIXED_RANDOM
    static const u64 FIXED_RANDOM = 0x9e3779b97f4a7c15;
    return mix64(x + FIXED_RANDOM);
  }
};

其亦可用于生成随机数:首先定义 splitmix64 rng64(seed),其中 seed 为种子,不填时使用默认种子;之后调用 rng64() 即可获得随机数.此时较为显然的,该函数具有 264 的周期(因为所有的大数都是奇数,所以 x 具有 264 的周期,且后面对 z 的所有操作都是完全可逆的).这也意味着即使你记不住如上的大参数,你也可以直接随便写几个足够大的奇数结合异或右移来达到相近的效果.

时滞斐波那契随机数生成器

利用下式来生成随机数序列 {Ri}(其中 0<j<k):

RiRijRikmodP

这里的 P 通常取 2 的幂(常用 232264), 表示二元运算符,可以使用加法,减法,乘法,异或.

该方法较传统的线性同余随机数生成器而言,拥有更长的周期,但随机性受初始条件影响较大.

参考实现
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
#include <iostream>
#include <random>
#include <vector>

using namespace std;

struct myrand {
  vector<unsigned int> vec;
  int l, j, k, cur;

  template <class URBG>
  myrand(int l, int j, int k, URBG &&rng) {
    this->l = l;
    this->j = j;
    this->k = k;
    cur = 0;
    for (int i = 0; i < l; i++) {
      vec.push_back(rng());  // 先用其他方法生成随机序列中的前几个元素
    }
  }

  unsigned int next() {
    vec[cur] = vec[(cur - j + l) % l] * vec[(cur - k + l) % l];
    // 这里用 unsigned 类型是为了实现自动对 2^32 取模
    return vec[cur++];
  }
};

// 最后一个可以换成其它随机数生成器,自行设置种子
myrand rnd(11, 4, 7, mt19937(SEED));

int main() {
  unsigned int x = rnd.next();
  cout << x << endl;
  return 0;
}

值得一提的是,在 C++11 中,std::subtract_with_carry_engine 是这一类随机数生成器的其中一种实现,而 std::ranlux24_basestd::ranlux48_base 是该模板类的特化,std::ranlux24std::ranlux48 则是在前者的基础上再套上一个适配器 std::discard_block_engine 得到.不过现在通常不再推荐使用它们.

随机算法

下面介绍在 C++ 标准中定义的一些依赖随机数生成的随机算法.

random_shuffle

std::random_shuffle 于 C++98 引入,用于随机打乱指定序列.使用时需要 #include<algorithm>,常用实现方法为 Fisher-Yates 算法

使用时传入指定区间的首尾指针或迭代器(左闭右开)即可:std::random_shuffle(first, last)std::random_shuffle(first, last, myrand)

内部使用的随机数生成器默认为 rand().当然也可以传入自定义的随机数生成器.

关于 random_shuffle 的随机性:

  • C++ 标准中要求 random_shuffle 在所有可能的排列中 等概率 随机选取,但 GCC2的默认标准库 libstdc++并未 严格执行.
  • GCC 中 random_shuffle 随机性上的缺陷的原因之一,是它使用了 rand()%n 这样的写法.如先前所述,这样生成的不是均匀随机的整数.
  • 原因之二,是 rand() 的值域有限.如果所传入的区间长度超过 RAND_MAX,将存在某些排列 不可能 被产生4
Warning

random_shuffle 已于 C++14 标准中被弃用,于 C++17 标准中被移除.

shuffle

std::shuffle 于 C++11 引入,效用同 random_shuffle.使用时需要 #include<algorithm>

区别在于必须使用自定义的随机数生成器:std::shuffle(first, last, myrand)

GCC2实现的 shuffle 符合 C++ 标准的要求,即在所有可能的排列中等概率随机选取.

下面是用 rand()random_shuffle() 编写的一个数据生成器.生成数据为 「ZJOI2012」灾难 的随机小数据.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
#include <algorithm>
#include <cstdlib>
#include <iostream>
using std::cout;

using std::random_shuffle;

int a[100];

int main() {
  srand(SEED);
  int n = rand() % 99 + 1;
  for (int i = 1; i <= n; i++) a[i] = i;
  cout << n << '\n';
  for (int i = 1; i <= n; i++) {
    random_shuffle(a + 1, a + i);
    int cnt = rand() % i;
    for (int j = 1; j <= cnt; j++) cout << a[j] << ' ';
    cout << 0 << '\n';
  }
}

下面是用 mt19937shuffle() 编写的同一个数据生成器.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
#include <algorithm>
#include <iostream>
#include <random>
using std::cout;
using std::mt19937;

using std::shuffle;

int a[100];

int main() {
  mt19937 rng(SEED);
  int n = rng() % 99 + 1;
  for (int i = 1; i <= n; i++) a[i] = i;
  cout << n << '\n';
  for (int i = 1; i <= n; i++) {
    shuffle(a + 1, a + i, rng);
    int cnt = rng() % i;
    for (int j = 1; j <= cnt; j++) cout << a[j] << ' ';
    cout << 0 << '\n';
  }
}

下面是随机排列前十个正整数的一个实现.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#include <algorithm>
#include <iostream>
#include <iterator>
#include <random>
#include <vector>
using std::copy;
using std::cout;
using std::mt19937;
using std::ostream_iterator;
using std::vector;

using std::shuffle;

int main() {
  vector<int> v = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10};

  mt19937 g(SEED);

  shuffle(v.begin(), v.end(), g);

  copy(v.begin(), v.end(), ostream_iterator<int>(cout, " "));
  cout << "\n";
}

sample

std::sample 于 C++17 引入,用于从序列里随机选择其中 n 个元素,与 shuffle 后提取前 n 个元素的主要区别在于选出来的元素输出后仍保持相对顺序.常用实现方法为 蓄水池抽样法

使用方法为 std::sample(first, last, dest, n, myrand),其中 [first, last) 表示要采样的范围,[dest, dest + n) 表示要输出到的地方,myrand 是指定的随机数生成器.

下面是从 11 个大写字母里随机抽取 4 个的实现.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
#include <algorithm>
#include <iostream>
#include <iterator>
#include <random>
#include <string>
using std::back_inserter;
using std::copy;
using std::cout;
using std::mt19937;
using std::string;

using std::sample;

int main() {
  string in{"ABCDEFGHIJK"}, out;
  mt19937 rng{SEED};
  sample(in.begin(), in.end(), back_inserter(out), 4, rng);
  cout << "Four random letters out of " << in << ": " << out << '\n';
}

参考资料与注释