A study for MT19937

寒假前被一个伪随机数卡了，终于有时间来看看

前言

如前文，不会就得学（悲

参考博客

MT19937

​ MT19937是许多语言的默认伪随机数生成器，其优点为：周期非常长（达到 $2^{19937}-1$ ）、在 $1\leq k\leq 623$ 的维度之间都可以均等分布、速度较快。

​ MT19937算法主要分为三个部分：

1
2
3
1、基础的梅森旋转链
2、输出处理算法（Tamper）
3、旋转链的旋转算法（Twist）

python中Random类采用MT19937算法，getrandbits(32)返回一个32位随机数，详细代码如下：

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
def _int32(x):
    return int(0xFFFFFFFF & x)

class MT19937:
    def __init__(self, seed):
        self.mt = [0] * 624
        self.mt[0] = seed
        self.mti = 0
        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):
        if self.mti == 0:
            self.twist()
        y = self.mt[self.mti]
        y = y ^ y >> 11
        y = y ^ y << 7 & 2636928640
        y = y ^ y << 15 & 4022730752
        y = y ^ y >> 18
        self.mti = (self.mti + 1) % 624
        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 >> 1) ^ self.mt[(i + 397) % 624]

            if y % 2 != 0:
                self.mt[i] = self.mt[i] ^ 0x9908b0df

针对题目

1.逆向extract_number

我们发现这其中有两种操作是较难的，我们先来看这一种

1
y = y ^ y >> 11

显而易见，这种操作并不会对y的高18位有任何影响，所以y'=y ^ y >> 11的高18位就是y的高18位，而y'接下来的18位正是y的高18位和y接下来的18位异或得来的，所以我们可以使用已知的高18位得到接下来的18位，同上，每操作一次我们就可得到y的18位，所以我们就能利用y'通过有限的操作后复原之前的y，实现代码如下：

1
2
3
4
5
6
y1 = 
y = y1^y1>>11
for i in range(32//11):
    y = y1^(y>>11)
print(y)
print(y1)

分析完这一种，我们来看下一种

1
y = y ^ y << 7 & 2636928640

类比上一种，y'的低7位实际上就是y^2636928640的结果，可以轻松的获得y的低7位，我们可以将y<<7&2636928640当做一个整体（也就是说掩码其实并没有实质性的作用），这样模仿前一种操作，复原代码如下：

1
2
3
4
5
6
7
y1 = 
y = y1^ y1 << 7 & 2636928640
tmp = y
for i in range(32 // 7):
    y = tmp ^ y << 7 & 2636928640
print(y)
print(y1)

我们将左移和右移（带掩码与不带掩码）的复原代码整理如下：

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
class unBitShift:
    def __init__(self):
        pass

    def RightXor(self,value,shift,bits=32):
        tmp = value
        for i in range(bits // shift):
            tmp = value ^ tmp >> shift
        return tmp

    def RightXorMasked(self,value,shift,mask,bits=32):
        tmp = value
        for i in range(bits // shift):
            tmp = value ^ tmp >> shift & mask
        return tmp

    def LeftXor(self,value, shift,bits=32):
        tmp = value
        for i in range(bits // shift):
            tmp = value ^ tmp << shift
        return tmp

    def LeftXorMasked(self,value, shift, mask,bits=32):
        tmp = value
        for i in range(bits // shift):
            tmp = value ^ tmp << shift & mask
        return tmp

2.逆向Twist

观察twist函数，可知新的state[i]只和原来的state[i]、state[i+1]与state[i+397]有关，首先异或掉state[i+397]，之后我们讨论是否将结果与0x9908b0df异或，因为0x9908b0df=0b10011001000010001011000011011111，又因为y>>1操作得知之前的数的高位一定是0，所以我们之间判断最高位即可知道是否需要异或0x9908b0df，也就是能得知原来state[i]的最高位，那么再对state[i-1]进行还原，即可得到完整的state[i]了，具体代码如下：

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
def backtrace(cur):
    high = 0x80000000
    low = 0x7fffffff
    mask = 0x9908b0df
    state = cur
    for i in range(623,-1,-1):
        tmp = state[i]^state[(i+397)%624]
        # recover Y(tmp)
        if tmp & high == high:
            tmp ^= mask
            tmp <<= 1
            tmp |= 1
        else:
            tmp <<=1
        # recover highest bit
        res = tmp&high
        # recover other 31 bits
        tmp = state[i-1]^state[(i+396)%624]
        # recover Y(tmp)
        if tmp & high == high:
            tmp ^= mask
            tmp <<= 1
            tmp |= 1
        else:
            tmp <<=1
        res |= (tmp)&low
        state[i] = res    
    return state