Поведение сброса для форсированных генераторов случайных чисел

Я обертываю импульсные генераторы случайных чисел классом адаптера для реализации процедуры Монте-Карло. При написании модульных тестов для функций-членов класса я предполагал, что поведение .discard (unsigned int N) состоит в том, чтобы рисовать N случайных чисел без их сохранения, тем самым улучшая состояние rng. Код повышения:

void discard(boost::uintmax_t z)
{
if(z > BOOST_RANDOM_MERSENNE_TWISTER_DISCARD_THRESHOLD) {
discard_many(z);
} else {
for(boost::uintmax_t j = 0; j < z; ++j) {
(*this)();
}
}
}

что подтверждает мое предположение. Тем не менее, я обнаружил, что последовательность, полученная из .discard (1), не отличается от числа на одну и ту же последовательность без сброса. Код:

#include <iostream>
#include <iomanip>
#include <random>
#include <boost/random.hpp>

int main()
{
boost::mt19937 uGenOne(1);
boost::variate_generator<boost::mt19937&, boost::normal_distribution<> > distOne(uGenOne, boost::normal_distribution<>());

boost::mt19937 uGenTwo(1);
boost::variate_generator<boost::mt19937&, boost::normal_distribution<> > distTwo(uGenTwo, boost::normal_distribution<>());
distTwo.engine().discard(1);

unsigned int M = 10;
std::vector<double> variatesOne(M);
std::vector<double> variatesTwo(M);

for (unsigned int m = 0; m < M; ++m) {
variatesOne[m] = distOne();
variatesTwo[m] = distTwo();
}

for (unsigned int m = 0; m < M; ++m)
std::cout << std::left << std::setw(15) << variatesOne[m] << variatesTwo[m] << std::endl;

return 0;
}

выходы

2.28493        0.538758
-0.668627      -0.0017866
0.00680682     0.619191
0.26211        0.26211
-0.806832      -0.806832
0.751338       0.751338
1.50612        1.50612
-0.0631903     -0.0631903
0.785654       0.785654
-0.923125      -0.923125

Моя интерпретация того, как .discard работает неправильно? Почему две последовательности отличаются в первых трех выходах, а затем идентичны?

(Этот код был скомпилирован в msvc 19.00.23918 и g ++ 4.9.2 в cygwin с одинаковыми результатами).

0

Решение

Проблема здесь заключается в том, что движок не изменяется должным образом или дистрибутив добавляет дополнительную работу. Если мы используем двигатель напрямую, как

int main()
{
boost::mt19937 uGenOne(1);

boost::mt19937 uGenTwo(1);
uGenTwo.discard(1);

unsigned int M = 10;
std::vector<double> variatesOne(M);
std::vector<double> variatesTwo(M);

for (unsigned int m = 0; m < M; ++m) {
variatesOne[m] = uGenOne();
variatesTwo[m] = uGenTwo();
}

for (unsigned int m = 0; m < M; ++m)
std::cout << std::left << std::setw(15) << variatesOne[m] << variatesTwo[m] << std::endl;

return 0;
}

Производит

1.7911e+09     4.28288e+09
4.28288e+09    3.09377e+09
3.09377e+09    4.0053e+09
4.0053e+09     491263
491263         5.5029e+08
5.5029e+08     1.29851e+09
1.29851e+09    4.29085e+09
4.29085e+09    6.30312e+08
6.30312e+08    1.01399e+09
1.01399e+09    3.96591e+08

Это последовательность с 1 сдвигом, как мы и ожидали, поскольку мы отбросили первый вывод.

Таким образом, вы правы в том, как сбрасывать работает. Я не совсем уверен, почему есть расхождение, когда вы делаете это через boost::variate_generator, Я не понимаю, почему первые три числа отличаются, но все остальные выходные данные совпадают.

0

Другие решения

Просто чтобы добавить важную деталь, которая была в комментариях к предыдущему ответу. Как упоминалось в @NathanOliver, .discard увеличивает генератор, который отправляет униформу в нормальное распределение, которое преобразует униформу в нормальное распределение. boost :: normal_distribution использует Алгоритм зиккурата который является алгоритмом типа «принятие / отклонение». Он рисует случайную форму, манипулирует ею, а затем проверяет, находится ли она в нужном распределении. Если нет, то отклоняется и новая случайная форма.

for(;;) {
std::pair<RealType, int> vals = generate_int_float_pair<RealType, 8>(eng);
int i = vals.second;
int sign = (i & 1) * 2 - 1;
i = i >> 1;
RealType x = vals.first * RealType(table_x[i]);
if(x < table_x[i + 1]) return x * sign;
if(i == 0) return generate_tail(eng) * sign;
RealType y = RealType(table_y[i]) + uniform_01<RealType>()(eng) * RealType(table_y[i + 1] - table_y[i]);
if (y < f(x)) return x * sign;
}

Важным моментом является то, что если последний if терпит неудачу тогда for цикл запускается снова и вызов generate_int_float_pair будет запущен снова. Это означает, что количество приращений базового генератора неизвестно.

Следовательно, нормальные последовательности будут иметь разные номера до тех пор, пока сумма, в которой сумма отклонений в подпоследовательности не будет одинаковой, в этот момент оставшаяся однородная последовательность будет идентичной. Это произошло на третьей позиции в примере, размещенном в вопросе. (На самом деле это немного более тонко, потому что базовый генератор может быть вызван один или два раза в алгоритме Зиккурата, но суть одна и та же — после того, как последовательности синхронизируются, они никогда не смогут производить разные переменные).

0

По вопросам рекламы ammmcru@yandex.ru
Adblock
detector