Генератор псевдослучайных чисел «умножение с переносом»

В информатике . умножение с переносом (MWC) — это метод, изобретенный Джорджем Марсальей ^{[ 1 ]} для генерации последовательностей случайных целых чисел на основе исходного набора от двух до многих тысяч случайно выбранных начальных значений. Основные преимущества метода MWC заключаются в том, что он использует простую компьютерную целочисленную арифметику и приводит к очень быстрой генерации последовательностей случайных чисел с огромными периодами, варьирующимися от примерно $2^{60}$ к $2^{2000000}$ .

Как и во всех генераторах псевдослучайных чисел , результирующие последовательности являются функциями предоставленных начальных значений.

Общая теория

Генератор MWC — это особая форма генератора случайных чисел Лемера. $x_{n}=bx_{n-1}{\bmod {p}}$ что позволяет эффективно реализовать простой модуль $p$ намного больше, чем размер машинного слова.

Обычные реализации генератора Лемера выбирают модуль, близкий к размеру машинного слова. Вместо этого генератор MWC сохраняет свое состояние в базе. $b$ , поэтому умножив на $b$ делается неявно путем сдвига одного слова. База $b$ обычно выбирается равным размеру слова компьютера, поскольку это делает арифметику по модулю $b$ тривиально. Это может варьироваться от $b=2^{8}$ чтобы микроконтроллер $b=2^{64}$ . (В этой статье используется $b=2^{32}$ для примеров.)

Значения начального состояния («начальные») являются произвольными, за исключением того, что они не должны быть полностью равными нулю или максимально допустимым значениям ( $x_{0}=b-1$ и $c_{0}=a-1$ ). (Обычно это делается путем выбора $c_{0}$ между 1 и $a-2$ .). Последовательность MWC тогда представляет собой последовательность пар. $x_{n},c_{n}$ определяется

x_{n}=(ax_{n-1}+c_{n-1})\,{\bmod {\,}}b,\ c_{n}=\left\lfloor {\frac {ax_{n-1}+c_{n-1}}{b}}\right\rfloor

Это называется последовательностью MWC с задержкой 1. Иногда предпочтительнее использовать нечетное основание, и в этом случае $b=2^{k}-1$ можно использовать, что почти так же просто реализовать. Отставание- $r$ последовательность представляет собой обобщение последовательности лаг-1, допускающую более длинные периоды ^{[ 2 ]}. Отставание- $r$ Последовательность MWC тогда представляет собой последовательность пар. $x_{n},c_{n}$ (для $n>r$ ) определяется

x_{n}=(ax_{n-r}+c_{n-1})\,{\bmod {\,}}b,\ c_{n}=\left\lfloor {\frac {ax_{n-r}+c_{n-1}}{b}}\right\rfloor

а выходной сигнал генератора MWC представляет собой последовательность $x$ х,

x_{r},x_{r+1},x_{r+2},...

В этом случае значения начального состояния («начальные») не должны быть ни равными нулю, ни равными нулю. $x_{0}=b-1$ и $c_{r-1}=a-1$ .

Множитель MWC $a$ и отставать $r$ определить модуль $p=ab^{r}-1$ . На практике, $a$ выбирается таким образом, чтобы модуль был простым, а последовательность имела большой период. Если модуль простой, период запаздывания $r$ Генератор MWC - это заказ $b$ в мультипликативной группе чисел по модулю $p$ . Хотя теоретически возможно выбрать непростой модуль, простой модуль исключает возможность того, что исходное начальное число будет иметь общий делитель с модулем, что уменьшит период генератора.

Поскольку 2 — квадратичный вычет чисел вида $8k\pm 1$ , $b=2^{k}$ не может быть примитивным корнем $p=ab^{r}-1$ . Поэтому генераторы MWC с базой $2^{k}$ их параметры выбраны так, чтобы их период был ( ab ^р−1)/2. Это одна из трудностей, связанных с использованием b = 2. ^к − 1 преодолевает.

Базовая форма генератора MWC имеет параметры a , b и r , а также r +1 слово состояния. Состояние состоит из r остатков по модулю b

0 ≤ x ₀ , x ₁ , x ₂ ,…, x _{r −1} < b,

и a перенос c _{r −1} < a .

