搜索
File 2026.10.28物理

正则系综与巨正则系综

Date / 日期
Updated / 更新
2026.10.28
Reading / 阅读
135 Min
Words / 字数
26917
Category / 分类
物理热力学

正则系综与巨正则系综

· 更新于 2026年10月28日· 约 135 分钟· 26917 字
物理热力学
#物理#热力学#正则系综#Ising模型#巨正则系综#化学势

点击跳转到英文版（English ver.)

上一章的微正则系综描述的是孤立系统，能量一次给定，每个平均值都要在一条很薄的能量壳层内计算。这种做法在概念上干净，在计算上却很麻烦，因为壳层限制会进入每一个积分，而且必须先算出熵 S(E,V)S(E,V)，别的量才算得下去。

本章放松对能量的约束。系统与热库接触，能量可以涨落，系综平均变成整个相空间上的无约束积分。随后再放松对粒子数的约束，这就引出巨正则系综与化学势。

正则系综

正则系综

考虑一个由两个子系统组成的孤立系统

1：N1,V1,H1N_1, V_1, H_1
2：N2,V2,H2N_2, V_2, H_2

其中

N=N1+N2V=V1+V2H=H1(p1,q1)+H2(p2,q2)\begin{align} N &= N_1 + N_2\\ V &= V_1 + V_2\\ H &= H_1(p_1,q_1) + H_2(p_2,q_2) \end{align}

并假设两个子系统之间的相互作用可以忽略不计。

如上一章所讨论，能量的分解中以最概然的分解 E=Eˉ1+Eˉ2E = \bar{E}_1 + \bar{E}_2 为主。我们关心 N2≫N1N_2 \gg N_1 且 Eˉ2≫Eˉ1\bar{E}_2 \gg \bar{E}_1 的特殊情形，此时子系统 2 是子系统 1 的一个大环境。

设 E=E1+E2E = E_1 + E_2 是能量在最概然分解附近涨落的任一分解。系统 1 处于 (p1,q1)(p_1,q_1) 且系统 2 处于 (p2,q2)(p_2,q_2) 的概率为 ∝dp1dq1dp2dq2\propto \dd p_1 \dd q_1 \dd p_2 \dd q_2。

如果只看系统 1，就应当把系统 2 在 E2=E−E1E_2 = E - E_1 下的所有可能微观态求和，

∝dp1dq1∫E2=E−E1dp2dq2=dp1dq1 Γ2(E−E1)\begin{align} \propto \dd p_1 \dd q_1 \int_{E_2 = E-E_1} \dd p_2 \dd q_2 = \dd p_1 \dd q_1\, \Gamma_2(E-E_1) \end{align}

也就是说，系统 1 处于 (p1,q1)(p_1,q_1) 处体积元 dp1dq1\dd p_1 \dd q_1 内某个微观态的概率为 ∝dp1dq1 Γ2(E−E1)\propto \dd p_1 \dd q_1\, \Gamma_2(E-E_1)，因此密度分布为

ρ(p1,q1)∝Γ2(E−E1)\begin{align} \rho(p_1,q_1) \propto \Gamma_2(E-E_1) \end{align}

由于 E1E_1 只在 Eˉ1\bar{E}_1 附近涨落，而 Eˉ1≪E\bar{E}_1 \ll E，

klog⁡Γ2(E−E1)=klog⁡Γ2(E)−E1∂ klog⁡Γ2(E2)∂E2∣E2=E+⋯=S2(E)−E1∂S2(E2)∂E2∣E2=E=S2(E)−E1T\begin{align} k \log \Gamma_2(E-E_1) &= k \log \Gamma_2(E) - E_1 \left.\frac{\partial\, k\log \Gamma_2(E_2)}{\partial E_2}\right|_{E_2 = E} + \cdots\\ &= S_2(E) - E_1 \left.\frac{\partial S_2(E_2)}{\partial E_2}\right|_{E_2 = E}\\ &= S_2(E) - \frac{E_1}{T} \end{align}

这里 TT 是整个系统的温度，因此也是系统 1 的温度，而 S2(E)S_2(E) 与系统 1 无关。

问题：为什么不直接展开 Γ2(E−E1)\Gamma_2(E-E_1) 本身？

答案：Γ2\Gamma_2 是相空间体积，随 N2N_2 指数增长。把这样的函数按 E1E_1 的幂次展开，得到的系数会随 N2N_2 增大，级数在第一项之后就没有用了。可加且随能量缓慢变化的量是 Γ2\Gamma_2 的对数，即熵，log⁡Γ2\log \Gamma_2 才是广延量。此外，Γ2\Gamma_2 在 Δ\Delta 的尺度上甚至不是光滑的，它只是一个粗粒化的体积，只在壳层宽度以内有定义。对 klog⁡Γ2k\log\Gamma_2 做展开是合法的，因为二阶项比一阶项小一个 E1/EE_1/E 量级的因子，而我们已假设 E1≪EE_1 \ll E。

略去下标 1，得到

ρ(p,q)∝e−H(p,q)/kT=e−βH(p,q),β=1kT\begin{align} \rho(p,q) \propto e^{-H(p,q)/kT} = e^{-\beta H(p,q)},\qquad \beta = \frac{1}{kT} \end{align}

这就是所谓的正则系综的概率分布。与微正则系综类似，正则系综的 ρ(p,q)\rho(p,q) 对 (p,q)(p,q) 的依赖同样只通过 Hamiltonian H(p,q)H(p,q) 实现。

正则系综的宏观条件：

  • NN 固定；
  • 系统与温度为 TT 的热库接触，TT 固定。

换句话说，能量不再是系统的宏观参数，它已被热库的温度取代。图 3.1(a) 画出了这个安排。

TikZ 源码
\begin{tikzpicture}[>=stealth, thick, font=\small]
  % ---- (a) canonical ensemble ----
  \draw (0,0) rectangle (2,2);
  \node at (1,1.35) {System 1};
  \node at (1,0.75) {$N_1,V_1$};


  \draw (2.9,0) rectangle (6.7,2.8);
  \node at (4.8,2.2) {Heat reservoir 2};
  \node at (4.8,1.6) {$N_2\gg N_1$};
  \node at (4.8,1.05) {$\bar{E}_2\gg\bar{E}_1$};
  \node at (4.8,0.5) {$T$ fixed};


  \draw[<->,red] (2.05,1.0) -- (2.85,1.0);
  \node[red,above] at (2.45,1.05) {$E$};
  \node at (1,-0.5) {(a)};


  % ---- (b) grand canonical ensemble ----
  \draw (0,-5.2) rectangle (6.7,-1.9);
  \node[anchor=north west] at (0.15,-2.0) {Reservoir: $T$, $\mu$};


  \draw (2.5,-4.3) rectangle (4.9,-2.35);
  \fill (2.9,-2.7) circle (2pt);
  \fill (4.4,-2.75) circle (2pt);
  \fill (3.0,-3.9) circle (2pt);
  \fill (4.5,-3.95) circle (2pt);
  \node at (3.7,-4.6) {open system: $V$ fixed, $N$ varies};


  \fill (0.6,-2.5) circle (2pt);
  \fill (1.1,-3.4) circle (2pt);
  \fill (1.7,-4.3) circle (2pt);
  \fill (0.6,-4.7) circle (2pt);
  \fill (5.9,-2.4) circle (2pt);
  \fill (6.2,-3.4) circle (2pt);
  \fill (5.5,-4.3) circle (2pt);
  \fill (6.2,-4.7) circle (2pt);


  \draw[->,red] (2.0,-3.25) -- (2.4,-3.25);
  \node[red,above left] at (2.3,-3.3) {$\Delta N$};
  \draw[->,red] (5.0,-3.25) -- (5.4,-3.25);
  \node[red,above right] at (5.1,-3.3) {$\Delta N$};


  \node at (1,-5.65) {(b)};
\end{tikzpicture}

把配分函数定义为 Γ\Gamma 空间中的有效体积，即微观态的有效数目，

QN(V,T)≡1N!∫d3Npd3Nq e−H/kT\begin{align} Q_N(V,T) \equiv \frac{1}{N!} \int \dd^{3N}p \dd^{3N}q\, e^{-H/kT} \end{align}

热力学量则由下式给出

QN(V,T)≡e−A(V,T)/kT\begin{align} Q_N(V,T) \equiv e^{-A(V,T)/kT} \end{align}

其中 A(V,T)A(V,T) 是自由能。要说明这一点，需要证明

