4. 单轴扭曲模型仿真

单轴扭曲（One-Axis Twisting，OAT）模型是自旋压缩理论中最经典、也最容易在仿真中复现的模型。它描述一个由 \(N\) 个等效二能级粒子组成的集体自旋，在非线性相互作用下绕某一固定轴发生“扭曲”的过程。理解这一模型的关键，不在于记住公式，而在于看清三件事：哈密顿量如何把相互作用强度写进相位、时间演化如何把初态“剪切”成压缩态、以及粒子数与相互作用强度分别从哪个方向影响压缩深度。

一、哈密顿量的建立与参照约定

本节统一采用如下约定：集体自旋算符 \(\hat{J}_x,\hat{J}_y,\hat{J}_z\) 满足 \([\hat{J}_i,\hat{J}_j]=i\varepsilon_{ijk}\hat{J}_k\)（取 \(\hbar=1\)），初始相干态沿 \(+x\) 方向极化，即 \(\langle\hat{J}_x\rangle=N/2\)，在 Bloch 球上位于赤道。单轴扭曲哈密顿量取为：

H_OAT = χ · Ĵ_z²

其中 χ 为非线性相互作用强度（单位：角频率），
Ĵ_z = (1/2) Σ σ_z^(k) 为集体 z 方向自旋算符。

这一形式只保留二次项，是因为在双模玻色—约瑟夫森结或原子内态干涉仪中，两模粒子数差对应的算符正是 \(\hat{J}_z\)，而相互作用能正比于粒子数差的平方。\(\hat{J}_z^2\) 的对角性意味着：在 \(z\) 基下，本征态 \(|m\rangle\)（\(m=-N/2,\dots,N/2\)）的能量为 \(\chi m^2\)。能量随 \(|m|\) 增大而升高，这正是“扭曲”的来源——不同 \(m\) 分量积累不同相位，相位差随 \(m\) 二次增长。

需要强调符号约定：若取 \(H=-\chi\hat{J}_z^2\)，则相位积累方向相反，压缩出现在另一正交分量上。仿真中必须固定一种约定并在全流程保持一致，否则会出现“压缩分量与预期相反”的假象。本节统一使用 \(H=+\chi\hat{J}_z^2\)，并约定在 \(t>0\) 后观察 \(\hat{J}_y\) 分量的方差压缩。

二、时间演化：从相干态到剪切椭圆

演化算符为 \(U(t)=e^{-iHt}=e^{-i\chi t\hat{J}_z^2}\)。在 \(z\) 基下它是对角算符，因此数值仿真最直接的做法是：把初态展开到 \(|m\rangle\) 基，乘以相位因子 \(e^{-i\chi t m^2}\)，再变换回所需表象。对初始相干态，其 \(z\) 基系数为二项式分布：

c_m(0) = sqrt( C(N, N/2+m) / 2^N )
c_m(t) = c_m(0) · exp(-i · χ · t · m²)

这一“相位二次调制”是 OAT 动力学的全部核心。直观上，Bloch 球赤道上的圆形不确定斑在演化中被剪切：原本各方位角分量相位一致，演化后相位随 \(m\) 二次变化，圆形斑被拉成倾斜椭圆。当剪切足够强时，椭圆短轴方向（即某一横向自旋分量）的方差降到散粒噪声极限以下，形成自旋压缩。

压缩程度用 Wineland 参数衡量：

ξ² = N · Var(Ĵ_⊥) / |⟨Ĵ_x⟩|²

其中 Ĵ_⊥ 取垂直于平均自旋方向、方差最小的横向分量。
ξ² < 1 表示压缩（相对散粒噪声极限），
以 dB 表示为 10·log10(ξ²)，负值即压缩量。

仿真时需注意：平均自旋方向在演化中会逐渐偏离 \(x\) 轴，因此不能固定用 \(\hat{J}_y\) 计算方差，而应先求 \(\langle\vec{J}\rangle\) 方向，再在垂直于它的平面内找最小方差分量。这是初学者最容易出错的一步。

三、粒子数 N 与相互作用强度 χ 的影响

在 OAT 模型中，最佳压缩时刻与压缩深度都有明确的标度关系。对初始相干态，短时近似下最小方差出现在 \(\chi t \approx N^{-2/3}\) 量级，对应最佳压缩参数 \(\xi^2_{\min}\approx N^{-2/3}\)，即压缩量约 \(10\log_{10}(N^{-2/3})\) dB。这意味着：

• 粒子数 N 增大：压缩极限按 \(N^{-2/3}\) 改善，但最佳压缩时刻 \(t_{\text{opt}}\propto N^{-2/3}\) 提前，对时间分辨率要求更高。例如 N=100 时理论极限约 −13.3 dB，N=1000 时约 −20 dB，但后者最佳时刻缩短约 4.6 倍。
• 相互作用强度 χ 增大：只改变时间尺度，不改变压缩极限。χ 翻倍，最佳时刻减半，\(\xi^2_{\min}\) 基本不变。因此 χ 是“快慢旋钮”，N 是“深度旋钮”。
• 过压缩：超过最佳时刻后，椭圆继续旋转并越过短轴方向，方差回升，甚至出现“压缩—反压缩”振荡。仿真中必须扫描时间找极小值，而不是取固定终点。

四、可执行的仿真流程与判读

下面给出一个可直接实现的数值流程（以 Python/NumPy 为例，示例参数 N=100、χ=1、扫描 200 个时间点）：

import numpy as np
N = 100
m = np.arange(-N//2, N//2 + 1)
coef0 = np.sqrt(np.array([np.math.comb(N, int(k + N/2)) for k in m]) / 2**N)
chi = 1.0
ts = np.linspace(0, 0.05, 200)
xi2 = []
for t in ts:
    coef = coef0 * np.exp(-1j * chi * t * m**2)
    # 构造 Jx, Jy, Jz 在 m 基下的矩阵并求期望与方差
    # 取平均自旋方向，在垂直平面内求最小方差
    ...
    xi2.append(N * var_min / mean_x**2)
t_opt = ts[np.argmin(xi2)]
print(t_opt, 10*np.log10(min(xi2)))

判读要点有三：其一，检查 \(\xi^2(t)\) 曲线是否先降后升，若单调下降说明时间窗口太短；其二，确认最小方差分量确实垂直于平均自旋，而非固定取 \(\hat{J}_y\)；其三，用 \(N^{-2/3}\) 标度做交叉验证——改变 N 后，\(t_{\text{opt}}\cdot N^{2/3}\) 应近似为常数（在 \(\chi\) 固定时）。若这一关系明显偏离，通常意味着基矢截断、归一化或算符构造有误。

把这三层关系串起来看：哈密顿量决定“相位如何随 m 二次积累”，时间演化决定“圆斑如何被剪成椭圆”，而 N 与 χ 分别控制“椭圆能压多扁”和“多快压到最扁”。掌握这三点，后续加入损耗、探测噪声或双轴扭曲项时，就能清楚判断每一处修改作用在链条的哪一环。