Хотя теория генераторов MWC допускает a > b , a почти всегда выбирается меньшим для удобства реализации.

Функция преобразования состояния генератора MWC представляет собой один шаг сокращения Монтгомери по модулю p . Состояние представляет собой большое целое число со старшим значащим словом c _{n −1} и наименее значимым словом x _{n − r} . На каждом шаге x _{n − r} ·( ab ^р−1) добавляется к этому целому числу. Это делается в два этапа: -1 · x _{n - r} добавляется к x _{n - r} , в результате чего наименее значимое слово равно нулю. Во-вторых, a · x _{n − r} к переносу добавляется . Это увеличивает целое число на одно слово, создавая два новых наиболее значимых слова x _n и c _n .

До сих пор это просто добавляло к состоянию число, кратное p , в результате чего получался другой представитель того же класса остатков по модулю p . Но, в конце концов, состояние сдвигается на одно слово вниз, деляся на b . При этом отбрасывается наименее значимое слово нуля (которое на практике вообще никогда не вычисляется) и эффективно умножается состояние на b. ⁻¹ (против п ).

Таким образом, генератор умножения с переносом — это генератор Лемера с модулем p и множителем b. ⁻¹ (мод р ). Это то же самое, что генератор с множителем b , но выдающий результат в обратном порядке, что не влияет на качество получаемых псевдослучайных чисел.

Кутюр и Л'Экуйер ^{[ 3 ]} доказали удивительный результат: решетка, связанная с генератором умножения с переносом, очень близка к решетке, связанной с генератором Лемера, который он моделирует. Таким образом, математические методы, разработанные для генераторов Лемера (такие как спектральный тест ), могут быть применены к генераторам умножения с переносом.

Эффективность

Линейный конгруэнтный генератор с базой b = 2 ³² реализуется как

x_{n+1}=(ax_{n}+c)\ {\bmod {\,}}2^{32},

где с — константа. Если a ≡ 1 (по модулю 4) и c нечетно, результирующая база 2 ³² конгруэнтная последовательность будет иметь период 2 ³². ^{[ 4 ]}

Это можно вычислить, используя только младшие 32 бита произведения a и текущего x . Однако многие микропроцессоры могут вычислить полный 64-битный продукт почти за то же время, что и младшие 32 бита. Действительно, многие вычисляют 64-битное произведение, а затем игнорируют старшую половину.

Генератор умножения с переносом с задержкой 1 позволяет нам сделать период почти 2 ⁶³ используя почти те же компьютерные операции, за исключением того, что верхняя половина 64-битного произведения не игнорируется после вычисления произведения. Вместо этого вычисляется 64-битная сумма, а верхняя половина используется в качестве нового значения переноса c, а не фиксированной аддитивной константы стандартной конгруэнтной последовательности: Вычислите ax + c в 64 битах, затем используйте верхнюю половину как новый c , а нижнюю половину как новый x .

Выбор множителя

С указанным множителем , каждая пара входных значений x , c преобразуется в новую пару,

x\leftarrow (ax+c)\,{\bmod {\,}}2^{32},\ \ c\leftarrow \left\lfloor {\frac {ax+c}{2^{32}}}\right\rfloor .

Если x и c не равны нулю, то период результирующей последовательности умножения с переносом будет порядка b = 2. ³² в мультипликативной группе остатков по модулю аб ^р − 1 , то есть наименьшее n такое, что б ^н ≡ 1 (против ab ^р − 1).

Если р = аб ^р − 1 простое число, то маленькая теорема Ферма гарантирует, что порядок любого элемента должен делить p − 1 = ab ^р − 2, поэтому один из способов обеспечить большой порядок — это выбрать такое , что p является « безопасным простым числом », то есть и p , и ( p − 1)/2 = ab ^р/2 − 1 простые. В таком случае при b = 2 ³² и r = 1, период будет ab ^р/2 − 1, приближаясь к 2 ⁶³, что на практике может быть приемлемо большим подмножеством возможных 32-битных пар ( x , c ).

Точнее, в таком случае порядок любого элемента делит p − 1, и возможных делителей всего четыре: 1, 2, ab. ^р/2 − 1, или ab ^р − 2. Первые два применимы только к элементам 1 и −1, а аргументы квадратичной взаимности показывают, что четвертый вариант не может применяться к b , поэтому остается только третий вариант.