(a) AA 是广延量；
(b) AA 与 UU 通过 A=U−TSA = U - TS 联系，其中 S≡−(∂A/∂T)VS \equiv -(\partial A/\partial T)_V。

证明 (b)：对下面的恒等式关于 β=1/kT\beta = 1/kT 求导，

1N!∫d3Npd3Nq eβ[A(V,T)−H(p,q)]=1\begin{align} \frac{1}{N!} \int \dd^{3N}p \dd^{3N}q\, e^{\beta[A(V,T)-H(p,q)]} = 1 \end{align}

注意到 H(p,q)H(p,q) 不依赖于 β\beta，得到

1N!∫d3Npd3Nq eβ[A−H][A−H+β∂A∂β]=0\begin{align} \frac{1}{N!} \int \dd^{3N}p \dd^{3N}q\, e^{\beta[A-H]} \left[ A - H + \beta \frac{\partial A}{\partial \beta} \right] = 0 \end{align}

由于

β∂A∂β=β∂A∂T∂T∂β=−T(∂A∂T)V=TS\begin{align} \beta \frac{\partial A}{\partial \beta} = \beta \frac{\partial A}{\partial T}\frac{\partial T}{\partial \beta} = -T \left(\frac{\partial A}{\partial T}\right)_V = TS \end{align}

这里用了 T=1/(kβ)T = 1/(k\beta)，因此 ∂T/∂β=−T/β\partial T/\partial \beta = -T/\beta。于是方括号内就是 A−H+TSA - H + TS，取平均给出 A−U+TS=0A - U + TS = 0，即 A=U−TSA = U - TS，其中 U=⟨H⟩U = \langle H \rangle 是内能。

证明 (a)：由两个相互独立的部分组成的系统，其配分函数因子化，QN1+N2=QN1QN2Q_{N_1+N_2} = Q_{N_1}Q_{N_2}，差别只在与上一章证明 SS 的广延性时相同的二项式计数和相同的 O(log⁡N)O(\log N) 项忽略。因此 −A/kT=log⁡QN-A/kT = \log Q_N 是可加量，AA 是广延量。

请注意，现在温度从一开始就已经定义好了，它就是热库的温度，是一个自变量。于是自由能 A(V,T)A(V,T) 是特性函数，可以直接导出

P=−(∂A∂V)TS=−(∂A∂T)VU=A+TSG=A+PV\begin{align} P &= -\left(\frac{\partial A}{\partial V}\right)_T\\ S &= -\left(\frac{\partial A}{\partial T}\right)_V\\ U &= A + TS\\ G &= A + PV \end{align}

正则系综中的自由能 A(V,T)A(V,T) 起到微正则系综中熵 S(E,V)S(E,V) 的作用。

原则上 dpdq\dd p \dd q 的积分不能扩展到整个 Γ\Gamma 空间，因为从一开始 ρ(p1,q1)\rho(p_1,q_1) 就要求 E1≤EE_1 \le E。但下一小节会说明，H≥EH \ge E 区域的贡献小到可以忽略。事实上，只有能量取某一个值时对配分函数的贡献是主要的。因此实际计算中积分可以遍及所有能量。这给正则系综中的计算带来很大方便。

此外，正则系综中的自由能 A(V,T)A(V,T) 可以显式计算。相比之下，微正则系综中的内能 U(S,V)U(S,V) 只能通过计算 S(E,V)S(E,V) 才有望得到，这种求解通常困难得多。

正则系综中配分函数 QN(V,T)Q_N(V,T) 的定义与微正则系综中微观态的总体积 Γ(E)\Gamma(E) 类似。事实上

klog⁡QN(V,T)=−A(V,T)/T=S−U/T\begin{align} k \log Q_N(V,T) = -A(V,T)/T = S - U/T \end{align}

问题：为什么可以这样做，又如何做到？

提示：只有能量为最概然值 Eˉ\bar{E} 的微观态，即能量为平均值 U=EˉU = \bar{E} 的微观态，对 QN(V,T)Q_N(V,T) 有贡献。因此 QN(V,T)=e−Eˉ/kTΓ(Eˉ)Q_N(V,T) = e^{-\bar{E}/kT}\Gamma(\bar{E})，取对数就直接得到 klog⁡QN=−U/T+Sk\log Q_N = -U/T + S。

正则系综的一个简单论证：设一个大孤立系统由许多相同的准独立子系统组成，每个子系统由 Hamiltonian H(p,q)H(p,q) 描述。按照上一章的结果，即 Maxwell-Boltzmann 分布，在 (p,q)(p,q) 找到某个子系统的概率分布为

ρ(p,q)∝e−βH(p,q)\begin{align} \rho(p,q) \propto e^{-\beta H(p,q)} \end{align}

因此不必重复热库的论证也能得到同样的结果，前提是子系统之间的相互作用可以忽略。

能量涨落

问题：

  • 微正则系综与正则系综在什么意义上等价？
  • 为什么需要两个系综？

任何宏观可观测量的计算都通过

⟨O⟩=1QNN!∫dpdq e−βHO\begin{align} \langle O \rangle = \frac{1}{Q_N N!} \int \dd p \dd q\, e^{-\beta H} O \end{align}

而能量为

U=⟨H⟩=1QNN!∫dpdq e−βHH\begin{align} U = \langle H \rangle = \frac{1}{Q_N N!} \int \dd p \dd q\, e^{-\beta H} H \end{align}

可以写成

1N!∫dpdq eβ(A−H)(U−H)=0\begin{align} \frac{1}{N!}\int \dd p \dd q\, e^{\beta(A-H)} (U - H) = 0 \end{align}

两边对 β\beta 求导，注意到 ∂U/∂β\partial U/\partial\beta 是唯一让导数作用在指数函数以外的东西上的项，得到

∂U∂β+1N!∫dpdq eβ(A−H)(U−H)(A−H+β∂A∂β)=0\begin{align} \frac{\partial U}{\partial \beta} + \frac{1}{N!}\int \dd p \dd q\, e^{\beta(A-H)} (U-H)\left( A - H + \beta \frac{\partial A}{\partial \beta} \right) = 0 \end{align}

进一步，由于

β∂A∂β=−T(∂A∂T)V=TS,U=A+TS\begin{align} \beta \frac{\partial A}{\partial \beta} = -T\left(\frac{\partial A}{\partial T}\right)_V = TS,\qquad U = A + TS \end{align}

方括号内是 A−H+TS=U−HA - H + TS = U - H，因此 ∂U/∂β+⟨(U−H)2⟩=0\partial U/\partial\beta + \langle (U-H)^2 \rangle = 0，即

⟨H2⟩−⟨H⟩2=−∂U∂β=kT2∂U∂T=kT2CV\begin{align} \langle H^2 \rangle - \langle H \rangle^2 = -\frac{\partial U}{\partial \beta} = kT^2 \frac{\partial U}{\partial T} = kT^2 C_V \end{align}

或者用偏差 ΔE=H−⟨H⟩\Delta E = H - \langle H \rangle 表示，

⟨(ΔE)2⟩=kT2CV\begin{align} \langle (\Delta E)^2 \rangle = kT^2 C_V \end{align}

事实上，这是著名的涨落耗散定理的一个特例，该定理把能量的涨落与对外部改变的响应联系起来。由于左边是方差，还可以得出平衡态下 CV≥0C_V \ge 0。

现在，由于 ⟨H⟩∼N\langle H \rangle \sim N 且 CV∼NC_V \sim N，当 N→∞N \to \infty 时

⟨H2⟩−⟨H⟩2⟨H⟩2∼1N→0\begin{align} \frac{\langle H^2 \rangle - \langle H \rangle^2}{\langle H \rangle^2} \sim \frac{1}{N} \to 0 \end{align}

如图 3.2 所示。因此，微正则系综与正则系综是等价的，意思是两种情形下能量都只在一个固定值附近涨落。不过在正则系综中能量并不给定，需要算出来。此外，这个定理还说明能量很大的那些状态的贡献可以忽略，权重 e−βHe^{-\beta H} 按指数衰减，涨落很小。

