產生偽隨機數兩種常用演算法
阿新 • • 發佈：2019-02-10
我們講的隨機數其實暗指偽隨機數。不少朋友可能想到C語言的rand()，可惜這個函式產生的隨機數隨機性非常差，而且速度很慢，相信幾乎不能勝任一般的應用。
古老的LCG(linear congruential generator)代表了最好的偽隨機數產生器演算法。主要原因是容易理解，容易實現，而且速度快。這種演算法數學上基於X(n+1) = (a * X(n) + c) % m這樣的公式，其中：
模m, m > 0
係數a, 0 < a < m
增量c, 0 <= c < m
原始值(種子) 0 <= X(0) < m
其中引數c, m, a比較敏感，或者說直接影響了偽隨機數產生的質量。
一般而言，高LCG的m是2的指數次冪(一般2^32或者2^64)，因為這樣取模操作截斷最右的32或64位就可以了。多數編譯器的庫中使用了該理論實現其偽隨機數發生器rand()。下面是部分編譯器使用的各個引數值：
Source m a c rand() / Random(L)的種子位
Numerical Recipes
2^32 1664525 1013904223
Borland C/C++
2^32 22695477 1 位30..16 in rand(), 30..0 in lrand()
glibc (used by GCC)
2^32 1103515245 12345 位30..0
ANSI C: Watcom, Digital Mars, CodeWarrior, IBM VisualAge C/C++
2^32 1103515245 12345 位30..16
Borland Delphi, Virtual Pascal
2^32 134775813 1 位63..32 of (seed * L)
Microsoft Visual/Quick C/C++
2^32 214013 2531011 位30..16
Apple CarbonLib
2^31-1 16807 0 見Park–Miller隨機數發生器
LCG不能用於隨機數要求高的場合，例如不能用於Monte Carlo模擬，不能用於加密應用。
LCG有一些嚴重的缺陷，例如如果LCG用做N維空間的點座標，這些點最多位於m1/n超平面上(Marsaglia定理)，這是由於產生的相繼X(n)值的關聯所致。
另外一個問題就是如果m設定為2的指數，產生的低位序列週期遠遠小於整體。
一般而言，輸出序列的基數b中最低n位，bk = m (k是某個整數)，最大週期bn.
有些場合LCG有很好的應用，例如記憶體很緊張的嵌入式中，電子遊戲控制檯用的小整數，使用高位可以勝任。
LCG的一種C++實現版本如下：
//************************************************************************
// Copyright (C) 2008 - 2009 Chipset
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU Affero General Public License as
// published by the Free Software Foundation, either version 3 of the
// License, or (at your option) any later version.
//
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Affero General Public License for more details.
//
// You should have received a copy of the GNU Affero General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.
//************************************************************************
#ifndef LCRANDOM_HPP_
#define LCRANDOM_HPP_
#include <ctime>
class lcrandom
{
public:
explicit lcrandom(size_t s = 0) : seed(s)
{
if (0 == seed) seed = std::time(0);
randomize();
}
void reset(size_t s)
{
seed = s;
if (0 == seed) seed = std::time(0);
randomize();
}
size_t rand()
{
//returns a random integer in the range [0, -1UL)
randomize();
return seed;
}
double real()
{
//returns a random real number in the range [0.0, 1.0)
randomize();
return (double)(seed) / -1UL;
}
private:
size_t seed;
void randomize() { seed = 1103515245UL * seed + 12345UL; }
};
class lcrand_help
{
static lcrandom r;
public:
lcrand_help() {}
void operator()(size_t s) { r.reset(s); }
size_t operator()() const { return r.rand(); }
double operator()(double) { return r.real(); }
};
lcrandom lcrand_help:: r;
extern void lcsrand(size_t s) { lcrand_help()(s); }
extern size_t lcirand() { return lcrand_help()(); }
extern double lcdrand() { return lcrand_help()(1.0); }
#endif // LCRANDOM_HPP_
如果需要高質量的偽隨機數，記憶體充足(約2kb)，Mersenne twister演算法是個不錯的選擇。Mersenne twister產生隨機數的質量幾乎超過任何LCG。不過一般Mersenne twister的實現使用LCG產生種子。
Mersenne twister(梅森旋轉演算法)是Makoto Matsumoto (松本)和Takuji Nishimura (西村)於1997年開發的偽隨機數產生器，基於有限二進位制欄位上的矩陣線性再生。可以快速產生高質量的偽隨機數，修正了古老隨機數產生演算法的很多缺陷。 Mersenne twister這個名字來自週期長度通常取Mersenne質數這樣一個事實。常見的有兩個變種Mersenne Twister MT19937和Mersenne Twister MT19937-64。
Mersenne Twister有很多長處，例如：週期2^19937 - 1對於一般的應用來說，足夠大了，序列關聯比較小，能通過很多隨機性測試。
關於Mersenne Twister比較詳細的論述請參閱http://www.cppblog.com/Chipset/archive/2009/01/19/72330.html
用Mersenne twister演算法實現的偽隨機數版本非常多。例如boost庫中的高質量快速隨機數產生器就是用Mersenne twister演算法原理編寫的。
古老的LCG(linear congruential generator)代表了最好的偽隨機數產生器演算法。主要原因是容易理解，容易實現，而且速度快。這種演算法數學上基於X(n+1) = (a * X(n) + c) % m這樣的公式，其中：
模m, m > 0
係數a, 0 < a < m
增量c, 0 <= c < m
原始值(種子) 0 <= X(0) < m
其中引數c, m, a比較敏感，或者說直接影響了偽隨機數產生的質量。
一般而言，高LCG的m是2的指數次冪(一般2^32或者2^64)，因為這樣取模操作截斷最右的32或64位就可以了。多數編譯器的庫中使用了該理論實現其偽隨機數發生器rand()。下面是部分編譯器使用的各個引數值：
Source m a c rand() / Random(L)的種子位
Numerical Recipes
2^32 1664525 1013904223
Borland C/C++
2^32 22695477 1 位30..16 in rand(), 30..0 in lrand()
glibc (used by GCC)
2^32 1103515245 12345 位30..0
ANSI C: Watcom, Digital Mars, CodeWarrior, IBM VisualAge C/C++
2^32 1103515245 12345 位30..16
Borland Delphi, Virtual Pascal
2^32 134775813 1 位63..32 of (seed * L)
Microsoft Visual/Quick C/C++
2^32 214013 2531011 位30..16
Apple CarbonLib
2^31-1 16807 0 見Park–Miller隨機數發生器
LCG不能用於隨機數要求高的場合，例如不能用於Monte Carlo模擬，不能用於加密應用。
LCG有一些嚴重的缺陷，例如如果LCG用做N維空間的點座標，這些點最多位於m1/n超平面上(Marsaglia定理)，這是由於產生的相繼X(n)值的關聯所致。
另外一個問題就是如果m設定為2的指數，產生的低位序列週期遠遠小於整體。
一般而言，輸出序列的基數b中最低n位，bk = m (k是某個整數)，最大週期bn.
有些場合LCG有很好的應用，例如記憶體很緊張的嵌入式中，電子遊戲控制檯用的小整數，使用高位可以勝任。
LCG的一種C++實現版本如下：
//************************************************************************
// Copyright (C) 2008 - 2009 Chipset
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU Affero General Public License as
// published by the Free Software Foundation, either version 3 of the
// License, or (at your option) any later version.
//
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Affero General Public License for more details.
//
// You should have received a copy of the GNU Affero General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.
//************************************************************************
#ifndef LCRANDOM_HPP_
#define LCRANDOM_HPP_
#include <ctime>
class lcrandom
{
public:
explicit lcrandom(size_t s = 0) : seed(s)
{
if (0 == seed) seed = std::time(0);
randomize();
}
void reset(size_t s)
{
seed = s;
if (0 == seed) seed = std::time(0);
randomize();
}
size_t rand()
{
//returns a random integer in the range [0, -1UL)
randomize();
return seed;
}
double real()
{
//returns a random real number in the range [0.0, 1.0)
randomize();
return (double)(seed) / -1UL;
}
private:
size_t seed;
void randomize() { seed = 1103515245UL * seed + 12345UL; }
};
class lcrand_help
{
static lcrandom r;
public:
lcrand_help() {}
void operator()(size_t s) { r.reset(s); }
size_t operator()() const { return r.rand(); }
double operator()(double) { return r.real(); }
};
lcrandom lcrand_help:: r;
extern void lcsrand(size_t s) { lcrand_help()(s); }
extern size_t lcirand() { return lcrand_help()(); }
extern double lcdrand() { return lcrand_help()(1.0); }
#endif // LCRANDOM_HPP_
如果需要高質量的偽隨機數，記憶體充足(約2kb)，Mersenne twister演算法是個不錯的選擇。Mersenne twister產生隨機數的質量幾乎超過任何LCG。不過一般Mersenne twister的實現使用LCG產生種子。
Mersenne twister(梅森旋轉演算法)是Makoto Matsumoto (松本)和Takuji Nishimura (西村)於1997年開發的偽隨機數產生器，基於有限二進位制欄位上的矩陣線性再生。可以快速產生高質量的偽隨機數，修正了古老隨機數產生演算法的很多缺陷。 Mersenne twister這個名字來自週期長度通常取Mersenne質數這樣一個事實。常見的有兩個變種Mersenne Twister MT19937和Mersenne Twister MT19937-64。
Mersenne Twister有很多長處，例如：週期2^19937 - 1對於一般的應用來說，足夠大了，序列關聯比較小，能通過很多隨機性測試。
關於Mersenne Twister比較詳細的論述請參閱http://www.cppblog.com/Chipset/archive/2009/01/19/72330.html
用Mersenne twister演算法實現的偽隨機數版本非常多。例如boost庫中的高質量快速隨機數產生器就是用Mersenne twister演算法原理編寫的。