Ниже приведены некоторые максимальные значения a для компьютерных приложений, которые удовлетворяют вышеуказанному условию безопасного простого числа для генераторов с задержкой 1:

в Биты	б	Максимум a такой, что ab −1 является безопасным простым числом.	Период
15	2 ¹⁶	2 ¹⁵−50 = 32,718	1,072,103,423
16	2 ¹⁶	2 ¹⁶−352 = 65,184	2,135,949,311
31	2 ³²	2 ³¹−563 = 2,147,483,085	4,611,684,809,394,094,079
32	2 ³²	2 ³²−178 = 4,294,967,118	9,223,371,654,602,686,463
64	2 ⁶⁴	2 ⁶⁴−742 = 18,446,744,073,709,550,874	2 ⁶³(2 ⁶⁴−742)−1 ≈ 1.701 × 10 ^ ³⁸
128	2 ¹²⁸	2 ¹²⁸−10,408	2 ¹²⁷(2 ¹²⁸−10,408)−1 ≈ 5.790 × 10 ^ ⁷⁶
256	2 ²⁵⁶	2 ²⁵⁶−9166	2 ²⁵⁵(2 ²⁵⁶−9166)−1 ≈ 6.704 × 10 ^ ¹⁵³
512	2 ⁵¹²	2 ⁵¹²−150,736	2 ⁵¹¹(2 ⁵¹²−150,736)−1 ≈ 8.988 × 10 ^ ³⁰⁷

В то время как безопасное простое число гарантирует, что почти любой элемент группы имеет большой порядок, период генератора конкретно равен порядку b . Для малых модулей можно использовать более ресурсоемкие методы для нахождения множителей a , где период равен ab /2 - 1. Ниже снова приведены максимальные значения a различных размеров.

в Биты	б ^р	Максимум такое , что b имеет порядок ab ^р/2−1	Период
8	2 ⁸	2 ⁸−7 = 249	31,871
8	(2 ⁸) ² = 2 ¹⁶	2 ⁸−32 = 224	7,340,031
15	2 ¹⁶	2 ¹⁵−29 = 32,739	1,072,791,551
16	2 ¹⁶	2 ¹⁶−22 = 65,514	2,146,762,751
8	(2 ⁸) ⁴ = 2 ³²	2 ⁸−64 = 192	412,316,860,415
15	(2 ¹⁶) ² = 2 ³²	2 ¹⁵−26 = 32,742	70,312,909,602,815
16	(2 ¹⁶) ² = 2 ³²	2 ¹⁶−2 = 65,534	140,733,193,388,031
31	2 ³²	2 ³¹−68 = 2,147,483,580	4,611,685,872,398,499,839
32	2 ³²	2 ³²−76 = 4,294,967,220	9,223,371,873,646,018,559
8	(2 ⁸) ⁸ = 2 ⁶⁴	2 ⁸−41 = 215	2 ⁶³(2 ⁸−41)−1 ≈ 1.983 × 10 ^ ²¹
15	(2 ¹⁶) ⁴ = 2 ⁶⁴	2 ¹⁵−50 = 32,718	2 ⁶³(2 ¹⁵−50)−1 ≈ 3.018 × 10 ^ ²³
16	(2 ¹⁶) ⁴ = 2 ⁶⁴	2 ¹⁶−56 = 65,480	2 ⁶³(2 ¹⁶−56)−1 ≈ 6.039 × 10 ^ ²³
31	(2 ³²) ² = 2 ⁶⁴	2 ³¹−38 = 2,147,483,610	2 ⁶³(2 ³¹−38)−1 ≈ 1.981 × 10 ^ ²⁸
32	(2 ³²) ² = 2 ⁶⁴	2 ³²−43 = 4,294,967,253	2 ⁶³(2 ³²−43)−1 ≈ 3.961 × 10 ^ ²⁸
63	2 ⁶⁴	2 ⁶³−140 = 9,223,372,036,854,775,668	2 ⁶³(2 ⁶³−140)−1 ≈ 8.507 × 10 ^ ³⁷
64	2 ⁶⁴	2 ⁶⁴−116 = 18,446,744,073,709,551,500	2 ⁶³(2 ⁶⁴−116)−1 ≈ 1.701 × 10 ^ ³⁸