TikZ 源码
\begin{tikzpicture}[>=stealth, thick, font=\small, x=0.8cm, y=2.8cm]
  \draw[->] (0,0) -- (10.5,0) node[below] {$N$};
  \draw[->] (0,0) -- (0,1.05) node[above] {$\frac{\sqrt{\langle(\Delta E)^2\rangle}}{U}$};
  \node[below left] at (0,0) {$O$};


  \draw[blue] plot[domain=0.95:10, samples=120] (\x, {0.85/sqrt(\x)});


  \draw[dashed] (1,0) -- (1,0.85);
  \draw[dashed] (4,0) -- (4,0.425);
  \draw[dashed] (4,0.425) -- (0,0.425);
  \draw[dashed] (1,0.85) -- (0,0.85);


  \fill (1,0.85) circle (2pt);
  \fill (4,0.425) circle (2pt);
  \fill (9,0.2833) circle (2pt);


  \node[right, xshift=4pt] at (1,0.85) {$N$};
  \node[right, xshift=4pt] at (4,0.425) {$4N$};
  \node[above] at (9,0.32) {$O(N^{-1/2})$};
\end{tikzpicture}

图中画出能量的相对涨落随 NN 的变化，曲线为 O(N−1/2)O(N^{-1/2})，所以系统增大到四倍，相对涨落减半。

Boltzmann 统计理论

考虑 NN 个相同的准独立粒子，

H=∑iϵ(pi,qi)\begin{align} H = \sum_i \epsilon(p_i,q_i) \end{align}

这里下标 ii 标记的是粒子。于是

QN(V,T)=1N!∑p,qe−βH=1N!∑p,q∏ie−βϵ(pi,qi)=1N!(∑pi,qie−βϵ(pi,qi))N=1N![Z(β,V)]N=e−βA(V,T)\begin{align} Q_N(V,T) &= \frac{1}{N!}\sum_{p,q} e^{-\beta H}\\ &= \frac{1}{N!}\sum_{p,q}\prod_i e^{-\beta\epsilon(p_i,q_i)} = \frac{1}{N!}\left(\sum_{p_i,q_i} e^{-\beta\epsilon(p_i,q_i)}\right)^N\\ &= \frac{1}{N!}\left[Z(\beta,V)\right]^N = e^{-\beta A(V,T)} \end{align}

其中 Z(β,V)=∑p,qe−βϵ(p,q)Z(\beta,V) = \sum_{p,q} e^{-\beta\epsilon(p,q)} 是单个粒子的配分函数，就是上一章引入的 Z(β,y)Z(\beta,y)，只是把 yy 换成了 VV。

用 log⁡N!≈Nlog⁡N−N\log N! \approx N\log N - N，得到 log⁡QN≈Nlog⁡Z−Nlog⁡N+N=Nlog⁡(Z/N)+N\log Q_N \approx N\log Z - N\log N + N = N\log(Z/N) + N，即

A(V,T)=−NkTlog⁡Z(β,V)N−NkT=−NkT(log⁡Z(β,V)N+1)\begin{align} A(V,T) = -NkT\log\frac{Z(\beta,V)}{N} - NkT = -NkT\left(\log\frac{Z(\beta,V)}{N} + 1\right) \end{align}

常数项 11 来自 Stirling 近似，在计算熵时不能丢掉。然后

S=−(∂A∂T)V=Nklog⁡ZN+NkT∂log⁡Z∂T+Nk=Nklog⁡ZN+NkT∂log⁡Z∂β∂β∂T+Nk=Nk(log⁡Z(β,V)N+1−β∂log⁡Z(β,V)∂β)\begin{align} S &= -\left(\frac{\partial A}{\partial T}\right)_V\\ &= Nk\log\frac{Z}{N} + NkT\frac{\partial \log Z}{\partial T} + Nk\\ &= Nk\log\frac{Z}{N} + NkT\frac{\partial \log Z}{\partial \beta}\frac{\partial \beta}{\partial T} + Nk\\ &= Nk\left(\log\frac{Z(\beta,V)}{N} + 1 - \beta\frac{\partial \log Z(\beta,V)}{\partial \beta}\right) \end{align}

这里用了 ∂β/∂T=−β/T\partial\beta/\partial T = -\beta/T。

把这里的 Z(β,V)Z(\beta,V) 与上一章引入的 Z(β,y)Z(\beta,y) 比较，唯一的差别是这里计算 Z(β,V)Z(\beta,V) 时能量不受总能量的限制。但这个差别小到可以忽略，因为 ϵ≪E∝N\epsilon \ll E \propto N。

正则系综的几个应用

能量均分定理

设 Hamiltonian 为

H=12∑iaipi2+U({qi})\begin{align} H = \frac{1}{2}\sum_i a_i p_i^2 + U(\{q_i\}) \end{align}

则 ⟨12aipi2⟩=12kT\left\langle \frac{1}{2} a_i p_i^2 \right\rangle = \frac{1}{2}kT。

证明：

⟨12aipi2⟩=∫−∞+∞12aipi2 e−βH d3Npd3Nq∫−∞+∞e−βH d3Npd3Nq=∫−∞+∞12aipi2 e−β12aipi2dpi∫−∞+∞e−β12aipi2dpi\begin{align} \left\langle \frac{1}{2} a_i p_i^2 \right\rangle &= \frac{\displaystyle\int_{-\infty}^{+\infty} \frac{1}{2} a_i p_i^2\, e^{-\beta H}\, \dd^{3N}p \dd^{3N}q}{\displaystyle\int_{-\infty}^{+\infty} e^{-\beta H}\, \dd^{3N}p \dd^{3N}q}\\ &= \frac{\displaystyle\int_{-\infty}^{+\infty} \frac{1}{2}a_i p_i^2\, e^{-\beta \frac{1}{2}a_i p_i^2}\dd p_i}{\displaystyle\int_{-\infty}^{+\infty} e^{-\beta \frac{1}{2}a_i p_i^2}\dd p_i} \end{align}

因为对其余所有变量的积分在分子与分母之间相消。对分子做分部积分，

numerator=∫12aipi2 e−β12aipi2dpi=∫12aipi⋅pie−β12aipi2dpi=−12β[pie−β12aipi2]−∞+∞+12β∫e−β12aipi2dpi=12kT⋅denominator\begin{align} \text{numerator} &= \int \frac{1}{2}a_i p_i^2\, e^{-\beta \frac{1}{2}a_i p_i^2}\dd p_i = \int \frac{1}{2}a_i p_i \cdot p_i e^{-\beta\frac{1}{2}a_i p_i^2}\dd p_i\\ &= -\frac{1}{2\beta}\left[ p_i e^{-\beta\frac{1}{2}a_i p_i^2} \right]_{-\infty}^{+\infty} + \frac{1}{2\beta}\int e^{-\beta\frac{1}{2}a_i p_i^2}\dd p_i\\ &= \frac{1}{2}kT \cdot \text{denominator} \end{align}

其中边界项为零，因为 pi→±∞p_i \to \pm\infty 时 H→∞H \to \infty。因此 ⟨12aipi2⟩=12kT\left\langle \frac{1}{2} a_i p_i^2 \right\rangle = \frac{1}{2}kT，命题得证。

证明所需的条件是 Boltzmann 因子在边界处衰减，即二次项系数为正，ai>0a_i > 0，并且对 pip_i 的积分可以延伸到 ±∞\pm\infty。同样的手法适用于 HH 中以二次形式出现的任何变量，无论它是坐标还是动量，势能 U({qi})U(\{q_i\}) 在平均值中不起作用。

练习：若 H(ξ→±∞)→∞H(\xi \to \pm\infty) \to \infty，则 ⟨ξ∂H∂ξ⟩=kT\left\langle \xi \frac{\partial H}{\partial \xi} \right\rangle = kT。

与上一章的推导比较。那里同一个恒等式需要仔细处理能量壳层，要在曲面 H=EH = E 上做分部积分，还要忽略一个含 (E−H)(E-H) 的项。这里壳层没有了，平均值就是普通的 Gaussian 积分，证明只用三行。这就是正则系综在实际计算上的好处。

对于由二次项之和构成的 Hamiltonian，H=∑iAiPi2+∑iBiQi2H = \sum_i A_i P_i^2 + \sum_i B_i Q_i^2，定理给出 ⟨H⟩=f2kT\langle H \rangle = \frac{f}{2}kT，其中 ff 是二次自由度的数目，与上一章的推导一致。

理想气体

对理想气体，H=12m∑ipi2H = \frac{1}{2m}\sum_i p_i^2，于是

