Вихор Мерсенна

Вихор Мерсенна (англ. Mersenne twister, MT) — генератор псевдовипадкових чисел (ГПВЧ), розроблений 1997 року японськими вченими Макото Мацумото (яп. 松本 眞) і Такудзі Нісімурою (яп. 西村 拓士). Вихор Мерсенна ґрунтується на властивостях простих чисел Мерсенна (звідси й назва) і забезпечує швидке генерування високоякісних за критерієм випадковості псевдовипадкових чисел.

Вихор Мерсенна позбавлений багатьох недоліків, властивих іншим ГПВЧ, таких як малий період, передбачуваність, легко відшукувані статистичні закономірності.

Проте, цей генератор не є криптостійким, що обмежує його використання в криптографії.

Існують принаймні два загальних варіанти алгоритму, що відрізняються тільки величиною використовуваного простого числа Мерсенна, найпоширенішим з яких є алгоритм MT19937, період якого становить 2¹⁹⁹³⁷ — 1 (приблизно 4,3·10⁶⁰⁰¹).

Властивості

Вихор Мерсенна є витковим регістром зсуву з узагальненою віддачею (TGFSR)^[1] (англ. twisted generalised feedback shift register). «Вихор» — це перетворення, яке забезпечує рівномірний розподіл генерованих псевдовипадкових чисел у 623 вимірах (для лінійних конгруентних генераторів його обмежено 5 вимірами). Тому функція кореляції між двома послідовностями вибірок у вихідній послідовності вихору Мерсенна мізерно мала.

Псевдовипадкова послідовність, породжена вихором Мерсенна, має дуже великий період, рівний числу Мерсенна, чого більш ніж достатньо для багатьох практичних застосувань.

Існують ефективні реалізації вихору Мерсенна, що перевершують за швидкістю багато стандартних ГПВЧ (зокрема, в 2-3 рази швидше від лінійних конгруентних генераторів). Алгоритм вихору Мерсенна реалізовано в програмній бібліотеці Boost^[2], Glib^[3] і стандартних бібліотеках для C++, Python^[4][5][6], Ruby^[7], R^[8], PHP^[9], MATLAB^[10], Autoit^[11].

Видавані вихором Мерсенна послідовності псевдовипадкових чисел успішно проходять статистичні тести Diehard, що підтверджує їхні хороші статистичні властивості.

Генератор не призначений для отримання криптографічно стійких випадкових послідовностей чисел. Для забезпечення криптостійкості вихідну послідовність генератора слід обробити одним із криптографічних алгоритмів хешування^[12].

K-розподіл

Запропоновано багато генераторів можливо «високої якості», але тільки деякі можна використати для серйозного моделювання через відсутність чіткого поняття «хаотичності» для генераторів псевдовипадкових чисел, оскільки кожен дослідник концентрується на конкретних параметрах хаотичності. Серед багатьох відомих мір, тести, засновані на вищому рівномірному розподілі, такі як спектральний тест і тест на k-розподілі, описаний нижче, вважаються найсильнішими.

Визначення

Кажуть, що псевдовипадкова послідовність x_i з періодом P, що складається з w-бітових цілих, має k-розподіл з v-бітовою точністю, якщо вона задовольняє такій умові: Нехай trunc_v(x) — число, утворене першими v бітами послідовності x_i, розглянемо P з k v-бітових векторів вигляду^[13] ${\displaystyle ({\text{trunc}}_{v}(x_{i}),\,{\text{trunc}}_{v}(x_{i+1}),\,...,\,{\text{trunc}}_{v}(x_{i+k-1}))\quad (0\leq i<P)}$. Тоді кожна з 2^kv можливих комбінацій бітів зустрічається рівне число разів, за винятком комбінації, що складається повністю з нулів, яка зустрічається на один раз менше.

Геометрична інтерпретація

Для кожного v = 1, 2,., w, нехай k(v) — найбільше число, таке, що послідовність є k-розподіленою з v-бітовою точністю. Поділимо кожне ціле x_i на 2^w для нормалізації в псевдовипадкове дійсне число x_i з інтервалу [0, 1]. Помістимо P точок в k-вимірний куб з координатами (x_i, x_i+1…, x_i+k-1)(i = 0, 1,…, P-1) для всього періоду. Кожна з осей даного k-вимірного куба поділена на 2^v інтервалів. Таким чином, ми поділили куб на 2^kv малих кубів. Послідовність є k-розподіленою з v-бітовою точністю, якщо кожен малий куб містить рівне число точок, крім куба в початку координат, який містить на одну точку менше. Отже, чим вище k(v) для кожного v, тим більше вимірів матиме розподіл з v-бітовою точністю. Під тестом на k-розподіл будемо розуміти отримання значень k(v).

