前一篇文章分析了Python中随机算法的实现细节,本文就来对其进行逆向。

由前文所述,MT19937提取随机数可分为两部分:twisttempering

def extract_number(self):
    if self._mti >= self.N:
        for i in range(self.N):
            self.twist(i)
        self._mti = 0

    y = self._mt[self._mti]
    y = self.tempering(y)

    self._mti += 1

    return _int32(y)

那么,其逆向过程就先从termpering操作开始。


逆向 tempering

def tempering(y):
    y ^= (y >> 11)
    y ^= (y << 7) & 0x9d2c5680
    y ^= (y << 15) & 0xefc60000
    y ^= (y >> 18)
    return y

我们倒着一步一步分析,约定记号如下:

  • 「异或」运算记为$\oplus$,「与」运算记为$\wedge$
  • 每一步运算前的变量为$y$,得到的结果为$z$​
  • 记变量最高位的下标为0,第二高位的下标为1,以此类推
  • 变量从高位到低位的连续一段切片,以上下标标记,下标为起点,上标为终点。例如$y$的高18字节记为$y_0^{17}$

先看最后一步:y ^= (y >> 18)

我们知道$z$是32位整数,根据这个公式,显而易见的结论有:

于是:

注意到这个$z\to y$的公式与前面$y\to z$的在形式上一模一样,故这一步的逆向我们只需照抄正向:

y ^= (y >> 18)

再看倒数第二步:y ^= (y << 15) & 0xefc60000

0xefc60000为$c$。

注意到bin(c) == 0b11101111110001100000000000000000,这个二进制数的低17位全为0。

故我们可以写出这一步的正向公式:

同理,容易写出逆公式:

发现形式也相同,故这一步也直接抄正向:

y ^= (y << 15) & 0xefc60000

接下来分析倒数第三步:y ^= (y << 7) & 0x9d2c5680

0x9d2c5680为$d_1$。

这里不容易像前面那样直接写出逆公式,不过我们可以用类似于递归的方法来求解。

首先我们有:

因此:

将此表达式直接代入右边的$y$,得到:

记上式为

我们来计算$X$:

这里$d_2=(d_1\ll7)\wedge d_1=\text{0x94284000}$​

同理,我们可以不断将下式

代入到等号右侧的$y$并展开,我们会得到:

我们记右侧的异或项序列为${X_i}$,即

其中,

计算得:

由此可知,我们在展开到第五项时,彻底消去了等号右侧的$y$,因此:

至此,我们已经可以写出这一步的逆向代码:

y ^= ((y << 7) & 0x9d2c5680) ^ ((y << 14) & 0x94284000) ^ ((y << 21) & 0x14200000) ^ ((y << 28) & 0x10000000)

最后,逆向第一步:y ^= (y >> 11)

类似于上一步,我们可以不断右移再异或,直到右侧的$y$变成0:

注意,$y$是32位整数,右移33位就归零了,因此,第一步的逆向如下:

y ^= (y >> 11) ^ (y >> 22)

综合上述内容,我们可以写出untempering

def untempering(y):
    y ^= (y >> 18)
    y ^= (y << 15) & 0xefc60000
    y ^= ((y << 7) & 0x9d2c5680) ^ ((y << 14) & 0x94284000) ^ ((y << 21) & 0x14200000) ^ ((y << 28) & 0x10000000)
    y ^= (y >> 11) ^ (y >> 22)
    return y

逆向 twist

def twist(self, i):
    y = (self._mt[i] & self.UPPER_MASK) | (self._mt[(i + 1) % self.N] & self.LOWER_MASK)
    self._mt[i] = y >> 1
    if y & 0x1 == 1:
        self._mt[i] ^= self.MATRIX_A
    self._mt[i] ^= self._mt[(i + self.M) % self.N]

逆向twist其实相当于恢复_mt[i]

我们首先写出最后一步的逆向:

tmp = self._mt[i] ^ self._mt[(i + self.M) % self.N]

这里tmp的值应为上面twist函数中经过3、4、5三行之后的_mt[i]的值,那么如何判断在正向过程中是否进入了这个if分支?其实我们只要关注tmp的最高位(32位整数的意义下)即可。

如果未曾进入if分支,那么tmp的值为y>>1,是某个32位整数右移得来的,故此时tmp最高位必为0;反之,若进入了if分支,其还会异或一个MATRIX_A,而我们知道MATRIX_A的最高位为1,因此这时tmp最高位也一定会变成1。

再考虑到进入if语句的条件是y的最低位为1,因此我们根据tmp的最高位的值,其实已经可以推导出twist函数里的变量y的值了。

接下来的几行代码:

if tmp & self.UPPER_MASK == self.UPPER_MASK:
    tmp ^= self.MATRIX_A
    tmp <<= 1
    tmp |= 1
else:
    tmp <<= 1

此时tmp的值相当于twist函数里的变量y,我们看看它包含了哪些信息:

y = (self._mt[i] & self.UPPER_MASK) | (self._mt[(i + 1) % self.N] & self.LOWER_MASK)

这说明y的最高位是_mt[i]的最高位,y的后31位则是_mt[i+1]的后31位。

因此我们已经恢复出_mt[i]的最高位了,接下来只要恢复其后31位即可。

显然,想要恢复_mt[i]的后31位,只需将前面所有操作的下标减去1,即可在tmp的后31位得到_mt[i]的后31位。这一步非常巧妙。

于是,untwist的完整代码就呼之欲出了:

def untwist(self, i):
    tmp = self._mt[i] ^ self._mt[(i + self.M) % self.N]
    if tmp & self.UPPER_MASK == self.UPPER_MASK:
        tmp ^= self.MATRIX_A
        tmp <<= 1
        tmp |= 1
    else:
        tmp <<= 1
    res = tmp & self.UPPER_MASK
	
    # 进行与前面一模一样的操作,不过将下标减去了1
    tmp = self._mt[i - 1] ^ self._mt[(i + self.M - 1) % self.N]
    if tmp & self.UPPER_MASK == self.UPPER_MASK:
        tmp ^= self.MATRIX_A
        tmp <<= 1
        tmp |= 1
    else:
        tmp <<= 1
    res |= tmp & self.LOWER_MASK
    self._mt[i] = res

逆向 extract_number

def extract_number(self):
    if self._mti >= self.N:
        for i in range(self.N):
            self.twist(i)
        self._mti = 0

    y = self._mt[self._mti]
    y = self.tempering(y)

    self._mti += 1

    return _int32(y)

这就很容易了,可以直接写出:

def unextract_number(self):
    self._mti = (self._mti - 1) % self.N
    if self._mti == 0:
        for i in range(self.N - 1, -1, -1):
            self.untwist(i)
        self._mti = self.N

调用此函数即可将当前的随机数内部状态倒回去一次迭代。


以上,我们已基本实现了对MT19937算法的逆向,以此为基础,我们便有了预测Python随机数的能力,具体内容见下一篇文章