QN(V,T)=1N!∫d3Npd3Nq e−β∑ipi2/2m=VNN!(∫d3p e−βp2/2m)N=VNN![q(β)]N\begin{align} Q_N(V,T) &= \frac{1}{N!}\int \dd^{3N}p \dd^{3N}q\, e^{-\beta \sum_i p_i^2/2m}\\ &= \frac{V^N}{N!}\left(\int \dd^3 p\, e^{-\beta p^2/2m}\right)^N\\ &= \frac{V^N}{N!}\left[q(\beta)\right]^N \end{align}

其中

q(β)=∫d3p e−βp2/2m=∫0∞dp p2e−p2/2mkT∫dΩ=(2πmkT)3/2\begin{align} q(\beta) = \int \dd^3 p\, e^{-\beta p^2/2m} = \int_0^\infty \dd p\, p^2 e^{-p^2/2mkT}\int \dd\Omega = (2\pi mkT)^{3/2} \end{align}

这里用了

∫0∞dp p2e−p2/2mkT=π4(2mkT)3/2,∫dΩ=4π\begin{align} \int_0^\infty \dd p\, p^2 e^{-p^2/2mkT} = \frac{\sqrt{\pi}}{4}(2mkT)^{3/2},\qquad \int \dd\Omega = 4\pi \end{align}

于是

A(V,T)=−kTlog⁡QN(V,T)=−NkT(log⁡VN+32log⁡(2πmkT)+1)\begin{align} A(V,T) &= -kT\log Q_N(V,T)\\ &= -NkT\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT) + 1\right) \end{align}

以及

S=−(∂A∂T)V=Nk(log⁡VN+32log⁡(2πmkT)+1)+NkT⋅32⋅1T=Nk(log⁡VN+32log⁡(2πmkT)+52)\begin{align} S &= -\left(\frac{\partial A}{\partial T}\right)_V\\ &= Nk\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT) + 1\right) + NkT\cdot\frac{3}{2}\cdot\frac{1}{T}\\ &= Nk\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT) + \frac{5}{2}\right) \end{align} P=−(∂A∂V)T=NkTV\begin{align} P = -\left(\frac{\partial A}{\partial V}\right)_T = \frac{NkT}{V} \end{align}

这就是状态方程。内能为 U=A+TS=32NkTU = A + TS = \frac{3}{2}NkT。用 UU 消去 TT，即 2πmkT=4πmU/(3N)2\pi mkT = 4\pi mU/(3N)，

S=Nk(log⁡VN+32log⁡4πmU3N+52)\begin{align} S = Nk\left(\log\frac{V}{N} + \frac{3}{2}\log\frac{4\pi mU}{3N} + \frac{5}{2}\right) \end{align}

这些结果与上一章微正则系综给出的完全一致，包括状态方程、U=32NkTU = \frac{3}{2}NkT，以及作为 VV 和 UU 的函数的熵，其中对数里的 V/NV/N 解决了 Gibbs 佯谬。

有两点需要说明。第一，因子 1/N!1/N! 是必需的；没有它，熵中会出现 log⁡V\log V 而不是 log⁡(V/N)\log(V/N)，就不再是广延量。第二，测度 d3pd3q\dd^3 p \dd^3 q 中没有出现 h3h^3 因子，相空间单位体积的选择只会把 SS 移动一个正比于 NN 的项，即每个粒子一个常数，不影响 PP、UU、CVC_V 以及 AA 的任何导数。这个常数在经典统计力学内部无法确定。

因此整个计算都是经典的，粒子被当作可以分辨的，只差交换，对状态的求和化为对相空间的积分。这描述的是高温下的稀薄气体，其中粒子之间的平均距离远大于其量子弥散。在低温或高密度下，态的计数必须用量子统计重做，Gibbs 因子中的 1/N!1/N! 则由多粒子波函数的正确对称性取代。

Ising 模型

Ising 模型是磁性材料最简单的模型，它能表现出相变。

图 3.3 画出二维的正方格子，用来描述一块晶体。每个格点上有一个磁矩，用 Si=±1S_i = \pm 1 表示。

TikZ 源码
\begin{tikzpicture}[>=stealth, thick, font=\small, x=1.35cm, y=1.35cm]
  % lattice
  \foreach \i in {0,1,2,3} {
    \draw[gray!55] (0,\i) -- (3,\i);
    \draw[gray!55] (\i,0) -- (\i,3);
  }


  % spins
  \foreach \x/\y/\s in {
    0/0/1, 1/0/-1, 2/0/1, 3/0/1,
    0/1/-1, 1/1/-1, 2/1/1, 3/1/-1,
    0/2/1, 1/2/1, 2/2/-1, 3/2/-1,
    0/3/-1, 1/3/1, 2/3/1, 3/3/-1} {
      \fill (\x,\y) circle (1.6pt);
      \draw[->, blue] (\x,\y-0.26*\s) -- (\x,\y+0.26*\s);
  }


  % nearest neighbour pair
  \draw[red, dashed, thick] (0.15,1) -- (0.85,1);
  \node[red, above] at (0.5,1.12) {$\langle ij\rangle$};


  % labels
  \node at (1.5,3.7) {square lattice};
  \draw[<->] (2.55,3.05) -- (3.25,3.05);
  \node[above] at (2.9,3.05) {$a$};
  \node at (1.5,-0.75) {$S_i=\pm 1$ on every site};
\end{tikzpicture}

Hamiltonian 写成

−1kTH=K∑⟨ij⟩SiSj+h∑iSi,Si=±1\begin{align} -\frac{1}{kT}H = K \sum_{\langle ij \rangle} S_i S_j + h \sum_i S_i, \qquad S_i = \pm 1 \end{align}

其中

⟨ij⟩:nearest neighbourJ=−kTK:couplinghˉ=−kTh:external magnetic field\begin{align} \langle ij \rangle &: \text{nearest neighbour}\\ J &= -kT K : \text{coupling}\\ \bar h &= -kT h : \text{external magnetic field} \end{align}

把整个 Hamiltonian 用 kTkT 作单位是一种约定，因为 Boltzmann 因子中只出现组合 H/kTH/kT。采用这个约定后，一个位形的权重是 eK∑SiSj+h∑Sie^{K\sum S_iS_j + h\sum S_i}，因此 K>0K > 0 使相邻自旋倾向于平行，即铁磁情形，而 hh 使 SiS_i 倾向于取某一个符号。一条键的能量是 JSiSj=−kTK SiSjJ S_i S_j = -kT K\, S_i S_j，一个自旋在外场中的能量是 hˉSi=−kTh Si\bar h S_i = -kT h\, S_i。

正则系综由 ρ∝e−H/kT\rho \propto e^{-H/kT} 描述。

略去 Gibbs 因子，配分函数写成

Z=∑{Si}e−H/kT\begin{align} Z = \sum_{\{S_i\}} e^{-H/kT} \end{align}

注意这里没有出现 1/N!1/N!，格点是标好的，各位形之间可以区分，求和遍及 NN 个自旋的 2N2^N 个位形。

平均磁化强度为

M=1Z∑{Si}(1N∑iSi)e−H/kT\begin{align} M = \frac{1}{Z}\sum_{\{S_i\}} \left( \frac{1}{N}\sum_i S_i \right) e^{-H/kT} \end{align}

可以求解的情形，这里 dd 表示空间维数：

  • d=1d = 1：没有真正的相变；
  • d=2d = 2：太困难，这里不处理。

h=0h = 0 时的结果如下。存在一个临界温度 TCT_C，使得

T>TC:M=0T<TC:M≠0T=TC:M=0 continuous,∂M/∂T discontinuous\begin{align} T > T_C&: \quad M = 0\\ T < T_C&: \quad M \neq 0\\ T = T_C&: \quad M = 0\ \text{continuous},\quad \partial M/\partial T\ \text{discontinuous} \end{align}

在 TCT_C 附近，磁化强度按如下方式变化

M∝(−τ)β,τ=T−TCTC\begin{align} M \propto (-\tau)^{\beta}, \qquad \tau = \frac{T - T_C}{T_C} \end{align}

所以 MM 在相变处连续而它的导数不连续，这是连续相变，∂M/∂T\partial M/\partial T 有跃变。这里的指数 β\beta 是序参量的常用临界指数；不幸的是它与 β=1/kT\beta = 1/kT 用了同一个字母，只能从上下文区分。

这样的模型如何求解，或者至少如何近似？平均场近似把一个自旋与近邻的相互作用换成与平均磁化强度的相互作用。

首先，一个好的近似应当满足什么要求？

  • 近似之后能够把计算做出来。
  • 有清楚的物理意义，并且在极端条件下近似应当正确或大致正确。

