Monte Carlo 核心原理¶

Monte Carlo 程序看起来是在不断随机更新构型，但它的数学骨架其实很简单。给定构型空间中的状态 \(s\)，一次更新由转移概率

\[ P(s\to s') \]

描述。整个算法要满足两件事：

\[ \sum_{s'}P(s\to s')=1 \]

和

\[ W(s)=\sum_{s'}W(s')P(s'\to s). \]

第一行保证这是一个合法随机过程，第二行保证 \(W(s)\) 是这个随机过程的不变分布。

合法的随机过程¶

对固定出发点 \(s\)，下一步一定会跳到某个构型 \(s'\)，包括可能留在原地。因此所有去向的概率必须加起来等于 \(1\)：

\[ \sum_{s'}P(s\to s')=1. \]

在 Ising 模型中，一个构型可以写为

\[ s=\{s_1,s_2,\cdots,s_N\},\qquad s_i=\pm 1. \]

一次单 spin flip 更新会尝试翻转某个 \(s_i\)。如果新构型被接受，就从 \(s\) 跳到 \(s'\)；如果被拒绝，就留在原地。也就是说，\(P(s\to s)\) 不是可有可无的项，而是拒绝更新的概率。没有这项，转移概率通常无法归一化。

平稳分布¶

假设第 \(t\) 步系统处于构型 \(s\) 的概率是 \(W_t(s)\)。经过一次更新后，第 \(t+1\) 步处于 \(s\) 的概率为

\[ W_{t+1}(s)=\sum_{s'}W_t(s')P(s'\to s). \]

这句话的意思是：下一步到达 \(s\)，可以从任意旧构型 \(s'\) 跳过来。于是要把“原来在 \(s'\) 的概率”和“从 \(s'\) 跳到 \(s\) 的概率”相乘，再对所有 \(s'\) 求和。

如果目标权重 \(W(s)\) 满足

\[ W(s)=\sum_{s'}W(s')P(s'\to s), \]

那么一旦系统分布达到 \(W(s)\)，经过一次更新后仍然是 \(W(s)\)。这就是平稳分布条件，也常叫 global balance。

严格来说，采样分布应该是

\[ \pi(s)=\frac{W(s)}{Z}, \qquad Z=\sum_s W(s). \]

但把它代入平稳条件：

\[ \frac{W(s)}{Z} = \sum_{s'}\frac{W(s')}{Z}P(s'\to s). \]

两边乘以 \(Z\)，就回到

\[ W(s)=\sum_{s'}W(s')P(s'\to s). \]

因此 Monte Carlo 不需要知道归一化常数 \(Z\)。这正是它能够处理配分函数、后验分布和复杂能量模型的关键。

Detailed Balance¶

直接构造满足 global balance 的 \(P(s\to s')\) 通常很难。实际中常用更强但更容易检查的条件：

\[ W(s)P(s\to s')=W(s')P(s'\to s). \]

这叫 detailed balance。它的含义是：在平衡状态下，从 \(s\) 流向 \(s'\) 的概率流，等于从 \(s'\) 流向 \(s\) 的概率流。

如果 detailed balance 对任意 \(s,s'\) 成立，那么对 \(s'\) 求和：

\[ \sum_{s'}W(s)P(s\to s') = \sum_{s'}W(s')P(s'\to s). \]

利用概率归一化

\[ \sum_{s'}P(s\to s')=1, \]

左边变成 \(W(s)\)，于是自动得到

\[ W(s)=\sum_{s'}W(s')P(s'\to s). \]

所以 detailed balance 是构造正确 Monte Carlo 更新的常用充分条件。

Proposal 与 Acceptance¶

实际算法通常把一次转移拆成两步：

\[ P(s\to s')=A(s\to s')P_{\rm acc}(s\to s'). \]

其中 \(A(s\to s')\) 是 proposal probability，即“尝试从 \(s\) 走到 \(s'\)”的概率；\(P_{\rm acc}(s\to s')\) 是接受率。

detailed balance 要求

\[ W(s)A(s\to s')P_{\rm acc}(s\to s') = W(s')A(s'\to s)P_{\rm acc}(s'\to s). \]

一种标准选择是 Metropolis-Hastings 接受率：

\[ P_{\rm acc}(s\to s') = \min\left[ 1, \frac{W(s')A(s'\to s)}{W(s)A(s\to s')} \right]. \]

如果 proposal 对称，

\[ A(s\to s')=A(s'\to s), \]

则接受率简化为

\[ P_{\rm acc}(s\to s') = \min\left[ 1, \frac{W(s')}{W(s)} \right]. \]

这就是最常见的 Metropolis 公式。

Ising 单 Spin Flip¶

考虑经典 Ising 模型的 reduced Hamiltonian：

\[ \mathcal H(s)=-K\sum_{\langle ij\rangle}s_is_j. \]

构型权重为

\[ W(s)=e^{-\mathcal H(s)} = \exp\left(K\sum_{\langle ij\rangle}s_is_j\right). \]

如果只翻转一个 spin

\[ s_i\to -s_i, \]

只有它和最近邻之间的局域能量会改变。翻转前局域能量为

\[ \mathcal H_{\rm old} = -Ks_i\sum_{j\in\partial i}s_j, \]

翻转后为

\[ \mathcal H_{\rm new} = Ks_i\sum_{j\in\partial i}s_j. \]

因此能量差是

\[ \Delta = \mathcal H_{\rm new}-\mathcal H_{\rm old} = 2K s_i\sum_{j\in\partial i}s_j. \]

权重比为

\[ \frac{W(s')}{W(s)} = e^{-\Delta}. \]

若单 spin proposal 是对称的，Metropolis 接受率就是

\[ P_{\rm acc} = \min(1,e^{-\Delta}) = \min\left(1,e^{-2K s_i\sum_{j\in\partial i}s_j}\right). \]

这条公式背后的完整逻辑是：

概率归一化保证 \(P(s\to s')\) 是合法转移概率。
平稳分布条件保证采样目标是 \(W(s)\)。
detailed balance 给出构造转移概率的实用方法。
proposal 与 acceptance 分解让算法可以只用权重比。
对 Ising 单点更新，权重比只依赖局域能量差。

程序中的对象¶

从代码角度看，一个 Monte Carlo 程序通常由下面几类对象组成：

构型：保存当前状态，以及局域能量、邻接关系、粒子数、世界线片段等辅助信息。
更新器：给出 proposal，计算权重比或能量差，并按接受率更新构型。
观测量：测量能量、磁化率、关联函数、绕数、Binder ratio 等。
随机数产生器：提供均匀随机数、离散抽样、高斯抽样等。
数据输出与分析：记录参数、热化步数、采样数、误差、关联时间和随机种子。

程序设计的关键不是把所有功能堆在一起，而是让“构型如何存储”“更新如何改变构型”“观测量如何从构型读出”三件事保持清晰对应。