Мерсенн Твистер

Mersenne Twister общего назначения — генератор псевдослучайных чисел (ГПСЧ), разработанный в 1997 году Макото Мацумото ( Макото Мацумото ) и Такудзи Нисимура ( Takuji Nishimura ) . Его название происходит от выбора простого числа Мерсенна в качестве длины периода.

Mersenne Twister был разработан специально для исправления большинства недостатков, обнаруженных в старых ГПСЧ.

Наиболее часто используемая версия алгоритма Твистера Мерсенна основана на простом числе Мерсенна. $2^{19937}-1$. Стандартная реализация этого, MT19937, использует длину слова 32 бита . Есть еще одна реализация (пять вариантов). ), использующий 64-битную длину слова, MT19937-64; он генерирует другую последовательность.

k -распределение

Псевдослучайная последовательность $x_{i}$ - битных W целых чисел периода P называется k-распределенным с точностью до v -бита, если выполняется следующее.

Пусть trunc _v ( x ) обозначает число, образованное старшими v битами x , и рассмотрим P k v . -битных векторов

$(\operatorname {trunc} _{v}(x_{i}),\operatorname {trunc} _{v}(x_{i+1}),\,\ldots ,\operatorname {trunc} _{v}(x_{i+k-1}))\quad (0\leq i<P)$.

Затем каждый из $2^{kv}$ возможные комбинации битов встречаются одинаковое количество раз за период, за исключением комбинации всех нулей, которая встречается реже.

Алгоритмическая деталь

Визуализация генерации псевдослучайных 32-битных целых чисел с использованием Mersenne Twister. В разделе «Извлечение числа» показан пример, где целое число 0 уже выведено, а индекс равен целому числу 1. «Создать числа» запускается, когда все целые числа были выведены.

Для длины слова w -бит Мерсенна Твистер генерирует целые числа в диапазоне $[0,2^{w}-1]$.

Алгоритм Mersenne Twister основан на матричной линейной рекуррентности над конечным двоичным полем. ${\textbf {F}}_{2}$. Алгоритм представляет собой регистр сдвига со скрученной обобщенной обратной связью. (twisted GFSR, или TGFSR) рациональной нормальной формы (TGFSR(R)), с отражением и смягчением битов состояния. Основная идея состоит в том, чтобы определить ряд $x_{i}$ через простое рекуррентное соотношение, а затем выводить числа в виде $x_{i}^{T}$, где T — обратимый ${\textbf {F}}_{2}$-матрица называется матрицей закалки .

Общий алгоритм характеризуется следующими величинами:

• w : размер слова (в битах)
• n : степень рецидива
• m : среднее слово, смещение, используемое в рекуррентном отношении, определяющем серию. $x$, $1\leq m<n$
• r : точка разделения одного слова или количество бит младшей битовой маски, $0\leq r\leq w-1$
• a : коэффициенты матрицы скручивания рациональной нормальной формы
• b , c : Битовые маски смягчения TGFSR(R)
• s , t : Сдвиг битов отпуска TGFSR(R)
• u , d , l : дополнительные смены/маски закалочных долот Mersenne Twister

с тем ограничением, что $2^{nw-r}-1$ является простым числом Мерсенна. Этот выбор упрощает тест на примитивность и тест на k -распределение , которые необходимы при поиске параметров.

Серия $x$ определяется как серия w -битных величин с рекуррентным соотношением:

$x_{k+n}:=x_{k+m}\oplus \left(({x_{k}}^{u}\mid {x_{k+1}}^{l})A\right)\qquad k=0,1,2,\ldots$

где $\mid$ обозначает конкатенацию битовых векторов (со старшими битами слева), $\oplus$ побитовое исключающее или (XOR), $x_{k}^{u}$ означает старшие w - r биты $x_{k}$, и $x_{k+1}^{l}$ означает младшие r биты $x_{k+1}$.

Все индексы могут быть смещены на -n.

$x_{k}:=x_{k-(n-m)}\oplus \left(({x_{k-n}}^{u}\mid {x_{k-(n-1)}}^{l})A\right)\qquad k=n,n+1,n+2,\ldots$

