Generator pseudolosowy Mersenne Twister


Tematy pokrewne
Przedziały liczbowe i liczby
Liczby parzyste i nieparzyste
Liczby podzielne lub niepodzielne przez zadane podzielniki
Ciągi arytmetyczne
NWD – algorytm Euklidesa
Liczby względnie pierwsze
Najmniejsza wspólna wielokrotność
Odwrotność modulo – rozszerzony algorytm Euklidesa
Liczby pierwsze – generacja przez sprawdzanie podzielności
Liczby pierwsze – generacja sitem Eratostenesa
Liczby pierwsze – generacja sitem liniowym
Liczby pierwsze – generacja sitem Atkina-Bernsteina
Czynniki pierwsze – metoda próbnych dzieleń
Czynniki pierwsze – metoda Fermata
Pierwszość liczby naturalnej – algorytmy naiwne
Pierwszość liczby naturalnej – Chiński Test Pierwszości
Pierwszość liczby naturalnej – Małe Twierdzenie Fermata
Pierwszość liczby naturalnej – test Millera-Rabina
Liniowe generatory liczb pseudolosowych
Generator pseudolosowy Park-Miller
Generator pseudolosowy Mersenne Twister
Wbudowane generatory liczb pseudolosowych
Generowanie liczb pseudolosowych
Liczby Fibonacciego
System liczbowy Fibonacciego
Całkowity pierwiastek kwadratowy

 

Opisane w poprzednim rozdziale liniowe generatory kongruencyjne liczb pseudolosowych posiadają wiele różnych wad eliminujących je z niektórych zastosowań. W celu usunięcia tych niedogodności, w roku 1997 dwaj japończycy, Makoto Matsumoto i Takuji Nishimura, opracowali nowy rodzaj generatora liczb pseudolosowych o nazwie Mersenne Twister – w skrócie MT. Generator oparty został na liniowej rekurencji macierzowej (ang. matrix linear recursion) w skończonej dziedzinie binarnej (ang. finite binary field). Umożliwia on generację wysokiej jakości liczb pseudolosowych o bardzo dużym okresie powtarzania, będącym liczbą pierwszą Mersenne'a.

Liczby Mersenne'a są potęgami liczby 2 pomniejszonymi o 1:

 

Mn = 2n - 1

 

Oto kilka z nich:

 

1, 3, 7, 15, 31, 63, 127, 255, 511, 1023, 2047, 4095, 8191, 16383, 32767, 65535, ...

 

Liczbą pierwszą Mersenne'a (ang. Mersenne Prime) nazywamy liczbę Mersenne'a, która jest pierwsza. Nie wszystkie liczby Mersenne'a są liczbami pierwszymi. Jedno z założeń mówi, iż liczba Mersenne'a może być liczbą pierwszą tylko wtedy, gdy wykładnik n sam jest liczbą pierwszą. Jeśli wykładnik jest liczbą złożoną, to 2n - 1 nie jest pierwsze. Pozwala to szybko eliminować z ciągu liczb Mersenne'a liczby złożone. Jednakże pozostałe liczby należy sprawdzać testem pierwszości, gdyż nie ma pewności, iż są one rzeczywiście pierwsze. Na przykład dla n = 11 otrzymujemy:

 

211 - 1 = 2047 = 23 × 89 – liczba złożona

 

Zadanie to jest trudne obliczeniowo, ponieważ liczby Mersenne'a bardzo szybko rosną. Na dzień dzisiejszy największą znaną liczbą Mersenne'a jest:

 

232582657 - 1

 

Jest to 44-ta liczba pierwsza Mersenne'a i posiada 9808358 cyfr dziesiętnych. Znaleziono ją 4 września 2006 roku dzięki projektowi GIMPS (ang. Great Internet Mersenne Prime Search), który polegał na rozproszonym (tj. wykorzystującym moc obliczeniową komputerów podłączonych do sieci Internet) sprawdzaniu pierwszości liczb Mersenne'a.

Omawiany tutaj generator MT wykorzystuje 24-tą liczbę Mersenne'a o wartości 19937 – stąd bierze się jego oznaczenie MT19937. Generator posiada olbrzymi okres powtarzania równy 219937 - 1. Ponieważ teoria będąca podstawą działania generatora MT jest bardzo skomplikowana matematycznie (wykracza poza poziom liceum), ograniczymy się jedynie do podania samego algorytmu, bez wnikania w podstawy matematyczne jego działania. Wszystkie opisane tutaj fakty pochodzą z oryginalnej publikacji:

 