Опис

x позначатимемо w-вимірні вектори над полем ${\displaystyle F_{2}}$ = {0,1}, відповідні машинним словам розміру w. Вихор Мерсенна генерує послідовність векторів, які є псевдовипадковими цілими з діапазону від 0 до 2^w — 1. Діленням на 2^w — 1 можна також отримати псевдовипадкове дійсне число з діапазону [0,1].

Алгоритм ґрунтується на такій рекурентній формулі:

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

де:

• n — ціле, яке позначає порядок рекурентності,
• m — ціле, 1 ≤ m < n,
• A — матриця розміру w×w, з елементами з ${\displaystyle F_{2}}$,
• ${\displaystyle \mid }$ — побітове АБО (OR),
• ${\displaystyle \oplus }$ — додавання за модулем два (XOR).

У правій частині x_k^u позначає «старші w-r біт» x_k і x_k+1^l «молодші r біт» x_k+1. Вектор (x_k^u | x_k+1^l) є конкатенацією старших w-r біт x_k і молодших r біт x_k+1. Взявши (x₀, x₁,…, x_n-1) для початкового заповнення. Тоді генератор обчислить x_n за рекурентним виразом при k= 0. Для k = 1,2,…, генератор обчислить x_n+1, x_n+2,… Форму матриці A обрано з розрахунку швидкості виконання множення на A:

${\displaystyle A={\begin{pmatrix}0&1&&&&\\0&0&1&&&\\0&\cdots &\cdots &\ddots &&\\&&&&\ddots &\\&&&&&1\\a_{w-1}&a_{w-2}&\cdots &\cdots &\cdots &a_{0}\end{pmatrix}}}$

І обчислення xA зводиться до побітових операцій:

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

де

${\displaystyle {\boldsymbol {x}}:=({x_{k}}^{u}\mid {x_{k+1}}^{l})\qquad \qquad k=0,1,\ldots }$

${\displaystyle {\boldsymbol {a}}:=({a_{w-1}},{a_{w-2}},\ldots ,{a_{0}})}$

${\displaystyle {\boldsymbol {x}}:=({x_{w-1}},{x_{w-2}},\ldots ,{x_{0}})}$

Неповні масиви

Неповні масиви

Послідовність МТ (x_k+n, x_k+n-1,…, x_k+1^u) утворює (n × w — r)-масив. Оскільки r бітів відкидаються, масив називають неповним масивом. Кожен елемент отримано з рекурентного співвідношення (1). Зміна стану MT також відбувається за лінійним співвідношенням і визначається за допомогою лінійного перетворення B.

Відповідно до загальної теорії лінійних рекурентних послідовностей, кожне значення в (n × w − r)-масиві є лінійною рекурентною послідовністю, що відповідає характеристичному многочлену φ_B(t) перетворення B. Матриця B має розміри (n × w − r) × (n × w − r) і таку форму:

${\displaystyle 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{\hbox{-th row}}\\\\\\\\\end{matrix}}}$

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

де ${\displaystyle I_{r}}$ — одинична матриця розміру r × r.

Для ${\displaystyle l=0,1,\ldots }$ виконується ${\displaystyle (x_{l+n},x_{l+n-1},\ldots ,x_{l+1}):=(x_{l+n-1},x_{l+n-2},\ldots ,x_{l})B}$.Характеристичний многочлен φ_B(t) перетворення B має вигляд:

${\displaystyle \Phi _{B}(t)=(t^{n}+t^{m})^{w-r}\cdot (t^{n-1}+t^{m-1})^{r}+a_{0}(t^{n}+t^{m})^{w-r}\cdot (t^{n-1}+t^{m-1})^{r-1}+...+a_{r-2}(t^{n}+t^{m})^{w-r}\cdot (t^{n-1}+t^{m-1})+a_{r-1}(t^{n}+t^{m})^{w-r}+a_{r}(t^{n}+t^{m})^{w-r-1}+...+a_{w-2}(t^{n}+t^{m})+a_{w-1}}$

Послідовність досягає максимального періоду 2^(nw-r)−1, тоді і тільки тоді коли φ_B(t) — примітивний.

Загартування

Необроблені послідовності, що генеруються рекурсією (1), мають поганий рівномірний розподіл на великих розмірностях. Щоб це виправити, використовується метод загартування (англ. tempering), на виході якого отримують кінцеву псевдовипадкову послідовність. Метод полягає в тому, що кожне згенероване слово множиться праворуч на спеціальну оборотну матрицю T розміром w×w. Для матриці T: x → z = x T, вибрано такі послідовні перетворення:

y := x ⊕ (x >> u)

y := :y ⊕ ((y << s) & b)

y := :y ⊕ ((y << t) & c)

z := y ⊕ (y >> l)

де l, s, t і u — цілі, а b і c — спеціально підібрані бітові маски розміру слова, і (x≫u) позначає побітову операцію зсуву вправо на u бітів.

Алгоритм

Вихор Мерсенна алгоритмічно реалізується двома основними частинами: рекурсивною і загартування. Рекурсивна частина являє собою регістр зсуву з лінійним зворотним зв'язком, у якому всі біти в його слові визначаються рекурсивно; потік вихідних бітів визначається також рекурсивно функцією бітів стану.

Блок-схема.

Регістр зсуву складається з 624 елементів, і, в цілому, з 19937 клітинок. Кожен елемент має довжину 32 біти за винятком першого елемента, який має тільки 1 біт завдяки відкиданню біта.

Процес генерування починається з логічного множення на бітову маску, що відкидає 31 біт (крім найбільш значущих).

Наступним кроком виконується ініціалізація (x₀, x₁,…, x₆₂₃) будь-якими беззнаковими 32-розрядними цілими числами. Наступні кроки включають об'єднання та перехідні стани.

Зміна стану МТ.

Простір станів має 19937 біт (624·32 — 31). Наступний стан генерується зсувом одного слова вертикально вгору і вставленням нового слова в кінець. Нове слово обчислюється гамуванням середньої частини з виключеною^[14]. Вихідна послідовність починається з x₆₂₄, x₆₂₅,…

Алгоритм має шість етапів:

 Крок 0. Попередньо ініціалізується значення сталих u, h, a
 u ← 10…0   бітова маска старших w-r бітів,
 h ← 01…1   бітова маска молодших r бітів,
 a ← aw-1aw-2…a0  останній рядок матриці A.

 Крок 1. x[0], x[1],…,x[n-1] ←  початкове заповнення

 Крок 2. Обчислення (xiu | xi+1l)
 y ← (x[i] AND u1) OR (x[(i + 1) mod n] AND h1)

 Крок 3. Обчислюється значення наступного елемента послідовності за рекурентним виразом (1)
 x[i] ← x[(i + m) mod n] XOR (y>>1) XOR a   якщо молодший біт y = 1

 Або
 x[i] ← x[(i + m) mod n] XOR (y>>1) XOR 0   якщо молодший біт y = 0
 Крок 4. Обчислення x[i]T
 y ← x[i]
 y ← y XOR (y>>u)
 y ← y XOR ((y<<s) AND b)
 y ← y XOR ((y<<t) AND c)
 z ← y XOR (y>>l)
 виведення z
 Крок 5. i ← (i + 1) mod n 

 Крок 6. Перейти до кроку 2.

Параметри 32-бітового генератора MT

Параметри MT ретельно підібрано для досягнення згаданих вище властивостей. Параметри n і r вибрано так, що характеристичний многочлен примітивний або nw — r дорівнює числу Мерсенна 19937. Зверніть увагу, що значення w еквівалентне розміру слова комп'ютера. У цьому випадку це 32 біти. Тоді як значення n, m, r і w фіксуються, значення останнього рядка матриці A вибирається випадковим чином. Для чисел Мерсенна тест на простоту цілих значно простіший. Так, відомо багато простих чисел Мерсенна (тобто простих виду 2^p−1), до p=43112609 [Архівовано 13 березня 2019 у Wayback Machine.] . Параметри загартування (англ. tempering) підібрано так, що ми можемо отримати хороший рівномірний розподіл. Всі параметри МТ наведено в таблиці 1.

Таблиця 1

n	624
w	32
r	31
m	397
a	9908B0DF16
u	11
s	7
t	15
l	18
b	9D2C568016
c	EFC6000016

Інші варіанти реалізації

У зв'язку зі змінами в технології, зокрема, зростанням продуктивності процесорів, винайдено 64-бітові і 128-бітові версії МТ. 64-розрядну версію запропонував 2000 року Такудзі Нісімура,^[15] 128-розрядну — 2006 року^[16][17] теж він. Обидві версії мають той самий період, що й оригінальний MT.

64-бітовий MT має дві версії. Перша версія використовує те саме рекурсивне співвідношення, в другу версію додано ще два вектори, що збільшило кількість членів характеристичного многочлена:

${\displaystyle x_{k+n}:=x_{k+m2}\oplus x_{k+m1}\oplus x_{k+m0}\oplus ({x_{k}}^{u}\mid {x_{k+1}}^{l})A\qquad (2)}$

Порівняно з 32-розрядною версією MT, 64-розрядна версія працює швидше. Всі параметри для 64-бітової версії наведено в таблиці 2.

Таблиця 2

ID	Для рекурсивного співвідношення (1)		Для рекурсивного співвідношення (2)
n	312	312	312	312	312
w	64	64	64	64	64
r	31	31	31	31	31
m	156	156
m0			63	63	63
m1			151	151	151
m2			224	224	224
a	B5026F5AA96619E916	F6A3F020F058B7A716	B3815B624FC82E2F16	8EBD4AD46CB39A1E16	CACB98F78EBCD4ED16
b	D66B5EF5B4DA000016	28AAF6CDBDB4000016	599CFCBFCA66000016	656BEDFFD9A4000016	A51DBEFFDA6C000016
c	FDED6BE00000000016	FDEDEAE00000000016	FFFAAFFE0000000016	FDFECE7E0000000016	FFEE9BF60000000016
u	29	29	26	26	26
s	17	17	17	17	17
t	37	37	33	33	33
l	41	41	39	39	39

Примітки

1. Twisted GFSR generators, 1992.
2. boost/random/mersenne_twister.hpp. Boost C++ Libraries. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
3. Changes to GLib. GLib Reference Manual. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
4. 9.6 random — Generate pseudo-random numbers. Python v2.6.8 documentation. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
5. 8.6 random — Generate pseudo-random numbers. Python v3.2 documentation. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
6. random — Generate pseudo-random numbers — Python 3.8.3 documentation. Python 3.8.3 documentation. Архів оригіналу за 28 липня 2021. Процитовано 23 червня 2020.
7. "Random" class documentation. Ruby 1.9.3 documentation. Архів оригіналу за 26 червня 2021. Процитовано 29 травня 2012.
8. Random Number Generators. CRAN Task View: Probability Distributions. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
9. mt_srand. php documentation. Архів оригіналу за 19 листопада 2012. Процитовано 29 травня 2012.
10. Control random number generation (RNG). Matlab Documentation. Архів оригіналу за 12 вересня 2010. Процитовано 23 листопада 2015.
11. Function Random. Архів оригіналу за 6 квітня 2021. Процитовано 28 липня 2021.
12. CryptMT Stream Cipher.
13. Matsumoto, Nishimura, 2017.
14. Cryptographic Mersenne Twister and Fubuki Stream/Block Cipher, 2005.
15. Nishimura, 2000.
16. SIMD-oriented Fast Mersenne Twister (SFMT). Архів оригіналу за 9 листопада 2020. Процитовано 28 липня 2021.
17. SFMT: Comparison of speed. Архів оригіналу за 8 січня 2020. Процитовано 28 липня 2021.

Література

• M. Matsumoto, T. Nishimura. Mersenne twister: A 623-dimensionally equidistributed uniform pseudorandom number generator // ACM Trans. on Modeling and Computer Simulations : journal. — 2017. — Vol. 8, no. 1 (25 December). — P. 3—30. — DOI:10.1145/272991.272995.
• Matsumoto, M.; Kurita, Y. Twisted GFSR generators // ACM Trans. on Modeling and Computer Simulations. — 1992. — Т. 2, № 3 (25 December). — С. 179—194. — DOI:10.1145/146382.146383.
• Matsumoto, Makoto; Nishimura, Takuji; Hagita, Mariko; Saito, Mutsuo. Cryptographic Mersenne Twister and Fubuki Stream/Block Cipher. — 2005. — 25 December. Архівовано з джерела 27 липня 2021. Процитовано 28 липня 2021.
• T. Nishimura. Tables of 64-bit Mersenne twisters // ACM Trans. on Modeling and Computer Simulations. — 2000. — Т. 10, № 4 (25 December). — С. 248—357. — DOI:10.1145/369534.369540.
• Matsumoto M., Saito M., Nishimura T., Hagita M. CryptMT Stream Cipher Ver. 3.The eSTREAM Project. Архівовано з джерела 20 липня 2020. Процитовано 28 липня 2021.

Посилання

• Makoto Matsumoto, Academic publications on MT [Архівовано 3 серпня 2020 у Wayback Machine.](англ.)
• Mersenne Twister home page, with C code(англ.)
• Cuba — a library for multidimensional numerical integration [Архівовано 7 квітня 2014 у Wayback Machine.]
• The Mersenne Twister [Архівовано 26 червня 2021 у Wayback Machine.](англ.)