Генераторы MWC как повторяющиеся десятичные дроби

Выходной сигнал генератора умножения с переносом эквивалентен по основанию - b разложению дроби со знаменателем p = ab. ^р − 1. Вот пример для простого случая b = 10 и r = 1, поэтому результатом является повторяющаяся десятичная дробь .

Начиная с ${\textstyle c_{0}=1,x_{0}=0}$ , последовательность MWC

x_{n}=(7x_{n-1}+c_{n-1})\,{\bmod {\,}}10,\ c_{n}=\left\lfloor {\frac {7x_{n-1}+c_{n-1}}{10}}\right\rfloor ,

создает следующую последовательность состояний:

10,01,07,49,67,55,40,04,28,58,61,13,22,16,43,25,37,52,19,64,34,31, 10,01,07,...

с периодом 22. Рассмотрим только последовательность x _i :

0,1,7,9,7,5,0,4,8,8,1,3,2,6,3,5,7,2,9,4,4,1, 0,1,7,9,7,5,0,...

Обратите внимание: если эти повторяющиеся сегменты значений x расположить в обратном порядке :

1449275\cdots 9710\,1449275\cdots 9710\,144\cdots

мы получаем расширение j /( ab −1) с a =7, b =10, j=10 :

{\frac {10}{69}}=0.1449275362318840579710\,14492753623\ldots

В целом это верно: последовательность x s, создаваемая генератором MWC с задержкой r :

x_{n}=(ax_{n-r}+c_{n-1}){\bmod {\,}}b\,,\ \ c_{n}=\left\lfloor {\frac {ax_{n-r}+c_{n-1}}{b}}\right\rfloor ,

в обратном порядке будет разложением по основанию b дроби j /( ab ^р − 1) для некоторого 0 < j < ab ^р.

Эквивалентность линейному конгруэнтному генератору

Продолжая приведенный выше пример, если мы начнем с ${\textstyle x_{0}=34}$ и сгенерируем обычную конгруэнтную последовательность

x_{n}=7x_{n-1}\,{\bmod {\,}}69

,

мы получаем последовательность периода 22

31,10,1,7,49,67,55,40,4,28,58,61,13,22,16,43,25,37,52,19,64,34, 31,10,1,7,...

и эта последовательность, уменьшенная по модулю 10, равна

1,0,1,7,9,7,5,0,4,8,8,1,3,2,6,3,5,7,2,9,4,4, 1,0,1,7,9,7,5,0,...

та же самая последовательность x , полученная из последовательности MWC.

В целом это верно (но, по-видимому, только для последовательностей MWC с задержкой 1). ^{[ нужна ссылка ]}):

Учитывая начальные значения $x_{0},c_{0}$ , последовательность $x_{1},x_{2},\ldots$ в результате последовательности MWC задержки-1

x_{n}=(ax_{n-1}+c_{n-1})\,{\bmod {b}}\,,\ \ c_{n}=\left\lfloor {\frac {ax_{n-1}+c_{n-1}}{b}}\right\rfloor

это в точности генератора случайных чисел Лемера выходная последовательность y _n = ay _{n - 1} mod ( ab - 1), уменьшено по модулю b .

Выбор другого начального значения y ₀ просто вращает цикл x' s.

Генераторы дополнительного умножения с переносом

Установление периода лаг- r генератора MWC обычно влечет за собой выбор множителя a так, чтобы p = ab ^р − 1 — простое число. Тогда p − 1 нужно будет разложить на множители, чтобы найти порядок b mod p . Если p — безопасное простое число , то это просто, и порядок b будет либо p − 1, либо ( p − 1)/2. В других случаях p - 1 может быть трудно факторизовать.

Однако алгоритм также допускает отрицательный множитель. Это приводит к небольшой модификации процедуры MWC и дает модуль p = | − аб ^р − 1 | = аб ^р + 1. Это делает p − 1 = ab ^р легко учитывать, что делает практическим определение периода очень больших генераторов.

Модифицированная процедура называется дополнительным умножением с переносом (CMWC), и ее настройка такая же, как и для lag- r MWC: множитель a , основание b и начальные числа r + 1,

