複製鏈接
請複製以下鏈接發送給好友

梅森旋轉算法

鎖定
梅森旋轉算法Mersenne twister)是一個偽隨機數發生算法。由松本真和西村拓士在1997年開發,基於有限二進制字段上的矩陣線性遞歸。可以快速產生高質量的偽隨機數,修正了古典隨機數發生算法的很多缺陷。
中文名
梅森旋轉算法
外文名
Mersenne twister

梅森旋轉算法應用[編輯]

梅森旋轉算法是R、PythonRubyIDLFree PascalPHP、Maple、Matlab、GNU多重精度運算庫和GSL的默認偽隨機數產生器。從C++11開始,C++也可以使用這種算法。在Boost C++,Glib和NAG數值庫中,作為插件提供。
在SPSS中,梅森旋轉算法是兩個PRNG中的一個:另一個是產生器僅僅為保證舊程序的兼容性 [1]  ,梅森旋轉被描述為“更加可靠”。梅森旋轉在SAS中同樣是PRNG中的一個,另一個產生器是舊時的且已經被棄用。

梅森旋轉算法優點

最為廣泛使用Mersenne Twister的一種變體是MT19937,可以產生32位整數序列。具有以下的優點:
  1. 週期非常長,達到219937−1。 [2]  儘管如此長的週期並不必然意味着高質量的偽隨機數,但短週期(比如許多舊版本軟件包提供的2)確實會帶來許多問題。
  2. 在1 ≤k≤ 623的維度之間都可以均等分佈(參見定義)。
  3. 除了在統計學意義上的不正確的隨機數生成器以外,在所有偽隨機數生成器法中是最快的(當時)

梅森旋轉算法缺點

為了性能,這個算法付出了巨大的空間成本(當時而言):需要 2.5KiB的緩存空間。2011年,松本真和西村拓士針對這一問題提出了一個更小的版本,僅佔127 bits的 TinyMT (Tiny Mersenne Twister)。

梅森旋轉算法算法詳細

整個算法主要分為三個階段:
第一階段:獲得基礎的梅森旋轉鏈;
第二階段:對於旋轉鏈進行旋轉算法;
第三階段:對於旋轉算法所得的結果進行處理;
算法實現的過程中,參數的選取取決於梅森素數,故此得名。

梅森旋轉算法相關代碼

下面的一段偽代碼使用MT19937算法生成範圍在[0, 2−1]的均勻分佈的32位整數:

梅森旋轉算法偽代碼

//創建一個長度為624的數組來存儲發生器的狀態 int[0..623] MT int index = 0  //用一個種子初始化發生器 function initialize_generator(int seed) {     i := 0     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) // 1812433253 == 0x6c078965     } }  // Extract a tempered pseudorandom number based on the index-th value, // calling generate_numbers() every 624 numbers function extract_number() {     if index == 0 {         generate_numbers()     }      int y := MT[index]     y := y xor (right shift by 11 bits(y))     y := y xor (left shift by 7 bits(y) and (2636928640)) // 2636928640 == 0x9d2c5680     y := y xor (left shift by 15 bits(y) and (4022730752)) // 4022730752 == 0xefc60000     y := y xor (right shift by 18 bits(y))     index := (index + 1) mod 624     return y }  // Generate an array of 624 untempered numbers function generate_numbers() {     for i from 0 to 623 {         int y := (MT[i] & 0x80000000)                       // bit 31 (32nd bit) of MT[i]                        + (MT[(i+1) mod 624] & 0x7fffffff)   // bits 0-30 (first 31 bits) of MT[...]         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) // 2567483615 == 0x9908b0df         }     } }

梅森旋轉算法Python 代碼

def _int32(x):    return int(0xFFFFFFFF & x)class MT19937:    def __init__(self, seed):        self.mt = [0] * 624        self.mt[0] = seed        for i in range(1, 624):            self.mt[i] = _int32(1812433253 * (self.mt[i - 1] ^ self.mt[i - 1] >> 30) + i)    def extract_number(self):        self.twist()        y = self.mt[0]        y = y ^ y >> 11        y = y ^ y << 7 & 2636928640        y = y ^ y << 15 & 4022730752        y = y ^ y >> 18        return _int32(y)    def twist(self):        for i in range(0, 624):            y = _int32((self.mt[i] & 0x80000000) + (self.mt[(i + 1) % 624] & 0x7fffffff))            self.mt[i] = y ^ self.mt[(i + 397) % 624] >> 1            if y % 2 != 0:                self.mt[i] = self.mt[i] ^ 0x9908b0df
調用函數MT19937(seed).extract_number()將會返回隨機數,其中seed是已確定的種子。
參考資料