Mersenne Twister: A 623-dimensionally equidistributed uniform pseudorandom number generator

Makoto Matsumoto and Takuji Nishimura

Keio University/Max Planck Institut für Mathematik
Keio University

Cały dokument możesz pobrać stąd.

 

Generator MT jest jednym, dużym rejestrem przesuwnym ze sprzężeniami zwrotnymi powodującymi rotację oraz mieszanie się bitów. Rejestr ma długość 19937 bitów (liczba Mersenne), które w pamięci zajmują 624 słowa 32-bitowe. Generowanie liczb pseudolosowych składa się z trzech etapów.

 

  1. Pierwszy etap polega na zainicjowaniu 624 słów 32-bitowych rejestru wartościami pseudolosowymi. Do tego celu można wykorzystać zwykły liniowy generator LCG (my wykorzystujemy generator Super-Duper LCG(232,69069,0,X0)).
  2. W drugim etapie generujemy liczbę pseudolosową na podstawie zawartości rejestru. Rejestr jest odpowiednio modyfikowany przy tworzeniu każdej nowej liczby pseudolosowej.
  3. W trzecim etapie wygenerowaną liczbę odpowiednio dopasowujemy, aby otrzymać równomierny rozkład bitów.

Algorytm inicjowania rejestru MT

Wejście
MT    tablica 624 liczb 32-bitowych przechowująca 19937 bitów rejestru przesuwnego
X0   ziarno dla generatora LCG(4294967296,69069,0), który posłuży do generacji zawartości rejestru przesuwnego. Wartość tę można pozyskiwać np z zegara czasu rzeczywistego.
Wyjście:

Tablica MT wypełniona liczbami pseudolosowymi

Zmienna pomocnicza

i – indeksuje kolejne liczby tablicy MT, i N

Lista kroków:
K01: MT[0] ← X0 and $FFFFFFFF ; ziarno ustawiamy w pierwszym elemencie
K02: Dla i = 1,2,...,623 wykonuj K03  
K03:     MT[i] ← (69069 × MT[i-1]) and $FFFFFFFF ; kolejne elementy wypełniamy liczbami pseudolosowymi
K04: Zakończ  

 

Algorytm generatora pseudolosowego Mersenne Twister

Wejście
MT    tablica 624 liczb 32-bitowych przechowująca 19937 bitów rejestru przesuwnego
mti  – indeks, wg którego algorytm wybiera elementy MT do modyfikacji oraz do generacji wartości wyjściowej. Jest to zmienna statyczna, która zachowuje swoją wartość pomiędzy kolejnymi wywołaniami funkcji generującej liczby pseudolosowe. Powinna być zainicjowana wartością 0.
Wyjście:

Liczba pseudolosowa z zakresu od 0 do 232 - 1.

Zmienne pomocnicze
MA    dwuelementowa tablica, która służy do mnożenia elementów MT przez wybraną macierz.
A[0] = 0
A[1] = $9908B0DF
mti  – indeks, wg którego algorytm wybiera elementy MT do modyfikacji oraz do generacji wartości wyjściowej. Jest to zmienna statyczna, która zachowuje swoją wartość pomiędzy kolejnymi wywołaniami funkcji generującej liczby pseudolosowe. Powinna być zainicjowana wartością 0.
y  – generowana liczba pseudolosowa
Lista kroków:
K01: y ← (MT[mti] and $80000000) or (MT[(mti + 1) mod 624] and $7FFFFFFF) ; łączymy dwa wyrazy z MT
K02: MT[mti] ← MT[(mti + 397) mod 624] xor (y shr 1) xor MA[y and 1] ; tworzymy sprzężenie zwrotne
K03: yMT[mti] ; wyznaczamy liczbę pseudolosową
K04: yy xor (y shr 11)  
K05: yy xor ((y shl 7) and $9D2C5680)  
K07: yy xor ((y shl 15) and $EFC60000)  
K08: yy xor (y shr 18)  
K09: mti ← (mti + 1) mod 624 ; obliczamy następny indeks
K10: Zakończ z wynikiem y  

 