Икс ₀ , Икс ₁ , Икс ₂ , ..., Икс _{р -1} и c _{р -1} .

Модификация заключается в создании новой пары ( x , c ). Перестраивая вычисления, чтобы избежать отрицательных чисел, новое значение x дополняется вычитанием его из b - 1:

x_{n}=(b-1)-(ax_{n-r}+c_{n-1})\,{\bmod {\,}}b,\ c_{n}=\left\lfloor {\frac {ax_{n-r}+c_{n-1}}{b}}\right\rfloor .

Результирующая последовательность x , полученная с помощью ГСЧ CMWC, будет иметь период порядка b в мультипликативной группе остатков по модулю ab. ^р+1, а выходные значения в x обратном порядке будут формировать по основанию b разложение j /( ab ^р+1) для некоторого 0 < j < ab ^р.

Использование lag -r CMWC значительно облегчает поиск периодов для r , таких больших, как 512, 1024, 2048 и т. д. (Присвоение r степени двойки немного упрощает доступ к элементам массива, содержащим самые последние x r .)

Еще одним преимуществом этой модифицированной процедуры является то, что период кратен b , поэтому выходные данные точно распределяются по модулю b . ^{[ 3 ]} (Обычный MWC в течение всего своего периода выдает каждый возможный результат одинаковое количество раз, за исключением того, что ноль выдается на один раз меньше, смещение, которым можно пренебречь, если период достаточно длинный.)

Одним из недостатков конструкции CMWC является то, что при базе степени двойки максимально достижимый период меньше, чем для генератора MWC аналогичного размера; вы потеряете несколько битов. Таким образом, генератор MWC обычно предпочтительнее для небольших лагов. Это можно исправить, используя b = 2. ^к−1 или выбрать задержку на одно слово длиннее для компенсации.

Некоторые примеры: При б = 2 ³², и a = 109111 или 108798 или 108517, период отставания - 1024 CMWC

x_{n}=(b-1)-(ax_{n-1024}+c_{n-1})\,{\bmod {\,}}b,\ c_{n}=\left\lfloor {\frac {ax_{n-1024}+c_{n-1}}{b}}\right\rfloor .

будет · 2 ³²⁷⁶² = аб ¹⁰²⁴/64, около 10 ⁹⁸⁶⁷.

При б = 2 ³² и а = 3636507990, р = аб. ¹³⁵⁹ − 1 — безопасное простое число, поэтому последовательность MWC, основанная на этом a, имеет период 3636507990·2. ⁴³⁴⁸⁷ ≈ 10 ¹³¹⁰¹.

При б = 2 ³², ГСЧ CMWC с почти рекордным периодом может быть основан на простом числе p = 15455296 b ⁴²⁶⁵⁸ + 1. Порядок b для этого простого числа равен 241489·2. ^1365056 ≈ 10 ⁴¹⁰⁹²⁸.

Более общие модули

Модуль MWC ab ^р−1 выбран для того, чтобы сделать вычисления особенно простыми, но имеет некоторые недостатки, в частности, то, что период составляет не более половины модуля. Есть несколько способов обобщить это за счет большего количества умножений на итерацию.

Во-первых, к произведению можно добавить дополнительные члены, получив модуль вида a _r b ^р+ а _с б ^с−1. Для этого необходимо вычислить c _n b + x _n = a _r x _{n - r} + a _s x _{n - s} . (Перенос ограничен одним словом, если a _r + a _s ≤ b .)

Однако это не решает проблему периода, которая зависит от младших битов модуля. К счастью, алгоритм приведения Монтгомери допускает другие модули, если они взаимно просты по отношению к основанию b , и это можно применить, чтобы разрешить модуль вида a _r b ^р− a ₀ , для широкого диапазона значений a ₀ . Горески и Клэппер ^{[ 5 ]} разработал теорию этих обобщенных генераторов умножения с переносом, доказав, в частности, что при выборе отрицательных значений a ₀ и a _r – a ₀ < b значение переноса всегда меньше b , что делает реализацию эффективной. Более общая форма модуля также улучшает качество генератора, хотя не всегда можно получить полный период.