平均场的替换是

∑⟨ij⟩SiSj→∑⟨ij⟩Si⟨S⟩=d⟨S⟩∑iSi\begin{align} \sum_{\langle ij \rangle} S_i S_j \to \sum_{\langle ij \rangle} S_i \langle S \rangle = d \langle S \rangle \sum_i S_i \end{align}

其中 dd 是配位数，即一个格点的最近邻数目；对超立方格子它等于空间维数的两倍，所以一维链的 d=2d = 2，正方格子的 d=4d = 4。在上面最后一个等号处，对 ⟨ij⟩\langle ij\rangle 的求和是这样计数的，每个格点贡献 dd 条键，也就是每条键算两次，两端各算一次；正是这个约定导致下面的 KC=1/dK_C = 1/d。

替换之后 Hamiltonian 不再耦合不同格点，

−1kTH=−1kT∑iH(Si),−1kTH(Si)=Kd⟨S⟩Si+hSi\begin{align} -\frac{1}{kT}H = -\frac{1}{kT}\sum_i H(S_i),\qquad -\frac{1}{kT}H(S_i) = Kd \langle S \rangle S_i + h S_i \end{align}

因此

Z=ZiN,Zi=∑Si=±1e−1kTH(Si)=∑Si=±1eKd⟨S⟩Si+hSi=2cosh⁡(Kd⟨S⟩+h)\begin{align} Z = Z_i^N, \qquad Z_i = \sum_{S_i = \pm 1} e^{-\frac{1}{kT}H(S_i)} = \sum_{S_i=\pm 1} e^{Kd\langle S\rangle S_i + hS_i} = 2\cosh\left( Kd\langle S\rangle + h \right) \end{align} M=⟨Si⟩=1Zi∑Si=±1Sie−1kTH(Si)=tanh⁡(Kd⟨S⟩+h)\begin{align} M = \langle S_i \rangle = \frac{1}{Z_i}\sum_{S_i = \pm 1} S_i e^{-\frac{1}{kT}H(S_i)} = \tanh\left( Kd\langle S \rangle + h \right) \end{align}

自洽条件是 M=⟨S⟩M = \langle S \rangle。对 h=0h = 0，

M=tanh⁡(KdM)\begin{align} M = \tanh(KdM) \end{align}

由于 tanh⁡x=x−x3/3+⋯\tanh x = x - x^3/3 + \cdots，有

M=KdM−13(KdM)3+⋯\begin{align} M = KdM - \frac{1}{3}(KdM)^3 + \cdots \end{align}

图解结果如图 3.4 所示，曲线 tanh⁡(KdM)\tanh(KdM) 与直线 MM 本身相交。

TikZ 源码
\begin{tikzpicture}[>=stealth, thick, font=\small, x=1.15cm, y=1.15cm]


  % ---------- (a) Kd < 1 ----------
  \begin{scope}[shift={(0,0)}]
    \draw[->] (-1.65,0) -- (1.65,0) node[below] {$M$};
    \draw[->] (0,-1.45) -- (0,1.45);
    \node[above] at (0,1.45) {$\tanh(KdM)$};
    \draw[gray] (-1.3,-1.3) -- (1.3,1.3);
    \draw[blue] plot[domain=-1.5:1.5, samples=150] (\x, {1-2/(exp(0.5*\x)+1)});
    \fill (0,0) circle (2pt);
    \node at (0,-1.85) {(a) $Kd<1$: only $M=0$};
  \end{scope}


  % ---------- (b) Kd > 1 ----------
  \begin{scope}[shift={(4.8,0)}]
    \draw[->] (-1.65,0) -- (1.65,0) node[below] {$M$};
    \draw[->] (0,-1.45) -- (0,1.45);
    \node[above] at (0,1.45) {$\tanh(KdM)$};
    \draw[gray] (-1.3,-1.3) -- (1.3,1.3);
    \draw[blue] plot[domain=-1.5:1.5, samples=150] (\x, {1-2/(exp(4*\x)+1)});
    \fill (0,0) circle (2pt);
    \fill (0.96,0.96) circle (2pt);
    \fill (-0.96,-0.96) circle (2pt);
    \node[above left] at (0.96,0.96) {$M_0$};
    \node[below, yshift=-3pt] at (-0.96,-1.02) {$-M_0$};
    \node at (0,-1.85) {(b) $Kd>1$: $M=0$ and $M=\pm M_0$};
  \end{scope}
\end{tikzpicture}

当 Kd<1Kd < 1 时只有一个解 M=0M = 0。当 Kd>1Kd > 1 时有三个根，M=0M = 0 和 M=±M0M = \pm M_0，也就是两个分支。

平衡态下自由能应当取极小，可以证明 M=±M0M = \pm M_0 是物理的解，而 M=0M = 0 只是驻点。因此

KCd=1,KC=1d\begin{align} K_C d = 1, \qquad K_C = \frac{1}{d} \end{align}

相变就发生了。

系统在 KCK_C 附近有什么特征？令

K=KC+ΔK,ΔK small⇒M small\begin{align} K = K_C + \Delta K, \qquad \Delta K\ \text{small} \Rightarrow M\ \text{small} \end{align}

于是

M=(KC+ΔK)dM−13(KC+ΔK)3d3M3≈M+ΔKdM−13M3−ΔKdM3≈M+ΔKdM−13M3\begin{align} M &= (K_C+\Delta K)dM - \frac{1}{3}(K_C+\Delta K)^3 d^3 M^3\\ &\approx M + \Delta K d M - \frac{1}{3}M^3 - \Delta K d M^3\\ &\approx M + \Delta K d M - \frac{1}{3}M^3 \end{align}

这里在得到最后一行时用了 KCd=1K_C d = 1，并因为 ΔKM3\Delta K M^3 是高阶项而把它略去。因此

M=±M0=±3ΔKKC\begin{align} M = \pm M_0 = \pm\sqrt{\frac{3\Delta K}{K_C}} \end{align}

由于 JJ 和 kk 固定时 K∝1/TK \propto 1/T，在 TCT_C 附近领头阶有 ΔK/KC≈−(T−TC)/TC=−τ\Delta K/K_C \approx -(T-T_C)/T_C = -\tau，于是 M∝(−τ)1/2M \propto (-\tau)^{1/2}，即平均场给出的指数 β\beta 的值是 1/21/2。

练习：求 KCK_C 处及其附近的 ∂M/∂ΔK\partial M/\partial \Delta K、UU 和 CVC_V，并讨论它们在 KCK_C 处及其附近的行为。

在什么条件下平均场近似几乎正确？

  • 空间维数较高，每个自旋有很多近邻，用平均值代替近邻的作法不至于太粗暴。
  • 长程相互作用，这时平均场算不上近似，它几乎就是相互作用本身。

要记住近似必须被论证，不能直接假设。d=1d = 1 时严格模型在任何有限温度下都没有相变，而上面用一维链 d=2d = 2 得到的平均场条件却预言 KC=1/2K_C = 1/2。

统计物理中的高等方法：

近似方法

  • 平均场方法，以及改进的平均场方法
  • 变分方法
  • 级数展开（高温展开与低温展开、累积量展开等）
  • 重整化群方法（实空间、动量空间）

严格方法

  • 数值求解
  • Monte Carlo 模拟
  • 分子动力学

巨正则系综

现在我们允许系统的粒子数不守恒。思路与推导正则系综时类似。

Γ\Gamma 空间：由 (p,q,N)(p,q,N) 张成，N=0,1,2,…N = 0,1,2,\dots

巨正则系综：一组系统，组成和宏观条件相同，但处于不同的粒子数和不同的微观态，彼此之间没有相互作用。

系综中的一个成员，即系统的一个微观态，用 Γ\Gamma 空间中的一个点表示。这些点的密度分布记为 ρ(p,q,N)\rho(p,q,N)。图 3.1(b) 画出了这个安排。

为求 ρ(p,q,N)\rho(p,q,N)，考虑一个 NN 个粒子、体积 VV、温度 TT 的正则系综。把系统分成两个子系统，

1：N1,V1N_1, V_1
2：N2=N−N1,V2=V−V1N_2 = N - N_1, V_2 = V - V_1

并考虑 V2≫V1V_2 \gg V_1 且 N2≫N1N_2 \gg N_1 的特殊情形。