Program

Ważne:

Zanim uruchomisz program, przeczytaj wstęp do tego artykułu, w którym wyjaśniamy funkcje tych programów oraz sposób korzystania z nich.

 

Program odczytuje z pierwszego wiersza trzy liczby a,b oraz n. Liczby a i b określają przedział całkowity, w którym mają być wygenerowane liczby pseudolosowe. Liczba n określa ile liczb pseudolosowych należy wygenerować. Liczby pseudolosowe są generowane w kolejnych wierszach.

 

Lazarus
// Generator Mersenne Twister
// Data   : 16.04.2008
// (C)2012 mgr Jerzy Wałaszek
//---------------------------

program prg;

uses SysUtils,DateUtils;

var

  MT  : array[0..623] of longword;
  mti : integer;

// Inicjuje MT[]
//--------------
procedure InicjujMT(x0 : longword);

var
  i : integer;
  x : qword;
begin
  MT[0] := x0;
  for i := 1 to 623 do
  begin
    x := MT[i-1];
    x := (23023 * x) and $00000000ffffffff;
    x := (    3 * x) and $00000000ffffffff;
    MT[i] := x;
  end;
  mti := 0;
end;

// Inicjuje MT wartościami losowymi
//---------------------------------
procedure UprzypadkowijMT;

var
  t  : TDateTime;
  x0 : longword;
begin
  t  := Now;
  x0 := HourOf(t);
  x0 := x0 *   60 + MinuteOf(t);
  x0 := x0 *   60 + SecondOf(t);
  x0 := x0 * 1000 + MillisecondOf(t);
  InicjujMT(x0);
end;

// Generator Mersenne Twister
//--------------------------
function MersenneTwister : longword;
const
  MA : array[0..1] of longword = (0,$9908b0df);
var
  y : longword;
  i1,i397 : integer;
begin
  i1      := mti +   1; if   i1 > 623 then i1 := 0;
  i397    := mti + 397; if i397 > 623 then dec(i397,624);
  y       := (MT[mti] and $80000000) or (MT[i1] and $7fffffff);
  MT[mti] := MT[i397] xor (y shr 1) xor MA[y and 1];
  y       := MT[mti];
  y       := y xor ( y shr 11);
  y       := y xor ((y shl  7) and $9d2c5680);
  y       := y xor ((y shl 15) and $efc60000);
  y       := y xor ( y shr 18);
  mti     := i1;
  MersenneTwister := y;
end;

var
  a,b,i,n : longint;
begin
  UprzypadkowijMT;
  readln(a,b,n);
  for i := 1 to n do write(a + (MersenneTwister mod (b - a + 1)),' ');
  writeln;
end.
Code::Blocks
// Generator Mersenne Twister
// Data   : 16.04.2008
// (C)2012 mgr Jerzy Wałaszek
//---------------------------

#include <iostream>

using namespace std;

typedef unsigned long long ulong;

unsigned int MT[624];
int mti = 0;

// Inicjuje MT[]
//--------------
void InicjujMT(unsigned int x0)
{
  ulong x;

  MT[0] = x0;
  for(int i = 1; i < 623; i++)
  {
    x = MT[i-1];
    x = (23023 * x) & 0xffffffffull;
    x = (    3 * x) & 0xffffffffull;
    MT[i] = x;
  }
}

// Inicjuje MT wartościami losowymi
//---------------------------------
void UprzypadkowijMT()
{
  InicjujMT((unsigned int)time(NULL));
}

// Generator Mersenne Twister
//--------------------------
unsigned int MersenneTwister()
{
  const unsigned int MA[] = {0,0x9908b0df};
  long int y;
  int i1,i397;

  i1      = mti +   1; if(  i1 > 623) i1 = 0;
  i397    = mti + 397; if(i397 > 623) i397 -= 624;
  y       = (MT[mti] & 0x80000000) | (MT[i1] & 0x7fffffff);
  MT[mti] = MT[i397] ^ (y >> 1) ^ MA[y & 1];
  y       = MT[mti];
  y       ^=  y >> 11;
  y       ^= (y <<  7) & 0x9d2c5680;
  y       ^= (y << 15) & 0xefc60000;
  y       ^=  y >> 18;
  mti     = i1;
  return y;
}