Чтобы реализовать генератор Горески-Клэппера, необходимо предварительно вычислить ⁻¹
₀ (mod b ) и меняет итерацию следующим образом: ^{[ 6 ]}

t=c_{n-1}+\sum _{1}^{r}a_{i}x_{n-i},\ x_{n}=a_{0}^{-1}t{\bmod {b}},\ c_{n}=\left\lfloor {\frac {t-a_{0}x_{n}}{b}}\right\rfloor .

В общем случае b = 2 ^к, 0 _{должен быть нечетным ,} чтобы существовало обратное.

Выполнение

Ниже представлена реализация алгоритма CMWC на языке программирования C. Также в программу включен пример функции инициализации. В этой реализации основание равно 2 ³²−1 и лаг r = 4096. Период результирующего генератора составляет около $2^{131104}$ .

// C99 Complementary Multiply With Carry generator
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>

// How many bits in rand()?
// https://stackoverflow.com/a/27593398
#define LOG_1(n) (((n) >= 2) ? 1 : 0)
#define LOG_2(n) (((n) >= 1<<2) ? (2 + LOG_1((n)>>2)) : LOG_1(n))
#define LOG_4(n) (((n) >= 1<<4) ? (4 + LOG_2((n)>>4)) : LOG_2(n))
#define LOG_8(n) (((n) >= 1<<8) ? (8 + LOG_4((n)>>8)) : LOG_4(n))
#define LOG(n)   (((n) >= 1<<16) ? (16 + LOG_8((n)>>16)) : LOG_8(n))
#define BITS_TO_REPRESENT(n) (LOG(n) + !!((n) & ((n) - 1)))
#if ((RAND_MAX | (RAND_MAX >> 1)) != RAND_MAX) 
#error "expected a RAND_MAX that is 2^n - 1!"
#endif
#define RAND_BITS BITS_TO_REPRESENT(RAND_MAX)

// CMWC working parts
#define CMWC_CYCLE 4096         // as Marsaglia recommends
#define CMWC_C_MAX 809430660    // as Marsaglia recommends
struct cmwc_state {
	uint32_t Q[CMWC_CYCLE];
	uint32_t c;	// must be limited with CMWC_C_MAX
	unsigned i;
};

// Collect 32 bits of rand(). You are encouraged to use a better source instead.
uint32_t rand32(void)
{
	uint32_t result = rand();
    for (int bits = RAND_BITS; bits < 32; bits += RAND_BITS)
        result = result << RAND_BITS | rand();
	return result;
}

// Init the state with seed
void initCMWC(struct cmwc_state *state, unsigned int seed)
{
	srand(seed);        
	for (int i = 0; i < CMWC_CYCLE; i++)
		state->Q[i] = rand32();
	do
		state->c = rand32();
	while (state->c >= CMWC_C_MAX);
	state->i = CMWC_CYCLE - 1;
}

// CMWC engine
uint32_t randCMWC(struct cmwc_state *state)  //EDITED parameter *state was missing
{
	uint64_t const a = 18782;	// as Marsaglia recommends
	uint32_t const m = 0xfffffffe;	// as Marsaglia recommends
	uint64_t t;
	uint32_t x;

	state->i = (state->i + 1) & (CMWC_CYCLE - 1);
	t = a * state->Q[state->i] + state->c;
	/* Let c = t / 0xffffffff, x = t mod 0xffffffff */
	state->c = t >> 32;
	x = t + state->c;
	if (x < state->c) {
		x++;
		state->c++;
	}
	return state->Q[state->i] = m - x;
}

int main()
{
	struct cmwc_state cmwc;
	unsigned int seed = time(NULL);

	initCMWC(&cmwc, seed);
	printf("Random CMWC: %u\n", randCMWC(&cmwc));
}

Ниже приведены реализации генераторов MWC с малым состоянием с 64-битным выходом с использованием 128-битных умножений.

// C99 + __uint128_t MWC, 128 bits of state, period approx. 2^127

/* The state must be neither all zero, nor x = 2^64 - 1, c = MWC_A1 -
   1. The condition 0 < c < MWC_A1 - 1 is thus sufficient. */

uint64_t x, c = 1;

#define MWC_A1 0xff3a275c007b8ee6

