コンテンツにスキップ
チュートリアル
モンテカルロシミュレーション
MC-02 帯磁率

MC-02 帯磁率

帯磁率 χ=∂⟨M⟩/∂h∣h=0\chi = \partial\langle M \rangle / \partial h \big|_{h=0} は、小さな外部磁場に対して系がどれだけ強く応答するかを表します。 その温度依存性は背後にあるスピン相関を調べる鋭敏なプローブです。相互作用の弱いスピン系はキュリー則 χ∝1/T\chi \propto 1/T に従いますが、相互作用や量子効果は特徴的なずれを生み出します。

このチュートリアルでは 4 つの系——1 次元鎖および 2 本足梯子上の古典および量子ハイゼンベルク模型——について χ(T)\chi(T) を計算し、その結果を 1 枚の図に重ねて示します。 この比較によって 2 つの重要な対比が浮かび上がります。量子ゆらぎが古典的描像をどのように変えるか、そして格子形状を鎖から梯子に変えることが低温での振る舞いにどう影響するかです。

符号の規約。 古典の spinmc コードも量子の looper コードも、ハミルトニアン H=J∑⟨i,j⟩S⃗i⋅S⃗jH = J \sum \langle i,j \rangle \vec{S}_i \cdot \vec{S}_j を用います。ここで J<0J < 0 は強磁性的な配列を、J>0J > 0 は反強磁性的な配列を好みます。以下の古典シミュレーションでは J=−1J = -1（強磁性体）を、量子シミュレーションでは J=+1J = +1（反強磁性体）を用います。どちらの選択でも低温で増大する正の帯磁率が得られるため、古典的な振る舞いと量子的な振る舞いを素直に比較できます。

古典ハイゼンベルク模型

1 次元鎖

コマンドラインでの設定と実行

パラメータファイル parm2a は、さまざまな温度における 60 サイトの鎖上の古典強磁性ハイゼンベルク模型のシミュレーションを設定します。

LATTICE="chain lattice"
L=60
J=-1
THERMALIZATION=15000
SWEEPS=500000
UPDATE="cluster"
MODEL="Heisenberg"
{T=0.05;}
{T=0.1;}
{T=0.2;}
{T=0.3;}
{T=0.4;}
{T=0.5;}
{T=0.6;}
{T=0.7;}
{T=0.8;}
{T=0.9;}
{T=1.0;}
{T=1.25;}
{T=1.5;}
{T=1.75;}
{T=2.0;}

標準的なコマンド列でシミュレーションを実行します。

parameter2xml parm2a
spinmc --Tmin 10 --write-xml parm2a.in.xml

Python での設定と実行

スクリプト tutorial2a.py は同じシミュレーションを設定して実行します。parm2a と同じフォルダに置いてください。

import pyalps
import matplotlib.pyplot as plt
import pyalps.plot

parms = []
for t in [0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.25, 1.5, 1.75, 2.0]:
    parms.append(
        {
            'LATTICE'        : "chain lattice",
            'T'              : t,
            'J'              : -1,
            'THERMALIZATION' : 10000,
            'SWEEPS'         : 500000,
            'UPDATE'         : "cluster",
            'MODEL'          : "Heisenberg",
            'L'              : 60
        }
    )

input_file = pyalps.writeInputFiles('parm2a', parms)
pyalps.runApplication('spinmc', input_file, Tmin=5, writexml=True)

評価と作図

出力ファイルから帯磁率を読み込み、温度に対してプロットします。

import pyalps
import matplotlib.pyplot as plt
import pyalps.plot

data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm2a'), 'Susceptibility')
susceptibility = pyalps.collectXY(data, x='T', y='Susceptibility')

plt.figure()
pyalps.plot.plot(susceptibility)
plt.xlabel('Temperature $T/J$')
plt.ylabel('Susceptibility $\chi J$')
plt.ylim(0, 0.22)
plt.title('Classical Heisenberg chain')
plt.show()

2 本足梯子

梯子形状では、レッグ方向の結合 J0J_0 に加えて、ラング方向の 2 つ目の結合 J1J_1 が導入されます。 格子の変更と 2 つの結合を除けば、シミュレーションの設定は鎖の場合とまったく同じです。

コマンドラインでの設定と実行

parm2b をダウンロードし、同じフォルダに置いてください。

LATTICE="ladder"
L=60
J0=-1
J1=-1
THERMALIZATION=15000
SWEEPS=500000
UPDATE="cluster"
MODEL="Heisenberg"
{T=0.05;}
{T=0.1;}
{T=0.2;}
{T=0.3;}
{T=0.4;}
{T=0.5;}
{T=0.6;}
{T=0.7;}
{T=0.8;}
{T=0.9;}
{T=1.0;}
{T=1.25;}
{T=1.5;}
{T=1.75;}
{T=2.0;}

シミュレーションを実行します。

parameter2xml parm2b
spinmc --Tmin 10 --write-xml parm2b.in.xml

Python での設定と実行

スクリプト tutorial2b.py は tutorial2a.py のコピーで、変更点は 3 つだけです。接頭辞を parm2b に変え、LATTICE を "ladder" に設定し、J を J0 と J1（どちらも -1）に置き換えます。

量子ハイゼンベルク模型

量子シミュレーションでは、ALPS モデルライブラリを用いて S=1/2S = 1/2 の量子スピン模型を指定し、looper QMC コードでシミュレーションを実行します。 古典の場合と比べた主なパラメータの変更点は次のとおりです。

  • MODEL="Heisenberg" の代わりに MODEL="spin" と local_S=1/2
  • looper コードを選ぶための ALGORITHM="loop"
  • 反強磁性結合としての J=+1（上記の符号の規約を参照）

