Monte Carlo módszerek a statisztikus fizikában. Az Ising modell.
8. előadás
Démon algoritmus az ideális gázra időátlag Å
fizikai mennyiségek átlagértéke Æ sokaságátlag E, V, N pl. molekuláris dinamika Monte Carlo módszerek Módszer: Creutz dolgozta ki a háló mértékelmélet alapján egy extra szabadságfokot adunk a rendszerhez Æ démon a démon energiát ad/vesz el, hogy megváltoztassa a rendszer dinamikai változóit ha a rendszer energiája csökken, a vátozás energiáját a démon kapja ha a rendszer energiája nő, a változást a démon adja, ha van neki elegendő megj. a démonnak nem lehet az energiája negatív Algoritmus: 1. Választunk véletlenszerűen egy részecskét és próbából megváltoztatjuk sebességét. 2. Kiszámoljuk a változás miatt bekövetkező ΔE energiaváltozást. 3. Ha ΔE ≤ 0, a rendszer a | ΔE | energiát a démonnak adja, akinek az energiája Ed = Ed – ΔE lesz és a változást elfogadjuk. 4. Ha ΔE > 0, és a démonnak elegendő energiája van (Ed ≥ ΔE), akkor a démon adja a szükséges energiát (Ed = Ed – ΔE) és a változást elfogadjuk. Ha a démonnak nincselég energiája, a próba konfigurációt elvetjük és a rendszer állapota változatlan marad.
Megj. Egy Monte Carlo lépés annyi próbálkozást jelent, ahány részecske van a rendszerben.
MC a statisztikus fizikában vizsgált mennyiség eloszlásfüggvény
állapottér eleme
teljes állapottér
rendszer Hamiltonfüggvénye
Fontossági mintavételezést használunk az integrálok számítására.
Jól megválasztott P(x)-el:
Metropolis MC szimuláció kanonikus sokaságban (Boltzmann statisztika) Egy mennyiség átlagértéke T hőmérsékleten:
f (H [x ]) ∝ e
−
H [x ] k BT
fontossági mintavételezés KIINDULÁS: a rendszer kezdetben a 0. állapotban van ε0, A0 lenullázzuk a következőket: N – MC próbálkozások száma Aössz – a tanulmányozott mennyiség összege ALGORITMUS: véletlenszerűen megváltoztatjuk a rendszer állapotát i. -ből i+1. -be ει+1, Ai+1 növeljük N-et összehasonlítjuk az ει+1 és az ει energiákat HA ει+1< ει
Aössz = Aössz + Ai+1
HA ει+1> ει generálunk egy r véletlen számot az egyenletes eloszlású [0, 1]-ből HA exp{-(εi+1 – εi)/(kBT)} > r
Aössz = Aössz + Ai+1
KÜLÖNBEN
Aössz = Aössz + Ai
Figyelem! A szimuláció elején a tranziens időszakban a mennyiséget nem összegezzük (nem vesszük bele az átlagérték számításba).
Ising modell anyagok mágnesességének tanulmányozása a mágnesezettség eredete: z töltött részecskék mozgása zárt pályákon z töltött részecskék saját tengely körüli forgása egyszerű klasszikus közelítés 2D mágnes
Ampere törvény mágneses momentum Ising spin
egy spin csakis a négy legközelebbi szomszédjával hat kölcsön (~r-3) kölcsönhatási energia:
J > 0 ferromágneses rendszer J < 0 antiferromágneses rendszer mágnesezettség:
Ferromágnesség, paramágnesség és a Curie hőmérséklet H = 0 két fázis lehetséges T függvényében
Curie hőmérséklet
ferromágneses paramágneses ferromágnes-paramágnes fázisátalakulás d>2 (magasabb dimenziókban) másodrendű fázisátalakulás, rend-rendezetlenség fázisátalakulás 2D 3D
Tc= 2.269 J/kB egzaktul megoldható Tc= 4.44 J/kB szimulációs eredmény négyzetrácson
Monte Carlo szimuláció FELADAT: véletlenszerűen generálni egy konfiguráció halmazt (minta) a lehetséges állapotok teréből (konfigurációs tér) egy bizonyos valószínűségi eloszlás szerint, és kiszámítani a megfigyelt mennyiség átlagértékét (pl. mágnesezettség) a minták alapján
kofiguráció, vagy minta: mágnesezettség átlagolása Boltzmann statisztika alapján
konfigurációs tér: 20x20-as rács:
óriási még kevés spin esetében is Ns = 400
lehetséges konfigurációk száma: 2400 = 2.58 x 10120 1millió konfiguráció/s 8.8 x 10103 év !!!
Monte Carlo szimuláció valószínűség-eloszlás, súlyfüggvény: fontossági mintavételezés, a Boltzmann-faktor alapján
Monte Carlo átlagérték: generálunk N független konfigurációt a Boltzmann-faktor alapján
a mágnesezettség és az energia átlagértékei
fajhő és mágneses szuszceptibilitás átlagértékei
Monte Carlo szimuláció - átlagértékek
Monte Carlo szimuláció - átlagértékek
Monte Carlo szimuláció - átlagértékek
0
Monte Carlo szimuláció Hogyan generáljunk mintákat a Boltzmann faktor alapján? pl. METROPOLIS algoritmus kiindulunk egy adott konfigurációból kiválasztunk egy spint:
si
próbából megfordítjuk:
si,próba= -si
kiszámítjuk a rendszer energiaváltozását: ΔE = E(s1,...,si,próba,...sN ) – E(s1,...,si,...sN ) generálunk egy r egyenletes eloszlású véletlenszámot HA
w = e
−
ΔE k BT
>r
si = si,próba
ismételjük ezeket a lépéseket MEGJEGYZÉSEK: egy MC lépés
N próbálkozás
akármilyen kezdeti állapotból kiindulva a rendszert néhány MC lépésen át hagyjuk termalizálódni, az átlagértékeket csak ezután kezdjük számolni erős végesméret effektus
periodikus határfeltételek használata
Monte Carlo szimuláció
ALGORITMUS: rögzítünk egy adott T hőmérsékletet rögzítünk egy kezdeti spinkonfigurációt (pl. véletlenszerűen) több MC lépést végzünk a rendszer termalizálására nagyszámú MC lépésre kiszámoljuk az E, E2, M, M2 átlagértékeket rögzítjük az <m(T)>, <E(T)>,
és <χ(T)> átlagértékeket megváltoztatjuk a T hőmérsékletet és megismételjük az algoritmust ábrázoljuk az átlagértékeket a hőmérséklet függvényében
Monte Carlo szimuláció
w=e
a w exponenciális faktor számítása:
− ΔE
k BT
exp függvény számolás minden lépésben lassú
számítsuk ki előre a H és T ismeretében!!! 4 szomszéd összege
s jj j szomszédok
x
, jy
sj
x
1, j y
megszorozva a spin értékével
sj
x
1, j y
si
sj
x
sj
, jy 1
+4 +2 s jx, jy 1 0 -2 -4 ugyanazok az értékek
j szomszédok
ha H≠0
+2H vagy -2H 10 lehetséges érték w-re
double w[17][3]; ... for (int i = -8; i <= 8; i += 4){ w[i+8][0] = exp( - (i * J + 2 * H) / T ); w[i+8][2] = exp( - (i * J - 2 * H) / T ); }
w 8 2 si
s j 1 si j szomszédok
Monte Carlo szimuláció Periodikus határfeltételek: int Lx, Ly; // spinek száma x és y irányban ... // véletlenszerűen kiválasztjuk az i,j -vel jellemzett spint int iPrev = i == 0 ? Lx-1 : i-1; int iNext = i == Lx-1 ? 0 : i+1; int jPrev = j == 0 ? Ly-1 : j-1; int jNext = j == Ly-1 ? 0 : j+1; int sumNeighbors = s[iPrev][j] + s[iNext][j] + s[i][jPrev] + s[i][jNext]; ...