把 V1V_1 中 N1N_1 个粒子的坐标记作 {p1,q1}\{p_1,q_1\}，V2V_2 中的记作 {p2,q2}\{p_2,q_2\}，

H(p,q,N)=H(p1,q1,N1)+H(p2,q2,N2)\begin{align} H(p,q,N) = H(p_1,q_1,N_1) + H(p_2,q_2,N_2) \end{align}

注意这里的 H(p1,q1,N1)H(p_1,q_1,N_1) 与 H(p2,q2,N2)H(p_2,q_2,N_2) 是同一个函数在不同自变量取值下的结果。这与讨论正则系综时不同，因为现在允许粒子在两个子系统之间交换，如果一个粒子从 V2V_2 移到 V1V_1，就把它的坐标从 {p2,q2}\{p_2,q_2\} 改成 {p1,q1}\{p_1,q_1\}。

整个正则系统的配分函数为

QN(V,T)=1N!∫Vdpdq e−βH(p,q)=1N!∫dp1dp2∑N1=0NN!N1!N2!∫V1dq1∫V2dq2 exp⁡{−β[H(p1,q1,N1)+H(p2,q2,N2)]}=∑N1=0N1N1!∫V1dp1dq1 e−βH(p1,q1,N1)⋅1N2!∫V2dp2dq2 e−βH(p2,q2,N2)=∑N1=0NQN2(V2,T)N1!∫V1dp1dq1 e−βH(p1,q1,N1)\begin{align} Q_N(V,T) &= \frac{1}{N!}\int_V \dd p \dd q\, e^{-\beta H(p,q)}\\ &= \frac{1}{N!}\int \dd p_1 \dd p_2 \sum_{N_1=0}^{N} \frac{N!}{N_1! N_2!} \int_{V_1}\dd q_1 \int_{V_2}\dd q_2\, \exp\{-\beta[H(p_1,q_1,N_1) + H(p_2,q_2,N_2)]\}\\ &= \sum_{N_1=0}^{N} \frac{1}{N_1!}\int_{V_1}\dd p_1 \dd q_1\, e^{-\beta H(p_1,q_1,N_1)} \cdot \frac{1}{N_2!}\int_{V_2}\dd p_2 \dd q_2\, e^{-\beta H(p_2,q_2,N_2)}\\ &= \sum_{N_1=0}^{N} \frac{Q_{N_2}(V_2,T)}{N_1!}\int_{V_1}\dd p_1 \dd q_1\, e^{-\beta H(p_1,q_1,N_1)} \end{align}

现在对 dq1\dd q_1 的积分限于 V1V_1，对 dq2\dd q_2 的积分限于 V2V_2，还需要对 N1N_1 求和。这个论证与推导正则系综时对能量的分解类似。组合数 N!/(N1!N2!)N!/(N_1!N_2!) 是把 N1N_1 个粒子分到 V1V_1、N2N_2 个粒子分到 V2V_2 的可能方式数。

例如 N=2N = 2 时，

∫Vdq1dq2=(∫V1+∫V2)(∫V1+∫V2)=∫V2∫V2+2∫V2∫V1+∫V1∫V1\begin{align} \int_V \dd q_1 \dd q_2 &= \left( \int_{V_1} + \int_{V_2} \right)\left( \int_{V_1} + \int_{V_2} \right)\\ &= \int_{V_2}\int_{V_2} + 2\int_{V_2}\int_{V_1} + \int_{V_1}\int_{V_1} \end{align}

其中二项式系数 2!/(N1!N2!)2!/(N_1!N_2!) 给出各项的重数。

引入分母 N!N! 是为了解决 Gibbs 佯谬。只要粒子数在物理现象中起作用，就需要这样的因子。这里 Gibbs 因子是必需的，对 QNQ_N、QN1Q_{N_1} 和 QN2Q_{N_2} 都是如此。

因此可以把在 V1V_1 中找到 N1N_1 个粒子的概率分布定义为

ρ(p1,q1,N1)=QN2(V2,T)QN(V,T)e−βH(p1,q1,N1)N1!\begin{align} \rho(p_1,q_1,N_1) = \frac{Q_{N_2}(V_2,T)}{Q_N(V,T)}\frac{e^{-\beta H(p_1,q_1,N_1)}}{N_1!} \end{align}

理解：计入 Gibbs 因子，

ρ(p1,q1,N1,p2,q2,N2)∝e−βH(p1,q1,N1)N1!e−βH(p2,q2,N2)N2!\begin{align} \rho(p_1,q_1,N_1,p_2,q_2,N_2) \propto \frac{e^{-\beta H(p_1,q_1,N_1)}}{N_1!}\frac{e^{-\beta H(p_2,q_2,N_2)}}{N_2!} \end{align}

如果只看系统 1，就应当对 p2p_2 和 q2q_2 积分，

ρ(p1,q1,N1)∝QN2(V2,T)e−βH(p1,q1,N1)N1!\begin{align} \rho(p_1,q_1,N_1) \propto Q_{N_2}(V_2,T)\frac{e^{-\beta H(p_1,q_1,N_1)}}{N_1!} \end{align}

因子

QN2(V2,T)QN(V,T)=exp⁡{−β[A(N−N1,V−V1,T)−A(N,V,T)]}\begin{align} \frac{Q_{N_2}(V_2,T)}{Q_N(V,T)} = \exp\{-\beta[A(N-N_1,V-V_1,T) - A(N,V,T)]\} \end{align}

保证归一化

∑N1=0N∫V1dp1dq1 ρ(p1,q1,N1)=1\begin{align} \sum_{N_1=0}^{N}\int_{V_1}\dd p_1 \dd q_1\, \rho(p_1,q_1,N_1) = 1 \end{align}

由于 N≫N1N \gg N_1 且 V≫V1V \gg V_1，

A(N−N1,V−V1,T)−A(N,V,T)=−N1μ+V1P\begin{align} A(N-N_1,V-V_1,T) - A(N,V,T) = -N_1 \mu + V_1 P \end{align}

其中 μ\mu 和 PP 分别是整个系统的化学势和压强，

μ=(∂A(N2,V,T)∂N2)N2=NP=−(∂A(N,V2,T)∂V2)V2=V\begin{align} \mu &= \left(\frac{\partial A(N_2,V,T)}{\partial N_2}\right)_{N_2 = N}\\ P &= -\left(\frac{\partial A(N,V_2,T)}{\partial V_2}\right)_{V_2 = V} \end{align}

注意 μ\mu 是一个宏观参数：

  • 它是强度量，μ=μ(N,V,T)=μ(V/N,T)\mu = \mu(N,V,T) = \mu(V/N,T)；
  • 平衡时 μ=μ(V1/N1,T)=μ(V2/N2,T)\mu = \mu(V_1/N_1,T) = \mu(V_2/N_2,T)。

略去下标 1，

ρ(p,q,N)=e−βPVeβμNN!e−βH(p,q,N),β=1kT\begin{align} \rho(p,q,N) = \frac{e^{-\beta PV} e^{\beta\mu N}}{N!}e^{-\beta H(p,q,N)}, \qquad \beta = \frac{1}{kT} \end{align}

在热力学极限下，允许 NN 的取值范围是 0≤N<∞0 \le N < \infty。这里 eβμN/N!e^{\beta\mu N}/N! 是加在正则分布 e−βH(p,q,N)e^{-\beta H(p,q,N)} 上的额外因子，这一点与微正则系综和正则系综都不同。换句话说，即使能量 H(p,q,N)H(p,q,N) 相同，概率 ρ(p,q,N)\rho(p,q,N) 也确实依赖于 NN。但非常重要的是，无论 μ\mu 为正还是为负，Gibbs 因子 1/N!1/N! 都会压低 NN 很大的概率。

除了 VV 和 TT，现在还多了一个独立宏观参数 Nˉ\bar N，即平均粒子数。如果 μ=0\mu = 0，Nˉ\bar N 就不是独立的宏观参数，它由 VV 和 TT 固定。

问题：正则系综中能量不守恒，这是否会带来一个额外的宏观参数？

答案：不会。正则系综中 NN 固定，TT 由热库给定，涨落的能量已经由 TT 和 Hamiltonian 决定，不会出现新的自变量。巨正则系综多出的参数 Nˉ\bar N 来自允许 NN 改变，与能量的涨落无关。

为了得到热力学，注意现在多了一个微观变量 NN。照例定义巨配分函数

