c - 避免蒙特卡罗模拟中的基本 rand() 偏差？

Question

我正在用 Objective C 重写蒙特卡罗模拟，以便在 VBA/Excel 的 dll 中使用。计算中的“引擎”是创建 0 到 10001 之间的随机数，并将其与 5000-7000 附近的变量进行比较。每次迭代使用 4-800 次，我使用了 100000 次迭代。因此，每次运行大约生成 50,000,000 代随机数。

虽然在 Objective C 中测试没有显示偏差，但我对 C 代码有很大的问题。 Objective C 是 C 的超集，所以 95% 的代码都是复制粘贴的，很难搞砸。昨天和今天一整天休息了很多次，没有发现任何问题。

使用 srand() 时，我留下了 arc4random_uniform() 和 rand() 之间的差异，特别是因为对 0 到 10000 的较低数字的偏差。我进行的测试与这种偏差一致.5 到 2 % 对于低于大约 5000 的数字。任何其他解释是如果我的代码避免重复，我猜它不会这样做。

代码非常简单(“spiller1evne”和“spiller2evne”是5500到6500之间的数字):

srand((unsigned)time(NULL));
for (j=0;j<antala;++j){
[..]
        for (i=1;i<450;i++){
            chance = (rand() % 10001);

[..]

             if (grey==1) {


                 if (chance < spiller1evnea) vinder = 1;
                 else vinder = 2;
            }
            else{
                if (chance < spiller2evnea) vinder = 2;
                else vinder = 1;
            }

现在我不需要真正的随机性，伪随机性就很好了。我只需要它在累积的基础上近似均匀分布(就像如果 5555 的可能性是 5556 的两倍，这并不重要。如果 5500-5599 的可能性比 5600-5699 的可能性高 5%，那么这确实很重要，并且如果 0-4000 比 6000-9999 存在明显 0.5-2% 的偏差。

首先，rand() 是我的问题听起来合理吗？是否有一个简单的实现可以满足我的低需求？

编辑:如果我的怀疑是合理的，我可以对此使用任何:

http://www.azillionmonkeys.com/qed/random.html

我可以将其复制粘贴作为替换吗(我正在用 C 语言编写并使用 Visual Studio，真的是新手)？:

#include <stdlib.h>

#define RS_SCALE (1.0 / (1.0 + RAND_MAX))

double drand (void) {
    double d;
    do {
       d = (((rand () * RS_SCALE) + rand ()) * RS_SCALE + rand ()) * RS_SCALE;
    } while (d >= 1); /* Round off */
    return d;
}

#define irand(x) ((unsigned int) ((x) * drand ()))

EDIT2:很明显，上面的代码在没有相同偏见的情况下工作，所以我建议任何有与我上面描述的相同“中间道路”需求的人。它确实会带来惩罚，因为它调用了 rand() 三次。所以我仍在寻找更快的解决方案。

Answer 1

rand() 函数生成 [0, RAND_MAX] 范围内的 int。如果您像原始代码那样通过模运算符 (%) 将其转换为不同的范围，那么就会引入不均匀性，除非目标范围的大小碰巧能均匀划分 RAND_MAX + 1。这听起来就像你所看到的那样。

您有多种选择，但如果您想坚持使用基于 rand() 的东西，那么我建议您对原始方法进行这种变体:

/*
 * Returns a pseudo-random int selected from the uniform distribution
 * over the half-open interval [0, limit), provided that limit does not
 * exceed RAND_MAX.
 */
int range_rand(int limit) {
    int rand_bound = (RAND_MAX / limit) * limit;
    int r;
    while ((r = rand()) >= rand_bound) { /* empty */ }
    return r % limit;
}

虽然原则上每次调用该函数将生成的 rand() 调用次数是无限的，但实际上，对于相对较小的 限制，平均调用次数仅略大于 1 值，并且每个 limit 值的平均值小于 2。它通过从 [0, RAND_MAX] 的子集中选择初始随机数来消除我之前描述的不均匀性，该子集的大小除以 limit。