1

我正在尝试为范围最大为 4 位的 C 应用程序的 TRNG 输出文件实现范围映射器。由于鸽子偏差问题，我决定使用丢弃算法。

我对简约算法的想法是这样的：

-- 从文件中读取 16 个字节并存储为索引的 128 位无符号整数位桶，以便一次选择 n 位作为位掩码。
-- 尽可能多地预先确定每个输入所需的范围/存储桶并存储在一个数组中。
-- 对于 bitbucket 中的每 n 位，从数组中选择一个输入，如果存在则不会丢弃它。如果 2 位找不到输入，请尝试 3 位，如果找不到输入，请尝试 4 位。起初，当有很多输入时，不丢弃应该很容易，但随着输入的选择变得低丢弃将变得更加普遍。我不完全确定从更少的位开始并按我的方式工作是否更好，或者相反。

这个位啜饮范围映射器的缺点似乎是我需要假设随机输入数据的数量大约是偏置缩放方法所需的两倍。例如，来自 4 位 rand 输出的 9 桶输入将丢失大约 43% 的时间。

现有的实现/算法： 这似乎是一个更复杂、更有效的简约范围映射方法的例子，但我发现他的解释完全难以理解。任何人都可以用英语向我解释或建议我可能读的一本书或我可能参加的大学课程，这会给我一个理解它的背景吗？

还有arc4random似乎是运行时优化的无偏模丢弃实现。像几乎所有无偏范围映射器实现一样，我发现这似乎并不特别关心它使用了多少数据。然而，这并不意味着它必然会降低数据效率，因为它具有更少未命中的优势。

arc4random 的基本思想似乎是，只要鸽子的数量（max_randvalue_output）可以被孔的数量（rangeupperbound）整除，模函数本身就是一个优雅且无偏的范围映射器。然而，模数似乎仅在您不啜饮时才相关​​，即当随机源的输出超过 ceil(log2(buckets)) 位时。

在“浪费”的随机位的数量和丢弃的百分比之间似乎存在权衡。未命中的百分比与范围映射器输入中的多余位数成反比。似乎应该有一种数学方法来比较一个有点小范围映射器的数据效率和一个更饥饿的版本和更少的失误，但我不知道。

所以我的计划是只写两个实现：有点吝啬类型的范围映射器，可能有点像也可能不像数学论坛的例子（我不明白）和一个接受字节输入的不变字节输入模范围映射器来自 TRNG 并使用从最大的顶部丢弃的模数除偏方法将 (x)n 只鸽子与 n 孔匹配，这类似于 arc4random。完成后，我计划将它们发布在 codereview 上。

我基本上是在寻找任何这些问题的帮助或建议，这些问题可能会帮助我编写一个更简洁但仍然不偏不倚的范围映射器，特别是在我的简约算法方面。运行时效率不是优先事项。

4

2 回答 2

2

有一种更简单的方法可以从随机比特流中生成一定范围内的随机数，这种方法不仅效率最佳，而且精确。它被称为 J. Lumbroso 的“快速掷骰子”方法：

“硬币翻转的最佳离散均匀生成和应用”，2013 年。

另请参阅此问题。

于 2020-03-22T10:03:36.983 回答
2

我查看了@ Peter.O指向的“快速骰子滚轮”（FDR），这确实很简单（并且避免了划分）。但是每次生成一个随机数时，这都会消耗一些位并丢弃它不使用的那些位的一部分。

“批处理”/“池化”技术似乎比 FDR 做得更好，因为（至少部分）保留了未使用的比特部分。

但有趣的是，您引用的DrMath与 FDR 基本相同，但它返回的每个随机值都不是从头开始的。

所以返回0..n-1的 FDR是：

  random(n):
    m = 1 ; r = 0 
    while 1:
        # Have r random and evenly distributed in 0..m-1
        # Need m >= n -- can double m and double r adding random bit until
        #                we get that.  r remains evenly distributed in 0..m-1 
        while m < n: r = 2*r + next_bit() ; m = m*2
        # Now have r < m and n <= m < n*2
        if r < n: return r   # Hurrah !
        # Have overshot, so reduce m and r to m MOD n and r MOD m
        m -= n ; r -= n ;

DrMath 是这样的：

  # Initialisation once before first call of random(m)
  ms = 1 ; rs = 0
  N = ... # N >= maximum n and N*2 does not overflow 

  # The function -- using the "static"/"global" ms, rs and N 
  random(n):
    m = ms ; r = rs
    while 1:
        # Same as FDR -- except work up to N not n
        while m < N: r = 2*r + next_bit() ; m = m*2 ;
        # Now have r < m and m >= N
        # Set nq = largest multiple of n <= m
        # In FDR, at this point q = 1 and nq = n
        q  = m DIV n ;
        nq = n * q
        if r < nq:             # all set if r < nq
            # in FDR ms = 1, rs = 0 
            ms = q             # keep stuff not used this time
            rs = r DIV n       # ditto
            return r MOD n     # hurrah !
        # Overshot, so reduce MOD n*q -- remembering, for FDR q == 1
        m = m - nq 
        r = r - nq

如前所述，它与 FDR 基本相同，但会跟踪未使用的随机性。

测试时我发现：

  FDR:    for 100000 values range=3 used 266804 bits cost=1.6833
  DrMath: for 100000 values range=3 used 158526 bits cost=1.0002

其中注意到 log2(3) = (1.58496) cost。bits-used / (100000 * log2(3))（所以cost是使用的位数除以希望使用的位数）。

