Simulace MONTE CARLO
Historie metody Monte Carlo sahá do čtyřicátých let dvacátého století (Stanislaw Marcin Ulam, John von Neuman, Národní laboratoř Los Alamos, USA) a poprvé byla použita při studiu chování neutronů, konkrétně jejich schopnosti procházet různými materiály, neboť právě tento problém nelze řešit teoreticky ani prakticky.
Metoda Monte Carlo má široké využití právě (a hlavně) v oblastech, kde klasické metody nepracují příliš efektivně, případně nejdou použít vůbec. Tato metoda má velmi široké použití například v oblasti simulace experimentů, výpočtů určitých integrálů (zejména vícerozměrných) nebo například při řešení diferenciálních rovnic.
Celá metoda Monte Carlo pracuje na velmi jednoduché myšlence – pokud je třeba určit hodnotu veličiny, která je výsledkem náhodného děje, stačí vytvořit počítačový model tohoto děje a nechat proběhnout dostatečné množství simulací a jejich výsledky poté vyhodnotit klasickými statistickými metodami.
Ukázka použití metody MONTE CARLO
Jako ukázku použití metody Monte Carlo jsem vybral jednoduchý výpočet Ludolfova čísla ?. Jako základ k určení hodnoty poslouží čtverec a jemu vepsaná čtvrtkružnice. Celé řešení spočívá v náhodných „hodech drobným předmětem“ do tohoto čtverce, kdy hodnotě ? odpovídá právě výsledný poměr počtu všech hodů do prostoru čtvrtkružnice (S1) k celkovému počtu hodů do prostoru čtverce (S2).
Přesnost určení hodnoty ? odpovídá počtu „hodů“ – tedy absolutně přesná hodnota bude určena jen při nekonečném počtu pokusů. Ve skutečnosti však stačí pro určení hodnoty s námi požadovanou přesností použít výrazně menšího počtu „hodů“.
Výpočet čísla PI dle postupu z předchozí kapitoly
Nyní se dle předcházející kapitoly podíváme na výpočet čísla ? v reálu. Dle velmi jednoduchého algoritmu jsem spočítal číslo ? s různou přesností, která byla ovlivněna především počtem hodů a zároveň hustotou sítě, ve které byly sledovány dopady jednotlivých hodů.
Na ukázce je také vysvětlen problém s generováním náhodných čísel pomocí dostupné funkce Rnd a navrženy a vyzkoušeny možnosti jeho řešení.
Pod hustotou sítě si můžeme představit množství jednotlivých čtverečků, na které je základní čtverec rozdělen. V našem případě se jedná o hodnoty od 10×10 po 10E06x10E06 čtverců. Počet hodů byl zkušebně určen v počtech od 10 do 1E09.
Nejlepší výsledek byl 3.1416336 pro počet hodů 1E07 a rozlišení čtverce 1E06x1E06. Celý výpočet – respektive jeho přesnost je limitovaná generátorem náhodných čísel. Pro případ požití jednoduchého generátoru přes funkci Rnd může dojít od určitého počtu „hodů“ k opakování náhodných čísel.
Jako ověření jsem zkusil použít při každém hodu restart generátoru náhodných čísel pomocí příkazu Randomize. Výsledky jsou uvedeny v následující tabulce:
Jak je vidět, restart generátoru náhodných čísel po každém hodu nebyl příliš prospěšný. Nyní tedy zkusíme restartovat generátor náhodných čísel vždy před začátkem jeho „opakování“, tedy každých 32768 hodů.
Toto vylepšení bude mít samozřejmě vliv jen při vyšším počtu hodů,
proto údaje pro nižší počet bodů nejsou zpřesněny – rozdíly v jejich
hodnotách jsou způsobeny pouze náhodností celého děje.
Prozatím bylo
dosaženo nejlepšího výsledku 3.14168, při nižším počtu hodů než bez
restartování generátoru náhodných čísel. Celkově mi přijde toto řešení
výhodné, zdá se, že rozložení náhodných čísel je výrazně lepší.
Pro seriózní generování náhodných čísel ale přesto doporučuji použít nějakého profesionálního generátoru náhodných čísel.
Oblasti použití metody MONTE CARLO
Tato metoda nalézá využití v širokém spektru aplikací, od matematiky (řešení integrálů, soustav lineárních rovnic, řešení diferenciálních rovnic…), fyziky, počítačové grafiky přes hazardní hry, finance a pojišťovnictví, analýzu rizik až po například výzkum polovodičů, studium životního prostředí, předpovědi struktury bílkovin nebo šíření ropných skvrn při haváriích.
Využití metody MONTE CARLO pro simulaci molekul
Pro pochopení konformace velkých molekul lze využít právě simulace vybrané molekuly pomocí metody Monte Carlo. V prvním kroku se vytvoří strukturní mřížka pro budoucí makromolekulu (dvourozměrná mřížka, koordinační číslo atomů je rovno čtyřem) a poté je zvolena startovní pozice prvního atomu. Při následné simulaci je zvoleno náhodné pokračování molekuly do jednoho z volných směrů (čtyři pro první atom, poté tři pro atomy další).
Pokud dojde při simulaci molekuly k tomu, že by měl přidávaný atom obsadit již zaplněnou pozici v mřížce, je simulace ukončena. Výsledkem těchto simulací je poté soubor molekul, který lze statisticky vyhodnotit a získat tak informace o příslušném polymeru (délka řetězce, její rozptyl…).
Všechny simulace pomocí metody Monte Carlo lze různě vylepšit a rozšířit a lze tak získat například více informací o vlivu energie konformací na řetězec polymeru. Dále je také možné například simulovat kinetiku změny konformace zvažováním různých voleb pro rotaci kolem zvolené vazby v molekule.
Ukázka simulace
Jako ukázku jsem vybral jednoduchou simulaci růstu molekuly, který je ukončen pouze případem, kdy molekula narazí do své již existující části. Pro jednoduchost mají všechny směry růstu stejnou pravděpodobnost.
Zdrojový kód je v tomto případě o něco složitější, je třeba především ověřovat, zda se přidávaný atom nevrací na předchozí pozici nebo zda nenarazí do již existující části celé molekuly. První případ je vyřešen tak, že se výběr nové pozice zopakuje pokud je zvolen směr, ve kterém je předchozí atom, v druhém případě simulace dané molekuly končí a je zaznamenána její délka.
V první simulaci jsem umožnil z každého bodu pokračování do 4 směrů (tedy čtyřvazný atom), respektive do tří směrů, neboť jeden směr je zabrán původní – již existující – vazbou.
Druhá simulace, sloužící pro srovnání, obsahuje jedinou změnu: atomu je osmivazný (osm možných směrů bylo zvoleno pro jednoduchost simulace a také jednoduchou prostorovou mřížku).
Následující tabulka uvádí výsledky jednotlivých simulací. Pro každou vaznost atomů byla vždy provedena simulace 10 až 100000 molekul a zaznamenána maximální dosažená délka molekuly a jejich průměrná délka (počet atomů v řetězci).
Jak vyplývá z hodnot v uvedených tabulkách, s vyšším počtem simulovaných molekul se zvyšuje šance na výskyt velmi dlouhé molekuly (odpovídá předpokladům) a zároveň se stabilizuje průměrná délka molekuly, v případě čtyřvazného atomu k hodnotě 5.6, v případě osmivazného atomu k hodnotě 7.0.
Závěrem přidávám i histogramy pro vybrané simulace pro počet molekul rovný 10000 a 100000 a obě vaznosti atomů.
Je patrné, že histogramy dodržují předpokládané rozložení, u osmivazných atomů dochází k zajímavému výskytu maxima u vyšší délky řetězce.
Poměrně zajímavou alternativou by byla simulace molekuly, kde by byla vyšší pravděpodobnost, že bude molekula lineární. Především z důvodů omezeného rozsahu této práce jsem se však tímto problémem dále nezabýval.
Ještě blíže realitě by byla simulace v prostoru, kdy by se základní mřížka rozrostla i ve třetím rozměru.