Я написал следующую функцию
bool random_bool(double probability)
{
double p_scaled = probability * (RAND_MAX+1) - rand();
if ( p_scaled >= 1 ) return true;
if ( p_scaled <= 0 ) return false;
return random_bool( p_scaled );
}
При условии rand()
генерирует число из равномерного распределения на {0,1,...,RAND_MAX-1,RAND_MAX}
и номера от последующих вызовов могут рассматриваться как независимые для всех практических целей, кроме криптографии, это должно вернуть true
с вероятностью p
: два if
заявления возвращаются true
с вероятностью чуть ниже p
, а также false
с вероятностью чуть выше 1-p
в то время как рекурсивный вызов имеет дело со всем остальным.
Однако следующий тест не пройден:
long long N = 10000000000; //1e10
double p = 10000.0 / N;
int counter = 0;
for (long long i=0;i<N;i++) if (random_bool(p)) counter++;
assert(9672 < counter && counter <= 10330);
Утверждение assert не выполняется только в 0,1% случаев. Однако он терпит неудачу все время (с counter
между 10600 и 10700).
В чем дело?
P.S .: я видел этот вопрос, но это не помогает …
Одним из распространенных дефектов в генераторах случайных чисел является небольшое смещение в сторону меньших результатов (в основном небольшое смещение в сторону 0 в старших разрядах). Это часто случается, когда перенос внутреннего состояния ГСЧ в выходной диапазон выполняется с использованием простого мода, который смещен в сторону высоких значений, если только RAND_MAX не является делителем размера внутреннего состояния. Вот типичная реализация смещения отображения:
static unsigned int state;
int rand() {
state = nextState(); /* this actually moves the state from one random value to the next, eg., using a LCG */
return state % RAND_MAX; /* biased */
}
Смещение происходит из-за того, что более низкие значения выводят и имеют еще одно отображение в режиме mod из состояния. Например, если состояние может иметь значения 0-9 (10 значений), а RAND_MAX равно 3 (поэтому значения 0-2), то % 3
результаты операции в зависимости от состояния
Output State
0 0 3 6 9
1 1 4 7
2 2 5 8
Результат 0 чрезмерно представлен, потому что у него есть шанс 4/10 быть выбранным против 3/10 для других значений.
В качестве примера с более вероятными значениями, если внутреннее состояние RNG является 16-целым, и RAND_MAX
35767 (как вы упомянули, это на вашей платформе), тогда все значения [0,6000] будут выведены для 3 различных значений состояния, но оставшиеся ~ 30000 значений будут выведены только для 2 различных значений состояния — значительный смещение. Этот тип смещения может привести к тому, что значение вашего счетчика будет выше ожидаемого (поскольку меньший, чем равномерный доход от rand () благоприятствует p_scaled >= 1
состояние.
Было бы полезно, если бы вы могли опубликовать точную реализацию rand () на вашей платформе. Если это оказывается смещением в старших битах, вы можете устранить это, передав значения, полученные из rand (), через хорошую хеш-функцию, но лучшим подходом, вероятно, будет просто использование высококачественного источника случайных чисел. числа, например, Мерсенн Твистер
. Лучший генератор также будет иметь больший выходной диапазон (эффективный, более высокий RAND_MAX), что означает, что ваш алгоритм будет подвергаться меньшему количеству повторов / рекурсии.
Даже если реализация среды выполнения Visual Studio страдает от этого дефекта, стоит отметить, что, вероятно, это был хотя бы частично преднамеренный выбор дизайна — использование RAND_MAX, подобного 35767, который относительно прост относительно размера состояния (обычно степени 2), обеспечивает лучшая случайность младших битов, поскольку операция% эффективно смешивает биты старшего и младшего разрядов, а наличие смещенных / неслучайных битов младшего разряда на практике часто является более серьезной проблемой, чем небольшое смещение в битах старшего разряда из-за повсеместности вызывающего rand()
уменьшение диапазона с использованием%, которое эффективно использует только биты младшего разряда для модулей, которые имеют степени 2 (также очень распространенные).
Я попробовал ваш код в Linux, и результаты были довольно приличными. Тем не менее, кажется, что вы находитесь в Windows, где RAND_MAX
около 32768, вероятно. Я говорю, так как GCC жаловался в Linux, что RAND_MAX+1
приводит к целочисленному переполнению, поэтому мне пришлось добавить приведение.
Так что проблема скорее всего в том, что либо RAND_MAX
слишком мала или реализация rand()
в вашей системе это не очень хорошо.
Если источником проблемы является реализация rand()
ваш единственный вариант — перейти на другую функцию из лучшей библиотеки. Однако, если проблема является первой, вы можете решить ее следующим образом.
/* change `rand()` to return two concatenated rands */
typedef long long rand_type; /* this type depends on your actual system, you might get away with `int` */
#define BIGGER_RAND_MAX ((RAND_MAX + 2) * RAND_MAX)
rand_type bigger_rand(void)
{
return (rand_type)rand() * (RAND_MAX + 1) + rand();
}
А затем попробуйте свою программу с этим рандом, который имеет более высокий диапазон. Если проблема не устранена, скорее всего, это ваша rand()
функция, которая далеко не случайна.
Примечание: ваш random_bool
должен вернуться bool
не double
! Так как вы проверяете double
против нуля, это также может быть источником проблемы, где у вас есть ложные срабатывания, потому что двойной может быть не совсем ноль.
я думаю, что результат этой функции связан со значением RAND_MAX, в этом случае p = 1e-6, если RAND_MAX равен 9999, то это никогда не вернет true