Я портировал Sebastiano Vigna xorshift1024 * PRNG совместим со стандартным форматом С++ 11 случайного числа генератора и заметил некоторое странное поведение с функцией jump(), которую он предоставляет.
Согласно Vigna, вызов jump() должен быть эквивалентен 2 ^ 512 вызовам next(). Поэтому серия вызовов jump() и next() должна быть коммутативной. Например, предположив, что генератор запускается в известном состоянии,
jump();
next();
должен оставить генератор в том же состоянии, что и
next();
jump();
так как оба должны быть эквивалентны
for (bigint i = 0; i < (bigint(1) << 512) + 1; ++i)
next();
Предполагая, что bigint - это целочисленный тип с чрезвычайно большим максимальным значением (и предполагая, что вы очень, очень, очень терпеливый человек).
К сожалению, это не работает с эталонной реализацией Vigna (которая я буду включать в конце для потомков, в случае, если реализация, связанная выше, изменится или будет снята в будущем). При тестировании первых двух параметров с использованием следующего тестового кода:
memset(s, 0xFF, sizeof(s));
p = 0;
// jump() and/or next() calls...
std::cout << p << ';';
for (int i = 0; i < 16; ++i)
std::cout << ' ' << s[i];
вызов jump() перед выходами next():
1; 9726214034378009495 13187905351877324975 10033047168458208082 990371716258730972 965585206446988056 74622805968655940 11468976784638207029 3005795712504439672 6792676950637600526 9275830639065898170 6762742930827334073 16862800599087838815 13481924545051381634 16436948992084179560 6906520316916502096 12790717607058950780
при вызове next() сначала результат:
1; 13187905351877324975 10033047168458208082 990371716258730972 965585206446988056 74622805968655940 11468976784638207029 3005795712504439672 6792676950637600526 9275830639065898170 6762742930827334073 16862800599087838815 13481924545051381634 16436948992084179560 6906520316916502096 12790717607058950780 9726214034378009495
Ясно, что либо мое понимание того, что делает jump(), так и неверно, или ошибка в функции jump(), либо неверные данные полиномиального перехода. Вигна утверждает, что такая функция скачка может быть рассчитана для любого шага периода, но не уточняет, как ее вычислить (в том числе в его статье на xorshift * генераторы). Как я могу вычислить правильные данные скачка, чтобы убедиться, что там нет опечатки?
Xorshift1024 * эталонная реализация; http://xorshift.di.unimi.it/xorshift1024star.c
/* Written in 2014-2015 by Sebastiano Vigna ([email protected])
To the extent possible under law, the author has dedicated all copyright
and related and neighboring rights to this software to the public domain
worldwide. This software is distributed without any warranty.
See <http://creativecommons.org/publicdomain/zero/1.0/>. */
#include <stdint.h>
#include <string.h>
/* This is a fast, top-quality generator. If 1024 bits of state are too
much, try a xorshift128+ generator.
The state must be seeded so that it is not everywhere zero. If you have
a 64-bit seed, we suggest to seed a splitmix64 generator and use its
output to fill s. */
uint64_t s[16];
int p;
uint64_t next(void) {
const uint64_t s0 = s[p];
uint64_t s1 = s[p = (p + 1) & 15];
s1 ^= s1 << 31; // a
s[p] = s1 ^ s0 ^ (s1 >> 11) ^ (s0 >> 30); // b,c
return s[p] * UINT64_C(1181783497276652981);
}
/* This is the jump function for the generator. It is equivalent
to 2^512 calls to next(); it can be used to generate 2^512
non-overlapping subsequences for parallel computations. */
void jump() {
static const uint64_t JUMP[] = { 0x84242f96eca9c41dULL,
0xa3c65b8776f96855ULL, 0x5b34a39f070b5837ULL, 0x4489affce4f31a1eULL,
0x2ffeeb0a48316f40ULL, 0xdc2d9891fe68c022ULL, 0x3659132bb12fea70ULL,
0xaac17d8efa43cab8ULL, 0xc4cb815590989b13ULL, 0x5ee975283d71c93bULL,
0x691548c86c1bd540ULL, 0x7910c41d10a1e6a5ULL, 0x0b5fc64563b3e2a8ULL,
0x047f7684e9fc949dULL, 0xb99181f2d8f685caULL, 0x284600e3f30e38c3ULL
};
uint64_t t[16] = { 0 };
for(int i = 0; i < sizeof JUMP / sizeof *JUMP; i++)
for(int b = 0; b < 64; b++) {
if (JUMP[i] & 1ULL << b)
for(int j = 0; j < 16; j++)
t[j] ^= s[(j + p) & 15];
next();
}
memcpy(s, t, sizeof t);
}