1 次元鎖

コマンドラインでの設定と実行

parm2c をダウンロードします。

LATTICE="chain lattice"
MODEL="spin"
local_S=1/2
L=60
J=1
THERMALIZATION=5000
SWEEPS=50000
ALGORITHM="loop"
{T=0.05;}
{T=0.1;}
{T=0.2;}
{T=0.3;}
{T=0.4;}
{T=0.5;}
{T=0.6;}
{T=0.7;}
{T=0.75;}
{T=0.8;}
{T=0.9;}
{T=1.0;}
{T=1.25;}
{T=1.5;}
{T=1.75;}
{T=2.0;}

spinmc ではなく loop アプリケーションを使って変換・実行します。

parameter2xml parm2c
loop parm2c.in.xml

Python での設定と実行

スクリプト tutorial2c.py は tutorial2a.py を量子のパラメータに合わせて書き換え、spinmc の代わりに loop を呼び出します。

input_file = pyalps.writeInputFiles('parm2c', parms)
pyalps.runApplication('loop', input_file)

2 本足梯子

2 本足梯子は擬 1 次元系の簡単な例です。結合 J0 を持つ 2 本の平行なハイゼンベルク鎖（レッグ）が、結合 J1 のラングボンドで結ばれています。

     J0      J0
  o-------o-------o--- ...     <- レッグ（結合 J0）
  |  J1   |  J1   |             <- ラング（結合 J1）
  o-------o-------o--- ...     <- レッグ（結合 J0）

比 α=J1/J0\alpha = J_1/J_0 は 2 本のレッグがどれだけ強く結合しているかを決め、切り離された鎖（α=0\alpha=0）と、ラングが強く結合したほぼ独立な二量体の系（α≫1\alpha\gg 1）との間を連続的に補間します。 ギャップのない鎖とは異なり、2 本足の反強磁性ハイゼンベルク梯子はスピンギャップを持ちます。その基底状態はラング一重項の積であり、χ\chi はギャップエネルギー以下では指数関数的に抑制され、低温では急激にゼロへ向かいます。

コマンドラインでの設定と実行

parm2d をダウンロードします。

LATTICE="ladder"
MODEL="spin"
local_S=1/2
L=60
J0=1
J1=1
THERMALIZATION=5000
SWEEPS=50000
ALGORITHM="loop"
{T=0.1;}
{T=0.2;}
{T=0.3;}
{T=0.4;}
{T=0.5;}
{T=0.6;}
{T=0.7;}
{T=0.8;}
{T=1.0;}
{T=1.25;}
{T=1.5;}
{T=1.75;}
{T=2.0;}
parameter2xml parm2d
loop parm2d.in.xml

Python での設定と実行

スクリプト tutorial2d.py は tutorial2c.py を書き換えたものです。接頭辞を parm2d に変え、LATTICE を "ladder" に変更し、J を J0 と J1（どちらも 1）に置き換えます。

4 つのシミュレーションの統合

同じフォルダで 4 つのシミュレーションをすべて実行したあと、スクリプト tutorial2full.py がすべての結果をまとめて読み込み、1 枚の図に重ねて表示します。

import pyalps
import matplotlib.pyplot as plt
import pyalps.plot

data = pyalps.loadMeasurements(pyalps.getResultFiles(), 'Susceptibility')
data = pyalps.flatten(data)

引数なしの pyalps.getResultFiles() は現在のフォルダにあるすべての出力ファイルを見つけます。pyalps.flatten は各データ点のシミュレーションパラメータを保持したまま、別々のファイルからの結果を 1 つのリストにまとめます。

MODEL と LATTICE の組み合わせごとに別々の曲線を作るために、foreach パラメータを付けて pyalps.collectXY を使います。

susceptibility = pyalps.collectXY(data, x='T', y='Susceptibility', foreach=['MODEL', 'LATTICE'])

保存されているパラメータを使って各曲線にラベルを付けます。

for s in susceptibility:
    if s.props['LATTICE'] == 'chain lattice':
        s.props['label'] = "chain"
    elif s.props['LATTICE'] == 'ladder':
        s.props['label'] = "ladder"
    if s.props['MODEL'] == 'spin':
        s.props['label'] = "quantum " + s.props['label']
    elif s.props['MODEL'] == 'Heisenberg':
        s.props['label'] = "classical " + s.props['label']

4 本の曲線をすべてプロットします。

plt.figure()
pyalps.plot.plot(susceptibility)
plt.xlabel('Temperature $T/J$')
plt.ylabel('Susceptibility $\chi J$')
plt.ylim(0, 0.25)
plt.legend()
plt.show()

得られる図は以下のようになるはずです。 高温では 4 本の曲線すべてがキュリー則 χ∝1/T\chi \propto 1/T に近づきます。 低温では、量子梯子はスピンギャップのために急激に落ち込みますが、古典鎖・量子鎖および古典梯子は有限のままか増大し続けます。

問題

  • 高温では 4 本の曲線すべてが同じ χ∝1/T\chi \propto 1/T のキュリー的振る舞いに近づきますか。前係数を決めているのは何でしょうか。
  • 低温において、古典鎖と量子鎖の最も目立つ違いは何でしょうか。
  • 低温では量子梯子の帯磁率が古典梯子よりはるかに速く落ち込みます。これを引き起こす物理的機構は何でしょうか。
  • より大きなシステムサイズや異なる格子（"cubic lattice"、"triangular lattice"——lattices.xml を参照）を試してください。結果はどう変わりますか。
MC-01b 平衡化MC-03 磁化