易辛模型 (英語:Ising model , ,德语: [ˈiːzɪŋ] ),是一個以物理學家恩斯特·伊辛 為名的數學模型 ,用於描述物質的鐵磁性 。該模型中包含了可以用來描述單個原子磁矩 的參數
σ
i
{\displaystyle \sigma _{i}}
,其值只能為+1或-1,分別代表自旋 向上或向下,這些磁矩通常會按照某種規則排列,形成晶格 ,並且在模型中會引入特定交互作用 的參數,使得相鄰的自旋互相影響。雖然該模型相對於物理現實是一個相當簡化的模型,但它卻和鐵磁性物質一樣會產生相變 。事實上,一個二維的方晶格易辛模型 是已知最簡單而會產生相變的物理系統。[ 1]
易辛模型最早是由物理學家威廉·楞次 在1920年發明的,他把該模型當成是一個給他學生恩斯特·易辛的問題。易辛在他一篇1924年的論文中求得了一維易辛模型的解析解 ,並且證明它不會產生相變。[ 2]
二維方晶格易辛模型相對於一維的難出許多,因此其解析的描述在一段時間之後才在1943年由拉斯·昂薩格 給出。一般來說,二維易辛模型的解析解可由傳遞矩陣法 求得,不過也有幾個和量子場論 有關的解法。對於大於三維的易辛模型目前還沒有找到解析解,但其近似解可由諸多方法求得,例如平均場論 。
令Λ為所有晶格點的集合 ,其中每個晶格點都有一個所有和它相鄰的晶格點的集合(在數學上稱之為圖 )並使這些晶格點形成一個d維的晶格。對於每個晶格點 k ∈ Λ 都有一個離散 變數 σk ,其中 σk ∈ {+1, −1},代表一個晶格點的自旋。而所有變數的集合σ = (σk )k ∈Λ 則稱作 自旋組態 。
對於兩個相鄰的晶格點i , j ∈ Λ ,我們可以引入一個交互作用 參數J ij ,此外,我們可以假設每個自旋j ∈ Λ都和外加的磁場 hj 作用。則整個系統的哈密頓量 可寫成:
H
(
σ
)
=
−
∑
<
i
j
>
J
i
j
σ
i
σ
j
−
μ
∑
j
h
j
σ
j
{\displaystyle H(\sigma )=-\sum _{<i~j>}J_{ij}\sigma _{i}\sigma _{j}-\mu \sum _{j}h_{j}\sigma _{j}}
其中 <ij> 代表晶格點 i 和晶格點 j 是相鄰的晶格點。因此哈密頓量的第一項為對每一對相鄰晶格點的總和(每一對只算一次),代表所有自旋之間交互作用的能量 ,而第二項則是磁場和自旋交互作用的能量。µ是晶格點磁矩 的值,值得注意的是,電子 的磁矩 和他的自旋方向相反,所以哈密頓量的第二項應該要是正號比較合理,但在習慣上,還是會令第二項為負號。[ 3]
該系統的組態機率 P(σ)為在熱平衡 下某個特定自旋組態 σ 的機率 ,為波茲曼分布 :
P
β
(
σ
)
=
e
−
β
H
(
σ
)
Z
β
,
{\displaystyle P_{\beta }(\sigma )={e^{-\beta H(\sigma )} \over Z_{\beta }},}
其中 β = (kB T )−1 ,而:
Z
β
=
∑
σ
e
−
β
H
(
σ
)
{\displaystyle Z_{\beta }=\sum _{\sigma }e^{-\beta H(\sigma )}}
是該機率分布的歸一化常數 ,在統計力學 中又稱做配分函數 。對於有為自旋組態函數 的物理量 f(σ) ,其期望值可表示為:
⟨
f
⟩
β
=
∑
σ
f
(
σ
)
P
β
(
σ
)
{\displaystyle \langle f\rangle _{\beta }=\sum _{\sigma }f(\sigma )P_{\beta }(\sigma )\,}
H (σ) 中兩項前的負號是約定俗成的。因為第一項為負號,因此參數 J ij 的正負號決定了該系統的性質,對於每一對 i 和 j :
J
i
j
>
0
{\displaystyle J_{ij}>0}
,則該系統為鐵磁性 。
J
i
j
<
0
{\displaystyle J_{ij}<0}
,則該系統為反鐵磁性 。
J
i
j
=
0
{\displaystyle J_{ij}=0}
表示自旋間無交互作用。
除此之外該系統為 非鐵磁性 。
在鐵磁性的易辛模型中,相鄰自旋 同方向時能量較低,因此自旋會傾向於同向排列,反之,在反鐵磁性的易辛模型中,相鄰自旋 反向的能量較低,因此自旋會頃向於反向排列。
H (σ) 中的第二項為負號,表示自旋頃向於和外加磁場 同向,因此 hj 的正負也決定自旋頃向排列的方向。對於所有的j ,如果:
h
j
>
0
{\displaystyle h_{j}>0}
, 則晶格點 j 頃向於朝向正向。
h
j
<
0
{\displaystyle h_{j}<0}
, 則晶格點 j 頃向於朝向負向。
h
j
=
0
{\displaystyle h_{j}=0}
, 表示沒有外加力場 作用在自旋 上 。
一個常見的簡化是假設沒有外加的磁場 作用在易辛模型上,也就是說,對於所有的 j∈Λ,hj = 0 。利用這項簡化,其哈密頓量 可以寫成:
H
(
σ
)
=
−
∑
<
i
j
>
J
i
j
σ
i
σ
j
.
{\displaystyle H(\sigma )=-\sum _{<i~j>}J_{ij}\sigma _{i}\sigma _{j}.}
此時易辛模型在反轉所有自旋之下是對稱 的:一個外加的力場會破壞這種對稱。
另一個常見的簡化是假設所有相鄰晶格點的交互作用都是相等的,因此可以設 Jij = J 對於所有相鄰的 i , j ∈Λ,而其哈密頓量可以寫:
H
(
σ
)
=
−
J
∑
<
i
j
>
σ
i
σ
j
.
{\displaystyle H(\sigma )=-J\sum _{<i~j>}\sigma _{i}\sigma _{j}.}
一維平移不變鐵磁性在無外加場下的易辛模型的所有組態之機率密度 投射在二維平面上。
被研究的最透徹易辛模型是在 d 維 晶格上,平移不變、鐵磁性 並且無外加場的模型。也就是 :Λ = Z d , J ij = 1, h = 0。
易辛在他1924年的博士論文中,解決了在 d=1 時的情況 ,這個一維的模型可以想像成一排的自旋 ,而每個自旋都只和它左右兩邊的自旋交互作用。這個一模型不會產生相變 [ 4]
,換句話說,對於所有正值的β,任意兩自旋的相關係數 <σi σj > 都對 |i −j | 呈指數衰減 :
⟨
σ
i
σ
j
⟩
β
≤
C
exp
(
−
c
(
β
)
|
i
−
j
|
)
{\displaystyle \langle \sigma _{i}\sigma _{j}\rangle _{\beta }\leq C\exp(-c(\beta )|i-j|)}
其中c(β)是一個只和β有關的函數。
由此可知。這個系統是無序 的。根據一維模型的結論,易辛錯誤的認為任何維度的易辛模型都不會有相變 。但事實上,在二維或二維以上的模型中,該系統可以歷經從無序相 轉變成有序相 的相變 。基本上在β值小(高溫)時,系統處在無序相,而β值大(低溫)時,系統處在有序相中。換句話說,當系統在有序相時:
⟨
σ
i
σ
j
⟩
β
≥
c
′
(
β
)
>
0.
{\displaystyle \langle \sigma _{i}\sigma _{j}\rangle _{\beta }\geq c^{\prime }(\beta )>0.\,}
其中c'(β)也是一個只和β有關的函數。
這個性質是首先被魯道夫·佩爾斯 (Rudolf Peierls)在1936年證明的[ 5]
,他的證法後來被稱為佩爾斯論述 (Peierls argument)。
在零磁下的二維方晶格易辛模型的解析解 後來在1943年被昂薩格 解出,他證明了該模型的相關函數 和自由能 可由一個無交互作用的格點費米子場 (noninteracting lattice fermion)來界定。昂薩格在1949年發表了一個決定了自發磁化現象的公式,但卻沒有給出推導過程。後來是楊振寧 在1951年發表了第一個正式的推導過程,其中裡面用到了包括塞格極限定理 和弗雷德霍姆行列式 等數學工具[ 6] 。
易辛模型一般來說很難直接進行數值計算,因為他的自旋組態非常之多。考慮一個擁有 L 個晶格點的模型,每個晶格點 σj 有 ±1 兩種可能,因此所有的自旋組態共有 2 L 種可能[ 7] ,這個數字會隨著 L 的增加而進行指數 增長。這也是為什麼一般在做易辛模型的數值模擬時,都會採用蒙地卡羅方法 (Monte Carlo method)[ 7]
用蒙特卡洛方法來模擬易辛模型所用的哈密頓量 為:
H
(
σ
)
=
−
J
∑
<
i
j
>
σ
i
σ
j
−
h
∑
j
σ
j
.
{\displaystyle H(\sigma )=-J\sum _{<ij>}\sigma _{i}\sigma _{j}-h\sum _{j}\sigma _{j}.}
此外,可以假設外加磁場 (h ) 為零以簡化模型,因為大部分的問題都只需要用到零磁下的模型。因此,近一步簡化的哈密頓量 為
H
(
σ
)
=
−
J
∑
<
i
j
>
σ
i
σ
j
.
{\displaystyle H(\sigma )=-J\sum _{<ij>}\sigma _{i}\sigma _{j}.}
可以從這個模型計算出的物理量 包括該系統在定溫下的比熱 和磁化強度 等等。[ 7]
Metropolis–Hastings演算法 是在數值模擬易辛模型時最常用的一種蒙地卡羅方法。[ 7] 這種方法首先要選定一個選擇機率 g (μ, ν),代表系統在狀態μ下,在所有可能狀態中選到狀態ν的機率。另外還需定義出一個接受機率 A (μ, ν) ,也就是說當系統在狀態μ下,接受系統跳到態ν的機率。如此的設計是為了讓系統達到細緻平衡 。如果狀態ν被接受了,則整個系統就會跳到狀態ν,並且選擇和決定下一個要跳到的狀態。如果狀態ν被沒被接受,系統則留在狀態μ,一樣重新選擇下一個要跳到的狀態。這樣的步驟一直重複直到某些條件達成為止,譬如說整個易辛模型完全被磁化 ,也就是所有的自旋 都只到同一個方向。[ 7] 另外在實行這種演算法,有一點需要注意的是需要選到適當 g (μ, ν) 以保證整個過程的遍歷性 (ergodicity)。
在熱平衡 時,整個系統的能量 只會有小幅度的擾動[ 7] ,這點促成了在演算時採用單一自旋反轉法 進行計算,也就是說每次系統轉換其狀態時,只改變其中一個自旋 的方向。對自旋數很多的一易辛模型來說,系統在不同的狀態之間跳躍時,其能量 改變的幅度都很小。事實上,對於每個晶格點都和 c 個晶格點相鄰的模型來說,每次能量改變的最大幅度為 2cJ 。此外,採用這種單一自旋反轉法可以保證演算過程的遍歷性 ,因為任意一個狀態都可以藉由逐次的反轉相異的自旋,而變成任意其他狀態[ 7] 。
在下面的演算過程我們會採用周期性边界条件 使的每個晶格點的相鄰數 c 都相等。[ 7]
整個用於數值模擬易辛模型的演算過程,可由下列的方法建立。
因為共有 L 個晶格點,而且單一自旋反轉法是唯一系統可以從一個狀態跳到另一個狀態的方法,所對於任何一個狀態 μ,有 L 個新狀態 ν 它可以跳去。可以假設 μ 跳到這 L 個新狀態的機率是相等的,因此 g (μ, ν) = 1/L 。 為了滿足細緻平衡 ,所以下面的等式要成立:
P
(
μ
,
ν
)
P
(
ν
,
μ
)
=
g
(
μ
,
ν
)
A
(
μ
,
ν
)
g
(
ν
,
μ
)
A
(
ν
,
μ
)
=
A
(
μ
,
ν
)
A
(
ν
,
μ
)
=
P
β
(
ν
)
P
β
(
μ
)
=
1
Z
e
−
β
(
H
ν
)
1
Z
e
−
β
(
H
μ
)
=
e
−
β
(
H
ν
−
H
μ
)
.
{\displaystyle {\frac {P(\mu ,\nu )}{P(\nu ,\mu )}}={\frac {g(\mu ,\nu )A(\mu ,\nu )}{g(\nu ,\mu )A(\nu ,\mu )}}={\frac {A(\mu ,\nu )}{A(\nu ,\mu )}}={\frac {P_{\beta }(\nu )}{P_{\beta }(\mu )}}={\frac {{\frac {1}{Z}}e^{-\beta (H_{\nu })}}{{\frac {1}{Z}}e^{-\beta (H_{\mu })}}}=e^{-\beta (H_{\nu }-H_{\mu })}.}
因此我們希望接受機率滿足:
A
(
μ
,
ν
)
A
(
ν
,
μ
)
=
e
−
β
(
H
ν
−
H
μ
)
.
{\displaystyle {\frac {A(\mu ,\nu )}{A(\nu ,\mu )}}=e^{-\beta (H_{\nu }-H_{\mu })}.}
如果 H ν > H μ 則 A (ν, μ) > A (μ, ν),因此 Metropolis 另 A (μ, ν) 和 A (ν, μ) 中較大的為 1。按照這樣的定可以得到 A (μ, ν) 的值為: [ 7]
A
(
μ
,
ν
)
=
{
e
−
β
(
H
ν
−
H
μ
)
,
if
H
ν
−
H
μ
>
0
1
,
otherwise
.
{\displaystyle A(\mu ,\nu )={\begin{cases}e^{-\beta (H_{\nu }-H_{\mu })},&{\text{if }}H_{\nu }-H_{\mu }>0\\1,&{\text{otherwise}}.\end{cases}}}
最後,整個演算法最基本的形式為:
用 g (μ, ν) 選出一個自旋,並且計算所有和其自旋相關的能量貢獻。(也就是所有和它相鄰自旋交互作用的能量貢獻。)
反轉其自旋,再計算一次所有和它相關的能量貢獻。
如果其能量貢獻下降,保持自旋反轉。
如果貢獻能量上升,則令自旋有
e
−
β
(
H
ν
−
H
μ
)
{\displaystyle e^{-\beta (H_{\nu }-H_{\mu })}}
的機率保持反轉。
重複步驟一。
整個系統能量的改變量 H ν −H μ 只和反轉的自旋和和它相鄰的自旋有關,所以只要相鄰的自旋數不要太多,計算的速度算是相當快的。而整個系統會逐漸的趨近於某個平衡。
將易辛模型比擬為馬可夫鏈 是一件很容易的事情,因為下一刻狀態 ν 的轉移機率 P β (ν) 只和目前狀態 μ有關 。事實上,Metropolis 演算法就是馬可夫蒙地卡羅 (Markov Chain Monte Carlo)的一個版本。而且因為演算時採用的是單一自旋反轉法,每個狀態可視為和另外 L 個狀態相連,而每一個狀態轉移則是將其中一個自旋 反轉。[ 8] 此外,因為不同狀態之間 H σ 變化只和相鄰自旋的交互作用 J 有關,所以易辛模型也可視為選民模型 的一種。
在一維易辛模型系統中,假設每個带有自旋的原子分布在一維的圓圈中,且原子仅和鄰居發生交互作用,交互作用均為
J
{\displaystyle J}
,能量可表示為
H
(
σ
)
=
−
J
∑
i
=
1
N
σ
i
σ
i
+
1
+
H
∑
i
=
1
N
σ
i
{\displaystyle H(\sigma )=-J\sum _{i=1}^{N}\sigma _{i}\sigma _{i+1}+H\sum _{i=1}^{N}\sigma _{i}}
其中h 为外加磁场的强度,J 为相邻原子的耦合强度。该系统的自由能 为:
Z
(
β
,
h
)
=
−
lim
N
→
∞
1
β
L
ln
(
Z
(
β
)
)
=
−
1
β
ln
(
e
β
J
cosh
β
h
+
e
2
β
J
(
sinh
β
h
)
2
+
e
−
2
β
J
)
{\displaystyle Z(\beta ,h)=-\lim _{N\to \infty }{\frac {1}{\beta L}}\ln(Z(\beta ))=-{\frac {1}{\beta }}\ln \left(e^{\beta J}\cosh \beta h+{\sqrt {e^{2\beta J}(\sinh \beta h)^{2}+e^{-2\beta J}}}\right)}
相邻自旋的相关函数 为:
⟨
σ
i
σ
j
⟩
−
⟨
σ
i
⟩
⟨
σ
j
⟩
=
C
(
β
)
e
−
c
(
β
)
|
i
−
j
|
{\displaystyle \langle \sigma _{i}\sigma _{j}\rangle -\langle \sigma _{i}\rangle \langle \sigma _{j}\rangle =C(\beta )e^{-c(\beta )|i-j|}}
其中,C(β)和c (β)是关于温度的函数,当温度T > 0时取正值。当T → 0时,c (β)趋于零。
易辛自己的论文中已给出了一维易辛模型配分函数 精确解。[ 2] 无外加磁场,即h = 0时,计算是很简单的。
此时能量简化为:
H
(
σ
)
=
−
J
(
σ
1
σ
2
+
⋯
+
σ
N
−
1
σ
N
)
.
{\displaystyle H(\sigma )=-J(\sigma _{1}\sigma _{2}+\cdots +\sigma _{N-1}\sigma _{N}).}
对于相邻两原子自旋朝向的四种情况,之间的能量只有两种状态:同向和反向。应用变量替换,
σ
j
′
=
σ
j
σ
j
−
1
j
≥
2.
{\displaystyle \sigma '_{j}=\sigma _{j}\sigma _{j-1}\qquad j\geq 2.}
代入配分函数定义式,得:
Z
(
β
)
=
∑
σ
1
,
…
,
σ
N
e
β
J
σ
1
σ
2
e
β
J
σ
2
σ
3
⋯
e
β
J
σ
N
−
1
σ
N
=
2
∏
j
=
2
N
∑
σ
j
′
e
β
J
σ
j
′
=
2
[
e
β
J
+
e
−
β
J
]
N
−
1
.
{\displaystyle Z(\beta )=\sum _{\sigma _{1},\ldots ,\sigma _{N}}e^{\beta J\sigma _{1}\sigma _{2}}\;e^{\beta J\sigma _{2}\sigma _{3}}\;\cdots e^{\beta J\sigma _{N-1}\sigma _{N}}=2\prod _{j=2}^{N}\sum _{\sigma '_{j}}e^{\beta J\sigma '_{j}}=2\left[e^{\beta J}+e^{-\beta J}\right]^{N-1}.}
因此,系统的自由能为
f
(
β
,
0
)
=
−
1
β
ln
[
e
β
J
+
e
−
β
J
]
.
{\displaystyle f(\beta ,0)=-{\frac {1}{\beta }}\ln \left[e^{\beta J}+e^{-\beta J}\right].}
同样的变量替换法,得间隔N-1个自旋间的相关函数:
⟨
σ
j
σ
j
+
N
⟩
=
[
e
β
J
−
e
−
β
J
e
β
J
+
e
−
β
J
]
N
{\displaystyle \langle \sigma _{j}\sigma _{j+N}\rangle =\left[{\frac {e^{\beta J}-e^{-\beta J}}{e^{\beta J}+e^{-\beta J}}}\right]^{N}}
T ≠ 0时相关性呈指数衰减。而仅仅在绝对零度时T = 0下,换句话说,β → ∞ 时,系统保持长程相关性。
如果外磁场h ≠ 0,配分函数的计算需要引入传递矩阵法 。周期性边界条件 的近似下,配分函数为
Z
(
β
)
=
∑
σ
1
,
…
,
σ
N
e
β
h
σ
1
e
β
J
σ
1
σ
2
e
β
h
σ
2
e
β
J
σ
2
σ
3
⋯
e
β
h
σ
N
e
β
J
σ
N
σ
1
=
∑
σ
1
,
…
,
σ
L
V
σ
1
,
σ
2
V
σ
2
,
σ
3
⋯
V
σ
N
,
σ
1
.
{\displaystyle Z(\beta )=\sum _{\sigma _{1},\ldots ,\sigma _{N}}e^{\beta h\sigma _{1}}e^{\beta J\sigma _{1}\sigma _{2}}\;e^{\beta h\sigma _{2}}e^{\beta J\sigma _{2}\sigma _{3}}\;\cdots e^{\beta h\sigma _{N}}e^{\beta J\sigma _{N}\sigma _{1}}=\sum _{\sigma _{1},\ldots ,\sigma _{L}}V_{\sigma _{1},\sigma _{2}}V_{\sigma _{2},\sigma _{3}}\cdots V_{\sigma _{N},\sigma _{1}}.}
系数
V
σ
,
σ
′
{\displaystyle V_{\sigma ,\sigma '}}
可视为2x2矩阵的元素。相邻两自旋各有四种可能状态,其玻尔兹曼因子分别为:
V
σ
,
σ
′
=
e
β
h
2
σ
e
β
J
σ
σ
′
e
β
h
2
σ
′
{\displaystyle V_{\sigma ,\sigma '}=e^{{\frac {\beta h}{2}}\sigma }e^{\beta J\sigma \sigma '}e^{{\frac {\beta h}{2}}\sigma '}}
或写作传递矩阵的四个元素:
V
=
[
e
β
(
h
+
J
)
e
−
β
J
e
−
β
J
e
−
β
(
h
−
J
)
]
.
{\displaystyle V={\begin{bmatrix}e^{\beta (h+J)}&e^{-\beta J}\\e^{-\beta J}&e^{-\beta (h-J)}\end{bmatrix}}.}
配分函数的展开式恰好和传递矩阵自乘后对角元素之和(矩阵的迹 )一致。而矩阵的迹可通过求解传递矩阵的特征值得出,因而:
Z
(
β
)
=
T
r
V
N
=
λ
1
N
+
λ
2
N
=
λ
1
N
[
1
+
(
λ
2
λ
1
)
N
]
{\displaystyle Z(\beta )={\rm {Tr}}V^{N}=\lambda _{1}^{N}+\lambda _{2}^{N}=\lambda _{1}^{N}\left[1+\left({\frac {\lambda _{2}}{\lambda _{1}}}\right)^{N}\right]}
其中λ1 是V 最大的特征值,λ2 是该矩阵另一个特征值。因而
λ
1
=
e
β
J
cosh
β
h
+
e
2
β
J
(
sinh
β
h
)
2
+
e
−
2
β
J
{\displaystyle \lambda _{1}=e^{\beta J}\cosh \beta h+{\sqrt {e^{2\beta J}(\sinh \beta h)^{2}+e^{-2\beta J}}}}
由于|λ2 | < λ1 ,当N 很大时对配分函数的贡献可忽略。
Kerson Huang, Introduction to Statistical Physics.
I. A. Stepanov. Exact Solutions of the One-Dimensional, Two-Dimensional, and Three-Dimensional Ising Models. - Nano Science and Nano Technology: An Indian Journal. 2012. Vol. 6. No 3. 118 - 122.(這篇文章可在該期刊的網站免費閱讀)
Barry A. Cipra, "The Ising model is NP-complete", SIAM News , Vol. 33, No. 6; online edition (.pdf) (一篇文章闡述為何任意的易辛模型無法有一般性的精確解,因為非平面的易辛模型是NP完全 的。)
存档副本 (PDF) . [2014-09-07 ] . (原始内容 (PDF) 存档于2014-04-07).
存档副本 . [2014-09-07 ] . (原始内容存档 于2015-04-03).
Newman MEJ, Barkema GT, "Monte Carlo Methods in Statistical Physics, Clarendon Press, 1999
Teif V.B. General transfer matrix formalism to calculate DNA-protein-drug binding in gene regulation. Nucleic Acids Res. 2007, 35 : e80. doi:10.1093/nar/gkm268 .