uint64_t inline next() {
    const __uint128_t t = MWC_A1 * (__uint128_t)x + c;
    c = t >> 64;
    return x = t;
}

// C99 + __uint128_t MWC, 256 bits of state, period approx. 2^255

/* The state must be neither all zero, nor x = y = z = 2^64 - 1, c =
   MWC_A3 - 1. The condition 0 < c < MWC_A3 - 1 is thus sufficient. */

uint64_t x, y, z, c = 1;

#define MWC_A3 0xff377e26f82da74a

uint64_t inline next() {
    const __uint128_t t = MWC_A3 * (__uint128_t)x + c;
    x = y;
    y = z;
    c = t >> 64;
    return z = t;
}

Ниже приведены реализации обобщенных генераторов MWC Горески-Клэппера с малым состоянием с 64-битным выходом с использованием 128-битных умножений.

// C99 + __uint128_t Goresky-Klapper GMWC, 128 bits of state, period approx. 2^127

/* The state must be neither all zero, nor x = 2^64 - 1, c = GMWC_A1 +
   GMWC_MINUS_A0. The condition 0 < c < GMWC_A1 + GMWC_MINUS_A0 is thus
   sufficient. */

uint64_t x = 0, c = 1;

#define GMWC_MINUSA0 0x7d084a4d80885f
#define GMWC_A0INV 0x9b1eea3792a42c61
#define GMWC_A1 0xff002aae7d81a646

uint64_t inline next() {
    const __uint128_t t = GMWC_A1 * (__uint128_t)x + c;
    x = GMWC_A0INV * (uint64_t)t;
    c = (t + GMWC_MINUSA0 * (__uint128_t)x) >> 64;
    return x;
}

// C99 + __uint128_t Goresky-Klapper GMWC, 256 bits of state, period approx. 2^255

/* The state must be neither all zero, nor x = y = z = 2^64 - 1, c =
   GMWC_A3 + GMWC_MINUS_A0. The condition 0 < c < GMWC_A3 + GMWC_MINUS_A0
   is thus sufficient. */

uint64_t x, y, z, c = 1; /* The state can be seeded with any set of values, not all zeroes. */

#define GMWC_MINUSA0 0x54c3da46afb70f
#define GMWC_A0INV 0xbbf397e9a69da811
#define GMWC_A3 0xff963a86efd088a2

uint64_t inline next() {
    const __uint128_t t = GMWC_A3 * (__uint128_t)x + c;
    x = y;
    y = z;
    z = GMWC_A0INV * (uint64_t)t;
    c = (t + GMWC_MINUSA0 * (__uint128_t)z) >> 64;
    return z;
}

Использование

Известно, что из-за простоты и скорости CMWC используется при разработке игр, особенно в современных играх- рогаликах . Он неофициально известен как Мать всех ГПСЧ — имя, первоначально придуманное самим Марсальей. ^{[ 7 ]} В libtcod CMWC4096 заменил MT19937 в качестве PRNG по умолчанию. ^{[ 8 ]}

См. также

Ссылки

^ Марсалья, Джордж ; Заман, Ариф (август 1995 г.). «Компакт-диск со случайными числами Марсальи, включая Непреклонную батарею тестов на случайность» . {{cite journal}}: Для цитирования журнала требуется |journal= ( помощь )
^ Марсалья, Джордж (май 2003 г.). «Генератор случайных чисел» . Журнал современных прикладных статистических методов . 2 (1): 2–13. дои : 10.22237/jmasm/1051747320 .
^ Перейти обратно: ^а ^б Кутюр, Раймонд; Л'Экуйер, Пьер (апрель 1997 г.). «Свойства распределения генераторов случайных чисел умножения с переносом» (PDF) . Математика вычислений . 66 (218): 591–607. Бибкод : 1997MaCom..66..591C . CiteSeerX 10.1.1.154.331 . дои : 10.1090/S0025-5718-97-00827-2 . Мы увидим, что для дополнительного MWC каждый бит выходного значения является справедливым, то есть две двоичные цифры будут появляться одинаково часто в течение всего периода - свойство, не свойственное генераторам MWC.
^ Халл, TE; Добелл, Арканзас (июль 1962 г.). «Генератор случайных чисел» (PDF) . Обзор СИАМ . 4 (3): 230–254. Бибкод : 1962SIAMR...4..230H . дои : 10.1137/1004061 . hdl : 1828/3142 .
^ Перейти обратно: ^а ^б Горески, Марк ; Клэппер, Эндрю (октябрь 2003 г.). «Эффективные генераторы случайных чисел умножения с переносом с максимальным периодом» (PDF) . ACM-транзакции по моделированию и компьютерному моделированию . 13 (4): 310–321. CiteSeerX 10.1.1.4.9190 . дои : 10.1145/945511.945514 . S2CID 13334372 .
^ Обратите внимание, что статья Горески и Клэппера ^{[ 5 ]} содержит ошибку в уравнении (4): последнее равенство неверно — нельзя исключить второе слагаемое из вычисления переноса.
^ «Мать всех генераторов случайных чисел Марсальи» . home.sandiego.edu .
^ «Генератор случайных чисел — RogueBasin» . www.roguebasin.com . Проверено 30 ноября 2016 г. .

