A munk´ ara fogott v´ eletlen ( Monte-Carlo-m´ odszer) Cserti J´ozsef E¨ otv¨ os Lor´and Tudom´anyegyetem, Komplex Rendszerek Fizik´aja Tansz´ek, H-1117 Budapest, P´azm´any P´eter s´et´any 1/A. (2003. szeptember 7.)
A p´ecsi Zipernowsky K´ aroly szakk¨ oz´episkolai tan´ arom, Balog J´ ozsef tisztelet´ere. ´ I. BEVEZETES
Mindennapi ´elet¨ unket gyakran befoly´asolj´ ak a v´eletlen esem´enyek. J´ol tudjuk, hogy a j´at´ekkaszin´okban a v´eletlen alapvet˝o szerepet j´atszik. Sz´amos v´eletlen jelens´eget figyelhet¨ unk meg a k¨ornyezet¨ unkben is (p´eld´aul a hegy´ere ´all´ıtott ceruza d˝ol´esi ir´anya). A term´eszetben is sz´amtalan v´eletlen jelens´eget ismer¨ unk. Egy tart´alyban l´ev˝ o g´azatomok v´eletlenszer˝ uen mozognak. Az atommagok boml´asa is v´eletlen folyamat. A v´eletlen seg´ıts´eg´evel k¨ozel´ıt˝ oleg meghat´arozhatjuk a π ´ert´ek´et. Dobjunk rizsszemeket (v´eletlenszer˝ uen) egy a oldal´ u n´egyzetlapra, amelybe egy a ´ atm´er˝oj˝ u k¨ort is berajzoltunk! V´egezz¨ unk N sz´am´ u k´ıs´erletet (csak azokat a dob´asokat tekintj¨ uk, melyekn´el a rizsszemek nem esnek ki a n´egyzetb˝ol), ´es sz´amoljuk meg h´any esik a k¨orbe! Jel¨olj¨ uk ezek sz´am´ at Nk -val! Ekkor az Nk /N ar´any nagy sz´am´ u k´ıs´erlet (N 1) eset´en j´o 2 k¨ozel´ıt´essel megegyezik a k¨or ´es a n´egyzet ter¨ ulet´enek ar´any´aval, azaz Nk /N = (a/2) π/a2 = π/4. ´Igy a π ´ert´eket π ≈ 4Nk /N alapj´an sz´am´ıthatjuk ki. Term´eszetesen ez a m´odszer nem adja meg teljesen pontosan a π-t. Min´el t¨obb k´ıs´erletet v´egz¨ unk, ann´al pontosabb eredm´enyt kaphatunk, felt´eve, hogy a rizsszemek val´oban v´eletlen¨ ul esnek a n´egyzetre. A fenti k´ıs´erletet nem sz¨ uks´eges a val´ os´ agban elv´egezni. A sz´am´ıt´og´eppel egyszer˝ u programot ´ırhatunk. Sz¨ uks´eg¨ unk van egy j´o v´eletlensz´ am gener´atorra. Ma m´ar sz´amos program l´etezik, mely a [0, 1] intervallumon egyenletes eloszl´asban gener´al egy v´eletlen sz´amot. Gener´aljunk egym´as ut´an kett˝ot, ´es jel¨olj¨ uk ezeket x-szel ill. y-nal! E k´et sz´amp´ arhoz egy pontot rendelhet¨ unk a koordin´atarendszer els˝o negyed´eben (rizsszem helye a dob´as ut´an). Ha a pont t´avols´ ag´ ara igaz, hogy x2 + y 2 < 1, akkor a pont az egys´egsugar´ u k¨or¨on bel¨ ul van. Tegy¨ uk fel, hogy a fenti algoritmust N -szer elv´egezve Nk sz´am´ u esetben esik a pont a k¨or¨on bel¨ ul. Hasonl´oan a rizsszemek eset´eben most is a 4Nk /N ar´any k¨ozel´ıti π ´ert´ek´et. Az al´abbi t´abl´azatban egyre n¨ovekv˝ o sz´am´ u k´ıs´erlet sor´an kapott π ´ert´eket ´es a hib´at t¨ untett¨ uk fel (a pontos ´ert´ek 9 tizedesjegyre π = 3.141592654).
N 10 102 103 104 105 106 107
π relat´ıv hiba % -ban 3.6 14.6 3.16 0.6 3.108 1.1 3.127 0.5 3.135 0.2 3.141 0.02 3.14155 0.001
L´athat´o, hogy N n¨ovel´es´evel egyre pontosabb ´ert´eket kapunk π-re. M´ar n´eh´any sz´azezer k´ıs´erlet is el´eg π-nek k´et tizedesjegy pontos kisz´am´ıt´ as´ ara. A mai sz´am´ıt´og´epekkel ak´ar 107 sz´am´ u k´ıs´erlet is egy percen bel¨ ul ´ elv´egezhet˝o. Erdemes kipr´ob´ alni! Teljesen v´eletlen jelens´eget felhaszn´alva egy j´ol meghat´arozott mennyis´eg ´ert´ek´et siker¨ ult k¨ozel´ıteni. A fenti m´ odszert tov´abb lehet fejleszteni ´es ´ıgy rendk´ıv¨ ul bonyolult feladatok megold´as´ara haszn´alhatjuk a v´eletlent. Az elj´ar´ast a Monte-Carloban tal´alhat´ o nevezetes kaszin´okra utalva, Monte-Carlo-m´odszernek h´ıvj´ak, ´es kiterjedten alkalmazz´ak mind a matematik´aban, mind a fizik´aban. A Monte-Carlo elnevez´est Metropolis ´es Ulam haszn´alt´ak el˝osz¨or egy 1949-es cikk¨ ukben arra utalva, hogy a m´odszerhez sz¨ uks´eges v´eletlen sz´amokat ak´ar egy j´at´ekkaszin´o j´at´ekeredm´enyeib˝ ol is vehetn´enk. A gyakorlatban viszont a v´eletlen sz´amokat a sz´am´ıt´og´epek maguk ´all´ıtj´ak el˝o. A m´odszert m´ar a sz´azad elej´en is haszn´alta n´eh´any statisztikus, de a Monte-Carlo-m´odszer csak akkor indult igaz´an fejl˝od´esnek, amikor Neumann, Ulam ´es Fermi atommagreakci´okra vonatkoz´o bonyolult matematikai probl´em´ ak sz´am´ıt´ og´eppel t¨ort´en˝ o k¨ozel´ıt˝o megold´as´ara haszn´alt´ak. 1
Sok esetben a feladatokat csak k¨ozel´ıt´esek alkalmaz´as´aval lehet megoldani. Szerencs´ere legt¨obbsz¨or nincs is sz¨ uks´eg¨ unk nagyon pontos ´ert´ekekre. Ilyenkor gyakran igen hat´ekonynak bizonyul a Monte-Carlo-m´odszer. A tov´ abbiakban n´eh´any p´eld´ at szeretn´enk bemutatni a m´odszer alkalmaz´as´ara a matematik´aban ´es a fizik´aban. Igyekezt¨ unk olyan probl´em´ akat v´alasztani, melyeket a mai sz´am´ıt´og´epes lehet˝os´egek mellett k¨oz´episkolai szinten vizsg´alhatunk. ´ II. MATEMATIKAI ALKALMAZAS
Gyakori feladat egy g¨orbe alatti ter¨ ulet meghat´aroz´asa. Az 1. ´abr´an l´athat´o f (x) f¨ uggv´enynek az [a, b] intervallumon a g¨orbe alatti A ter¨ ulet´et u ´gy hat´arozhatjuk meg, hogy az [a, b] szakaszt felosztjuk N egyenl˝o PN ∆x = (b − a)/N intervallumra ´es a ter¨ uletet az u ´n. t´egl´any¨osszeggel k¨ozel´ıtj¨ uk: A ≈ i=1 f (xi )∆x, ahol xi az i-edik intervallum k¨ozepe. Az egyszer˝ us´eg kedv´e´ert feltessz¨ uk, hogy a f¨ uggv´eny pozit´ıv az [a, b] szakaszon, ´es legyen a f¨ uggv´eny legnagyobb ´ert´eke M ezen a szakaszon. y M
f(x)
a
x
b
FIG. 1. Egy f (x) f¨ uggv´enynek a g¨ orbe alatti ter¨ ulete az [a, b] intervallumon.
A t´egl´any¨osszeg ann´al pontosabb, min´el nagyobb N . Ez az egyik legismertebb (´es legyeszer˝ ubb) m´odszer egy f¨ uggv´eny g¨orbe alatti ter¨ ulet´enek meghat´aroz´as´ara. A 2. ´abr´an l´athat´o f (x) = sin2 x1 f¨ uggv´eny azonban nagyon gyorsan oszcill´al az orig´o k¨or¨ ul, ´es ´ıgy a ter¨ ulet kiel´eg´ıt˝ o pontoss´ag´ u kisz´am´ıt´as´ahoz nagyon nagyra kellene v´alasztani N ´ert´ek´et. 1
1
0.9
0.9
0.8
0.8
0.7
0.7
0.6
0.6
0.5
0.5
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0 -1
-0.5
FIG. 2. Az f (x) = sin
2
1 x
0 0
0.5
1
-0.4
-0.2
0
0.2
0.4
f¨ uggv´eny a [−1, 1] (bal oldali ´ abra) ill. a [−0.5, 0.5] intervallumon (jobb oldali ´ abra).
A gyorsan v´altoz´ o f¨ uggv´enyek g¨orbe alatti ter¨ ulet´et a Monte-Carlo-m´odszer seg´ıts´eg´evel hat´ekonyabban becs¨ ulhetj¨ uk meg. Ism´et ’rizsszemek’ dob´al´ as´ at alkalmazzuk. Gener´aljunk egy x v´eletlen sz´amot (´altal´aban az egyes programnyelvekben tal´alhat´ o olyan v´eletlensz´am gener´ator, mely a [0, 1] intervallumon egyenletes eloszl´asban ad egy v´eletlen sz´amot)! Ekkor az x → a + x(b − a) transzform´aci´oval a kapott v´eletlen sz´am az [a, b] intervallumba ker¨ ul. Gener´aljunk egy m´asik y v´eletlen sz´amot! Az y → yM transzform´aci´o ut´an az y v´eletlen sz´am az [0, M ] intervallumba esik. Tekints¨ uk a k´et sz´amot egy pont (x, y) koordin´at´ainak! Ez a pont az 1. ´abr´an l´athat´o (b − a) ´es M hossz´ us´ ag´ u t´eglalapon bel¨ ul tal´alhat´o. Ha viszont a pont (x, y) koordin´at´aira teljes¨ ul az y < f (x) egyenl˝ otlens´eg, akkor a pont (rizsszem) az f (x) g¨orbe al´a esik. Ism´etelj¨ uk meg a fenti algoritmust N¨ossz alkalommal ´es sz´amoljuk meg h´anyszor esik az adott pont a g¨orbe al´a. Jel¨olj¨ uk ezeknek az esem´enyeknek a sz´am´ at Nbe -vel! Elegend˝oen sok k´ıs´erlet ut´an azt v´arjuk, hogy a f¨ uggv´enynek a g¨orbe alatti A ter¨ ulet´enek ´es az M (b − a) ter¨ ulet˝ u t´eglalapnak az ar´anya j´o k¨ozel´ıt´essel megegyezik az Nbe /N¨ ar´ a nnyal ´ e s ´ ıgy ossz
2
A ≈ M (b − a)
Nbe . N¨ossz
Alkalmazzuk a fenti Monte-Carlo-m´ odszert a 2. ´abr´an l´athat´o f¨ uggv´enyre! Kisz´am´ıtottuk a g¨orbe alatti ter¨ uletet az [a, b] = [0, b] intervallumon u ´gy, hogy b ´ert´ek´et 0 ´es 1 k¨oz¨ott v´altoztatjuk. A k¨ ul¨onb¨oz˝o sz´am´ u k´ıs´erlet sor´an kapott eredm´enyek a 3. ´abr´ an l´athat´ok. A g¨orb´ek ter¨ ulet´et ´altal´aban fels˝obb matematikai m´ odszerrel, az u ´n. integr´ alsz´ am´ıt´ assal lehet kisz´am´ıtani. Ezt a m´odszert tekinthetj¨ uk az egzakt, azaz pontos sz´am´ıt´asnak. Az ´abr´ akon berajzoltuk az ´ıgy kapott eredm´enyeket is a Monte-Carlo-m´odszer hat´ekonys´ag´anak illusztr´al´asa ´erdek´eben. 0.19
0.6
0.6
0.18
0.17
0.4
0.4
0.2
0.2
0
0
0.2 2
1 x
0.4
0 0.6
0.8
1
0.16 0.4
0
0.41
0.42
0.2
0.43
0.4
0.44
0.45
0.6
0.8
1
f¨ uggv´enynek a g¨ orbe alatti ter¨ ulete N¨ abra) ´es N¨ FIG. 3. Az f (x) = sin ossz = 100 (bal oldali ´ ossz = 10000 (jobb oldali ´ abra) k´ıs´erlet sor´ an. Az egzakt eredm´enyt a szaggatott g¨ orbe mutatja.
L´athat´o, hogy N¨ossz = 100 k´ıs´erlet eset´en a Monte-Carlo-m´odszer m´eg el´eg nagy hib´aval k¨ozel´ıti a pontos eredm´enyt. Ugyanakkor N¨ossz = 10000 k´ıs´erlet sor´an m´ar nagyon j´o egyez´est kapunk az integr´alsz´am´ıt´asb´ol nyert pontos eredm´ennyel. A tov´abbiakban sz¨ uks´eg¨ unk lesz a 4. ´abr´ an l´athat´o f (x) = sin2 (sin (x)) ’kalapszer˝ u’ f¨ uggv´enynek a g¨orbe alatti ter¨ ulet´ere a [0, π] intervallumon. 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0
0.5
1
1.5
2
2.5
3
FIG. 4. ’A kalapf¨ uggv´eny’. A [0, π] intervallumon a f¨ uggv´eny maximuma: M = sin2 (1) = 0.7081.
Ez a ter¨ ulet a Monte-Carlo-m´ odszer alapj´an A ≈ 1.2194, melyet N¨ossz = 222 = 4194304 sz´am´ u k´ıs´erlet elv´egz´ese ut´an kaptunk, h´arom tizedesjegy pontosan egyezik az integr´alsz´am´ıt´assal kapott A = 1.2191 pontos eredm´ennyel. A Monte-Carlo-m´ odszert rendk´ıv¨ ul sokf´ele ter¨ uleten haszn´alj´ak a fizikusok is. Cikk¨ unk II. r´esz´eben a v´eletlenen alapul´o sz´am´ıt´ asi m´odszerek fizikai alkalmaz´asaib´ol mutatunk be n´eh´anyat.
3
´ III. FIZIKAI ALKALMAZASOK
Cikk¨ unk els˝o r´esz´eben bemutattuk a Monte-Carlo-m´odszer alapjait ´es n´eh´any egyszer˝ u matematikai alkalmaz´as´at. Ezt az elj´ar´ ast – t¨obbek k¨oz¨ ott – a fizikusok is nagyon sok ter¨ uleten haszn´alj´ak. A sz´amos fizikai alkalmaz´as k¨oz¨ ul a tov´ abbiakban k´et p´eld´ at fogunk bemutatni.
A. S´ ulypont meghat´ aroz´ asa
A k¨oz´episkol´as f¨ uggv´enyt´ abl´ azatban megtal´alhat´o a s´ ulypont meghat´aroz´as´anak m´odja, ´es n´eh´any speci´alis alakzatra z´art alakban megadhat´o k´epleteket is fel´ırhatunk. (P´eld´aul egy r sugar´ u, homog´en t¨omegeloszl´as´ u f´elk¨ or alak´ u lap s´ ulypontja a k¨oz´eppontt´ ol 4r/(3π) t´avol van.) Bonyolultabb alakzatokra m´ar igen neh´ez kisz´am´ıtani a s´ ulypont hely´et, az ´altal´ aban csak integr´alsz´am´ıt´assal hat´arozhat´o meg. Tekints¨ uk p´eld´aul az I. r´esz 4. ´ abr´ aj´ an l´ athat´ o, kalapf¨ uggv´enynek nevezett alakzat. Az egyszer˝ us´eg kedv´e´ert most csak olyan s´ıkbeli alakzatokat vizsg´alunk, amelyeknek van szimmetriatengelye. A kalapf¨ uggv´eny szimmetriatengelye az [a, b] = [0, π] intervallumon az x = π/2 f¨ ugg˝oleges egyenes, ´ıgy nyilv´anval´o, hogy a s´ ulypont v´ızszintes koordin´at´ aja π/2. Azonban a f¨ ugg˝ oleges koordin´at´at azonban m´ar nem lehet elemi u ´ton kisz´am´ıtani. Az ´altal´anos m´odszer szerint felosztjuk az alakzatot sok kicsi r´eszre, melyek t¨omege mi , ´es elegend˝oen finom feloszt´as eset´en a s´ ulypont f¨ ugg˝ oleges hely´et a k¨ovetkez˝o k´eplettel sz´am´ıthatjuk ki P yi mi ys = Pi . i mi A nevez˝o ar´anyos az alakzat ter¨ ulet´evel. Cikk¨ unk I. r´esz´eben Monte-Carlo-m´odszer seg´ıts´eg´evel m´ar kisz´am´ıtottuk a kalapf¨ uggv´eny ter¨ ulet´et: A ≈ 1,219. A fenti k´eplet sz´aml´al´oj´at is kisz´am´ıthatjuk a MonteCarlo-m´odszerrel. Gener´aljunk egy x ´es egy y v´eletlen sz´amot a [0, 1] intervallumban. Az x → x(b − a) = xπ ´es y → M y transzform´aci´ oval az (x, y) koordin´at´aj´ u pont a π ´es M oldalhossz´ us´ag´ u t´eglalapba fog esni (M a kalapf¨ uggv´eny maximuma: M = sin2 (1) = 0,7081). Ekkor a k¨ovetkez˝o igen egyszer˝ u algoritmussal sz´am´ıthatjuk ki a fenti k´eplet sz´aml´ al´ oj´ at: if (y < f(x) ) then S = S + y N_be = N_be + 1 endif Az algoritmust N¨ossz alkalommal v´egrehajtva a fenti k´eplet sz´aml´al´oj´at az SM (b − a)/N¨ossz ´ert´ekkel k¨ozel´ıthetj¨ uk. Itt Nbe azoknak az (x, y) pontoknak a sz´am´at jel¨oli, amelyek az f (x) al´a esnek. Cikk¨ unk I. r´esz´eben l´ attuk, hogy a g¨orbe alatti ter¨ uletet k¨ozel´ıt˝ oleg az A = M (b − a)Nbe /N¨ossz kifejez´es alapj´an hat´arozhatjuk meg. ´Igy a s´ ulypont f¨ ugg˝ oleges koordin´at´aj´ at az ys =
S Nbe
k´epletb˝ol sz´am´ıthatjuk ki. Ne felejts¨ uk el az S ´es Nbe v´altoz´okat null´azni a ciklus elej´en! Term´eszetesen a fenti ´ programr´eszletet a haszn´alt programnyelv szab´alyai szerint esetleg ´at kell ´ırni. Erdemes megjegyezni, hogy a s´ ulypontnak a Monte-Carlo-m´ odszerrel t¨ort´en˝ o kisz´am´ıt´as´ahoz nem volt sz¨ uks´eg a g¨orbe alatti ter¨ uletre. Integr´alsz´am´ıt´assal a s´ ulypontra a pontosnak mondhat´o ys = 0,274 975 ´ert´eket kapjuk. Az 4. a s´ ulypontnak k¨ ul¨ onb¨oz˝o sz´am´ u k´ıs´erletekb˝ ol kapott ´ert´ekeit ´abr´ azoltuk, ´es felt¨ untett¨ uk az egzakt ´ert´eket is (szaggatott vonal).
4
0.29
0.285
ys 0.28
0.275
0.27
0
2
4
6
8
lg N
FIG. 5. A kalapf¨ uggv´eny s´ ulypontja Monte-Carlo-m´ odszer alapj´ an. A szaggatott vonal az egzakt ´ert´eket jelzi.
L´athat´o, hogy N¨ossz = 106 k´ıs´erlet eset´en m´ar nagyon j´ol egyezik a ’m´er´es¨ unk’ a pontos ´ert´ekkel (0,02%-os az elt´er´es).
B. Ising-modell
Egyes anyagok m´agneses tulajdons´agokkal rendelkeznek. M´ar az ´okori g¨or¨og¨ok ismert´ek azokat a Magnezia k¨orny´ek´en tal´alhat´ o ’k¨oveket’, melyek vonzz´ak a vasat. Eg´eszen a XX. sz´azadig, a kvantummechanika t¨orv´enyeinek felismer´es´eig nem siker¨ ult kiel´eg´ıt˝ oen magyar´azni az anyagok m´agneses viselked´es´et. Az Ising ´altal kidolgozott modell fontos el˝orel´ep´est jelentett ebben a kutat´asban. Ezen egyszer˝ us´ıtett modell szerint az anyagban l´ev˝o kis m´agnesek (mai ismereteink szerint az atomok spinjei) k´et ´allapot valamelyik´eben lehetnek: vagy azonos, vagy ellent´etes ir´anyban ´allnak a k¨ uls˝o m´agneses t´erhez k´epest. K´epzelj¨ uk el, hogy a kis m´agnesek egy n´egyzetr´acs pontjaiban helyezkednek el. Feltessz¨ uk, hogy a n´egyzetr´ acs i-edik r´acspontj´aban egy kis m´agnes (spin) van, melynek ´ert´eke Si = ±1, att´ol f¨ ugg˝oen, hogy a spin azonos vagy ellent´etes ir´any´ u a k¨ uls˝ o H t´errel: 1, ha ↑ Si = −1, ha ↓ . A spinek egyik lehets´eges be´all´ asa l´athat´ o a 6. ´abr´an.
H
FIG. 6. A spinrendszerek Ising-modellje. A k¨ uls˝ o H m´ agneses t´erhez k´epest a spinek vagy azonos, vagy ellent´etes ir´ anyban ´ allnak.
Az egyes spinek k¨olcs¨ onhatnak egym´assal. A modellben csak a k¨ozvetlen szomsz´edok (´ un. els˝oszomsz´edok) k¨ ozti k¨olcs¨onhat´ast vessz¨ uk figyelembe, a t´avolabbi m´agnesek egym´asra hat´asait elhanyagoljuk. ´Igy minden r´acspontban a spin a 4 els˝oszomsz´edj´ aval hat csak k¨olcs¨on. A k¨ uls˝o H m´agneses t´er igyekszik a spineket a t´er ir´ any´aba ´all´ıtani. A fenti k¨olcs¨ onhat´ asok alapj´an a rendszer teljes E energi´aj´at a k¨ovetkez˝o alakban ´ırhatjuk: X X E = −J Si Sj − H Si ,
i
5
ahol J a szomsz´edos spinek k¨ozti k¨olcs¨ onhat´ as er˝oss´eg´ere jellemz˝o sz´am, melyet csatol´asi ´alland´onak szoktak nevezni. Az els˝o tagban az ¨osszegz´est csak az els˝oszomsz´edokra vessz¨ uk, erre utal a <, > jel¨ol´es. A rendszer M m´agnesezetts´ege az ¨osszes spinv´altoz´o ¨osszege (vagy ezzel az ¨osszeggel ar´anyos): X M = Si . i
H˝ uts¨ uk le a rendszert (gondolatban) abszol´ ut z´erus fokra, ´es fokozatosan kapcsoljuk ki a k¨ uls˝o m´agneses teret is! Ekkor a rendszer alap´allapotba (a legalacsonyabb energi´aj´ u, m´assz´oval energetikailag ’legkedvez˝obb’ ´allapotba) ker¨ ul. Ha a J csatol´ asi ´alland´ o pozit´ıv, akkor energetikailag az a kedvez˝o ´allapot, ha minden spin azonos ir´ anyba ´all (hiszen ebben az esetben a szomsz´edos spinek szorzata +1), teh´at egy rendezett ´allapot alakul ki. Alap´allapotban ´es J > 0 eset´en a rendszer m´agnezetts´ege maxim´alis. Ezt az anyag er˝osen m´agneses, ferrom´ agneses ´allapot´anak nevezik. A h˝om´ers´eklet n¨ovel´es´evel egyes spinek ’´atbillennek’ az ellent´etes ir´any´ u ´allapotba, ´es ´ıgy lecs¨okken a rendszer m´ agnesezetts´ege. Egy bizonyos Tc h˝om´ers´eklet´ert´ek ut´an – k¨ uls˝o m´agneses t´er n´elk¨ ul – ´atlagosan a spinek fele felfel´e, m´asik fele lefel´e ´all: a rendszer M m´agnezetts´ege z´erus lesz. Ezt a h˝om´ers´ekletet kritikus h˝om´ers´ekletnek nevezik, ´es Pierre Curie francia fizikus kutat´asai nyom´an Curie-h˝om´ers´ekletnek is h´ıvj´ak. A h˝om´ers´eklet teljesen ’sz´etzil´alja’ a rendezett ´allapotot m´eg z´erus H t´er mellett is. Tc alatt a rendszer m´agnesezetts´ege v´eges ´ert´ek˝ u, m´ıg Tc f¨ol¨ott z´erus. Ez a jelens´eg (a sok tekintetben hasonl´o halmaz´allapotv´altoz´asokkal egy¨ utt) a f´azis´atalakul´asok csal´adj´ ahoz tartozik. Ezen egyszer˝ u modell alapj´an is meg´erthetj¨ uk fizikailag, hogy mi´ert veszti el egy m´agnes a m´agnesess´eg´et, ha t˝ uzbe dobjuk. A m´agnesezetts´eg h˝om´ers´ekletf¨ ugg´es´et (v´azlatosan) a 7. ´abra mutatja. M
H>0
H=0
T H<0
FIG. 7. Ising-modell m´ agnesezetts´eg´enek h˝ om´ers´eklett˝ ol val´ o f¨ ugg´ese k¨ ul¨ onb¨ oz˝ o k¨ uls˝ o m´ agneses terek eset´en.
´ Erdemes megjegyezni, hogy ha a J csatol´asi ´alland´o negat´ıv, akkor alap´allapotban ´es H = 0 eset´en a szomsz´edos spinek ellent´etes ir´anyban ´allnak (ekkor k´et szomsz´edos spin szorzata −1), m´as sz´oval egy ’sakkt´ablaszer˝ u’ ´allapot alakul ki. Az ilyen (a term´eszetben is megfigyelhet˝o) rendszereket antiferrom´ agneseknek nevezik. Tov´abbiakban egy a Monte-Carlo-m´ odszeren alapul´o algoritmust mutatunk be, mellyel meghat´arozhatjuk a 7. ´abr´an l´athat´o m´agnesezetts´egnek h˝om´ers´eklett˝ol ´es a k¨ uls˝o m´agneses t´ert˝ol val´o f¨ ugg´es´et. Az algoritmust el˝osz¨or N. Metropolis g¨or¨ og sz´armaz´ as´ u amerikai matematikus alkalmazta spinrendszerekre a m´ ult sz´azad k¨ozep´en, ´es az´ota Metropolis-algoritmusnak is nevezik1 . A spineket gondolatban egy n´egyzetr´acsra helyezz¨ uk a 6. ´abr´anak megfelel˝oen, ´es peri´odikus hat´arfelt´etelt alkalmazunk. Ez azt jelenti, hogy a teljes n´egyzetr´acsot gondolatban megism´etelj¨ uk a n´egyzet oldal´eleinek ir´any´aban, vagy ami ezzel egyen´ert´ek˝ u: a n´egyzetr´acsot egy t´oruszra (aut´ogumi alak´ u fel¨ uletre) rajzoljuk. Ezzel a ’tr¨ ukkel’ azt ´erj¨ uk el, hogy az eredeti n´egyzet sz´elein l´ev˝o spineknek is n´egy els˝oszomsz´edja legyen, teh´at egyen´ert´ek˝ uek legyenek a r´acs belsej´eben tal´alhat´o t´arsaikkal. Induljunk ki a spinek valamilyen kezdeti elrendez˝od´es´eb˝ol (konfigur´aci´oj´ab´ol). V´alaszthatjuk pl. a z´erus h˝om´ers´ekletnek megfelel˝o rendezett ´allapotot, az alap´allapotot kiindul´asi konfigur´aci´ok´ent. Ezut´an a k¨ovetkez˝o l´ep´eseket ism´etelj¨ uk: uen kiv´alasztunk egy spint. Kisz´am´ıtjuk, hogy az ´atford´ıt´as´aval mennyit v´altozna meg a 1. V´eletlenszer˝ rendszer energi´aja. Legyen ez az energiav´ altoz´as ∆E. 6
2. Kisz´amoljuk annak val´ osz´ın˝ us´eg´et, hogy T h˝om´ers´ekleten a kiv´alasztott spin ´atfordul az ellent´etes ir´any´ u ´allapotba. Ez a val´ osz´ın˝ us´eg (az u ´n. ´atmeneti val´osz´ın˝ us´eg) az al´abbi k´eplettel adhat´o meg: ( W =
1,
ha ∆E < 0
e− kT , ha ∆E > 0, ∆E
ahol k a Boltzmann-´alland´ ot jel¨oli. 3. Gener´alunk egy r v´eletlen sz´amot a [0, 1] intervallumon! Ha azt tapasztaljuk, hogy r < W , akkor megford´ıtjuk a kiszemelt spint, k¨ ul¨ onben nem v´altoztatunk az ´all´as´an. 4. Visszat´er¨ unk a 1. ponthoz. Bel´ athat´o, hogy a fenti algoritmussal elegend˝o sz´am´ u ciklus ut´an a m´agnesezetts´eg be´all valamekkora hM i egyens´ ulyi ´ert´ekre, melyet a T h˝ om´ers´eklet ´es a k¨ uls˝o H m´agneses t´er hat´aroz meg. Bizonyos sz´am´ u iter´aci´ora (id˝ ore) mindig sz¨ uks´eg van, hogy el´erj¨ uk az egyens´ ulyi ´allapotot. Ez az ´allapot Trel id˝o (az u ´n. relax´aci´os id˝o) ut´an ´all be. Az id˝ot Monte-Carlo-l´ep´esekben (MCS) szokt´ak sz´amolni. Egy MCS az az id˝o, ami alatt ´atlagosan minden spint egyszer kiv´alasztunk az iter´aci´ ok sor´an. A 8. ´abr´an a m´agnesezetts´eg tipikus id˝obeli v´altoz´asa (a ciklusok sz´am´at´ol val´ o f¨ ugg´ese) l´athat´ o.
<M>
Tcorr
Trel
MCS
FIG. 8. A m´ agnesezetts´eg id˝ obeli v´ altoz´ asa az iter´ aci´ ok sor´ an MCS egys´egekben.
A relax´aci´os id˝ot tapasztalat u ´tj´ an ´allap´ıthatjuk meg. Ha az M m´agnesezetts´eg m´ar nem v´altozik jelent˝osen, csak kiss´e ingadozik, akkor el´ert¨ uk az egyens´ ulyi ´allapotot. Ekkor az ´atlagos hM i m´agnesezetts´eget u ´gy hat´arozhatjuk meg, hogy bizonyos gyakoris´ aggal kisz´amoljuk M ´ert´ek´et, ´es vessz¨ uk ezek ´atlag´at. A h˝om´ers´eklet v´altoztat´as´aval vizsg´alhatjuk az ´atlagos m´agnesezetts´eg viselked´es´et, ’m´er´eseinkkel’ meghat´arozhatjuk annak h˝om´ers´ekletf¨ ugg´es´et. A fenti k´etdimenzi´ os v´egtelen Ising-modell viselked´es´et els˝ok´ent L. Onsager 2 sz´amolta ki analitikusan (k¨ozel´ıt´esmentesen, z´art alakban) k¨ uls˝ o m´agneses t´er n´elk¨ uli esetben. A rendszer kritikus h˝om´ers´eklete kTc = 2,27 J. Az egzakt megold´as ismeret´eben lehet˝os´eg van a Monte-Carlo-m´odszer pontoss´ag´anak tanulm´anyoz´as´ara. Jelenleg nem ismeretes analitikus megold´as h´arom dimenzi´oban (vagyis t´erbeli r´acsokra), Monte-Carlo-m´odszerrel azonban m´ar alaposan tanulm´anyozt´ak ezeket a rendszereket is. A Monte-Carlo-m´ odszer alkalmaz´as´ an´ al figyelni kell a lehets´eges hibaforr´asokra. P´eld´aul a rendszer v´eges m´erete miatt a kritikus h˝om´ers´eklet s´ıkbeli r´acsokn´al nem egyezik meg az ismert elm´eleti ´ert´ekkel. Ilyenkor vizsg´alni kell, hogy hogyan f¨ ugg a v´eges m´eret˝ u mint´an kapott kritikus h˝om´ers´eklet a minta m´eret´et˝ol, ´es ebb˝ol tudunk k¨ovetkeztetni a v´egtelen rendszerre. A m´asik gyakori gond a ’j´o’ v´eletlensz´ amgener´ator. El˝ofordul, hogy a gener´ator egym´ast k¨ovet˝o sz´amai nem teljesen f¨ uggetlenek egym´ast´ ol; ekkor biztosan hib´as eredm´enyeket kapunk. Rendk´ıv¨ ul fontos teh´at a v´eletlensz´amgener´ator tesztel´ese. V´egezet¨ ul megeml´ıtj¨ uk, hogy viszonylag egyszer˝ uen programozhat´o probl´em´ak tal´alhat´ok3 cikkben.
K¨osz¨onetemet szeretn´em kifejezni koll´eg´ aimnak, Kert´esz J´ anosnak ´es Vicsek Tam´ asnak, akik megismertettek a Monte-Carlo-m´odszerrel, Geszti Tam´ asnak ´es Pollner P´eternek a hasznos tan´acsaik´ert, valamint egyetemi
7
hallgat´oimnak, Koltai J´ anosnak, Pafka Szil´ ardnak ´es Sz˝ ucs J´ ozsefnek, akik seg´ıts´egemre voltak ebben a munk´aban.
1
N. Metropolis, A.W. Rosenbluth, M.N. Rosenbluth, A.H. Teller, E. Teller: J. Chem. Phys. 21, 1087 (1953). L. Onsager: Phys. Rev. 65, 117 (1944). 3 J. Kert´esz, J. Cserti, J. Sz´ep: Monte Carlo simulation programs for microcomputer, Eur. J. Phys. 6, 232 (1985). 2
8