int main()
{
  int a,b,i,n;

  UprzypadkowijMT();
  cin >> a >> b >> n;
  for(i = 1; i <= n; i++) cout << (a + (MersenneTwister() % (b - a + 1))) << " ";
  cout << endl << endl;

  return 0;
}
Free Basic
' Generator Mersenne Twister
' Data   : 16.04.2008
' (C)2012 mgr Jerzy Wałaszek
'---------------------------

Dim Shared MT(623) As Uinteger
Dim Shared mti As Integer = 0

' Inicjuje MT[]
'--------------
Sub InicjujMT(Byval x0 As Uinteger)

  Dim i As Integer, x As Ulongint

  MT(0) = x0
  For i = 1 To 623
    x = MT(i-1)
    x = (23023 * x) And &Hffffffff
    x = (    3 * x) And &Hffffffff
    MT(i) = x
  Next  
End Sub

' Inicjuje MT wartościami losowymi
'---------------------------------
Sub UprzypadkowijMT()
  InicjujMT(Timer * 1000)
End Sub

' Generator Mersenne Twister
'--------------------------
Function MersenneTwister() As Uinteger

  Dim MA(1) As Uinteger => {0,&H9908b0df}
  Dim y As Uinteger
  Dim As Integer i1,i397

  i1      = mti +   1: If   i1 > 623 Then i1 = 0
  i397    = mti + 397: If i397 > 623 Then i397 -= 624
  y       = (MT(mti) And &H80000000) Or (MT(i1) And &H7fffffff)
  MT(mti) = MT(i397) Xor (y Shr 1) Xor MA(y And 1)
  y       = MT(mti)
  y       = y Xor ( y Shr 11)
  y       = y Xor ((y Shl  7) And &H9d2c5680)
  y       = y Xor ((y Shl 15) And &Hefc60000)
  y       = y Xor ( y Shr 18)
  mti     = i1
  MersenneTwister = y
End Function

Dim As Integer a,b,i,n 

UprzypadkowijMT()

Open Cons For Input As #1

Input #1,a,b,n

Close #1

For i = 1 To n
  Print a + (MersenneTwister() Mod (b - a + 1));
Next

Print

End
Wynik
1 6 200
6 2 3 6 2 5 1 5 4 3 2 3 2 4 1 1 3 6 4 2 4 5 5 4 4 5 3 4 4 3 3 6 1 3 5 1 1 4 4 3
4 3 2 5 1 2 2 1 6 3 1 1 1 3 4 6 6 4 1 5 1 3 4 2 5 1 6 3 2 5 5 5 5 6 5 6 3 2 5 3
5 5 4 6 5 1 5 4 4 6 5 6 1 5 3 1 6 3 4 5 3 6 5 3 5 6 5 5 6 2 2 1 1 4 6 1 5 3 4 5
2 1 1 3 6 1 1 2 4 5 4 5 4 4 1 1 1 1 4 2 3 3 1 1 5 1 2 3 2 1 4 4 6 5 3 6 6 4 2 3
6 5 1 6 4 1 5 2 2 6 2 1 6 2 6 2 4 3 4 1 6 3 2 6 2 3 3 4 4 1 6 3 2 5 5 3 3 4 2 6

 



List do administratora Serwisu Edukacyjnego Nauczycieli I LO

Twój email: (jeśli chcesz otrzymać odpowiedź)
Temat:
Uwaga: ← tutaj wpisz wyraz  ilo , inaczej list zostanie zignorowany

Poniżej wpisz swoje uwagi lub pytania dotyczące tego rozdziału (max. 2048 znaków).

Liczba znaków do wykorzystania: 2048

 

W związku z dużą liczbą listów do naszego serwisu edukacyjnego nie będziemy udzielać odpowiedzi na prośby rozwiązywania zadań, pisania programów zaliczeniowych, przesyłania materiałów czy też tłumaczenia zagadnień szeroko opisywanych w podręcznikach.



   I Liceum Ogólnokształcące   
im. Kazimierza Brodzińskiego
w Tarnowie

©2017 mgr Jerzy Wałaszek

Dokument ten rozpowszechniany jest zgodnie z zasadami licencji
GNU Free Documentation License.