ALPS 2.1 Tutorials:TEBD-02 kink/ja
|
Languages: |
English • 日本語 (ja) • 繁體中文 (zh-tw) • 简体中文 (zh) |
Contents
磁壁の展開
このチュートリアルでは、非平衡状態で準備したS=1/2スピン鎖の時間発展について学習します。選択した特殊な状態は、次のように、真ん中で上下スピンをもち、左右でそれぞれ上スピン、下スピンを持つような構造です。
この状態は初期状態としてINITIAL_STATEに'kink'と設定されています。いくつかの厳密解は一次元XXモデルに関して知られています。
XX modelの完全解
XXモデルでのkink初期状態の時間発展はJordan-WignerPhys. Rev. E 59, 4912 (1999)の自由フェルミオンへの相転移によって厳密解が解かれています。時間関数として、どんなサイトでも磁化の期待値はBessel関数の合計として表現できることが明らかとなりました。そして、磁化と初期磁壁から大きな距離は、変数
でスケーリングします。
は中心からの距離、
は時間です。
ここで、
は磁化で、は中心からの距離、
はオーダ
のBessel関数を表します。このチュートリアルの最初のパートでは、これらの2つの結果を演習します。
Pythonを用いた実行法
Pythonスクリプトを用いた実行方法です。tutorial2a.pyを参照してください。
import pyalps import matplotlib.pyplot as plt import pyalps.plot import numpy as np import copy import math import scipy.special
#prepare the input parameters
parms = [{
'L' : 50,
'MODEL' : 'spin',
'local_S' : 0.5,
'CONSERVED_QUANTUMNUMBERS' : 'Sz',
'Jxy' : 1,
'INITIAL_STATE' : 'kink',
'CHI_LIMIT' : 40,
'TRUNC_LIMIT' : 1E-12,
'NUM_THREADS' : 1,
'TAUS' : [20.0],
'POWS' : [0.0],
'GS' : ['H'],
'GIS' : [0.0],
'GFS' : [0.0],
'NUMSTEPS' : [500],
'STEPSFORSTORE' : [2]
}]
mathとscipy.specialモジュールは、この厳密解を比較するために必要な特別な関数を生成するために使用されます。POWSを0としていることに注意してください。まったくクエンチしないことに対応しています。GS、GIS、GFSの値は任意です。TAUS、NUMSTEPSはそれぞれ、全計算時間とタイムステップを与えます。入力ファイルを準備し、計算を実行し、出力を得て、最後にプロットをおこないます。
baseName='tutorial_2a' nmlname=pyalps.writeTEBDfiles(parms, baseName) res=pyalps.runTEBD(nmlname)
#Get the results of the simulation Data=pyalps.load.loadTimeEvolution(pyalps.getResultFiles(prefix='tutorial_2a'), measurements=['Local Magnetization'])
厳密解と比較するために、計算結果の後処理が必要です。 まず最初に、後処理データを保持するための空の配列を定義します。
#define a dataset numericalSolution to contain the numerical result numericalResult=[] #define a dataset exactSolution to contain the exact solution exactResult=[] #define a dataset scalingForm to contain the scaling form scalingForm=[]
時間データから厳密解を計算します。比較するために各サイトの磁化の計算された値を使用します。
#Compute the exact result M(n,t)=<S_n^z>=-(1/2)*sum_{i=1-n}^{n-1} j_i(t)^2, where
# j_i(t) is the Bessel function of order i and compare to the numerically obtained result
for q in Data:
syssize=q[0].props['L']
#Assign a label 'Distance' denoting the distance from the center n (only do the first two sites
#to avoid cluttering the plot)
for n in range(1,3):
#Create copies of the data for postprocessing
numericalCopy=copy.deepcopy(q)
exactCopy=copy.deepcopy(q)
numericalCopy[0].props['Distance']=n
numericalCopy[0].props['SIMID']='Numerical at n='+str(n)
exactCopy[0].props['Distance']=n
exactCopy[0].props['SIMID']='Exact at n='+str(n)
#compute the exact result of the manetization n sites from the center
loc=0.0
for i in range(1-n,n):
loc-=0.5*scipy.special.jn(i,q[0].props['Time'])*scipy.special.jn(i,q[0].props['Time'])
exactCopy[0].y=[loc]
#add to the the exact dataset
exactResult.extend(exactCopy)
#get the numerical result of the magnetization n sites from the center
numericalCopy[0].y=[q[0].y[syssize/2+n-1]]
#add to the the numerical dataset
numericalResult.extend(numericalCopy)
次に、正確なスケーリング関数を計算します。スケーリング変数の関数として磁化を計算し、厳密解と比較します。
#compute the scaling form
# \phi(n/t)=-(1/pi)*arcsin(n/t) that M(n,t) approaches as n->infinity and t->infinity
# and compare it with the numerically computed values of M(n/t)
for q in Data:
syssize=q[0].props['L']
#Assign a label 'Distance' denoting the distance from the center n (only do the first few sites
#to avoid cluttering the plot)
for n in range(0,5):
#Create a copy of the data for postprocessing
scalingCopy=copy.deepcopy(q)
scalingCopy[0].props['Distance']=n
#The first distance contains the exact scaling form \phi(n/t)=-(1/pi)*arcsin(n/t)
if n==0:
scalingCopy[0].props['Time']=1.0/scalingCopy[0].props['Time']
scalingCopy[0].y=[-(1.0/3.1415926)*math.asin(min(scalingCopy[0].props['Time'],1.0))]
scalingCopy[0].props['SIMID']='Exact'
#The other distances contain the numerical data as a function of the scaling variable M(n/t)
else:
scalingCopy[0].props['Time']=n/scalingCopy[0].props['Time']
scalingCopy[0].y=[scalingCopy[0].y[syssize/2+n-1] ]
scalingCopy[0].props['SIMID']='Numerical at n='+str(n)
#add to the scaling dataset
scalingForm.extend(scalingCopy)
最後に、厳密解、数値解を比較のためプロットします。
#Plot the numerical and exact magnetization for comparison
exactMag=pyalps.collectXY(exactResult, x='Time', y='Local Magnetization',foreach=['SIMID'])
for q in exactMag:
q.props['label']=q.props['SIMID']
numericalMag=pyalps.collectXY(numericalResult, x='Time', y='Local Magnetization',foreach=['SIMID'])
for q in numericalMag:
q.props['label']=q.props['SIMID']
plt.figure()
pyalps.plot.plot([exactMag, numericalMag])
plt.xlabel('Time $t$')
plt.ylabel('Magnetization')
plt.legend(loc='lower right')
plt.title('Magnetization vs. time')
#Plot the scaling form with the numerical data for comparison Scal=pyalps.collectXY(scalingForm, x='Time', y='Local Magnetization', foreach=['SIMID']) for q in Scal: q.props['label']=q.props['SIMID']
plt.figure()
pyalps.plot.plot(Scal)
plt.xlabel('Scaling variable $n/t$')
plt.ylabel('Magnetization$(n,t)$')
plt.legend()
plt.xlim(0,1.5)
plt.title('Magnetization scaling function; numerical and exact results')
plt.show()
視覚的には磁化が非常によく一致し、関連する制限されている正確なスケーリングフォームに近づくことがわかります。
Vistrailsを用いた実行法
Vistrailsを用いた実行は、tutorial1a.vtを開き、ワークフロー"tutorial2b"を選択してください。
TEBDの誤差解析：タイムステップ誤差
時間の関数としてTEBDシミュレーションの誤差を計算するために、厳密解を使用しています。最初に "無限小"時間ステップdtを変化の影響を調査します。
Pythonを用いた実行法
Pythonスクリプトを用いた実行方法です。tutorial2b.pyを参照してください。
import pyalps import matplotlib.pyplot as plt import pyalps.plot import numpy as np import math import scipy.special
#prepare the input parameters
parms=[]
count=0
for nsteps in [100, 250, 500, 750, 1000]:
count+=1
parms.append({
'L' : 50,
'MODEL' : 'spin',
'local_S' : 0.5,
'CONSERVED_QUANTUMNUMBERS' : 'Sz',
'Jxy' : 1,
'INITIAL_STATE' : 'kink',
'CHI_LIMIT' : 20,
'TRUNC_LIMIT' : 1E-12,
'NUM_THREADS' : 1,
'TAUS' : [20.0],
'POWS' : [0.0],
'GS' : ['H'],
'GIS' : [0.0],
'GFS' : [0.0],
'NUMSTEPS' : [nsteps],
'STEPSFORSTORE' : [int(math.floor(nsteps/100))],
'SIMID': count
})
NUMSTEPSパラメータを代えることで、暗黙的にタイムステップを変化させます。それはTAUが固定値なので、結果、タイムステップを変化させることになります。入力ファイルを準備し、計算を実行させ、データを収集します。
baseName='tutorial_2b_' nmlnameList=pyalps.writeTEBDfiles(parms, baseName) res=pyalps.runTEBD(nmlnameList)
#Get magnetization data Magdata=pyalps.load.loadTimeEvolution( pyalps.getResultFiles(prefix='tutorial_2b'), measurements=['Local Magnetization'])
時間データから正確な結果を計算します。そして、磁化の数値計算と厳密解を比較します。
#Postprocessing-get the exact result for comparison
for q in Magdata:
syssize=q[0].props['L']
#Get the exact result of M(1,t)=-(1/2)*(j_0(t)^2), where j_0(t) is the 0^{th} order
# bessel function and M(1,t) is the magnetization one site to the right of the chain center
loc=-0.5*scipy.special.jn(0,q[0].props['Time'])*scipy.special.jn(0,q[0].props['Time'])
#Get the difference between the computed and exact results
q[0].y=[abs(q[0].y[syssize/2+1-1]-loc)]
最後に磁化の誤差をプロットします。
#Plot the Error in the magnetization one site to the right of the chain center Mag=pyalps.collectXY(Magdata, x='Time', y='Local Magnetization', foreach=['SIMID']) for q in Mag: dt=round(q.props['TAUS']/q.props['NUMSTEPS'],3) q.props['label']='dt='+str(dt)
plt.figure()
pyalps.plot.plot(Mag)
plt.xlabel('Time $t$')
plt.yscale('log')
plt.ylabel('Magnetization Error')
plt.title('Error in the magnetization vs. time')
plt.legend(loc='lower left')
plt.show()
次のことが言えます。短い時間だと、指数関数のトロッター分解からエラーへの寄与を反映し、、誤差はおおまかにdt^2に比例しています。しかし長い時間では、最小のdtでのシミュレーションでは、大きなdtと比較して大きな誤差が発生しており、最終的に誤差は膨大なものになってしまします。次のセクションでさらに追求します。
Vistrailsを用いた実行法
Vistrailsを用いた実行は、tutorial1a.vtを開き、ワークフロー"tutorial2a"を選択してください。
TEBDの誤差解析：エンタングルメントカットオフ誤差
続いて、エンタングルメントカットオフ-パラメータ
の誤差への影響を調査します。
Pythonを用いた実行法
Pythonスクリプトを用いた実行方法です。tutorial2c.pyを参照してください。
import pyalps import matplotlib.pyplot as plt import pyalps.plot import math import scipy.special
#prepare the input parameters
parms=[]
count=0
for chi in [10, 20, 30, 40]:
count+=1
parms.append({
'L' : 50,
'MODEL' : 'spin',
'local_S' : 0.5,
'CONSERVED_QUANTUMNUMBERS' : 'Sz',
'Jxy' : 1,
'INITIAL_STATE' : 'kink',
'CHI_LIMIT' : chi,
'TRUNC_LIMIT' : 1E-12,
'NUM_THREADS' : 1,
'TAUS' : [20.0],
'POWS' : [0.0],
'GS' : ['H'],
'GIS' : [0.0],
'GFS' : [0.0],
'NUMSTEPS' : [500],
'STEPSFORSTORE' : [5],
'SIMID': count
})
入力ファイルを準備し、計算を実行し、出力を得て、最後にプロットをおこないます。
baseName='tutorial_2c_'
nmlnameList=pyalps.writeTEBDfiles(parms, baseName)
res=pyalps.runTEBD(nmlnameList)
#Get magnetization data
Magdata=pyalps.load.loadTimeEvolution( pyalps.getResultFiles(prefix='tutorial_2c'), measurements=['Local Magnetization'])
#Postprocessing-get the exact result for comparison
for q in Magdata:
syssize=q[0].props['L']
#Get the exact result of M(1,t)=-(1/2)*(j_0(t)^2), where j_0(t) is the 0^{th} order
# bessel function and M(1,t) is the magnetization one site to the right of the chain center
loc=-0.5*scipy.special.jn(0,q[0].props['Time'])*scipy.special.jn(0,q[0].props['Time'])
#Get the difference between the computed and exact results
q[0].y=[abs(q[0].y[syssize/2+1-1]-loc)]
最後に、磁化誤差のプロットをおこないます。
#Plot the Error in the magnetization one site to the right of the chain center
Mag=pyalps.collectXY(Magdata, x='Time', y='Local Magnetization', foreach=['SIMID'])
for q in Mag:
q.props['label']='$\chi$='+str(q.props['CHI_LIMIT'])
plt.figure()
pyalps.plot.plot(Mag)
plt.xlabel('Time $t$')
plt.yscale('log')
plt.ylabel('Magnetization Error')
plt.title('Error in the magnetization vs. time')
plt.legend(loc='lower left')
plt.show()
短い時間だと、指数関数のトロッター分解からエラーへの寄与を反映し、、誤差はおおまかにdt^2に比例しています。時間が増加すると、多くの発散誤差が発生します。最初のシミュレーション
では
あたりで、
では
で発散します。この時、非常に状態が複雑になっている時の時間発展状態を近似する行列積を求めるための手法が近似式を使用していることに起因します。この近似は、波動関数の繰り込みを含み、誤差が時間内にほぼ指数関数的に蓄積されていってしまいます。
この誤差の指数関数的爆発は、小さなdtのときでも発生しています。dtがより小さい時でも、対応できるような近似伝搬法を採用しなくてはなりません。それは打ち切り誤差の蓄積による指数関数的な爆発を意味します。よって、タイムステップの増加によって発生する誤差と、さらなるタイムステップを取ることによって発生する誤差との間のデリケートなバランスを調整する必要があります。すべての結果はとdtで収束のチェックを注意深くおこなうべきです。
Vistrailsを用いた実行法
Vistrailsを用いた実行は、tutorial1a.vtを開き、ワークフロー"tutorial2c"を選択してください。
XXZモデルの解法
磁化は、速度
で弾道敵に拡大していくことが厳密解からわかりました。XXモデルでは多くの特別な性質を持っているので、このような磁化の振る舞いが一般的なことかどうか興味が持たれます。このチュートリアルのパートでは、ハミルトニアンに
を加え、その影響を調査します。これはXXZモデルに相当します。この追加によって、平行構造のスピンは凍結スピンとなり、初期状態はハミルトニアンの厳密な固有状態となります。ハミルトニアンのXX項はスピンを反転しようとし、ピュアなXXモデルでみられたように磁化の波面を伝搬する役割をはたします。スピン輸送する系の特性の定量的な測定として、Phys. Rev. E 71 036102(2005)の
で定義された手法を元に磁化の総合的なフローを考慮します。
Pythonを用いた実行法
Pythonスクリプトを用いた実行方法です。tutorial2d.pyを参照してください。
import pyalps import matplotlib.pyplot as plt import pyalps.plot import math import scipy.special
#prepare the input parameters
parms=[]
count=0
for z in [0.0, 0.3, 0.9, 1.0, 1.1, 1.5]:
count+=1
parms.append({
'L' : 50,
'MODEL' : 'spin',
'local_S' : 0.5,
'CONSERVED_QUANTUMNUMBERS' : 'Sz',
'Jxy' : 1,
'Jz' : z,
'INITIAL_STATE' : 'kink',
'CHI_LIMIT' : 40,
'TRUNC_LIMIT' : 1E-12,
'NUM_THREADS' : 1,
'TAUS' : [20.0],
'POWS' : [0.0],
'GS' : ['H'],
'GIS' : [0.0],
'GFS' : [0.0],
'NUMSTEPS' : [500],
'STEPSFORSTORE' : [5],
'SIMID': count
})
Jz-結合のシミュレーションをおこないます。入力ファイルを準備し、計算を実行し、出力結果を得ます。
baseName='tutorial_2d' nmlnameList=pyalps.writeTEBDfiles(parms, baseName) res=pyalps.runTEBD(nmlnameList)
#Get magnetization data Magdata=pyalps.load.loadTimeEvolution( pyalps.getResultFiles(prefix='tutorial_2d'), measurements=['Local Magnetization'])
計算された磁化データから、上記で定義したように総合的な磁化を計算します。
#Compute the integrated magnetization across the center
for q in Magdata:
syssize=q[0].props['L']
#Compute the integrated flow of magnetization through the center \Delta M=\sum_{n>L/2}^{L} (<S_n^z(t)>+1/2)
#\Delta M= L/4
loc=0.5*(syssize/2)
#\Delta M-=<S_n^z(t)> from n=L/2 to L
q[0].y=[0.5*(syssize/2)+sum(q[0].y[syssize/2:syssize])]
最後に、Jzカップリング幅の総合磁化をプロットします。
#Plot the integrated magnetization Mag=pyalps.collectXY(Magdata, x='Time', y='Local Magnetization', foreach=['Jz'])
plt.figure()
pyalps.plot.plot(Mag)
plt.xlabel('Time $t$')
plt.ylabel('Integrated Magnetization $\Delta M(t)$')
plt.title('Integrated Magnetization vs. Time')
plt.legend(loc='upper left')
plt.show()
Jz>1では、総合磁化はおおよそリニアに増加していきます。その磁化の相転移は、XXの場合のように弾道的(ballistic)です。1のJzでは、統合磁化は最終的に飽和するといった定性的な振る舞いがみられます。
Vistrailsを用いた実行法
Vistrailsを用いた実行は、tutorial1a.vtを開き、ワークフロー"tutorial2d"を選択してください。
質疑
- Jz=1での総合磁化が異なる性質に変化する振る舞いをするポイントは、XXZモデルは臨界点から反強磁性相に転移しているということです。しかし、この相転移は、定常状態に影響与える低エネルギーでの先験的現象です。この低エネルギーでの変化がどのように、今回おこなった高エネルギーの初期状態のダイナミックスに影響を与えるのか、推測できますか？
