Part0:隨機數的性質

隨機數一般來說符合下面這幾個性質.

(馬爾科夫性)(1.)它產生時后面那個數與前面的毫無關系.

(不確定性)(2.)給定樣本的一部分和隨機算法,無法推出樣本的剩余部分.

(不可再現性)(3.)其隨機樣本不可重現.

另外還要說一下統計學偽隨機數概念.

統計學偽隨機性.統計學偽隨機性指的是在給定的隨機比特流樣本中,1的數量大致等于0的數量,同理,"10""01""00""11"四者數量大致相等.類似的標準被稱為統計學隨機性.滿足這類要求的數字在人類"一眼看上去”是隨機的.(摘自百度詞條)

實際上這也是在計算機中對偽隨機數優劣的概念.

Part1:偽隨機數

偽隨機數,也就是我們C++中常用的隨機數.為什么說它是"偽"隨機呢?其實只要稍微詳細的了解過C++rand函數的人應該都能懂得這一點.

因為計算機本身并不能夠產生和隨機數,只能通過產生一組循環節很長的數來"偽造"隨機.

C++的rand函數實際上只是根據你提供的種子計算出來的一組循環節很長的數.只要兩次提供的種子是一樣的,那么rand函數提供的隨機數也是一樣的.

那么,循環節到底長到什么程度呢?

事實上,rand()函數在LINUX系統下的實現如下：

static unsigned long next=1;

// RAND_MAX assumed to be 32767 

int rand()
{
    next=next*1103515245+12345;
    return ((unsigned)(next/65536)%32768);
}

void srand(unsigned seed)
{
    next=seed;
}

通過這個算法我們可以推知,rand函數的循環節應該是在(32768(2^{15}))以內.

因此,在計算機方面,目前來說,如果不借助外部幫助,是無法達到真隨機的.

Part2:隨機數的優劣判定

在講隨機數算法之前,應該先講講隨機數優劣的判定.畢竟只有清除了隨機數的優劣,我們才能說如何生成優質隨機數.

在這里我們就要用到前面說的統計學偽隨機性:

統計學偽隨機性.統計學偽隨機性指的是在給定的隨機比特流樣本中,1的數量大致等于0的數量,同理,"10""01""00""11"四者數量大致相等.類似的標準被稱為統計學隨機性.滿足這類要求的數字在人類"一眼看上去”是隨機的.(摘自百度詞條)

結合計算機隨機數的特性,我們能夠得出以下三項判斷隨機數優劣的性質:

(1.)隨機程度,即隨機算法是否足夠復雜,隨機數之間會不會存在明顯聯系.

(2.)分布范圍,即是否存在隨機數在分布區域內大量出現偏大偏小的現象,分布是否足夠平均.

(3.)循環長度,即是否會在大量調用時很快地出現循環的情況.

有了這些評判規則,我們就比較好學習優質隨機數的生成.

Part3:基于C++rand的優質隨機數的生成

我們先來講一下基于C++rand函數的隨機數生成.

1.來回擺動式

這種隨機數主要是針對退火算法之類的需要用隨機數來修正答案的.既然是修正答案,那么我們希望最好是來回擺動,一正一負.這種隨機數的特點便是通過一部分人工處理,將原本的rand函數產生的隨機數變成正負交替的.

static int f=3000;
static double del=0.999;// f和del是用來控制隨機數幅度不斷變小的 
static int con=-1;
static int g=1; // 控制正負交替 

int rand1()
{
	f*=del;
	g*=con;
	int tmp=f*g*rand();
	return tmp;
}

這種隨機數的產生引入了退火的思路,當然,你也可以直接使用算法中現成的溫度來控制.

2.平均式

這種主要是用于平衡樹treap的,特點就是在保證單個數隨機的情況下在整體上保證分布比較平均.

int p; // 希望的分布位置