Марсалья, Джордж (4 июля 2003 г.). «Xorshift RNG» . Журнал статистического программного обеспечения . 8 (14): 1–6. дои : 10.18637/jss.v008.i14 .
Марсалья, Джордж (октябрь 2005 г.). «О случайности числа Пи и других десятичных разложений» . Интерстат . CiteSeerX 10.1.1.694.4783 .
Пресс, Уильям Х .; Теукольский, Саул А. ; Веттерлинг, Уильям Т.; Фланнери, Брайан П. (2007). «Раздел 7.1.2.B Умножение с переносом (MWC)» . Численные рецепты: искусство научных вычислений (3-е изд.). Нью-Йорк: Издательство Кембриджского университета. ISBN 978-0-521-88068-8 .

[Marsaglia95-1] Марсалья, Джордж ; Заман, Ариф (август 1995 г.). «Компакт-диск со случайными числами Марсальи, включая Непреклонную батарею тестов на случайность» . {{cite journal}}: Для цитирования журнала требуется |journal= ( помощь )

[Marsaglia03-2] Марсалья, Джордж (май 2003 г.). «Генератор случайных чисел» . Журнал современных прикладных статистических методов . 2 (1): 2–13. дои : 10.22237/jmasm/1051747320 .

[Couture97-3] Перейти обратно: ^а ^б Кутюр, Раймонд; Л'Экуйер, Пьер (апрель 1997 г.). «Свойства распределения генераторов случайных чисел умножения с переносом» (PDF) . Математика вычислений . 66 (218): 591–607. Бибкод : 1997MaCom..66..591C . CiteSeerX 10.1.1.154.331 . дои : 10.1090/S0025-5718-97-00827-2 . Мы увидим, что для дополнительного MWC каждый бит выходного значения является справедливым, то есть две двоичные цифры будут появляться одинаково часто в течение всего периода - свойство, не свойственное генераторам MWC.

[4] Халл, TE; Добелл, Арканзас (июль 1962 г.). «Генератор случайных чисел» (PDF) . Обзор СИАМ . 4 (3): 230–254. Бибкод : 1962SIAMR...4..230H . дои : 10.1137/1004061 . hdl : 1828/3142 .

[Goresky03-5] Перейти обратно: ^а ^б Горески, Марк ; Клэппер, Эндрю (октябрь 2003 г.). «Эффективные генераторы случайных чисел умножения с переносом с максимальным периодом» (PDF) . ACM-транзакции по моделированию и компьютерному моделированию . 13 (4): 310–321. CiteSeerX 10.1.1.4.9190 . дои : 10.1145/945511.945514 . S2CID 13334372 .

[6] Обратите внимание, что статья Горески и Клэппера ^{[ 5 ]} содержит ошибку в уравнении (4): последнее равенство неверно — нельзя исключить второе слагаемое из вычисления переноса.

[7] «Мать всех генераторов случайных чисел Марсальи» . home.sandiego.edu .

[8] «Генератор случайных чисел — RogueBasin» . www.roguebasin.com . Проверено 30 ноября 2016 г. .

[ 1 ]

[ 2 ]

[ 3 ]

[ 4 ]

[ 5 ]

[ 6 ]

[ 7 ]

[ 8 ]