通过探索性种子反向工程线性 PRNG

• 下载源代码 - 4.6 KB

引言

线性伪随机数生成器（PRNG）在安全性不是问题的应用中非常受欢迎：它们速度快，而且，更重要的是，它们允许进行一些超越简单生成随机数流的技巧，例如以对数复杂度沿序列的双向跳跃。这在我们进行分布式和/或并行计算时尤其有用。

执行此类跳跃的算法基于通用理论，该理论已在我的先前文章中进行了说明，例如：System.Random and Infinite Monkey Theorem 或 XorShift Jump 101。但是，在应用此通用理论之前，我们必须进行一些准备步骤，这些步骤并不简单，并且取决于特定的PRNG类型。

自然，这会引发一个问题：我们能做得更好吗？幸运的是，答案是肯定的：我们确实可以进一步推广，并使整个过程与PRNG的依赖性更小。

背景

假定读者熟悉 System.Random and Infinite Monkey Theorem 中介绍的线性PRNG理论基础。与那篇文章不同，本文的演示代码是用C++ 11编写的，因此读者应在中级水平上精通此方言。

该代码已故意编写为尽可能独立于库，除了用于单元测试基础架构的 Google Test 1.8.1 之外。

构建演示需要G++及其配套的工具链，例如Cygwin（实际上已使用），但任何现代Linux环境也同样适用。

探索性加种：概念

让我们从理论的简要回顾开始。正如我们在之前的文章中所见，线性PRNG的转移函数 [2] 可以表示为

其中 \(T\) 是某个可以称为转移矩阵的矩阵，而 \(s^{(i+1)}\) 和 \(s^{(i)}\) 都是属于状态集 \(S\) 的 \(N\) 维向量。

这个转移矩阵因PRNG算法的不同而不同，我们必须弄清楚它。虽然对于滞后斐波那契生成器，这个矩阵相当明显，并且可以简单地写出，但对于xorshift，则需要编写一个程序来对指定xorshift族中特定PRNG的公式进行符号变换。

现在考虑一个特殊值 \(s=\left\langle 1, 0, \ldots, 0\right\rangle\)。在这种情况下

对 \(s=\left\langle 0, 1, 0, \ldots, 0\right\rangle\) 执行相同操作，我们可以看到模式

也就是说，如果我们设置PRNG状态为 \(\left\langle 1, 0, \ldots, 0\right\rangle\) 并将转移函数应用于该状态，则后续状态将等于相应转移矩阵的第一列。重复此过程 \(N\) 次，使单位值在状态向量中移动，我们可以收集整个矩阵 \(T\)。实际上，我们计算矩阵乘法 \(T I = T\) 的每一列，从而获得 \(T\) 的值。

这项技术如此明显，我绝不会声称是第一个发明的。它肯定以前在其他地方被使用过，但不幸的是，我找不到任何可以引用的学术论文。

回到技术术语，使用某个状态初始化PRNG称为“加种”，而且，由于我们这样做是为了获取有关PRNG内部的一些知识，因此也可以将其视为一种探索。这解释了将整个过程称为探索性加种的原因。

与我们早期获取转移矩阵的方式相比，这种方法所需的PRNG算法信息要少得多：事实上，这种PRNG仍然必须是线性的，并且能够设置和检索其内部状态向量，但不再需要知道转移函数的精确公式。实际上，我们甚至可以反向工程以“黑箱”形式呈现给我们的某些PRNG。然而，在现实生活中，进行这种反向工程的需求非常罕见：你总能找到一些描述清晰、值得信赖的PRNG算法。因此，探索性加种的主要价值在于它能够以更少的精力获取转移矩阵，而不是解开未知的奥秘。

软件实现

本文附带的源代码旨在提供一个简单的概念证明演示，而不是一个复杂的库。因此，在适当的时候，我将牺牲效率来换取简洁和清晰。

整个代码分为两部分：首先，我们将实现一些流行的PRNG算法，以展示探索性加种方法如何可以推广。第二部分是该方法的核心，将是加种过程的通用实现，适用于所有这些算法。