int rand2()
{
    int tmp=(p+rand())/2; // 通過取于分布位置的平均數,是產生的數更加靠近期望分布 
    return tmp;
}

3.多次調用不重復式

當然,如果有人真的需要非常接近真隨機的數,也就是多次運行程序也不會出現相同的情況,那就需要用到一定的外部干擾了.

首先是clock函數.上文已經說過,一個程序在不斷調用期間.每一次的運行時間都會有細小的變化.我們就可以利用好這個變化.每次調用完后都重置一次隨機數種子.

還有一個可能大家都會忽視的方法,計算機本身的誤差.眾所周知,計算機在做浮點運算時是會產生精度損失的,那么我們也可以利用這個特點輔助``clock調整種子(畢竟程序調用時間相同其實可能性也不小,畢竟clock```只精確到( ext{ms})).

int count;

int rand3()
{
	++count;
	int t=clock()+1; // 使用當前時間 
	
	for(int i=1;i<12121307;++i) // 降速
		t+=rand();
	
	t+=clock();  // 降速后擴大時間變化 
	t*=-1234;
	srand(t*count+rand()); // 重置隨機數種子 
	return rand();
}

經過大量實驗,筆者發現該函數前三個數出現重復幾率相對會比較大(7~9%)建議從第四個開始使用.

上面的代碼并沒用用精度損失來隨機化,因為同一個式子的進度損失值太小,以至于幾乎不會發生什么改變,所以并沒有使用.

優劣度分析

首先,隨機程度方面,雖然之前看過rand()函數代碼,可能清楚數字之間的關聯.但在實際運用中,這個數字之間的關聯還是基本可以忽略的.所以在隨機程度方面,rand()函數還是能夠勉強通過的.

在平均分布方面,單看代碼可能感覺不出來.

那么,筆者就做一個測試:

int data[10007];

int main()
{
    for(int i=1;i<=1000000;++i)
    {
        int tmp=rand()%100000; // 生成一個100000以內的隨機數 
        ++data[tmp/10]; // 統計出現次數 
    }

    for(int i=1;i<=1000;++i)
        printf("%d
",data[i]);
}

最后結果：

從中我們可以看到,這個分布還是非常平均的.

循環長度...

這個主要就是rand()函數的硬傷了,(32768)這個長度真的挺不夠用的.在需要大量調用rand()函數的算法中(比如退火),基本都會把rand()卡出循環.

那有沒有既優質又循環節長的算法呢?

Part4:梅森旋轉算法(MT19937)

梅森旋轉算法是目前產生優質偽隨機數的普遍算法,在C++11的標準中,也增加了MT19937的實現,在random庫中.

其實,這個算法是由松本真和和西村拓士在1997年開發的一種算法,和梅森沒有多大關系.它之所以有這個名字,是因為它有長達(2^{19937}-1)的循環節,這是一個梅森素數.況且,這種算法還能在如此長的循環節下產生均勻的隨機數.

那么,MT19937的原理究竟是什么呢?

這個旋轉算法實際上是對一個(19937)位的二進制序列作變換.我們知道,一個長度為(n)的二進制序列,它的排列長度最長為(2^n).但是,有時候會因為某些操作不當,導致循環節小于(2^n).而如何將這個序列的排列恰好達到(2^n)個,就是這個旋轉算法的精髓.

如果反饋函數的本身(+1)是一個既約多項式,那么它的循環節達到最長,即(2^n-1).

我們拿一個4位寄存器模擬一下:

我們這里使用的反饋函數是(y=x^4+x^2+x+1)(這個不是既約多項式,只是拿來好理解)

這個式子中(x^4),(x^2),(x)的意思就是我們每次對這個二進制序列的從后往前數第4位和第2位做異或運算,然后(x)的意思是我們再拿結果和最后一位做異或運算.把最后的結果放到序列的開頭,整個序列后移一位,最后一位舍棄.

1.初始數組({1,0,0,0}).

2.將它的第四位和第二位取出來做異或運算.

3.把剛剛的運算結果和最后一位再做一次運算

4.把最后的運算結果放到第一位,序列后移.最后一位被拋棄.

這就是一次運算,然后這個算法就是不斷循環這幾步,從而不斷偽隨機改變這個序列.

因為它所使用的反饋函數(y=x^4+x+1)是既約多項式,所以最后的循環節為(2^4-1=15),運算結果如下：

[egin{array}
{|c|c|c|c|}\
a_3&a_2&a_1&a_0\
1&0&0&0\
1&1&0&0\
1&1&1&0\
1&1&1&1\
0&1&1&1\
1&0&1&1\
0&1&0&1\
1&0&1&0\
1&1&0&1\
0&1&1&0\
0&0&1&1\
1&0&0&1\
0&1&0&0\
0&0&1&0\
0&0&0&1\
-&-&-&-\
1&0&0&0
end{array}
]

大家可以看到,這個運算結果包含了(1,2,...,2^4-1)中的所有整數,并且沒有循環,同時擁有很好的隨機性.

Part5:MT19937的偽代碼及C++實現

初始化隨機種子:

(indexleftarrow 0)
(MTleftarrow new array with size 624//624 imes 32-31=19937)
(// ext{Above are global variables})
( ext{MT19937_SRAND}(seed):)
(indexleftarrow 0)
(MT[0]leftarrow seed)
(mathbf{for} ileftarrow 1 mathbf{to} 623:)
(quad tleftarrow 1812433253cdot(MT[i-1]oplus(MT[i-1]gg 30))+i//oplus ext{ is the xor operation}, gg ext{ is the right-shift operation})
(quad MT[i]leftarrow t& ext{0xffffffff}// ext{get the last 32 bits})
(//& ext{ is the bit-and operation}, ext{0x means that the number next is a hex number})

梅森旋轉:

( ext{MT19937_GENERATE}():)
(mathbf{for} ileftarrow 0 mathbf{to} 623:)
(quad yleftarrow (MT[i]& ext{0x80000000})+(MT[(i+1)mod 624]& ext{0x7fffffff}))
(quad MT[i]leftarrow MT[(i+397)mod 624]oplus(ygg 1))
(quadmathbf{if} y&1=1:)
(quadquad MT[i]leftarrow MT[i]oplus 2567483615)

生成隨機數:

( ext{MT19937_RAND}():)
(mathbf{if} index=0:)
(quad ext{MT19937_GENERATE}())
(yleftarrow MT[index])
(yleftarrow yoplus (ygg 11))
(yleftarrow yoplus ((yll 7)&2636928640))
(yleftarrow yoplus ((yll 15)& 4022730752))
(yleftarrow yoplus (ygg 18))
(indexleftarrow (index+1)mod 624)
(mathbf{return} y)

C++實現:

int index;
int MT[624];

// 設置隨機數種子
inline void srand(int seed)
{
	index=0;
	MT[0]=seed;
	
	for(int i=1;i<=623;++i)
	{
		int t=1812433253*(MT[i-1]^(MT[i-1]>>30))+i;
		MT[i]=t&0xffffffff;
	}	
} 

// 梅森旋轉 
inline void generate()
{
	for(int i=0;i<=623;++i)
	{
		int y=(MT[i]&0x80000000)+(MT[(i+1)%624]&0x7fffffff);
		MT[i]=MT[(i+397)%624]^(y>>1);
		
		if(y&1)
			MT[i]^=2567483615;
	}
}

// 生成隨機數
inline int rand()
{
	if(index==0)
		generate();
	
	int y=MT[index];
	y=y^(y>>1);
	y=y^((y<<7)&2636928640);
	y=y^((y<<15)&4022730752);
	y=y^(y>>18);
	index=(index+1)%624;
	return y;
 } 

本文完

總結

以上是生活随笔為你收集整理的浅谈随机数的生成的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。