ALPS 2.1 Tutorials:MC-07 Phase Transition/ja
|
Languages: |
English • 日本語 (ja) • 繁體中文 (zh-tw) • 简体中文 (zh) |
Contents
最初に
このチュートリアルの目的は、有名な2次元イジングモデルを例にとり、有限サイズのシミュレーションから二次相転移を再現する方法を学ぶことです。本当の相転移現象は有限システムでは求められませんが、有限サイズスケーリング法を用いて、相転移現象を明確に予測することが可能です。
相転移現象は2次元イジングモデルでは正確に解けることが知られています。このチュートリアルでは、臨界点並びに臨界指数の計算手順を説明します。このチュートリアルの計算は、ある程度の時間を要しますので、効率的に勧めてください。
2番目のシミュレーションのパラメータファイルはparm7bです。次のコマンドで実行開始してください。
parameter2xml parm7b spinmc --Tmin 10 parm7b.in.xml &
Pythonでは、tutorial7b.py
import pyalps
import matplotlib.pyplot as plt
import pyalps.pyplot
import numpy as np
import pyalps.fit_wrapper as fw
#prepare the input parameters
parms = []
for l in [32,48,64]:
for t in [2.24, 2.25, 2.26, 2.27, 2.28, 2.29, 2.30, 2.31, 2.32, 2.33, 2.34, 2.35]:
parms.append(
{
'LATTICE' : "square lattice",
'T' : t,
'J' : 1 ,
'THERMALIZATION' : 5000,
'SWEEPS' : 150000,
'UPDATE' : "cluster",
'MODEL' : "Ising",
'L' : l
}
)
#write the input file and run the simulation
input_file = pyalps.writeInputFiles('parm7b',parms)
pyalps.runApplication('spinmc',input_file,Tmin=5)
Vistrailを用いる場合は、mc-07-tutorial.vtを参照してください。
相転移のおおよその特定
最初に、おおよその臨界範囲を特定するために、小さい系で温度による探索をおこないます。使用するパラメータファイルはparm7aです。
parameter2xml parm7a spinmc --Tmin 5 parm7a.in.xml
Pythonを用いる場合、スクリプトはtutorial7a.pyです。
import pyalps
import matplotlib.pyplot as plt
import pyalps.pyplot
#prepare the input parameters
parms = []
for l in [4,8,16]:
for t in [5.0,4.5,4.0,3.5,3.0,2.9,2.8,2.7]:
parms.append(
{
'LATTICE' : "square lattice",
'T' : t,
'J' : 1 ,
'THERMALIZATION' : 1000,
'SWEEPS' : 400000,
'UPDATE' : "cluster",
'MODEL' : "Ising",
'L' : l
}
)
for t in [2.6, 2.5, 2.4, 2.3, 2.2, 2.1, 2.0, 1.9, 1.8, 1.7, 1.6, 1.5, 1.2]:
parms.append(
{
'LATTICE' : "square lattice",
'T' : t,
'J' : 1 ,
'THERMALIZATION' : 1000,
'SWEEPS' : 40000,
'UPDATE' : "cluster",
'MODEL' : "Ising",
'L' : l
}
)
#write the input file and run the simulation
input_file = pyalps.writeInputFiles('parm7a',parms)
pyalps.runApplication('spinmc',input_file,Tmin=5)
磁化、帯磁率、比熱
低音度におけるイジングモデルを表現するパラメータのオーダーは、サイトmあたりの磁化です。有限システムでは、
は対称性により0ですが、磁化のサイト毎の平均絶対値に注目すべきです。
次の記述はPythonシェルによる、実行、データ解析、計算結果をロードする一連の手順です。
pyalps.evaluateSpinMC(pyalps.getResultFiles(prefix='parm7a')) #load the observables and collect them as function of temperature T data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm7a'),['|Magnetization|', 'Connected Susceptibility', 'Specific Heat', 'Binder Cumulant', 'Binder Cumulant U2']) magnetization_abs = pyalps.collectXY(data,x='T',y='|Magnetization|',foreach=['L']) connected_susc = pyalps.collectXY(data,x='T',y='Connected Susceptibility',foreach=['L']) spec_heat = pyalps.collectXY(data,x='T',y='Specific Heat',foreach=['L']) binder_u4 = pyalps.collectXY(data,x='T',y='Binder Cumulant',foreach=['L']) binder_u2 = pyalps.collectXY(data,x='T',y='Binder Cumulant U2',foreach=['L'])
に対する|Magnetization|のプロットをおこないます。
plt.figure()
pyalps.plot.plot(magnetization_abs)
plt.xlabel('Temperature $T$')
plt.ylabel('Magnetization $|m|$')
plt.title('2D Ising model')
磁化は０から低Tでの飽和値まで上昇することが明らかにわかります。システムサイズを大きくなるにつれて、上昇が鋭くなるにつれて、これらが起きる温度はそれほど明白ではありません。 次に、磁化の変動を見てみます。
ここで、磁化率
を見てみます。
は全スピン数です。
plt.figure()
pyalps.plot.plot(connected_susc)
plt.xlabel('Temperature $T$')
plt.ylabel('Connected Susceptibility $\chi_c$')
plt.title('2D Ising model')
周辺でピークが観測され、最も変動が強い箇所です。ピークはシステムサイズにともなって発散する傾向にあります。詳細は後述しますが、この発散は臨界指数によって特徴付けられます。
比熱の振る舞い、 エネルギーの摂動
（
はサイトの内部エネルギー）をみることで、おおまかな臨界温度を決定します。
plt.figure()
pyalps.plot.plot(spec_heat)
plt.xlabel('Temperature $T$')
plt.ylabel('Specific Heat $c_v$')
plt.title('2D Ising model')
Binderキュムラント
シミュレーションをおこない、常に容易に曲線の最大値（帯磁率や比熱のピーク）を求められるとは限りません。また、前途のプロットにあるように、ピーク温度はシステムサイズによって変動することもあります。これらは後述する有限サイズスケーリング法によって対処が可能です。
相転移点を特定する他の効率的な方法にキュムラントを用いる方法があります。ここでBinderキュムラントと呼ばれる比率
を取り上げます。低温度で、この比率は、低温度では一定ですが、高温度ではパラメータのガウシアンゆらぎ
は３となります。Binderキュムラントの特徴は、その値が、システムサイズから独立していることにあります。すなわち、温度に対するBinderキュムラントのプロットは、Tcの決定に非常に有益です。
plt.figure()
pyalps.plot.plot(binder_u4)
plt.xlabel('Temperature $T$')
plt.ylabel('Binder Cumulant U4 $g$')
plt.title('2D Ising model')
このプロットから、相転移は![T_c \in [2.2,2.3]](/mediawiki/images/math/2/7/2/272a421f40fbe231c52020c273135d66.png)
の狭い範囲であると結論づけられます。また、別の課題として、キュムラント
の場合も考えて見てください。
相転移場所の特定、臨界指数、collapse plotsの評価
密な温度グリッド点や大規模システムを用いることで、より正確な相転移の性質をシミュレーションすることができます。チュートリアルで用いるファイルはparm7bです。結果ファイルから観測量を取り出し解析します。
pyalps.evaluateSpinMC(pyalps.getResultFiles(prefix='parm7b')) #load the observables and collect them as function of temperature T data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm7b'),['|Magnetization|', 'Connected Susceptibility', 'Specific Heat', 'Binder Cumulant', 'Binder Cumulant U2']) magnetization_abs = pyalps.collectXY(data,x='T',y='|Magnetization|',foreach=['L']) connected_susc = pyalps.collectXY(data,x='T',y='Connected Susceptibility',foreach=['L']) spec_heat = pyalps.collectXY(data,x='T',y='Specific Heat',foreach=['L']) binder_u4 = pyalps.collectXY(data,x='T',y='Binder Cumulant',foreach=['L']) binder_u2 = pyalps.collectXY(data,x='T',y='Binder Cumulant U2',foreach=['L'])
Collapsing data
Bider cumulant交差の計算をおこないます。
plt.figure()
pyalps.plot.plot(binder_u4)
plt.xlabel('Temperature $T$')
plt.ylabel('Binder Cumulant U4 $g$')
plt.title('2D Ising model')
plt.show()
このプロットから
Hint : Try 
close to unity ...
#perform a data collapse of the Binder cumulant:
Tc=... #your estimate
a=... #your estimate
for d in binder_u4:
d.x -= Tc
d.x = d.x/Tc
l = d.props['L']g
d.x = d.x * pow(float(l),a)
plt.figure()
pyalps.plot.plot(binder_u4)
plt.xlabel('$L^a(T-T_c)/T_c$')
plt.ylabel('Binder Cumulant U4 $g$')
plt.title('2D Ising model')
plt.show()
現在、何をやっているかというと、臨界指数を得るため有名なテクニックで、スケーリング用のデータを調整します。この場合、相関距離の臨界指数を読み取ることが可能です。
です。
を決定するためにBinder cumulantを考慮することは有益です。
臨界指数
比熱と磁化率に関して、チュートリアルをおこないます。
#make a plot of the specific heat and connected susceptibility:
plt.figure()
pyalps.plot.plot(connected_susc)
plt.xlabel('Temperature $T$')
plt.ylabel('Connected Susceptibility $\chi_c$')
plt.title('2D Ising model')
plt.figure()
pyalps.plot.plot(spec_heat)
plt.xlabel('Temperature $T$')
plt.ylabel('Specific Heat $c_v$')
plt.title('2D Ising model')
plt.show()
バインダーキュムラント交差で得られた
と多少異なるTでの磁化率と比熱は関連性があります。この僅かな差異は小さな有限サイズ上で予想されます。ピークの位置で定義される臨界温度
はシステムサイズでドリフとすると予想され、次の用に定義されます。
。定数
はの値によって、異なる値をとります。
次に、臨界指数について考えます。、またはバインダーキュムラントから得られたのいずれかで、磁化率や比熱の値から読み取ることができます。
and 
のようなスケーリングを想定するでしょう。例えば、異なる曲線での磁化率の最大値を見つけ、
の関係を用いています。th>\eta</math>はどのように求めることができるか分かりますか？
#make a fit of the connected susceptibility as a function of L:
cs_mean=[]
for q in connected_susc:
cs_mean.append(np.array([d.mean for d in q.y]))
peak_cs = pyalps.DataSet()
peak_cs.props = pyalps.dict_intersect([q.props for q in connected_susc])
peak_cs.y = np.array([np.max(q) for q in cs_mean])
peak_cs.x = np.array([q.props['L'] for q in connected_susc])
sel = np.argsort(peak_cs.x)
peak_cs.y = peak_cs.y[sel]
peak_cs.x = peak_cs.x[sel]
pars = [fw.Parameter(1), fw.Parameter(1)]
f = lambda self, x, pars: pars[0]()*np.power(x,pars[1]())
fw.fit(None, f, pars, peak_cs.y, peak_cs.x)
prefactor = pars[0].get()
gamma_nu = pars[1].get()
plt.figure()
plt.plot(peak_cs.x, f(None, peak_cs.x, pars))
pyalps.plot.plot(peak_cs)
plt.xlabel('System Size $L$')
plt.ylabel('Connected Susceptibility $\chi_c(T_c)$')
plt.title('2D Ising model, $\gamma$ is %.4s' % gamma_nu)
plt.show()
比熱に関しても同様の手順で求めます。
は何ですか？ は、磁化率とは対照的に、一般に比熱が連続的な遷移で発散する必要がないことを意味し、正または負であることに注意してください。
#make a fit of the specific heat as a function of L:
sh_mean=[]
for q in spec_heat:
sh_mean.append(np.array([d.mean for d in q.y]))
peak_sh = pyalps.DataSet()
peak_sh.props = pyalps.dict_intersect([q.props for q in spec_heat])
peak_sh.y = np.array([np.max(q) for q in sh_mean])
peak_sh.x = np.array([q.props['L'] for q in spec_heat])
sel = np.argsort(peak_sh.x)
peak_sh.y = peak_sh.y[sel]
peak_sh.x = peak_sh.x[sel]
pars = [fw.Parameter(1), fw.Parameter(1)]
f = lambda self, x, pars: pars[0]()*np.power(x,pars[1]())
fw.fit(None, f, pars, peak_sh.y, peak_sh.x)
prefactor = pars[0].get()
alpha_nu = pars[1].get()
plt.figure()
plt.plot(peak_sh.x, f(None, peak_sh.x, pars))
pyalps.plot.plot(peak_cs)
plt.xlabel('System Size $L$')
plt.ylabel('Specific Heat $c_v(T_c)$')
plt.title(r'2D Ising model, $\alpha$ is %.4s' % alpha_nu)
plt.show()
練習としては、
としてのスケールする絶対磁化に対して同様の解析をおこなうことができます。また、関連する磁化率、磁化のデータ解析をおこなってみてください。
hint: 対応するスケーリングフォームは、
、 
です。
#make a data collapse of the connected susceptibility as a function of (T-Tc)/Tc:
for d in connected_susc:
d.x -= Tc
d.x = d.x/Tc
l = d.props['L']
d.x = d.x * pow(float(l),a)
two_minus_eta=... #your estimate
for d in connected_susc:
l = d.props['L']
d.y = d.y/pow(float(l),two_eta)
plt.figure()
pyalps.plot.plot(connected_susc)
plt.xlabel(' $L^a(T-T_c)/T_c$')
plt.ylabel(r'$L^{\gamma/\nu}\chi_c$')
plt.title('2D Ising model')
plt.show()
#make a data collapse of the |magnetization| as a function of (T-Tc)/Tc
for d in magnetization_abs:
d.x -= Tc
d.x = d.x/Tc
l = d.props['L']
d.x = d.x * pow(float(l),a)
beta_over_nu=... #your estimate
for d in magnetization_abs:
l = d.props['L']
d.y = d.y / pow(float(l),-beta_over_nu)
plt.figure()
pyalps.plot.plot(magnetization_abs)
plt.xlabel(' $L^a(T-T_c)/T_c$')
plt.ylabel(r'$L^{-\beta/\nu} |m|$')
plt.title('2D Ising model')
plt.show()
臨界指数の正確な推定値
二次元イジングモデルの正確な臨界指数は、
、 
です。比熱はべき乗ではありませんが、対数的に発散します。それゆえ、比熱指数は０次元です。厳密解とどのように比較しましたか？
このチュートリアルでは、計算時間が膨大なるため、おこないませんでしたが、大規模なシステムサイズで計算をおこなうと、得られる臨界指数や温度の推定値の信頼性が大幅に増加します。また、それらは臨界温度の推定値の品質に依存します。臨界指数で、より高精度な求めかたは、臨界温度
の正確な値を用い、大規模サイズで計算をおこなうことです。
パラメータとして一つの温度しか無いので、を決める前回の方法は、あまり役に立ちません。前に述べたバインダーキュムラントのスケーリング手法は\nuを決める別の方法です。バインダーキュムラントをTに関して微分すると、で
の微分は
でスケールすることが分かります。計算で得られたデータから数値微分して求めることも出来ますし、また、モンテカルロシミュレーションで熱力学的平均値をとっても求めることが可能です。これらの計算には十分な統計量が必要であることに注意してください。計算されたデータで、指数は、
の関数として、温度で、この微分のべき乗近似によって決定することができます。