还：

  FDR:    for 100000 values range=17: 576579 bits cost=1.4106
  DrMath: for 100000 values range=17: 408774 bits cost=1.0001

和：

  FDR:    for 100000 values ranges=5..60: 578397 bits cost=1.2102
  DrMath: for 100000 values ranges=5..60: 477953 bits cost=1.0001

其中构造了 100000 个值，并为每个值选择了一个范围5..60（含）。

在我看来 DrMath 有它！尽管对于更大的范围，它的优势较小。

请注意... DrMath 每个返回的随机值至少使用 2 个除法，这给了我运行时间方面的建议。但是您确实说过您对运行时效率不感兴趣。

它是如何工作的 ？

因此，我们希望一系列随机值r均匀分布在一个范围内0..n-1。不方便的是，我们只有一个随机源，它为我们提供了均匀分布在 中的随机值0..m-1。通常m是 2 的幂——让我们假设n < m（如果n == m问题是微不足道的，如果n > m问题是不可能的）。对于任何r我们可以r MOD n在所需范围内给出随机值的值。如果我们只使用rwhen r < nthen（很简单），我们就有了我们想要的均匀分布。如果我们只使用rwhenr < (n * q)并且(n * q) < m我们也有一个均匀分布。我们在这里“拒绝” r“太大”的东西。越少r我们拒绝，更好。所以我们应该选择q这样(n * q) <= m < (n * (q-1))——所以n * q是n小于或等于的最大倍数m。反过来，这告诉我们n“少得多”比m是首选。

当我们“拒绝”一个给定的r东西时，我们可以把它全部扔掉，但事实证明这并不是完全必要的。此外，m不一定是 2 的幂。但我们稍后会讨论。

这是一些工作的Python：

M = 1
R = 0
N = (2**63)    # N >= maximum range

REJECT_COUNT = 0

def random_drmath(n):
    global M, R, REJECT_COUNT

    # (1) load m and r "pool"
    m = M
    r = R
    while 1:
        # (2) want N <= m < N*2
        #     have 0 <= r < m, and that remains true.
        #     also r uniformly distributed in 0..m-1, and that remains true
        while m < N:
            r = 2*r + next_bit()
            m = m*2

        # (3) need r < nq where nq = largest multiple of n <= m
        q  = m // n
        nq = n * q
        if r < nq:
            # (4) update the m and r "pool" and return 0..n-1 
            M = q
            R = r // n
            return r % n       # hurrah !

        # (5) reject: so reduce both m and r by MOD n*q
        m = m - nq 
        r = r - nq
        REJECT_COUNT += 1

必须有N>= 最大范围，最好更大。 2**31或者2**63是显而易见的选择。

在第 (2) 步的第一次调用中，random_drmath()将读取随机位以“填充池”。对于,N = 2**63将以63 个随机位结束。显然是随机且均匀分布在 中的。[到现在为止还挺好。]m = 2**63rr0..m-1

现在（以及对 的所有进一步调用）我们希望从中random_drmath()均匀地提取随机值，如上所述。所以——步骤（3）——构造小于或等于的最大倍数。如果我们不能使用它，因为其中的值少于-- 这是通常的“拒绝”标准。0..n-1rnqnmr >= nqnnq..m-1

那么，哪里r < nq可以返回一个值——步骤（4）。这里的诀窍是将m和r视为数字“base-n”。的 ls “数字”r被提取 ( r % n) 并返回。然后m和r右移一位“数字”（q = m // n和r // n），并存储在“池”中。我认为很明显，此时r和m仍然是随机r < m且r均匀分布在0..m-1. 但是m不再是 2 的幂——但这没关系。

但是，如果r >= nq必须减少r并m一起——步骤（5）——再试一次。琐碎，可以设置m = 1 ; r = 0并重新开始。但是我们所做的是nq从两者中减去m和 ，r 这使得r均匀分布在 中0..m-1。这最后一步感觉很神奇，但我们知道rinnq..m-1并且每个可能的值具有相等的概率，所以r-nqis 在范围内0..m-nq-1并且每个可能的值仍然具有相等的概率！while[请记住，循环顶部的“不变量”r是随机且均匀分布在 中的0..m-1。]

对于小n的拒绝步骤将丢弃大部分r，但对于小的n（与 相比N）我们希望不要经常拒绝。相反，对于大n（与 相比N），我们可能希望更频繁地拒绝，但这至少保留了我们迄今为止吃过的一些随机位。我觉得可能有一种方法可以保留更多r……但是还没有想到一种简单的方法来做到这一点……如果读取一个随机位的成本很高，那么可能值得尝试找到一个不简单的方法！

FWIW：设置N = 128我得到：

  FDR:    for 100000 values ranges=3.. 15: 389026 bits cost=1.2881
  DrMath: for 100000 values ranges=3.. 15: 315815 bits cost=1.0457

  FDR:    for 100000 values ranges 3.. 31: 476428 bits cost=1.2371
  DrMath: for 100000 values ranges 3.. 31: 410195 bits cost=1.0651

  FDR:    for 100000 values ranges 3.. 63: 568687 bits cost=1.2003
  DrMath: for 100000 values ranges 3.. 63: 517674 bits cost=1.0927

  FDR:    for 100000 values ranges 3..127: 664333 bits cost=1.1727
  DrMath: for 100000 values ranges 3..127: 639269 bits cost=1.1284

所以随着n接近N，单位价值的成本上升。

于 2020-03-22T18:49:46.593 回答