易辛模型

易辛模型（英語：Ising model，/ˈaɪsɪŋ/，德語：[ˈiːzɪŋ]），是一個以物理學家恩斯特·易辛為名的數學模型，用於描述物質的鐵磁性。該模型中包含了可以用來描述單個原子磁矩的參數 $\sigma _{i}$ ，其值只能為+1或-1，分別代表自旋向上或向下，這些磁矩通常會按照某種規則排列，形成晶格，並且在模型中會引入特定交互作用的參數，使得相鄰的自旋互相影響。雖然該模型相對於物理現實是一個相當簡化的模型，但它卻和鐵磁性物質一樣會產生相變。事實上，一個二維的方晶格易辛模型（英語：square-lattice Ising model）是已知最簡單而會產生相變的物理系統。^[1]

易辛模型最早是由物理學家威廉·楞次在1920年發明的，他把該模型當成是一個給他學生恩斯特·易辛的問題。易辛在他一篇1924年的論文中求得了一維易辛模型的解析解，並且證明它不會產生相變。^[2] 二維方晶格易辛模型相對於一維的難出許多，因此其解析的描述在一段時間之後才在1943年由拉斯·昂薩格給出。一般來說，二維易辛模型的解析解可由傳遞矩陣法求得，不過也有幾個和量子場論有關的解法。對於大於三維的易辛模型目前還沒有找到解析解，但其近似解可由諸多方法求得，例如平均場論。

定義

令Λ為所有晶格點的集合，其中每個晶格點都有一個所有和它相鄰的晶格點的集合(在數學上稱之為圖)並使這些晶格點形成一個d維的晶格。對於每個晶格點 k ∈ Λ 都有一個離散變數 σ_k ，其中 σ_k ∈ {+1, −1}，代表一個晶格點的自旋。而所有變數的集合σ = (σ_k)_k∈Λ則稱作 自旋組態。

對於兩個相鄰的晶格點i, j ∈ Λ ，我們可以引入一個交互作用參數J_ij，此外，我們可以假設每個自旋j ∈ Λ都和外加的磁場 h_j 作用。則整個系統的哈密頓量可寫成:

H(\sigma )=-\sum _{<i~j>}J_{ij}\sigma _{i}\sigma _{j}-\mu \sum _{j}h_{j}\sigma _{j}

其中 <ij> 代表晶格點 i 和晶格點 j 是相鄰的晶格點。因此哈密頓量的第一項為對每一對相鄰晶格點的總和（每一對只算一次），代表所有自旋之間交互作用的能量，而第二項則是磁場和自旋交互作用的能量。µ是晶格點磁矩的值，值得注意的是，電子的磁矩和他的自旋方向相反，所以哈密頓量的第二項應該要是正號比較合理，但在習慣上，還是會令第二項為負號。^[3]

該系統的組態概率 P(σ)為在熱平衡下某個特定自旋組態 σ 的概率，為波茲曼分佈:

$P_{\beta }(\sigma )={e^{-\beta H(\sigma )} \over Z_{\beta }},$

其中 β = (k_B T)⁻¹，而:

$Z_{\beta }=\sum _{\sigma }e^{-\beta H(\sigma )}$

是該概率分佈的歸一化常數，在統計力學中又稱做配分函數。對於有為自旋組態函數的物理量 f(σ) ，其期望值可表示為:

$\langle f\rangle _{\beta }=\sum _{\sigma }f(\sigma )P_{\beta }(\sigma )\,$

參數

H(σ) 中兩項前的負號是約定俗成的。因為第一項為負號，因此參數 J_ij 的正負號決定了該系統的性質，對於每一對 i 和 j :

J_{ij}>0

，則該系統為鐵磁性。

J_{ij}<0

，則該系統為反鐵磁性。

J_{ij}=0

表示自旋間無交互作用。

除此之外該系統為 非鐵磁性。

在鐵磁性的易辛模型中，相鄰自旋同方向時能量較低，因此自旋會傾向於同向排列，反之，在反鐵磁性的易辛模型中，相鄰自旋反向的能量較低，因此自旋會頃向於反向排列。

H(σ) 中的第二項為負號，表示自旋頃向於和外加磁場同向，因此 h_j 的正負也決定自旋頃向排列的方向。對於所有的j ，如果:

h_{j}>0

, 則晶格點 j 頃向於朝向正向。

h_{j}<0

, 則晶格點 j 頃向於朝向負向。

h_{j}=0

, 表示沒有外加力場作用在自旋上。

簡化

一個常見的簡化是假設沒有外加的磁場作用在易辛模型上，也就是說，對於所有的 j∈Λ，h_j = 0 。利用這項簡化，其哈密頓量可以寫成:

$H(\sigma )=-\sum _{<i~j>}J_{ij}\sigma _{i}\sigma _{j}.$

此時易辛模型在反轉所有自旋之下是對稱的:一個外加的力場會破壞這種對稱。

另一個常見的簡化是假設所有相鄰晶格點的交互作用都是相等的，因此可以設 J_ij = J 對於所有相鄰的 i, j ∈Λ，而其哈密頓量可以寫: $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| 呈指數衰減:

$\langle \sigma _{i}\sigma _{j}\rangle _{\beta }\leq C\exp(-c(\beta )|i-j|)$

其中c(β)是一個只和β有關的函數。

由此可知。這個系統是無序的。根據一維模型的結論，易辛錯誤的認為任何維度的易辛模型都不會有相變。但事實上，在二維或二維以上的模型中，該系統可以歷經從無序相轉變成有序相的相變。基本上在β值小(高溫)時，系統處在無序相，而β值大(低溫)時，系統處在有序相中。換句話說，當系統在有序相時：

$\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年發表了第一個正式的推導過程，其中裏面用到了包括塞格極限定理（英語：Szegő limit theorems）和弗雷德霍姆行列式（英語：Fredholm determinant）等數學工具^[6]。

楊李定理

數值模擬

蒙地卡羅方法

易辛模型一般來說很難直接進行數值計算，因為他的自旋組態非常之多。考慮一個擁有 L 個晶格點的模型，每個晶格點 σ_j 有 ±1 兩種可能，因此所有的自旋組態共有 2^L 種可能^[7]，這個數字會隨着 L 的增加而進行指數增長。這也是為什麼一般在做易辛模型的數值模擬時，都會採用蒙地卡羅方法（Monte Carlo method）^[7]

用蒙地卡羅方法來模擬易辛模型所用的哈密頓量為:

H(\sigma )=-J\sum _{<ij>}\sigma _{i}\sigma _{j}-h\sum _{j}\sigma _{j}.

此外，可以假設外加磁場 (h) 為零以簡化模型，因為大部分的問題都只需要用到零磁下的模型。因此，近一步簡化的哈密頓量為

H(\sigma )=-J\sum _{<ij>}\sigma _{i}\sigma _{j}.

可以從這個模型計算出的物理量包括該系統在定溫下的比熱和磁化強度等等。^[7]

Metropolis 演算法

Metropolis–Hastings演算法是在數值模擬易辛模型時最常用的一種蒙地卡羅方法。^[7]這種方法首先要選定一個選擇概率 g(μ, ν)，代表系統在狀態μ下，在所有可能狀態中選到狀態ν的概率。另外還需定義出一個接受概率 A(μ, ν) ，也就是說當系統在狀態μ下，接受系統跳到態ν的概率。如此的設計是為了讓系統達到細緻平衡。如果狀態ν被接受了，則整個系統就會跳到狀態ν，並且選擇和決定下一個要跳到的狀態。如果狀態ν被沒被接受，系統則留在狀態μ，一樣重新選擇下一個要跳到的狀態。這樣的步驟一直重複直到某些條件達成為止，譬如說整個易辛模型完全被磁化，也就是所有的自旋都只到同一個方向。^[7]另外在實行這種演算法，有一點需要注意的是需要選到適當 g(μ, ν) 以保證整個過程的遍歷性(ergodicity)。

在熱平衡時，整個系統的能量只會有小幅度的擾動^[7]，這點促成了在演算時採用單一自旋反轉法進行計算，也就是說每次系統轉換其狀態時，只改變其中一個自旋的方向。對自旋數很多的一易辛模型來說，系統在不同的狀態之間跳躍時，其能量改變的幅度都很小。事實上，對於每個晶格點都和 c 個晶格點相鄰的模型來說，每次能量改變的最大幅度為 2cJ。此外，採用這種單一自旋反轉法可以保證演算過程的遍歷性，因為任意一個狀態都可以藉由逐次的反轉相異的自旋，而變成任意其他狀態^[7]。

在下面的演算過程我們會採用周期性邊界條件使的每個晶格點的相鄰數 c 都相等。^[7]

演算過程

整個用於數值模擬易辛模型的演算過程，可由下列的方法建立。

因為共有 L 個晶格點，而且單一自旋反轉法是唯一系統可以從一個狀態跳到另一個狀態的方法，所對於任何一個狀態 μ，有 L 個新狀態 ν 它可以跳去。可以假設 μ 跳到這 L 個新狀態的概率是相等的，因此 g(μ, ν) = 1/L。為了滿足細緻平衡（英語：detailed balance），所以下面的等式要成立:

{\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 })}.

因此我們希望接受概率滿足:

{\frac {A(\mu ,\nu )}{A(\nu ,\mu )}}=e^{-\beta (H_{\nu }-H_{\mu })}.

如果 H_ν > H_μ 則 A(ν, μ) > A(μ, ν)，因此 Metropolis 另 A(μ, ν) 和 A(ν, μ) 中較大的為 1。按照這樣的定可以得到 A(μ, ν) 的值為: ^[7]

$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^{-\beta (H_{\nu }-H_{\mu })}$ 的概率保持反轉。
重複步驟一。

整個系統能量的改變量 H_ν−H_μ 只和反轉的自旋和和它相鄰的自旋有關，所以只要相鄰的自旋數不要太多，計算的速度算是相當快的。而整個系統會逐漸的趨近於某個平衡。

將易辛模型視為馬可夫鏈

將易辛模型比擬為馬可夫鏈是一件很容易的事情，因為下一刻狀態 ν 的轉移概率 P_β(ν) 只和目前狀態 μ有關。事實上，Metropolis 演算法就是馬可夫蒙地卡羅(Markov Chain Monte Carlo)的一個版本。而且因為演算時採用的是單一自旋反轉法，每個狀態可視為和另外 L 個狀態相連，而每一個狀態轉移則是將其中一個自旋反轉。^[8] 此外，因為不同狀態之間 H_σ 變化只和相鄰自旋的交互作用 J 有關，所以易辛模型也可視為選民模型（英語：Voter model）的一種。

一維易辛模型

在一維易辛模型系統中，假設每個帶有自旋的原子分佈在一維的圓圈中，且原子僅和鄰居發生交互作用，交互作用均為 $J$ ，能量可表示為

H(\sigma )=-J\sum _{i=1}^{N}\sigma _{i}\sigma _{i+1}+H\sum _{i=1}^{N}\sigma _{i}

其中h為外加磁場的強度，J為相鄰原子的耦合強度。該系統的自由能為：

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)

相鄰自旋的相關函數（英語：Correlation function (statistical mechanics)）為：

\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(\sigma )=-J(\sigma _{1}\sigma _{2}+\cdots +\sigma _{N-1}\sigma _{N}).

對於相鄰兩原子自旋朝向的四種情況，之間的能量只有兩種狀態：同向和反向。應用變量替換，

\sigma '_{j}=\sigma _{j}\sigma _{j-1}\qquad j\geq 2.

代入配分函數定義式，得：

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(\beta ,0)=-{\frac {1}{\beta }}\ln \left[e^{\beta J}+e^{-\beta J}\right].

同樣的變量替換法，得間隔N-1個自旋間的相關函數：