Ξ(μ,V,T)≡∑N=0∞eβμNQN(V,T)\begin{align} \Xi(\mu,V,T) \equiv \sum_{N=0}^{\infty} e^{\beta\mu N} Q_N(V,T) \end{align}

则

βPV=log⁡Ξ(μ,V,T)\begin{align} \beta P V = \log \Xi(\mu,V,T) \end{align}

其中巨热力学势为 Φ(μ,V,T)≡−kTlog⁡Ξ(μ,V,T)=−PV\Phi(\mu,V,T) \equiv -kT \log \Xi(\mu,V,T) = -PV。

这看起来已经是状态方程，但 μ\mu 不便于实验测量。于是再计算平均粒子数，

Nˉ=⟨N⟩=1Ξ(μ,V,T)∑N=0∞NeβμNQN(V,T)=kT∂∂μlog⁡Ξ(μ,V,T)\begin{align} \bar N = \langle N \rangle = \frac{1}{\Xi(\mu,V,T)}\sum_{N=0}^{\infty} N e^{\beta\mu N} Q_N(V,T) = kT \frac{\partial}{\partial \mu}\log \Xi(\mu,V,T) \end{align}

因此 μ=μ(Nˉ,V,T)\mu = \mu(\bar N,V,T)。事实上，这个关系也可以通过在 μ\mu 的定义中把 NN 换成 Nˉ\bar N 得到。

问题：这一点可以证明吗？

提示：一个量的平均值等于它的最概然值。

在 βPV=log⁡Ξ\beta PV = \log\Xi 与 Nˉ=kT ∂μlog⁡Ξ\bar N = kT\,\partial_\mu\log\Xi 之间消去 μ\mu，就得到真正的状态方程 P=P(Nˉ,V,T)P = P(\bar N,V,T)。

这个系综里粒子数会涨落，涨落的大小再由一次求导给出。把 Nˉ=∂log⁡Ξ/∂(βμ)\bar N = \partial\log\Xi/\partial(\beta\mu) 再对 βμ\beta\mu 求导一次，

⟨N2⟩−⟨N⟩2=kT∂⟨N⟩∂μ\begin{align} \langle N^2 \rangle - \langle N \rangle^2 = kT \frac{\partial \langle N \rangle}{\partial \mu} \end{align}

因此

⟨N2⟩−⟨N⟩2⟨N⟩2∝1⟨N⟩→0\begin{align} \frac{\langle N^2 \rangle - \langle N \rangle^2}{\langle N \rangle^2} \propto \frac{1}{\langle N \rangle} \to 0 \end{align}

所以正则系综与巨正则系综也是等价的。对理想气体这给出 ⟨N2⟩−⟨N⟩2=Nˉ\langle N^2 \rangle - \langle N \rangle^2 = \bar N，用下面算出的显式关系 μ(Nˉ,V,T)\mu(\bar N,V,T) 可以验证。

内能仍定义为 Hamiltonian 的平均值。为计算内能，记 y=eβμy = e^{\beta\mu}，

U=−(∂∂βlog⁡Ξ(y,V,T))y,V,U=U(Nˉ,V,T)\begin{align} U = -\left(\frac{\partial}{\partial \beta}\log \Xi(y,V,T)\right)_{y,V}, \qquad U = U(\bar N,V,T) \end{align}

这里求导时固定 yy，所以 μ\mu 对 β\beta 的依赖没有被求导。于是

CV=(∂U∂T)V,NˉS=∫0TCVTdTA=U−TS\begin{align} C_V &= \left(\frac{\partial U}{\partial T}\right)_{V,\bar N}\\ S &= \int_0^T \frac{C_V}{T}\dd T\\ A &= U - TS \end{align}

与正则系综相比，多了一个宏观参数 μ\mu，或者说是 Nˉ\bar N。Nˉ\bar N 比 μ\mu 更方便实验测量。一旦固定 Nˉ\bar N，热力学量的计算与正则系综相同，例如可以由自由能 A(Nˉ,V,T)A(\bar N,V,T) 导出所有常见的量。

请注意，内能 UU 对 Nˉ\bar N 的依赖只通过 μ\mu。Hamiltonian 本身确实依赖于粒子数 NN。但 NN 是微观变量，在计算平均值时已经求和掉了。

练习：计算粒子数的相对涨落。

化学势

热力学

化学势 μ\mu 的定义是，在 TT 和 VV 固定时，粒子数改变 dNˉ\dd \bar N，自由能改变 μdNˉ\mu \dd \bar N。也就是说，μ\mu 是单个粒子携带的自由能。因此

dA=−PdV−SdT+μdNˉ\begin{align} \dd A = -P \dd V - S \dd T + \mu \dd \bar N \end{align}

再由 U=A+TSU = A + TS 可以推出第一定律的推广形式，

dU=−PdV+TdS+μdNˉ\begin{align} \dd U = -P \dd V + T \dd S + \mu \dd \bar N \end{align}

换句话说，如果 SS 和 VV 保持不变，μ\mu 就是一个粒子携带的能量，能量守恒必须计入粒子交换所转移的能量。

当 μ>0\mu > 0 时，在 TT、VV 固定下增加 Nˉ\bar N 会使 AA 增加 μdNˉ>0\mu \dd\bar N > 0，因此 μ\mu 倾向于把 Nˉ\bar N 推向更小的值，以便在 TT、VV 固定下降低自由能，反之亦然。

例 1：理想气体。

QN=VN(2πmkT)3N/2N!\begin{align} Q_N = \frac{V^N (2\pi mkT)^{3N/2}}{N!} \end{align} A(N,V,T)=−NkT(log⁡V+32log⁡(2πmkT)−log⁡N+1)=−NkT(log⁡VN+32log⁡(2πmkT)+1)\begin{align} A(N,V,T) &= -NkT\left(\log V + \frac{3}{2}\log(2\pi mkT) - \log N + 1\right)\\ &= -NkT\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT) + 1\right) \end{align}

这里用了 log⁡N!=Nlog⁡N−N\log N! = N\log N - N。于是

μ=(∂A∂N)V,T=−kT(log⁡VN+32log⁡(2πmkT)+1)−NkT(−1N)=−kT(log⁡VN+32log⁡(2πmkT))\begin{align} \mu &= \left(\frac{\partial A}{\partial N}\right)_{V,T}\\ &= -kT\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT) + 1\right) - NkT\left(-\frac{1}{N}\right)\\ &= -kT\left(\log\frac{V}{N} + \frac{3}{2}\log(2\pi mkT)\right) \end{align}

这里不能忘记配分函数 QNQ_N 中的因子 1/N!1/N!，它解决了 Gibbs 佯谬。否则例如会得到 μ∼log⁡V\mu \sim \log V 这样一个有问题的关系，而合理的应是 μ∼log⁡(V/N)\mu \sim \log(V/N)。没有这个因子也无法导出巨正则系综。

可以写成

eβμ=NV(2πmkT)3/2\begin{align} e^{\beta\mu} = \frac{N}{V(2\pi mkT)^{3/2}} \end{align}

因此 μ\mu 的正负取决于温度和气体的密度；当 N/V=(2πmkT)3/2N/V = (2\pi mkT)^{3/2} 时 μ\mu 为零。

在这个简单情形下可以直接验证，把 NN 换成 Nˉ\bar N 就得到 μ=μ(Nˉ,V,T)\mu = \mu(\bar N,V,T)。

例 2：NN 个相同的准独立粒子，QN=1N!Z(β,V)NQ_N = \frac{1}{N!}Z(\beta,V)^N，其中 Z(β,V)Z(\beta,V) 是单个粒子的配分函数，1/N!1/N! 是 Gibbs 因子。因此

−βA=log⁡QN≈−N(log⁡N−1)+Nlog⁡Z(β,V)\begin{align} -\beta A = \log Q_N \approx -N(\log N - 1) + N\log Z(\beta,V) \end{align}

于是

−βμ=−β∂A∂N=log⁡Z(β,V)N\begin{align} -\beta\mu = -\beta\frac{\partial A}{\partial N} = \log\frac{Z(\beta,V)}{N} \end{align}

另一方面，由最概然分布的理论，N=e−αZ(β,V)N = e^{-\alpha}Z(\beta,V)，由此认定 α=−βμ\alpha = -\beta\mu。

巨正则系综的一个简单论证：设一个大系统由许多相同的准独立子系统组成，每个子系统有 NN 个粒子，由 Hamiltonian H(p,q)H(p,q) 描述。按照 Maxwell-Boltzmann 分布，在 (p,q,N)(p,q,N) 找到某个子系统的概率分布为

