计算科学 - 生成边界外的随机数： - 吾爱随笔录

生成边界外的随机数：

计算科学算法随机数生成

2021-12-05 03:18:58

假设我有一个随机数生成器，它在[0, RAND_MAX]和内生成一个数字RAND_MAX < UINT_MAX。

如何在和内生成一个随机数，[0, i]同时保持均匀分布，并且不超过任何计算？i>RAND_MAXi<UINT_MAXUINT_MAX

第一次尝试是将范围分成n相等的部分，以便i mod n == 0通过找到最小分母i（使用整数分解筛算法之一）。每个子列表具有相等的分布，可以随机选择任何一个。如果i是素数，则失败：

0 1 2 3 4 5 6
|___| |_____|

第二种尝试是将范围分成两部分a = [0, floor(i/2)]和b = [floor(i/2)+1, i]，然后将分布计算为两个子范围的长度之比：

g = gcd(length(a), length(b))
if random(0, length(a)/g + length(b)/g - 1) >= length(a)/g
    random(0, floor(i/2))
else
    random(floor(i/2)+1, i)

给定以下范围，选择右侧子范围的频率是左侧子范围的两倍：

0 1 2 3 4 5
|_| |_____|

不幸的是，如果第一次调用random()可能会永远递归gcd()。

最后的尝试是纠正通过拆分奇数长度范围引入的分布偏斜：

if random(1)
    if random(1)
        random(0, floor(i/2))
    else
        random(floor(i/2) + 1, i)
else
    if random(1)
        random(0, floor(i/2) - 1)
    else
        random(floor(i/2), i)

这可能行不通。

编辑：

对于任何有类似问题的人，两个统一变量的总和遵循三角分布（请参阅Irwin-Hall 分布）。这是无效的：

# triangular distribution:
random(a, floor(b/2)) + random(a, b - floor(b/2))

编辑（接受）：

我考虑的两个答案都是高质量的。此答案以每位或每位为基础考虑 RNG。接受的答案以数字方式考虑 RNG，而不考虑数字或位。在这两种情况下，阅读评论都很重要。

我还发布了一个后续问题维护跨子范围的均匀分布。

2个回答

为了您的示例，假设 RAND_MAX=12 和 i=17。然后执行以下程序：选择两个随机数 $r_1,r_2$ 并将它们组合成一个均匀分布的随机数 $[0,144)$ 通过计算 $r=r_1*12+r_2$ . 这当然是错误的间隔。你得到一个均匀分布的随机数 $[0,17)$ 通过重复这个过程直到你得到一个数字 $r$ 在这个区间；换句话说，每次你得到一个随机数 $r\ge i$ ，你丢弃它，然后再试一次。这保证了你最终得到一个统一的随机数。

与其修复提供的生成器，不如使用合理的现代选择。它将更快并且具有更好的统计质量。可能的示例包括 xorshift+、xorshift* 和 PCG 的各种变体。

直接回答提出的问题。您可以生成一个样本，屏蔽掉两个可用位的最大功率，生成另一个并移动那个位或前一个。对于两个“n”的非幂执行拒绝方法。

编辑 2：因此使用 64 位 xorshift+ 生成 64 位统一序列可能如下所示：

// multiple 'state' blocks to allow for friendly multi-threading
typedef struct {
  uint64_t s0;
  uint64_t s1;
} rng_state;

// get 64-bits
inline uint64_t rng_next(rng_state_t* s)
{
  uint64_t s1 = s->s0;
  uint64_t s0 = s->s1;

  s->s0 = s0;
  s1   ^= s1 << 23;
  s->s1 = s1 ^ s0 ^ (s1 >> 18) ^ (s0 >> 5);

  return s->s1 + s0; 
}

尝试“修补”随机，假设简化假设它返回一个统一的 N 位数字并将其扩展到 2N 位：

inline uint64_t rng_next() {
   uint64_t r0 = random();
   return (r0 << RANDOM_BITS) | random();
 }

快速浏览一下源代码，看起来（比如说）glibc 的当前版本将返回 31 位，并且默认情况下它使用二次幂 LCG（质量非常低）。修补的最终结果是质量显着降低，运行时成本显着增加。

注意：所有代码都是在帖子中输入的，甚至可能无法编译。

其它你可能感兴趣的问题

上一篇特征值数字的已知问题？下一篇在无通量边界条件的半无限域上数值求解耦合偏微分方程