伊辛模型 (英语: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 .