ρ(p,q,N)∝e−α−βH(p,q)\begin{align} \rho(p,q,N) \propto e^{-\alpha - \beta H(p,q)} \end{align}

其中 α=−βμb\alpha = -\beta\mu_b，μb\mu_b 是大系统的化学势。另一方面，μb\mu_b 正是子系统的自由能。由于自由能是广延量，可以写成 μb=Nμ\mu_b = N\mu，μ\mu 是子系统的化学势。计入 Gibbs 因子，

ρ(p,q,N)∝1N!eβμN−βH(p,q,N)\begin{align} \rho(p,q,N) \propto \frac{1}{N!}e^{\beta\mu N - \beta H(p,q,N)} \end{align}

练习：导出理想气体的热力学，即状态方程，以及 UU、CVC_V、SS、AA。提示：要导出理想气体的熵，只用等容过程不够，还需要一个等温过程。

化学平衡

设在巨正则系综中有两种粒子，H=H1+H2H = H_1 + H_2，两种粒子之间的相互作用可以忽略，

Ξ=∑N1=0∞∑N2=0∞QN1QN2eβ(μ1N1+μ2N2)\begin{align} \Xi = \sum_{N_1=0}^{\infty}\sum_{N_2=0}^{\infty} Q_{N_1}Q_{N_2} e^{\beta(\mu_1 N_1 + \mu_2 N_2)} \end{align}

其中

QNi=e−βAi(Ni,V,T),A=∑i=12Ai,Ai=−kTlog⁡QNi,μi=(∂Ai∂Ni)T,V,Nj≠i\begin{align} Q_{N_i} = e^{-\beta A_i(N_i,V,T)}, \qquad A = \sum_{i=1}^{2} A_i, \qquad A_i = -kT \log Q_{N_i}, \qquad \mu_i = \left(\frac{\partial A_i}{\partial N_i}\right)_{T,V,N_{j\neq i}} \end{align}

例如在水和冰的混合物中，一个水分子可以变成冰分子，反之亦然。

假设粒子总数守恒，即 N=N1+N2N = N_1 + N_2 为常量。平衡时自由能取极小，

0=δA=∑iδAi=∑i∂Ai∂NiδNi=∑iμiδNi\begin{align} 0 = \delta A = \sum_i \delta A_i = \sum_i \frac{\partial A_i}{\partial N_i}\delta N_i = \sum_i \mu_i \delta N_i \end{align}

由粒子数守恒 0=δN=∑iδNi0 = \delta N = \sum_i \delta N_i，得到平衡条件

μ1δN1+μ2δN2=μ1δN1−μ2δN1=0\begin{align} \mu_1 \delta N_1 + \mu_2 \delta N_2 = \mu_1 \delta N_1 - \mu_2 \delta N_1 = 0 \end{align}

即 μ1=μ2(V/N,T)\mu_1 = \mu_2(V/N,T)。

一般地，可以有多种粒子，以及一个反应

ν1X1+ν2X2+⋯⇌ν1′Y1+ν2′Y2+⋯\begin{align} \nu_1 X_1 + \nu_2 X_2 + \cdots \rightleftharpoons \nu_1' Y_1 + \nu_2' Y_2 + \cdots \end{align}

把 YY 改写成 XX，ν′\nu' 改写成 ν\nu，则

∑i=1KνiXi=0\begin{align} \sum_{i=1}^{K} \nu_i X_i = 0 \end{align}

例如 2H2+O2⇌2H2O2H_2 + O_2 \rightleftharpoons 2H_2O，其中 ν1=2\nu_1 = 2，ν2=1\nu_2 = 1，ν3=−2\nu_3 = -2。

现在粒子数守恒写成

δN1ν1=δN2ν2=⋯=δNKνK=const≡δN\begin{align} \frac{\delta N_1}{\nu_1} = \frac{\delta N_2}{\nu_2} = \cdots = \frac{\delta N_K}{\nu_K} = \text{const} \equiv \delta N \end{align}

例如 δNH22=δNO21=δNH2O−2\frac{\delta N_{H_2}}{2} = \frac{\delta N_{O_2}}{1} = \frac{\delta N_{H_2O}}{-2}。

注意现在粒子数守恒不止给出一个方程，任意一对粒子都给出一个方程。于是让自由能 AA 取极值，

0=δA=∑i∂Ai∂NiδNi=∑iμiνiδN\begin{align} 0 = \delta A = \sum_i \frac{\partial A_i}{\partial N_i}\delta N_i = \sum_i \mu_i \nu_i \delta N \end{align}

得到

∑i=1Kμiνi=0\begin{align} \sum_{i=1}^{K} \mu_i \nu_i = 0 \end{align}

这就是化学平衡的一般条件，把化学势按上面的符号约定用化学计量系数加权，其和为零。

总结

正则系综：系统 1 与温度为 TT 的热库接触，热库被积掉。系统的能量只通过 Hamiltonian 进入分布。

→ρ(p,q)∝e−βH(p,q)\rightarrow \rho(p,q) \propto e^{-\beta H(p,q)}, β=1/kT\beta = 1/kT

→QN(V,T)=1N!∫d3Npd3Nq e−βH=e−βA\rightarrow Q_N(V,T) = \frac{1}{N!}\int \dd^{3N}p \dd^{3N}q\, e^{-\beta H} = e^{-\beta A}

→A=U−TS\rightarrow A = U - TS 作为特性函数，PP、SS、UU、GG 由求导得到，积分中的能量壳层已经消失。

→\rightarrow 能量涨落 ⟨(ΔE)2⟩=kT2CV∝N\langle(\Delta E)^2\rangle = kT^2 C_V \propto N，相对涨落 O(N−1/2)O(N^{-1/2})，因此正则系综与微正则系综等价。

→\rightarrow 应用：每个二次项贡献 12kT\frac{1}{2}kT 的能量均分；经典理想气体，计入 Gibbs 因子后 S(V,U)S(V,U) 与微正则系综相同。

Ising 模型：格子上的正则系综，用平均场近似求解；KC=1/dK_C = 1/d，相变附近 M=±3ΔK/KCM = \pm\sqrt{3\Delta K/K_C}，这个近似只在高维或长程相互作用下才有依据。

巨正则系综：允许粒子数涨落，ρ∝1N!eβμN−βH\rho \propto \frac{1}{N!}e^{\beta\mu N - \beta H}。

→Ξ(μ,V,T)=∑NeβμNQN\rightarrow \Xi(\mu,V,T) = \sum_N e^{\beta\mu N}Q_N, βPV=log⁡Ξ\beta PV = \log\Xi, Φ=−PV\Phi = -PV

→Nˉ=kT ∂μlog⁡Ξ\rightarrow \bar N = kT\,\partial_\mu \log\Xi，⟨(ΔN)2⟩=kT ∂μ⟨N⟩\langle(\Delta N)^2\rangle = kT\,\partial_\mu\langle N\rangle，相对涨落同样是 O(Nˉ−1/2)O(\bar N^{-1/2})，因此这个系综与正则系综等价。

→\rightarrow 化学势：dA=−PdV−SdT+μdNˉ\dd A = -P\dd V - S\dd T + \mu\dd\bar N，dU=−PdV+TdS+μdNˉ\dd U = -P\dd V + T\dd S + \mu\dd\bar N；化学平衡 ∑iνiμi=0\sum_i \nu_i\mu_i = 0。

正则系综的缺点是 NN 必须固定，因此要过渡到巨正则系综。代价是多了一个参数 μ\mu，它不能直接测量，通常要用 Nˉ\bar N 消掉。

练习

  • 用正则系综证明，若 H(ξ→±∞)→∞H(\xi \to \pm\infty) \to \infty，则
⟨ξ∂H∂ξ⟩=kT\begin{align} \left\langle \xi \frac{\partial H}{\partial \xi} \right\rangle = kT \end{align}

  • 用 Ising 模型的平均场近似，求 KCK_C 处及其附近的 ∂M/∂ΔK\partial M/\partial \Delta K、UU 和 CVC_V，并讨论它们在 KCK_C 处及其附近的行为。

  • 计算巨正则系综中粒子数的相对涨落。

  • 用巨正则系综导出理想气体的热力学，即状态方程，以及 UU、CVC_V、SS、AA。

  • 教材中的习题 7.1、7.2 和 7.3。

  • 讨论：Gibbs 佯谬

Comment评论

添加新评论