где сейчас ЛХС, $x_{k}$, — это следующее сгенерированное значение в серии с точки зрения значений, сгенерированных в прошлом, которые находятся в правой части.

Твист-преобразование A определяется в рациональной нормальной форме как:

$$A={\begin{pmatrix}0&I_{w-1}\\a_{w-1}&(a_{w-2},\ldots ,a_{0})\end{pmatrix}}$$

с $I_{w-1}$ как $(w-1)(w-1)$ идентификационная матрица. Преимущество рациональной нормальной формы состоит в том, что умножение на A можно эффективно выразить как: (помните, что здесь умножение матриц выполняется в ${\textbf {F}}_{2}$, и поэтому место сложения занимает побитовое исключающее ИЛИ)

$${\boldsymbol {x}}A={\begin{cases}{\boldsymbol {x}}\gg 1&x_{0}=0\\({\boldsymbol {x}}\gg 1)\oplus {\boldsymbol {a}}&x_{0}=1\end{cases}}$$

где $x_{0}$ это бит младшего порядка $x$.

Как и в случае с TGFSR(R), «Твистер Мерсенна» каскадируется с умеренным преобразованием , чтобы компенсировать пониженную размерность равнораспределения (из-за того, что A находится в рациональной нормальной форме). Обратите внимание, что это эквивалентно использованию матрицы A , где $A=T^{-1}*AT$ для T обратимая матрица, и поэтому приведенный ниже анализ характеристического полинома остается в силе.

Как и в случае с A , мы выбираем смягчающее преобразование так, чтобы оно было легко вычислимым, и поэтому фактически не конструируем T само по себе. В случае Mersenne Twister этот отпуск определяется как

${\begin{aligned}y&\equiv x\oplus ((x\gg u)~\And ~d)\\y&\equiv y\oplus ((y\ll s)~\And ~b)\\y&\equiv y\oplus ((y\ll t)~\And ~c)\\z&\equiv y\oplus (y\gg l)\end{aligned}}$

где $x$ это следующее значение из ряда, $y$ является временным промежуточным значением, и $z$ значение, возвращаемое алгоритмом, с $\ll$и $\gg$ при побитовом сдвиге влево и вправо , и $\&$ побитовое И. как Первое и последнее преобразования добавлены для улучшения равнораспределения младших разрядов. Из собственности ТГФСР, $s+t\geq \left\lfloor {\frac {w}{2}}\right\rfloor -1$ требуется для достижения верхней границы равнораспределения старших битов.

Коэффициенты для MT19937:

${\begin{aligned}(w,n,m,r)&=(32,624,397,31)\\a&={\textrm {9908B0DF}}_{16}\\(u,d)&=(11,{\textrm {FFFFFFFF}}_{16})\\(s,b)&=(7,{\textrm {9D2C5680}}_{16})\\(t,c)&=(15,{\textrm {EFC60000}}_{16})\\l&=18\\\end{aligned}}$

Обратите внимание, что 32-битные реализации Mersenne Twister обычно имеют d = FFFFFFFF ₁₆ . В результате d иногда опускается в описании алгоритма, поскольку побитовое и с d в этом случае не имеет никакого эффекта.

Коэффициенты для MT19937-64:

${\begin{aligned}(w,n,m,r)=(64,312,156,31)\\a={\textrm {B5026F5AA96619E9}}_{16}\\(u,d)=(29,{\textrm {5555555555555555}}_{16})\\(s,b)=(17,{\textrm {71D67FFFEDA60000}}_{16})\\(t,c)=(37,{\textrm {FFF7EEE000000000}}_{16})\\l=43\\\end{aligned}}$

Инициализация

Состояние, необходимое для реализации Mersenne Twister, представляет собой массив из n значений по w бит каждое. Для инициализации массива w -битное начальное значение. используется $x_{0}$ через $x_{n-1}$ установив $x_{0}$ до начального значения и после этого установка

$x_{i}=f\times (x_{i-1}\oplus (x_{i-1}\gg (w-2)))+i$

