モンテカルロ法による2次元イジング模型のシミュレーション
理論系の研究室に配属されたので早速Fortran90で簡単なプログラムを書いてみた。
で表されるからJ>0として考えると、隣接するスピンが平行の場合に安定となる。
アルゴリズムを簡単に解説すると
①すべてのスピンをランダムに決める
②温度を決める
③スピンを1つランダムに選び、それをフリップさせる
④フリップによるエネルギー変化が0以下の場合そのまま、0より大の場合はボルツマン分布から求めた温度に依存した確率でそのままか元に戻すか決める
⑤必要な物理量の値(アンサンブル)を調べる
⑥③に戻ってこの操作を繰り返す
⑦アンサンブルから知りたい物理量(エネルギー、磁化、比熱、キュムラントなど)を求める
⑧②に戻ってこの操作を繰り返す
というような感じ。
スピン数がLxL個の系で周期的境界条件を使った。
物理量を求める際はアンサンブルの平均を取るんだけど、その時③~⑥の繰り返しのうち初めの数万回(系のスケールによる)は平衡状態になる前の状態なのでデータとしては捨てる。
簡単なプログラムなんだけど物性実験で初めて触って、それ以来やってなかったから結構時間かかった。
これからXY模型やHeisenberg模型への拡張はスラスラやれたかな。
結果を載せていくと、
磁化の転移が綺麗に出たのがスピンの数が60x60の場合で、そのグラフは
こんな感じ。
低温ではスピンか完全に一方向に揃っていて、高温では熱運動によりスピンが任意の方向を向けるから全体としては磁化0になる。
転移点ではスピンがどちらでも向けるようになるので磁化が不の領域にも振れてるっぽい。
比熱はエネルギーの分散から定義できるから
より
横軸が温度で縦軸が比熱。
右上の数字が系のスケール。
比熱は転移点で発散して、系のスケールにはよらないように見えるんだけど、ちょっとずつずれてるらしい。
磁化もそうなんだけど、サンプルを多く取りすぎたせいでやや見にくい(醜い)。
温度はそこまで細かく区切る必要は無いらしく、XYmodelを作るときはもっと見やすくする。
ところどころ飛んでるところは気にしない。
最後にキュムラントは
で定義されるから、グラフは
となり、無限系では各線が交わる温度で相転移が起こることがわかる。
この方法は複数の有限系のサンプルから無限系の転移点を求める方法らしく、正直理論はよく分からない(笑)
温度変化に対してちょっとずつ磁化が下がっていくのでなく、ある温度を境に急激に0になるのって不思議。
イジングモデルに関しては統計力学をやってれば式の上で転移点を理解することは出来るんだけど、こういう単純な繰り返しでそれを再現できるっていうのが何とも面白かった。
追々XY模型についても書いていきたい。
おしまい。