Mersenne Twister随机数发生器及随机性测试
一、实验目的
用MATLAB实现Mersenne Twister随机数发生器,并对其随机性进行测试。二、实验原理
伪随机数的产生，首先是选取种子，然后是在此种子基础上根据具体的生成算法计算得到一个伪随机数，然后利用此伪随机数再根据生成算法递归计算出下二个伪随机数，直到将所有不重复出现的伪随机数全部计算出来。这个伪随机数序列就是以后要用到的伪随机数序列。上面的计算过程可以一次性计算完毕，也可以使用一次递归计算一次，每次生成的伪随机数就是这个伪随机数序列中的一个，不过不管怎么样，只要确定了种子，确定了生成算法，这个序列就是确定的了。所谓种子，就是一个对伪随机数计算的初始值。
Mersenne Twister算法是一种随机数产生方法，它是移位寄存器法的变种。该算法的原理：Mersenne Twister算法是利用线性反馈移位寄存器(LFSR)产生随机数的，LFSR的反馈函数是寄存器中某些位的简单异或，这些位也称之为抽头序列。一个n位的LFSR能够在重复之前产生2^n-1位长的伪随机序列。只有具有一定抽头序列的LFSR才能通过所有2^n-1个内部状态，产生2^n - 1位长的伪随机序列，这个输出的序列就称之为m序列。为了使LFSR成为最大周期的LFSR，由抽头序列加上常数1形成的多项式必须是本原多项式。一个n阶本原多项式是不可约多项式，它能整除x^(2*n-1)+1而不能整除x^d+1，其中d能整除2^n-1。例如(32,7,5,3,2,1,0)是指本原多项式x^32+x^7+x^5+x^3+x^2+x+1，把它转化为最大周期LFSR就是在LFSR小邓第32，7，5，2，1位抽头。利用上述两种方法产生周期为m的伪随机序列后，只需要将产生的伪随机序列除以序列的周期，就可以得到(0，1)上均匀分布的伪随机序列了。
伪代码如下：
// 建立624位随机序列数组
int[0..623] MT
int index = 0
//初始化随机序列数组
function initializeGenerator(int seed) {
MT[0] := seed
for i from 1 to 623 {
MT[i] := last 32 bits of(1812433253 * (MT[i-1] xor(right shift by 30 bits(MT[i-1]))) + i) // 0x6c078965
}
}
//根据index的值提取数组中的某个数来生成随机数
// 每624个数调用一次generateNumbers() 函数
function extractNumber() {
if index == 0 {
generateNumbers()
}
int y := MT[index]
y := y xor (right shift by 11 bits(y))
y := y xor(left shift by 7 bits(y) and (2636928640)) // 0x9d2c5680 y := y xor(left shift by 15 bits(y) and(4022730752)) // 0xefc60000 y := y xor (right shift by 18 bits(y))
index := (index + 1) mod 624
return y
}
//产生随机数
function generateNumbers() {
for i from 0 to 623 {
int y := 32nd bit of(MT[i]) + last 31 bits of(MT[(i+1) mod624]) MT[i] := MT[(i + 397) mod 624] xor (right shift by 1 bit(y)) if (y mod 2) != 0 { // y is odd
MT[i] := MT[i] xor (2567483615) // 0x9908b0df
}
}
}
三、程序结构说明
函数的伪代码已经在实验原理中给出。
1，main.m文件：
输入随机序列的种子(默认为1)；输入随机码个数(默认为100)。
调用initializeGenerator函数，初始化随机序列数组。
调用extractNumber()函数生成M个随机序列，存储到result中。
进行随机性测试，画出随机数分布图。
2，initializeGenerator.m文件：
initializeGenerator函数：初始化随机序列数组的长度和值
3，extractNumber.m文件：
extractNumber()函数:根据index的值提取数组中的某个数，调用generateNumbers()函数来生成随机数
4，generateNumbers.m文件：
generateNumbers()函数:产生随机数
四、代码：
五、结果及分析
1.输入随机序列种子和随机码个数：
2.随机码的分布图如下
六、心得体会
MT算法的中文资料很罕见，很多教材、参考书上在讲述随机数发生器时只提到线性同余算法，所以在理解实验原理时遇到了很大困难。
后来终于查找到了关于MT算法的外文资料，认真阅读后消除了编写代码前的障碍。
本实验中用到了MT19937，这是MT的一种变体，可以产生32位整数序列。具有以下的优点:
有219937− 1的非常长的周期。在许多软件包中的短周期—232随机数产生器在一个长周期中不能保证生成随机数的质量。
在1 ≤k≤ 623的维度之间都可以均等分布(参见定义).
除了在统计学意义上的不正确的随机数生成器以外，在所有伪随机数生成器法中是最快的(当时)