для $i$ от $1$ к $n-1$.

• Первое значение, которое затем генерирует алгоритм, основано на $x_{n}$, не на $x_{0}$.
• Константа f образует еще один параметр генератора, но не является частью самого алгоритма.
• Значение f для MT19937 составляет 1812433253.
• Значение f для MT19937-64 составляет 6364136223846793005.

Код C

#include <stdint.h>

#define n 624
#define m 397
#define w 32
#define r 31
#define UMASK (0xffffffffUL << r)
#define LMASK (0xffffffffUL >> (w-r))
#define a 0x9908b0dfUL
#define u 11
#define s 7
#define t 15
#define l 18
#define b 0x9d2c5680UL
#define c 0xefc60000UL
#define f 1812433253UL

typedef struct
{
    uint32_t state_array[n];         // the array for the state vector 
    int state_index;                 // index into state vector array
} mt_state;


void initialize_state(mt_state* state, uint32_t seed) 
{
    uint32_t* state_array = &(state->state_array[0]);
    
    state_array[0] = seed;                                  // suggested initial seed = 19650218UL
    
    for (int i=1; i<n; i++)
    {
        seed = f * (seed ^ (seed >> (w-2))) + i;    // Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier.
        state_array[i] = seed; 
    }
    
    state->state_index = 0;
}


uint32_t random_uint32(mt_state* state)
{
    uint32_t* state_array = &(state->state_array[0]);
    
    int k = state->state_index;      // point to current state location
    
//  int k = k - n;                   // point to state n iterations before
//  if (k < 0) k += n;               // modulo n circular indexing
                                     // the previous 2 lines actually do nothing
                                     //  for illustrative purposes only
    
    int j = k - (n-1);               // point to state n-1 iterations before
    if (j < 0) j += n;               // modulo n circular indexing

    uint32_t x = (state_array[k] & UMASK) | (state_array[j] & LMASK);
    
    uint32_t xA = x >> 1;
    if (x & 0x00000001UL) xA ^= a;
    
    j = k - (n-m);                   // point to state n-m iterations before
    if (j < 0) j += n;               // modulo n circular indexing
    
    x = state_array[j] ^ xA;         // compute next value in the state
    state_array[k++] = x;            // update new state value
    
    if (k >= n) k = 0;               // modulo n circular indexing
    state->state_index = k;          // 0 <= state_index <= n-1   always
    
    uint32_t y = x ^ (x >> u);       // tempering 
    y = y ^ ((y << s) & b);
    y = y ^ ((y << t) & c);
    uint32_t z = y ^ (y >> l);
    
    return z; 
}

с

Чтобы достичь $2^{nw-r}-1$ теоретический верхний предел периода в Т GFSR , $\phi _{B}(t)$ должен быть примитивным полиномом , $\phi _{B}(t)$ являющийся полиномом характеристическим

$B={\begin{pmatrix}0&I_{w}&\cdots &0&0\\\vdots &&&&\\I_{w}&\vdots &\ddots &\vdots &\vdots \\\vdots &&&&\\0&0&\cdots &I_{w}&0\\0&0&\cdots &0&I_{w-r}\\S&0&\cdots &0&0\end{pmatrix}}{\begin{matrix}\\\\\leftarrow m{\text{-th row}}\\\\\\\\\end{matrix}}$

$S={\begin{pmatrix}0&I_{r}\\I_{w-r}&0\end{pmatrix}}A$

Преобразование твиста улучшает классический GFSR следующими ключевыми свойствами:

• Период достигает теоретического верхнего предела $2^{nw-r}-1$ (кроме случая, когда инициализировано значением 0)
• Равнораспределение в n измерениях (например, линейные конгруэнтные генераторы в лучшем случае могут обеспечить разумное распределение в пяти измерениях)

Варианты

CryptMT — это поточный шифр и криптографически безопасный генератор псевдослучайных чисел , который внутри использует Mersenne Twister. Он был разработан Мацумото и Нисимурой вместе с Марико Хагитой и Муцуо Сайто. Он был представлен проекту eSTREAM сети eCRYPT . В отличие от Mersenne Twister или других его производных, CryptMT запатентован .

