1. MOLEKULÁRIS SZIMULÁCIÓK
1.7. A dinamikus Monte Carlo módszer
A molekuláris dinamika (MD) a rendszert alkotó részecskék mozgásegyenleteinek megoldásával időben követi a rendszer változásait, ezért egyensúlyi és nemegyensúlyi rendszerek vizsgálatára egyaránt alkalmas. Annak ellenére, hogy az MC módszert főként egyensúlyi rendszerekre alkalmazzák, nemegyensúlyi vagy relaxációs folyamatok szimulációjára is használható [36-373839404142]. A módszert ilyenkor gyakran dinamikus MC módszerként (DMC) emlegetik. A technika azon alapul, hogy a sztochasztikus folyamat során hagyományos MC lépésekkel generált, a Markov-lánc tulajdonságával rendelkező állapotok sorozata megfeleltethető a rendszer időbeli változásának [43, 44]. Ez alapján a ξi
konfigurációk ξ1(t1)→ξ2(t2)→...ξi(ti)→ξi+1(ti+1)→... Markov-lánca úgy értelmezhető, hogy a rendszer ti ideig ξi konfigurációban marad, azaz ξi és ξi+1 konfigurációk között
i
i t
t t= +1−
∆ idő telik el. A fentiek leginkább a következő gondolatmenettel szemléltethetők. Vizsgáljuk azt az esetet, mikor egy MD vagy DMC szimuláció kezdeti konfigurációját adott gáz- vagy folyadékfázis részecskéinek valamilyen szabályos rácsszerű elrendeződése reprezentálja. Az MD szimuláció során az idő előrehaladtával a rendszer egy
rendezetlenebb, a kezdeti állapothoz képest valószínűbb állapotba jut. Az MD szimulációban a részecskék a rájuk ható erők irányában mozdulnak el, tehát a rendszer idővel egy energetikailag kedvezőbb állapotba jut. Az ugyanebből a kezdeti konfigurációból indított DMC szimuláció a rendszer állapotait szintén az energiaminimum felé (tehát mindig valószínűbb állapotok elérésének irányába) hajtja. Ekkor ugyan egyetlen részecskeelmozdítási lépés nem feleltethető meg egy molekuláris dinamikai időlépéshez tartozó részecskeelmozdításnak, de mivel a vizsgált rendszer evolúciója a két számításban átlagosan hasonló, a számítandó időfüggő tulajdonságok átlagai (pl. diffúziós állandó) összehasonlíthatók.
A DMC módszerrel, az MD módszertől eltérően, az eltelt idő mértékegysége nem a megszokott időegységekben (pl.: ps), hanem MC szimulációs lépésekben (MCS) van kifejezve. Ha az MC szimulációban csak a részecskék transzlációs léptetéseivel generáljuk a konfigurációkat, akkor az eltelt idő a részecskeelmozdítási kísérletek számával adható meg. Ilyen esetben a DMC módszer algoritmusa megegyezik a hagyományos kanonikus sokaságon végzett MC szimuláció algoritmusával (ld. 1.3. fejezet). Eszerint véletlenszerűen kiválasztunk egy részecskét, majd eredeti pozíciójából véletlen távolságba elmozdítjuk. Az elmozdítás mértékét az rmax parameter határozza meg (emlékezzünk, hogy rmax annak a kockának az oldalhossza, amelyen belül az aktuálisan kiválasztott i-edik részecskét véletlenszerűen elmozgatjuk). Az így kialakult konfigurációt - a mikroszkopikus reverzibilitás elégséges kritériuma alapján alkotott szabály szerint - elfogadjuk, ha a megkísérelt mozdítás után a rendszer potenciális energiája csökkent, és exp(–∆U/(kBT)) valószínűséggel fogadjuk el, ha nőtt.
A DMC módszer alkalmazásakor lehetőség van a molekuláris dinamikába nehezen implementálható speciális technikák (részecskeeltávolítás és -beillesztés; identitáscsere-eljárás; „parallel tempering”, stb.) beépítésére, melyekkel nagyobb időlépések ölelhetőek fel. Hátrányuk azonban, hogy minden egyes szimulációs lépést illetve eljárást ugyanarra az MCS egységre kell vonatkoztatni (a különböző technikák lépéseinek más-más az „abszolút időigénye”). Noha a helyes időszükségletek megtalálása összetett feladat, az irodalomban léteznek módszerek, melyeket már sikeresen alkalmaztak DMC számítások teljesítményének növelésére [45, 46].
A DMC szimulációban eltelt idő és ennek megfelelően a számított időfüggő (dinamikus) tulajdonságok nem függetlenek a különböző szimulációs lépések paramétereinek értékeitől (pl.: a részecskeelmozdítás nagyságának maximuma, a komponensek kiválasztásának gyakorisága). Egykomponensű rendszerekben, ha csak a részecskék transzlációja az egyetlen szimulációs lépés, akkor az rmax szimulációs paraméter határozza meg a transzlációs lépések átlagos nagyságát adott időegységre vonatkoztatva.
Többkomponensű rendszerek esetén rmax értéke illetve a komponensek kiválasztásának gyakorisága az egyes komponensek mért időfüggő tulajdonságait eltérően befolyásolja: két különböző komponens részecskéjének egységnyi időre eső átlagos elmozdulásának aránya változhat az alkalmazott rmax függvényében, azaz az egyik részecske mozgékonyabb lesz a másikhoz képest. Elvben a szimulációs paraméterek megválasztásának garantálnia kell, hogy adott rendszerre a számított dinamikus tulajdonságok DMC szimulációban és MD szimulációban megegyezzenek. Ha az időskálák különbözősége miatt az összehasonlítás direkt módon nem lehetséges, akkor az egyezés szükségszerűen az időfüggő tulajdonságok arányának egyezését kell, hogy jelentse. Két különböző DMC számítás eredményeként kapott időfüggő tulajdonságok abszolút értékei, továbbá egy DMC számításon belül az egyes komponensek dinamikus tulajdonságai akkor hasonlíthatóak össze, ha e tulajdonságok aránya megegyezik az MD referenciaszimulációk megfelelő tulajdonságainak arányával (ilyen például a mozgékonyságok aránya). A helyes paraméterértékek megtalálásának egy lehetséges módja az MD szimulációk segítségével történő kalibráció [37], melyet minden egyes esetben el kell végezni, ha a rendszer valamely tulajdonsága (sűrűsége, összetétele, nyomása) megváltozik.
A dinamikus Monte Carlo módszer néhány előnye és hátránya a molekuláris dinamikával szemben az alábbi pontokban foglalható össze.
A DMC módszer előnyei:
• merevgömbi potenciálok és merev falak kezelhetősége egyszerűbb;
• a molekuláris dinamikába nehezen beépíthető speciális technikákkal nagyobb időlépések ölelhetőek fel;
• a szimulációs paraméterek megválasztásával olyan hatások is figyelembe vehetők, amelyekre az alkalmazott modellek keretein belül esetleg nincs lehetőség (rmax
értékét oly módon is megválaszthatjuk, hogy például implicit módon kezelt oldószer használatakor is figyelembe vehessük az oldószer molekuláinak direkt hatását).
A DMC módszer hátrányai:
• az eltelt idő mértékegysége MC szimulációs lépés (MCS);
• minden különböző típusú szimulációs lépést és eljárást ugyanarra az MCS egységre kell vonatkoztatni;
• az előzőekből következően a szimulációs paramétereket hangolni (kalibrálni) kell.
1.7.1. Ellenőrzés és kalibráció elegyekre
Többkomponensű rendszer esetén a részecskék egységnyi idő alatt megtett elmozdulásainak arányára a következő összefüggésnek kell teljesülnie [37]:
MD
ahol MSD az átlagos négyzetes elmozdulás, és az α illetve β két különböző komponenst jelöl. ∆t idő alatt az átlagos négyzetes elmozdulás a következő összefüggéssel számítható [37]: részecskéjének pozíciója t0 időpillanatban. A megfelelő MSD-arány beállítható a részecsketranszlációt befolyásoló szimulációs paraméterek beállításával. Huitema és van der Eerden többkomponensű rendszerekben az MSD-arányok finomhangolását a transzlációs lépésekre vonatkozó komponenskiválasztási gyakoriságok, illetve rmax
értékének változtatásával valósította meg [37]. A módszer, noha egzakt, mindig járulékos MD számítást (referenciaszámítást) igényel.
1.7.2. Saját szerkezeti függvény
Az irodalomban létezik egy független módszer, melyet elegyek esetén az rmax
elegy főkomponensének saját szerkezeti függvényéből állapította meg rmax optimális értékét [42]. Az α komponens részecskéinek relaxációs viselkedését jellemző saját szerkezeti függvény írja le a részecskepozíciók változása és az eltelt idő közti viszonyt („self-intermediate scattering function”) [47]:
( )
A fenti kifejezésben ∆t jelöli az eltelt időt, ri(t0) az α komponens i-edik részecskéjének pozíciója t0 időpillanatban, szimulációs átlagot jelent, q pedig a szórási változó, melynek nagysága meghatározható a parciális szerkezeti függvényből (más néven: parciális struktúrafaktor vagy Faber-Ziman-féle szerkezeti függvény) [47]: