伊辛模型

伊辛模型（英语：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]