MTGP — это вариант Mersenne Twister, оптимизированный для графических процессоров, опубликованный Муцуо Сайто и Макото Мацумото. Базовые операции линейной рекурсии расширены из MT, а параметры выбраны так, чтобы позволить множеству потоков вычислять рекурсию параллельно, одновременно разделяя свое пространство состояний для уменьшения нагрузки на память. В документе утверждается улучшенное равномерное распределение по сравнению с MT и производительность на старом графическом процессоре (2008 года выпуска) ( Nvidia GTX260 со 192 ядрами) на уровне 4,7 мс для 5×10. ⁷ случайные 32-битные целые числа.

SFMT ( SIMD -ориентированный быстрый вихрь Мерсенна) представляет собой вариант вихря Мерсенна, представленный в 2006 году. спроектирован так, чтобы быть быстрым при работе на 128-битном SIMD.

• Он примерно в два раза быстрее Mersenne Twister.
• Он имеет лучшее свойство равнораспределения точности v-бит, чем MT, но хуже, чем WELL («Well Equidistributed Long- period Linear») .
• Он быстрее восстанавливается из начального состояния с нулевым превышением, чем MT, но медленнее, чем WELL.
• Он поддерживает различные периоды от 2 ⁶⁰⁷ − 1 к 2 ²¹⁶⁰⁹¹  − 1.

Intel SSE2 и PowerPC AltiVec поддерживаются SFMT. Он также используется для игр с Cell BE на PlayStation 3 .

TinyMT — это вариант Mersenne Twister, предложенный Сайто и Мацумото в 2011 году. TinyMT использует всего 127 бит пространства состояний, что значительно меньше по сравнению с исходными 2,5 КиБ состояния. Однако он имеет период $2^{127}-1$, намного короче оригинала, поэтому авторы рекомендуют его только в тех случаях, когда память ограничена.

Характеристики

Преимущества:

• Лицензионно-лицензированный и непатентованный для всех вариантов, кроме CryptMT.
• Проходит многочисленные тесты на статистическую случайность, включая тесты Дайхарда и большинство, но не все, тестов TestU01 .
• Очень длительный период $2^{19937}-1$. Обратите внимание, что длительный период не является гарантией качества генератора случайных чисел, короткие периоды, такие как $2^{32}$ часто встречающийся во многих старых пакетах программного обеспечения, может быть проблематичным.
• k – распределяется с точностью до 32 бит для каждого $1\leq k\leq 623$ (определение k -распределенного см. ниже )
• Реализации обычно создают случайные числа быстрее, чем аппаратно реализованные методы. Исследование показало, что Mersenne Twister создает 64-битные случайные числа с плавающей запятой примерно в двадцать раз быстрее, чем аппаратно реализованный процессорный набор инструкций RDRAND .

Недостатки:

• Относительно большой буфер состояний, почти 2,5 КБ , если не используется вариант TinyMT.
• Посредственная пропускная способность по современным стандартам, если не используется вариант SFMT (обсуждаемый ниже).
• Демонстрирует две явные ошибки (линейная сложность) в Crush и BigCrush в пакете TestU01. Тест, как и Mersenne Twister, основан на ${\textbf {F}}_{2}$-алгебра.
• Множественные экземпляры, которые отличаются только начальным значением (но не другими параметрами), обычно не подходят для моделирования Монте-Карло , требующего независимых генераторов случайных чисел, хотя существует метод выбора нескольких наборов значений параметров.
• Плохая диффузия: может потребоваться много времени, чтобы начать генерировать выходные данные, которые проходят тесты на случайность , если начальное состояние сильно неслучайно, особенно если начальное состояние имеет много нулей. Следствием этого является то, что два экземпляра генератора, запущенные с почти одинаковыми начальными состояниями, обычно выдают почти одну и ту же последовательность в течение многих итераций, прежде чем в конечном итоге разойдутся. В обновлении алгоритма MT от 2002 года улучшена инициализация, поэтому начало с такого состояния маловероятно. Говорят, что версия с графическим процессором (MTGP) еще лучше.
• Содержит подпоследовательности, в которых нулей больше, чем единиц. Это добавляет к плохому свойству диффузии, что затрудняет восстановление из состояний со многими нулями.
• Не является криптографически безопасным , если не CryptMT используется вариант (обсуждаемый ниже). Причина в том, что соблюдение достаточного количества итераций (624 в случае MT19937, поскольку это размер вектора состояния, из которого производятся будущие итерации) позволяет предсказать все будущие итерации.