\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(\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_{\sigma ,\sigma '}$ 可視為2x2矩陣的元素。相鄰兩自旋各有四種可能狀態，其玻爾茲曼因子分別為：

V_{\sigma ,\sigma '}=e^{{\frac {\beta h}{2}}\sigma }e^{\beta J\sigma \sigma '}e^{{\frac {\beta h}{2}}\sigma '}

或寫作傳遞矩陣的四個元素：

V={\begin{bmatrix}e^{\beta (h+J)}&e^{-\beta J}\\e^{-\beta J}&e^{-\beta (h-J)}\end{bmatrix}}.

配分函數的展開式恰好和傳遞矩陣自乘後對角元素之和（矩陣的跡）一致。而矩陣的跡可通過求解傳遞矩陣的特徵值得出，因而：

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]

其中λ₁是V最大的特徵值，λ₂是該矩陣另一個特徵值。因而

\lambda _{1}=e^{\beta J}\cosh \beta h+{\sqrt {e^{2\beta J}(\sinh \beta h)^{2}+e^{-2\beta J}}}

由於|λ₂| < λ₁，當N很大時對配分函數的貢獻可忽略。

結論

系統最低能量為−N，此時所有自旋朝向相同。而其它構型相較基態的能量增量等於一維序列中自旋轉向的次數k。例如，能量次低的狀態相較基態，能量差為2k。由於能量與轉向數之間是線性關係，則兩個相鄰自旋方向相反的概率符合玻爾茲曼分佈：

{p \over 1-p}=e^{-2\beta J}.

藉由統計力學的配分函數可以計算再給定溫度下 $T$ ( $\beta =1/(kT)$ )的每個原子的磁矩期望值為

M(H,T)={\frac {\sinh(\beta H)}{\sqrt {\sinh ^{2}(\beta H)+e^{-4\beta J}}}}

.

所以一維易辛模型並沒有居里溫度、不會發生相變，即沒有自發磁化（英語：spontaneous magnetization）的現象。

M(0,T)=0

.

二維易辛模型

設二維方晶格易辛模型橫縱兩方向的交互作用能分別為 $J_{1}$ and $J_{2}$ 。拉斯·昂薩格求得無外磁矩，即 $h=0$ 時自由能的解析解：

-\beta f=\ln 2+{\frac {1}{8\pi ^{2}}}\int _{0}^{2\pi }d\theta _{1}\int _{0}^{2\pi }d\theta _{2}\ln[\cosh(2\beta J_{1})\cosh(2\beta J_{2})-\sinh(2\beta J_{1})\cos(\theta _{1})-\sinh(2\beta J_{2})\cos(\theta _{2})].

從自由能的偏導數可得到各種熱力學函數。

特別地，二維易辛模型有一個相變點，臨界溫度 $T_{c}$ 滿足以下方程：

\sinh \left({\frac {2J_{1}}{kT_{c}}}\right)\sinh \left({\frac {2J_{2}}{kT_{c}}}\right)=1

.

若 $J_{1}=J_{2}$ ，則 $T_{c}={\frac {2J}{k\ln(1+{\sqrt {2}})}}$ ，或 $\beta ={\frac {1}{kT_{c}}}\approx {\frac {0.44}{J}}$ .

參見

延伸閱讀

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完全的。)

註解

^ See Gallavotti (1999), Chapters VI-VII.
^ ^2.0 ^2.1 Ernst Ising, Contribution to the Theory of Ferromagnetism. [2014-08-14]. （原始內容存檔於2020-11-22）.
^ See Baierlein (1999), Chapter 16.
^ 存档副本 (PDF). [2014-09-07]. （原始內容 (PDF)存檔於2014-04-07）.
^ 存档副本. [2014-09-07]. （原始內容存檔於2015-04-03）.
^ Montroll, Potts & Ward 1963，第308–309頁
^ ^7.0 ^7.1 ^7.2 ^7.3 ^7.4 ^7.5 ^7.6 ^7.7 ^7.8 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.

[1] See Gallavotti (1999), Chapters VI-VII.

[Ising-2] 2.0 ^2.1 Ernst Ising, Contribution to the Theory of Ferromagnetism. [2014-08-14]. （原始內容存檔於2020-11-22）.

[3] See Baierlein (1999), Chapter 16.

[4] 存档副本 (PDF). [2014-09-07]. （原始內容 (PDF)存檔於2014-04-07）.

[5] 存档副本. [2014-09-07]. （原始內容存檔於2015-04-03）.

[Montroll_1963_pages=308-309-6] Montroll, Potts & Ward 1963，第308–309頁

[Newman-7] 7.0 ^7.1 ^7.2 ^7.3 ^7.4 ^7.5 ^7.6 ^7.7 ^7.8 Newman MEJ, Barkema GT, "Monte Carlo Methods in Statistical Physics, Clarendon Press, 1999

[8] 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.

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]