到目前为止，我们的示例PRNG包括：一种所谓的乘法同余PRNG（MCG），即Park和Miller提出的“最小标准”MCG，并在[4, § 7.1]中进行了描述；一个来自[1, § 3.2.2]的滞后斐波那契生成器家族的代表；最后，128位xorshift [3]，这是我先前文章中的主力。

请注意：这些生成器仅因其简单性和方法的不同而选入本文。它们不代表当前伪随机数生成技术的最高水平，存在已知问题，因此不建议在生产环境中使用。

项目结构非常简单：每种PRNG类型都在单独的头文件中定义，文件名为park_miller_mcg.h，依此类推。头文件tr_utils.h包含探索性加种本身的实现。由于广泛使用C++模板，所有类定义都存在于这些头文件中。最后，项目中唯一的.cc文件explseed.cc包含单元测试，用于检查上述每种PRNG系列的转移矩阵的正确性。

PRNG实现将在某种程度上是不完整的。事实上，由于我们更关注生成器状态的操作，因此我们不必实现输出函数[2]。因此，PRNG类的公共接口仅限于两个函数：init()，它根据PRNG特定的数据参数初始化PRNG状态；以及step_fwd()，它将生成器的内部状态向前推进一步。

然而，这个公共PRNG接口不足以进行探索性加种。事实上，我们需要一种方法来读取PRNG的内部状态并将其设置回。更重要的是，这种方式应该是统一的，以便我们能够开发通用代码。定义两个互补的函数：get_state()和set_state()是自然的。但是，我们如何做到这一点并满足我们的要求呢？一种可能的解决方案是将这些函数直接添加到表示每个PRNG的类中，但这种功能相当特定，并非旨在普通使用。另一种方法允许我们将代码按其受众分离，称为类型特征。

基本思想非常简单：我们将辅助功能移到一个单独的伴随类中。请注意，实现将这些类关联起来的能力很重要，这样它们就可以在模板中一起使用，其中PRNG类型本身被作为模板参数传递。解决此问题的标准方法是模板特化。首先，我们声明一个代表特定类型特征族的模板类

template <class R> class rng_traits;

然后，对于每个PRNG类，我们特化该模板

class xorshift128
{
  ...
};

template<>
class rng_traits<xorshift128>
{
  ...
};

最终的特化集稍后可以与固定的PRNG类型一起使用

void do_something()
{
  using traits = rng_traits<xorshift128>;

  xorshift128 rng;

  auto s = traits::get_state(rng);

  ...
}

甚至可以与当前未知的类型一起使用

template <class R>
void do_something()
{
  using traits = rng_traits<R>;

  R rng;

  auto s = traits::get_state(rng);

  ...
}

请注意，使用类型特征，我们可以多次以不同的方式扩展某个类的功能，并且可以非侵入性地做到这一点，只要我们可以通过其公共接口操作该类的对象。然而，在本例中，我决定允许类型特征访问PRNG类的内部，从而减小了它们的公共接口。此外，这种方法非常类似于Haskell中可用的类型类的概念。

本文实现的演示代码将PRNG类及其类型特征都放在同一个头文件中，但在实际应用中，可能最好使用单独的文件和命名空间约定来将这些内部信息与客户端代码隔离。

实现探索性加种的代码对PRNG类型特征做出以下假设

• PRNG类型特征定义了一个类似向量的容器类型state_t，它表示PRNG状态。该容器的大小是固定的，并且可以通过包含单个项目N的匿名枚举在编译时获得。

• 状态向量元素的类型可以是任意的，但行为类似于数字（即，整数常量0和1被隐式转换为该类型的相应元素，乘法和就地加法运算以自然的方式工作）。

  我们的演示使用了两种类型的状态向量：标准的uint32_t和用户定义的int_mod<M>，实现了整数模\(M\)有限环上的所需操作子集。

• 状态向量可以从PRNG实例使用get_state()读取，并使用set_state()写入。

• 不对state_t实例的初始值做任何假设。因此，需要一个附加函数zero_state()来创建在探索性加种过程中使用的零状态向量。

