Metody Monte Carlo są rozumiane na różne sposoby. Można powiedzieć, że stanowią one narzędzia numeryczne oparte na teorii rachunku prawdopodobieństwa służące do rozwiązywania problemów, w tym problemów o podłożu deterministycznym, najczęściej wtedy, gdy zawodzą metody analityczne czy klasyczne metody numeryczne. Jednym z twórców współczesnych metod Monte Carlo był polski matematyk Stanisław ULam pracujący w latach 40-tych XX wieku w Stanach Zjednoczonych nad konstrukcją broni jądrowej. Wraz z rozwojem technologii komputerowych metody Monte Carlo osiągnęły zawrotny rozwój i są wykorzystywane obecnie w wielu dziedzinach nauki i techniki.
Już przykład 10.21 oraz przykład 10.22 podane przez na przy okazji omawiania znaczenia mocnego prawa wielkich liczb ilustrowały podstawowe idee metod Monte Carlo. Sygnalizujemy teraz inny klasyczny przykład. Różni się on tym od tych dwóch, że używa narzędzia probabilistyczne do rozwiązania problemu o charakterze deterministycznym.
Przykład – 12.1 Należy podać przybliżoną wartość liczby \(\pi \).
Schemat rozwiązania: Losujemy dużo punktów, powiedzmy \(n\), \(n \ge 1000\), z kwadratu \([-1,1]^2\) i wyznaczamy liczbę, \(k\), punktów, które spełniają nierówność \(x^2 + y^2 \le 1\). Liczbę \(4\frac {k}{n}\) uważamy za przybliżenie \(\pi \).
Uzasadnimy później sensowność powyższego postępowania. Pokażemy istotne uogólnienie. Zaczniemy od wyjaśnienie kwestii losowania.
Metody Monte Carlo bazują na możliwości generowania przez komputer w bardzo krótkim czasie bardzo wielu liczb, które można uważać za wielkości losowe niezależnie od tego, że najczęściej są one generowane przez algorytmy deterministyczne. Takie liczby nazywamy liczbami pseudolosowymi. Po omówieniu sposobów pozyskiwania liczb pseudolosowych powiemy jak mogą być one użyte do przybliżonego wyznaczania całek, a także w problemach optymalizacyjnych.
Problem pozyskiwania liczb pseudolosowych o dobrych własnościach jest kluczowy dla wielu obliczeń metodami Monte Carlo oraz w symulacjach i nadal pozostaje w centrum uwagi specjalistów.
Jeszcze do niedawna programy komputerowe używały najczęściej następującego algorytmu: dla ustalonych liczb naturalnych \(a\), \(b\) i \(p\) wybieramy dowolną liczbę naturalną \(X_0\), zwaną ziarnem (ang. seed), a następnie określamy rekurencyjnie ciąg:
\[ X_{n+1} = a X_n + b \ (\mbox { mod }\, p). \]
Mówiąc inaczej, za każdym razem obliczamy \(X_{n+1}' = a X_n + b\), a jako \(X_{n+1}\) bierzemy resztę z dzielenia \(X_{n+1}'\) przez \(p\) – tak więc wszystkie wyrazy naszego ciągu są liczbami naturalnymi mniejszymi od \(p\).
Jeżeli parametry \(a\), \(b\), \(p\) i \(X_0\) są odpowiednio dobrane, to okazuje się, że liczby:
\[ U_n = X_n /p \]
mają własności niemal takie same, jak liczby wylosowane niezależnie jedna od drugiej z rozkładu jednostajnego na przedziale \((0,1)\).
Przykład – 12.2 \(a = 5\), \(b = 0\), \(p=7\). \(X_0 = 4\). Wtedy otrzymujemy kolejno:
\(4, 6, 2, 3, 1, 5, 4, 6, \dots \) i po podzieleniu przez 7 odpowiednio: \(.5714285716,\) \(.8571428574,\) \(.2857142858,\) \(.4285714287, \) \(.1428571429,\) \(.7142857145,\) \(.5714285716,\) \(.8571428574, \dots \).
Istnieją pewne zasady wybierania parametrów. W szczególności, \(p\) powinno być bardzo duże, aby jak najbardziej ograniczyć zjawisko okresowości. Z podobnych względów także \(a\) powinno być dużą liczbą, najlepiej względnie pierwszą z \(p\). Natomiast wybór \(b\) ma mniejsze znaczenie – często przyjmuje się \(b=0\).
Okazuje się, że przy odpowiednio wybranych parametrach oraz przy zastosowaniu dodatkowych procedur liczby pseudolosowe i ich zestawy mają bardzo dobre własności – potwierdzają to także odpowiednie testy statystyczne.
Przykładowo, program Maple (już w wersji 5) używał generatora liczb pseudolosowych z bardzo dużymi parametrami \(a\) oraz \(p\), mianowicie:
\[a = 427419669081, \;\; p = 999999999989.\]
Wartość ziarna \(X_0\) można w każdej chwili zadać zgodnie z aktualnymi potrzebami. Może nam na przykład zależeć, aby przy powtórzeniach danego losowania zawsze otrzymywać te same liczby pseudolosowe – zadajemy wtedy taką samą (stałą) wartość \(X_0\). Z drugiej strony, możemy żądać, aby w każdym losowaniu otrzymywać inne liczby pseudolosowe – można to na przykład osiągnąć, zadając wartość ziarna w zależności od upływu czasu zużytego przez procesor od rozpoczęcia aktualnej sesji.
Mając liczby pseudolosowe z rozkładu jednostajnego na odcinku \([0,1]\), można, stosując odpowiednią metodę, uzyskać liczby pseudolosowe z innego zadanego rozkładu.
Przykład – 12.3 (liczby pseudolosowe z rozkładu dyskretnego) Niech rozkład \(Q\) będzie zadany przez ciągi: \(x_1,x_2, \dots , \), \(p_1,p_2, \dots , \). Czyli \(Q(x_i) = p_i\). Definiujemy liczby \(F_0 = 0\), \(F_i = F_{i-1} + p_i\). \(i = 1,2,3, \dots \). Przypuśćmy, że liczby pseudo-losowe \(u_1, u_2, u_3, \dots \) są wylosowane według rozkładu \(U(0,1)\). Wtedy punkty \(y_1, y_2, y_3, \dots \) wybrane według reguły:
\(y_j = x_k\), gdzie \(k\) jest takie, że: \(u_j \in (F_{k-1}, F_k]\).
są wybrane z rozkładu \(Q\). Uzasadnić (ćwiczenie).
Przykład – 12.4 (liczby pseudolosowe z rozkładu ciągłego) Załóżmy, że zmienna losowa \(X\) ma rozkład \(Q\) oraz dystrybuantę ciągłą, \(F\). Dla uproszczenia załóżmy, że \(F\) jest silnie rosnąca. Zauważmy, że zmienna losowa \(U = F(X) = F\circ X\) ma rozkład \(U(0,1)\).
Rzeczywiście, dla każdego \(0 < u < 1\): \(F_U(u) = P(U \le u) = P(F\circ X \le u) = P(X \le F^{-1}(u)) = F(F^{-1}(u)) = u\), a to oznacza, że \(F_U\) jest dystrybuantą rozkładu \(U(0,1)\).
Aby więc wylosować liczbę \(x\) z rozkładu \(Q\) można wylosować liczbę z rozkładu \(U(0,1)\) i wziąć \(x = F^{-1}(u)\).
Ponieważ odwracanie dystrybuanty może być trudne (czasochłonne), istnieją również inne metody pozyskiwania liczb pseudo-losowych z określonych typów rozkładów.
Obecnie chyba najczęściej stosowaną metodę pozyskiwania liczb pseudo losowych z rozkłądu \(U(0,1)\) stanowią różne wersje: Mersenne Twister Algorithm.