Приложения

Mersenne Twister используется в качестве ГПСЧ по умолчанию в следующем программном обеспечении:

• Языки программирования: Dyalog APL , ИДЛ , Р , Руби , Свободный Паскаль , PHP , Python (также доступен в NumPy значение по умолчанию было изменено на PCG64 ). , однако в версии 1.17 ), CMU Common Lisp , Встраиваемый Common Lisp , Steel Bank Common Lisp , Julia (до версии Julia 1.6 LTS, будет доступна позже, но начиная с версии 1.7 по умолчанию используется лучший/более быстрый ГСЧ)
• Linux Библиотеки и программное обеспечение : GLib , Библиотека арифметики множественной точности GNU , ГНУ Октава , Научная библиотека ГНУ
• Другое: Microsoft Excel , ГАУСС , Гретель , Был , МудрецМатематика , Скилаб , Клен , МАТЛАБ

Он также доступен в Apache Commons , в стандартной библиотеке C++ (начиная с C++11 ), и в Математике . Дополнительные реализации предоставляются во многих программных библиотеках, включая библиотеки Boost C++ . библиотека CUDA , и Цифровая библиотека NAG .

Mersenne Twister — один из двух ГПСЧ в SPSS : другой генератор сохранен только для совместимости со старыми программами, а Mersenne Twister считается «более надежным». Mersenne Twister также является одним из PRNG в SAS : остальные генераторы устарели и устарели. Mersenne Twister — это PRNG по умолчанию в Stata , другой — KISS , для совместимости со старыми версиями Stata.

Альтернативы

Альтернативный генератор WELL («Well Equidistributed Long- period Linear») предлагает более быстрое восстановление, равную случайность и почти равную скорость.

Марсальи Генераторы и варианты xorshift являются самыми быстрыми в классе LFSR.

64-битные MELG («64-битные максимально равнораспределенные ${\textbf {F}}_{2}$-Линейные генераторы с простым периодом Мерсенна") полностью оптимизированы с точки зрения свойств k -распределения.

Семейство ACORN (опубликовано в 1989 г.) — это еще один PRNG с k -распределением, который демонстрирует скорость вычислений, аналогичную MT, и лучшие статистические свойства, поскольку удовлетворяет всем текущим (2019 г.) критериям TestU01; при использовании с соответствующим выбором параметров ACORN может иметь сколь угодно большой период и точность.

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

Дальнейшее чтение

• Харасе, С. (2014), «На $\mathbb {F} _{2}$-линейные отношения генераторов псевдослучайных чисел Мерсенна Твистера», Mathematics and Computers in Simulation , 100 : 103–113, arXiv : 1301.5435 , doi : 10.1016/j.matcom.2014.02.002 , S2CID   6984431 .
• Харасе, С. (2019), «Преобразование Мерсенна Твистера в числа с плавающей запятой двойной точности», Mathematics and Computers in Simulation , 161 : 76–83, arXiv : 1708.06018 , doi : 10.1016/j.matcom.2018.08.006 , S2CID   19777310 .

Внешние ссылки

• Научная статья по MT и связанные с ней статьи Макото Мацумото.
• Домашняя страница Mersenne Twister с кодами на C, Fortran, Java, Lisp и некоторых других языках.
• Примеры Mersenne Twister — коллекция реализаций Mersenne Twister на нескольких языках программирования — на GitHub.
• SFMT в действии: Часть I — Создание DLL, включая поддержку SSE2 — в Code Project