在完成准备工作后，我们就可以实现旅程的目标：计算任意线性PRNG转移矩阵的通用函数makeT()。请注意，我们的实现是完全抽象的，可以用于任何大小的PRNG状态，甚至任何类型的状态项。这种级别的灵活性可以通过将所有PRNG特定实现细节转移到相应的类型特征类来实现。

template<class R>
class tr_utils
{
  using traits = rng_traits<R>;
  
  enum { N = traits::N };
  
public:
  using state_t = typename traits::state_t ;
  using T_t = std::array<std::array<typename state_t::value_type, N>, N> ;
  
  static T_t make_T()
  {
    T_t T;

    R rng;

    for(size_t i = 0; i < N; i++)
    {
      state_t st = traits::zero_state();

      st[i] = 1;

      traits::set_state(rng, st);

      rng.step_fwd();

      st = traits::get_state(rng);

      for(size_t j = 0; j < N; j++)
        T[j][i] = st[j];
    }

    return T;
  }
  ...
};

请注意，makeT()已成为tr_utils模板类的static成员。该类还实现了其他一些联合使用的函数：matmul()，它简单地将转移矩阵与状态向量相乘，从而将PRNG向前推进一步；以及print_T()，它输出转移矩阵，可用于调试。将所有这些函数放入同一个外部类允许我们引入通用的类型别名，使代码更清晰。

最后，我们来考虑explseed.cc，这是演示代码中唯一的.cc文件。它包含三个单元测试，可以用来检查一切是否正常工作，并且另外提供了一个使用上述类和函数的示例。

TEST(ExplSeed, XorShift128)
{
  using traits = rng_traits<xorshift128>;
  using utils = tr_utils<xorshift128>;
  
  xorshift128 rng;

  rng.init({
    123456789,
    362436069,
    521288629,
    88675123
  });
  
  auto T = utils::make_T();
  auto s = traits::get_state(rng);

  for(int i = 0 ; i < num_steps; i++)
  {
    rng.step_fwd();

    s = utils::matmul(T, s);

    ASSERT_EQ(s, traits::get_state(rng))
       << "Failed at step " << i;
  }
}

除了使用的PRNG类型外，所有测试的结构都完全相同。我们能否添加一个模板函数来进一步减小源代码大小？我不会这样做，因为在这种情况下，错误消息中显示的所有行号对于所有测试都将相同，并且将指向该辅助函数而不是某个特定的失败测试。

测试过程本身非常简单：我们创建一个PRNG类型的实例，复制其状态，然后使用为该PRNG类型计算的转移矩阵及其成员函数step_fwd()执行num_steps次状态更新。在每一步之后，提取生成器的内部状态，并将其与计算出的状态向量进行比较，如果不匹配则报告错误。

结论

在本文中，我们讨论了一种创建某些线性PRNG转移矩阵的又一种方法。与我在先前文章中考虑过的其他方法不同，这种方法更简单，并且在某种程度上可以以“黑箱”方式使用，即，无需完全了解PRNG的内部。然而，仍然需要一些信息，例如状态向量的维度及其组件的性质（它们是实数还是整数，整数算术使用什么模数，等等）。

本文的实现部分也可以用作C++编程技术（类型特征）的简单用例，这种技术更适合学术目的，而不是引用生产级库（如Boost）的源代码。

参考文献

• [1] Donald E. Knuth. The Art of Computer Programming, Volume 2 (3rd Ed.): Seminumerical Algorithms. Addison-Wesley Longman Publishing Co., Inc., Boston, MA, USA, 1997.
• [2] P. L’Ecuyer. Uniform random number generation. Annals of Operations Research, 53:77–120, 1994.
• [3] George Marsaglia. Xorshift RNGs. Journal of Statistical Software, Articles, 8(14):1–6, 2003. (Available online.)
• [4] William H. Press, Saul A. Teukolsky, William T. Vetterling, and Brian P. Flannery. Numerical Recipes in C (2nd Ed.): The Art of Scientific Computing. Cambridge University Press, USA, 1992.