来源论文: https://arxiv.org/abs/2604.02769v1 生成时间: Apr 06, 2026 09:49

哈伯德模型的半经典表示：基于非常规相干态的深度解析

0. 执行摘要

哈伯德模型作为强关联电子系统的基石，因其复杂性而难以精确求解。传统的近似方法，如摄动论、变分法和基于哈伯德-斯特拉托诺维奇（Hubbard-Stratonovich, HS）变换的半经典方法，在捕捉某些关键物理现象（如非局部关联和纠缠序）方面存在局限性。本文所解析的论文提出了一种创新的半经典近似方案，其核心在于采用一种非常规的相干态表示，旨在最小化格拉斯曼变量的使用。该方法将局部希尔伯特空间中的费米子奇偶性编码为一个单一的格拉斯曼变量，而自旋和电荷自由度则由玻色子变量描述。这种构造策略带来了多重优势：它提供了一种非微扰、适用于有限温度、并能够自然地扩展到多轨道系统的框架，同时能够更好地描述同时涉及自旋和电荷对称性破缺的纠缠序。

论文通过与一格点和二格点哈伯德模型的精确解进行比较，对该方法的有效性进行了严格评估。结果表明，这种半经典方法在定性上成功复现了粒子数、双占据、跳跃振幅和自旋关联等关键物理量的行为。尽管在定量上存在一些偏差（主要源于底层态密度从离散变为连续的特性），但其在处理强关联系统方面展现了显著的潜力和独特性。此外，论文还深入探讨了这种相干态构造与非线性费米子变换之间的精确联系，揭示了哈伯德模型表示的新视角，并为将其推广到大M极限和映射到马约拉纳-近藤格点模型提供了理论基础。

1. 核心科学问题，理论基础，技术难点，方法细节

1.1 核心科学问题

哈伯德模型（Hubbard model）是凝聚态物理领域研究强关联电子系统最基础且最重要的模型之一。其哈密顿量，看似简洁，却能捕捉从费米液体行为、磁性、莫特金属-绝缘体转变到高温超导等一系列丰富而复杂的量子现象。然而，作为一个真正的量子多体问题，哈伯德模型的精确求解，无论是解析还是数值，都极具挑战性。其核心难点在于电子之间的局部库仑相互作用（U项）导致的多体纠缠效应，使得单粒子近似失效。因此，开发有效的近似方法来理解其物理性质变得至关重要。

现有的一些近似方法，如在相互作用强度U或跳跃振幅t上进行的摄动展开，以及基于试探波函数的变分法，都各有其适用范围和局限性。特别是在半经典近似中，通常采用哈伯德-斯特拉托诺维奇（Hubbard-Stratonovich, HS）变换引入辅助玻色子场，然后忽略这些场的虚时依赖性（即将其视为静态场）。这种标准的HS变换存在一定的任意性，并且由于辅助场是局域的，其静态处理无法捕捉非局部关联，例如d波超导所需的长程关联。这促使研究人员重新思考相干态的构造，以期开发出能够克服这些限制的更普适、更受控的近似方案。

1.2 理论基础：非常规相干态表示

论文的核心创新在于提出了一种基于“分级相干态”（Graded Coherent States）的非常规路径积分表述。其核心思想是最小化格拉斯曼（Grassmann）变量的使用，以简化路径积分的复杂性，并更好地处理局部希尔伯特空间的自旋和电荷自由度。

相干态的构造：

对于每个格点i，局部希尔伯特空间包含四个基态：空位$|0\rangle$、自旋向上$|\uparrow\rangle$、自旋向下$|\downarrow\rangle$和双占据$|\uparrow\downarrow\rangle$。传统的费米子相干态通常涉及多个格拉斯曼变量。而论文提出的分级相干态，表示为$|\Omega_i, \psi_i\rangle$，只引入了一个单一的格拉斯曼变量$\psi_i$来描述费米子的宇称（即奇偶性），而自旋和电荷扇区的方向则由一组玻色子参数$\Omega_i$来编码。具体形式如下：

$|\Omega_i, \psi_i\rangle = C_1(\Omega_i)|0\rangle - C_2(\Omega_i)\psi_i c_{i\downarrow}^\dagger|0\rangle - C_3(\Omega_i)\psi_i c_{i\uparrow}^\dagger|0\rangle + C_4(\Omega_i)\psi_i c_{i\downarrow}^\dagger c_{i\uparrow}^\dagger|0\rangle = |\alpha_i\rangle + \psi_i|\beta_i\rangle$

其中，$|\alpha_i\rangle$代表偶费米子宇称扇区（空位和双占据），$|\beta_i\rangle$代表奇费米子宇称扇区（单占据）。系数$C_a(\Omega_i)$是玻色子参数$\Omega_i$的函数。

宇称分离与SU(2)对称性：

这种构造使得格拉斯曼变量$\psi_i$能够充当偶宇称和奇宇称扇区之间的媒介。通过施加特定的归一化条件：

$\langle\Omega_i, \psi_i|\Omega_i, \psi_i\rangle = e^{\psi_i\psi_i}$，

这进一步限制了系数，得到$|C_1(\Omega_i)|^2 + |C_4(\Omega_i)|^2 = |C_2(\Omega_i)|^2 + |C_3(\Omega_i)|^2 = 1$。这种结构精确反映了半满单轨道哈伯德模型中电荷（Holon-Doublon）和自旋扇区的SU(2)η × SU(2)spin对称性。因此，玻色子自由度$\Omega_i$可以被参数化为两个独立的布洛赫球面（Bloch spheres），分别对应于电荷和自旋自由度。

费米子算符的分解：

为了在新的相干态基底下处理哈密顿量，原始的费米子产生/湮灭算符$c_{i\sigma}$被分解为$

c_{i\sigma} = \xi_{i\sigma} + \pi_{i\sigma}$

其中$\xi_{i\sigma}$作用于偶宇称态$|\alpha_i\rangle$，而$\pi_{i\sigma}$作用于奇宇称态$|\beta_i\rangle$。这种分解与Ref. [29]中引入的复合费米子算符紧密相关，其中$\xi_{i\sigma}$描述下哈伯德带（单占据位）上的电子，而$\pi_{i\sigma}$描述上哈伯德带（双占据位）上的电子。

路径积分表述：

通过插入完备性关系：

$1 = \mathcal{C} \int d\Omega_i d\psi_i e^{-\bar{\psi}_i\psi_i} |\Omega_i, \psi_i\rangle\langle\Omega_i, \psi_i|$

并在每个虚时切片上应用铃木-特罗特（Suzuki-Trotter）分解，可以构建路径积分。关键在于，这种相干态表示下，哈伯德模型的相互作用项$U n_{i\uparrow} n_{i\downarrow}$（其中$n_{i\sigma} = c_{i\sigma}^\dagger c_{i\sigma}$）通过特定的数学操作和格拉斯曼变量的代数性质，最终在有效作用量中变为格拉斯曼场（$\psi_i$）的双线性形式。这与标准HS变换引入辅助玻色子场来处理四费米子相互作用项形成了鲜明对比，极大地简化了格拉斯曼变量的积分。

1.3 技术难点与半经典近似细节

尽管非常规相干态的引入简化了路径积分中格拉斯曼变量的处理，但玻色子变量$\Omega_i$的虚时依赖性（即量子涨落）仍然使问题难以求解。为了获得可处理的框架，论文引入了半经典近似。

半经典近似：

核心思想是忽略玻色子变量$\Omega_i$的虚时依赖性，即假设$\Omega_i(\tau) \rightarrow \Omega_i$。这意味着玻色子场被视为静态场，其量子涨落被忽略。这种近似类似于量子自旋模型中的大S展开，或标准HS变换后的静态辅助场近似。

有效作用量与玻色子变量积分：

在半经典近似下，玻色子变量$\Omega_i$的虚时导数项消失，有效作用量$S_{eff}$简化为：

$S_{eff} = \int d\tau \sum_i \bar{\psi}_i \partial_\tau \psi_i + \sum_i (\dots) + \int d\tau \sum_{(ij)\sigma} (\dots)$

其中，玻色子变量$\Omega_i$的虚时导数项消失，只留下与玻色子场本身相关的项和与跳跃项相关的项。更重要的是，格拉斯曼场部分的有效作用量现在完全是二次的。这意味着对于给定的玻色子场配置$\Omega = \{\Omega_i\}$，格拉斯曼场可以被精确积分。积分后的结果是一个有效的单体费米子哈密顿量，其形式为：

$H_{eff}[\Omega, d] = \sum_i [\mu + (U - 2\mu)|C_4(\Omega_i)|^2] d_i^\dagger d_i + \sum_{(ij)} (T_{ij} d_j^\dagger d_i + \Delta_{ij} d_i^\dagger d_j^\dagger + \text{H.c.})$

其中$d_i$是湮灭一个无自旋费米子的算符。跳跃矩阵$T_{ij}$和配对势$\Delta_{ij}$都由玻色子变量$\Omega_i$的系数$C_a(\Omega_i)$决定。这些系数通过参数化$\Omega_i$（例如，利用两个布洛赫球的极角和方位角$\theta_{1i}, \phi_{1i}, \theta_{2i}, \phi_{2i}$）来确定。

求解过程：

  1. 对于固定的玻色子场配置$\Omega$： 求解上述单体费米子哈密顿量$H_{eff}[\Omega, d]$是一个标准问题。这可以通过对$H_{eff}[\Omega, d]$进行对角化来获得其本征值和本征态。一旦得到，费米子部分的配分函数$Tr[e^{-\beta H_{eff}[\Omega, d]}] $可以精确计算。

  2. 对玻色子场进行积分： 最终的配分函数$Z_{cl}$通过对所有可能的玻色子场配置$\Omega$进行积分来获得：

    $Z_{cl} = \int d\Omega e^{-\beta(U-2\mu) \sum_i |C_4(\Omega_i)|^2} Tr[e^{-\beta H_{eff}[\Omega, d]}]$

    由于玻色子场的维度可能很高，对于格点数量超过一或二的系统，这个积分通常需要通过蒙特卡罗（Monte Carlo）模拟来完成。

物理量的计算：

任何物理量$A$的期望值$\langle A \rangle$都可以通过对玻色子场和费米子场进行平均来计算。在半经典近似下，费米子场的期望值可以利用维克定理（Wick’s theorem）精确计算，而玻色子场的平均则通过蒙特卡罗采样完成。

1.4 大M构造与多轨道泛化

受控的半经典极限（大M构造）：

为了将上述半经典近似置于更严格的理论基础之上，论文提出了一种“大M构造”。该方法通过将局部希尔伯特空间复制M次，然后投影到完全对称扇区，从而构建一系列层次化的哈伯德模型。在这种大M极限下，玻色子场的量子涨落被1/M的因子抑制，使得半经典近似渐近精确。这类似于量子自旋模型中的大S展开，其中1/M扮演了与$\hbar_{eff}$相似的角色。

多轨道系统泛化：

该非常规相干态构造可以自然地推广到多轨道系统。在多轨道情况下，局部希尔伯特空间维度更高。论文提出，仍然可以通过引入L个格拉斯曼变量来描述费米子奇偶性，其中L可以从1（最小格拉斯曼变量数）到N（标准费米子相干态的格拉斯曼变量数）变化，N是自旋轨道味的数目。对于每个轨道块（block），可以关联一个格拉斯曼变量，从而保持费米子奇偶性的分级结构，同时用玻色子参数描述自旋和电荷自由度。

1.5 到马约拉纳-近藤格点模型的映射

论文还展示了该方法的有效作用量，在不进行半经典近似的情况下，可以精确地映射到一个包含无自旋费米子$d_i$和局部赝自旋算符$Q_i$的哈密顿量。通过引入马约拉纳费米子$y_{1i}, y_{2i}$来表示$d_i$和$d_i^\dagger$，以及利用马约拉纳表示来处理量子自旋$Q_i$，该模型最终可以被重写为一个马约拉纳-近藤格点模型。这个模型描述了巡游的马约拉纳费米子与局部量子自旋之间的耦合，其中甚至出现了三体相互作用项。更进一步的代数变换表明，这个马约拉纳-近藤格点模型最终可以变回原始的哈伯德模型。这种映射不仅为理解哈伯德模型中的强关联物理提供了新的视角，也揭示了不同物理模型之间的深层联系。

2. 关键 benchmark 体系，计算所得数据，性能数据

论文通过与哈伯德模型在一格点（原子极限）和二格点系统上的精确解进行比较，全面评估了所提出的半经典近似方法的性能。这些基准测试体系虽然简单，但足以揭示该方法在定性上捕获关键物理行为的能力以及在定量上存在的局限性。

2.1 一格点哈伯德原子

在一格点情况下，系统不存在跳跃项，因此哈密顿量只包含化学势和局域相互作用项。这是最简单的哈伯德模型，其精确解可以解析获得。

配分函数 $Z_{cl}(\beta; U', \mu)$：

  • 半经典结果： 经过对玻色子变量$\Omega$的积分和对d费米子的追踪，半经典配分函数被解析地计算出来，其表达式包含一系列指数和代数项。例如，在$U' = U - 2\mu$的半经典形式中，它为$(4\pi)^2 \frac{1-e^{-\beta U'}}{(\beta U')^2} + e^{\beta\mu}$。
  • 精确解： 精确配分函数为$1 + 2e^{\beta\mu} + e^{-\beta U'}$，对应于空位、单占据（自旋向上或向下）和双占据这三种态的玻尔兹曼因子之和。这三项明确反映了系统的离散能级结构。
  • 比较： 半经典结果在形式上与精确解相似，但在参数$U'$不为零时，其函数形式略有不同。最显著的区别在于，半经典解的渐近行为（例如当$\delta\mu \rightarrow +\infty$时）呈代数衰减，而精确解呈指数衰减。这直接反映了半经典近似中态密度的连续性与精确解中离散能级的根本差异。

态密度 $W_{cl}(E)$：

  • 半经典结果： 通过对半经典配分函数进行逆拉普拉斯变换，得到的态密度是连续且均匀分布的。具体地，其由两个Heaviside步函数（阶跃函数）和Dirac delta函数组成，覆盖了$E \in [0, U']$（对于$U' > 0$）或$E \in [U', 0]$（对于$U' < 0$）的能量范围。这直接反映了半经典变量$\Omega$的连续性质。
  • 精确解： 精确态密度由一系列Dirac delta函数组成，对应于系统的离散能级：$2\delta(E + \mu) + \delta(E) + \delta(E - U')$。
  • 比较： 尽管在形式上存在连续与离散的根本差异，但半经典和精确态密度的总谱权重是守恒的（除去一个归一化因子）。特别是在半满填充（$U' = 0$）时，半经典和精确态密度在形式上完全一致，都简化为两个$2\delta(E + \mu) + 2\delta(E)$峰。

电子数 $\langle n \rangle$：

  • 计算： 电子数通过对局部数算符$n_{i\sigma} = c_{i\sigma}^\dagger c_{i\sigma}$在半经典近似下进行求值获得。结果在定性上（例如随化学势$\delta\mu$的增加而增加，或随温度T的变化）与精确解一致。图2（上排）展示了填充数对化学势的依赖性，半经典和精确结果都呈现出类似的阶梯状行为。
  • 定量差异： 与配分函数和态密度类似，在定量上，半经典和精确结果之间存在差异。特别是在$\delta\mu \rightarrow +\infty$时，半经典电子数呈代数衰减到2（$2 - 1/(\beta\delta\mu')$），而精确解呈指数衰减到2（$2 - 2e^{-2\beta\delta\mu'}$）。这再次归因于半经典近似中连续态密度与精确解中离散能级的区别。

双占据 $\langle D \rangle$：

  • 计算： 双占据算符$D_i = n_{i\uparrow} n_{i\downarrow}$的期望值也在半经典近似下计算。结果在半满填充时（$\mu = U/2$ 或 $U' = 0$）与精确解完全一致。
  • 定性行为： 图2（下排）显示，在大多数情况下，半经典结果定性地捕捉了双占据随填充数的变化。然而，在高温且$U \neq 0$的情况下，半经典结果在$n < 1$时出现了非单调行为，这是一个由于连续态密度导致的人为产物。虽然这个伪影在低温下会消失且幅度较小，但它凸显了半经典近似的局限性。
  • U=0时的关联： 当$U=0$时，精确解满足$\langle D \rangle_{exact} = \langle n \rangle_{exact}^2 / 4$，反映了无相互作用体系中没有关联。但半经典近似下不满足此关系，存在一个人工关联，偏差范围在$0 \le \langle D \rangle - \langle n \rangle^2 / 4 < 0.07957$。

单粒子格林函数 $G_0(\tau)$：

  • 动态量： 作为动态量的代表，格林函数在半经典近似中计算。在$\tau=0$处的跳跃（$G_0(-0) - G_0(+0)$）在精确解中为1，但在半经典近似中小于1，存在“权重损失”。此外，半经典近似中$C_{1,2,3,4}$的$\tau$依赖性被忽略，导致格林函数的$\tau$依赖性无法完全反映量子力学行为。

2.2 二格点哈伯德模型

对于二格点系统，由于存在跳跃项$t \neq 0$，系统展现出更复杂的关联效应。此时，玻色子场$\Omega_i$的积分需要通过蒙特卡罗模拟来完成。

模拟算法：

蒙特卡罗模拟采用Metropolis算法更新玻色子场配置$\Omega_i$。在每个蒙特卡罗步骤中，对提议的新配置$\Omega_i'$，需要计算相应的有效哈密顿量$H_{eff}[\Omega, d]$，并精确求解d费米子问题以获得其配分函数$Tr[e^{-\beta H_{eff}[\Omega, d]}] $。此费米子迹项与玻色子场的有效作用量结合，共同决定Metropolis接受率。这种方法确保了即使玻色子场是半经典的，费米子部分也得到了量子处理。

自旋关联 $\langle S_1 \cdot S_2 \rangle$：

  • 定性行为： 图3（上排）展示了二格点系统在半满填充下自旋关联随温度T的变化。半经典结果和精确解在高温度区域表现出相似的趋势。
  • 定量差异与能量尺度： 在低温度区域，精确解渐近地呈指数衰减到零温极限，而半经典结果则呈幂律行为。在$U \rightarrow \infty$的极限下，哈伯德模型映射到海森堡自旋模型。精确解中，两个量子自旋$S=1/2$形成单线态，零温极限下$\langle S_1 \cdot S_2 \rangle_{exact} = -3/4$。然而，在半经典描述中，自旋被视为长度为$S=1/2$的经典矢量，导致$\langle S_1 \cdot S_2 \rangle_{cl} = -1/4$。这表明半经典近似在自旋关联的幅度上存在约1/3的损失，且特征能量尺度（$J = 4t^2/U$）也被削减（半经典$J_{cl} = t^2/U = J/4$）。图4中通过对半经典结果的轴进行缩放（垂直轴乘以3，水平轴乘以16）后，可以使半经典与精确解的特征能量尺度视觉上对齐，进一步证实了能量尺度的差异。
  • 填充依赖性： 图5（上排和中排）展示了自旋关联随填充数n或化学势$\delta\mu$的变化。半经典和精确结果都显示在半满填充$n=1$处有一个峰值，随后随着填充偏离半满而下降。然而，半经典结果的幅度仍然小于精确解。在$n=1$附近，半经典结果在$T=0.1$时表现出一种“奇异”行为，这再次与连续态密度的性质以及有限温度玻尔兹曼因子对这些连续态的采样有关。

跳跃振幅 $\sum_\sigma \langle c_{i\sigma}^\dagger c_{j\sigma} + \text{H.c.} \rangle$：

  • 定性行为： 图3（下排）展示了跳跃振幅随温度T的变化。与自旋关联类似，半经典和精确结果在高温度区域表现出相似的行为，而在低温度区域则存在幂律与指数行为的差异。
  • 定量差异： 在非相互作用极限$U \rightarrow 0$时，精确解的跳跃振幅渐近趋近于2，而半经典近似趋近于1。这种差异可归因于半经典近似中跳跃由无自旋费米子携带，而精确解中跳跃由自旋费米子携带，后者具有两个内部自由度。这意味着半经典近似低估了电子在格点间的离域程度。
  • 填充依赖性： 图5（下排）展示了跳跃振幅随填充数n的变化。在小U时，跳跃振幅在$n=1$处增强；在大U且半满时，强库仑斥力抑制跳跃，导致跳跃振幅较小。随着电子或空穴掺杂，跳跃振幅逐渐增加。半经典近似能够捕捉这些定性趋势，但在幅度上仍然小于精确解。

2.3 性能数据（未提供）

值得注意的是，论文中并未提供关于计算性能的详细数据，例如模拟所需的CPU/GPU时间、内存消耗或可扩展性基准。它主要关注理论框架的建立和与精确解的定性比较。对于更复杂的系统，例如二维格点，蒙特卡罗模拟的计算成本将是一个重要考虑因素。论文提到“更复杂的更新算法可以提高模拟效率”，暗示了当前实现可能在计算效率上面临挑战，但没有给出具体的性能提升数据。

论文主要聚焦于理论框架的构建和概念验证，并未提供具体的代码实现细节、所用软件包清单或任何开源代码仓库链接。因此，本节将侧重于根据论文描述，对该方法的概念性实现步骤进行详细阐述，以指导潜在的复现工作。

3.1 核心算法概述

该方法的计算流程主要分为两个阶段：一是玻色子场的蒙特卡罗采样，二是对每个玻色子场配置下费米子问题的精确求解。

  1. 参数初始化：

    • 定义系统参数：跳跃积分$t$、相互作用强度$U$、化学势$\mu$和温度$T$。确定格点数量$N_{site}$。
    • 为每个格点$i$初始化玻色子场$\Omega_i$。$\Omega_i$通常通过两个布洛赫球的极角$\theta_{1i}, \theta_{2i}$和方位角$\phi_{1i}, \phi_{2i}$进行参数化，这些角度可以在$[0, \pi]$和$[0, 2\pi)$范围内随机初始化。这些角度定义了相干态中的系数$C_1, C_2, C_3, C_4$。

  2. 蒙特卡罗（Metropolis）循环：

    • 迭代步骤： 重复进行一系列蒙特卡罗步骤，直至系统达到平衡并收集足够多的样本。
    • 玻色子场提议： 对于当前配置$\Omega = \{\Omega_i\}$，随机选择一个格点$i$，并对该格点的$\Omega_i$提出一个小的随机扰动，得到新的提议配置$\Omega_i'$。例如，可以对$\theta_{1i}, \theta_{2i}, \phi_{1i}, \phi_{2i}$施加小的随机增量。
    • 构建有效费米子哈密顿量：
      • 根据当前提议的玻色子场配置$\Omega'$（或旧配置$\Omega$），计算每个格点$i$的系数$C_a(\Omega_i)$。
      • 利用这些$C_a(\Omega_i)$系数，根据论文公式（48）和（49）计算有效的跳跃矩阵$T_{ij}$和配对势$\Delta_{ij}$。
      • 构造一个$2N_{site} \times 2N_{site}$的单体费米子哈密顿量矩阵$H_{eff}[\Omega', d]$。这个哈密顿量描述的是无自旋费米子$d_i$。
    • 求解费米子问题：
      • 对$H_{eff}[\Omega', d]$矩阵进行对角化，得到其本征值$\epsilon_k$（$k=1, \dots, 2N_{site}$）。
      • 计算费米子部分的配分函数$Z_d(\Omega') = Tr[e^{-\beta H_{eff}[\Omega', d]}] = \prod_k (1 + e^{-\beta\epsilon_k})$。
    • 计算总权重与接受/拒绝采样：
      • 计算玻色子场部分的有效作用量$S_{eff, \Omega}(\Omega') = \beta(U - 2\mu) \sum_i |C_4(\Omega_i')|^2$。
      • 计算总权重$W(\Omega') = Z_d(\Omega') e^{-S_{eff, \Omega}(\Omega')}$。
      • 使用Metropolis-Hastings准则决定是否接受新配置$\Omega'$。如果$W(\Omega') / W(\Omega) \ge 1$，则接受；否则以概率$W(\Omega') / W(\Omega)$接受。如果接受，则$\Omega \leftarrow \Omega'$。
    • 收集数据： 在达到平衡后，开始收集玻色子场配置$\Omega$的样本。对于每个接受的配置，计算所需的物理量（如粒子数、双占据、自旋关联、跳跃振幅），通过对$H_{eff}[\Omega, d]$的本征向量求和来获得费米子部分的期望值，然后进行累加平均。

3.2 复现指南与关键数学步骤

复现这项工作需要对论文中的数学推导有深刻理解，并细致实现以下几个关键部分：

  1. 相干态系数 $C_a(\Omega_i)$ 的参数化：

    • 论文使用$\theta_{1i}, \phi_{1i}, \theta_{2i}, \phi_{2i}$这四个角度来参数化$C_1, C_2, C_3, C_4$。例如，$C_1(\Omega_i) = \cos(\theta_{1i}/2)$，$C_4(\Omega_i) = e^{i\phi_{1i}}\sin(\theta_{1i}/2)$，以及$C_2(\Omega_i) = \cos(\theta_{2i}/2)$，$C_3(\Omega_i) = e^{i\phi_{2i}}\sin(\theta_{2i}/2)$（注意论文中$C_1$和$C_4$与$C_2$和$C_3$的下标不同，具体对应关系需参考公式(53)）。确保这些系数满足归一化条件。

  2. 构建跳跃矩阵 $T_{ij}$ 和配对势 $\Delta_{ij}$：

    • 根据公式(48)和(49)，精确计算$T_{ij}$和$\Delta_{ij}$。这些矩阵元素由$C_a(\Omega_i)$系数的组合以及玻色子场之间的乘积构成。例如， $T_{ij} = -t \sum_\sigma (C_{1i}^*(C_{2i})C_{3j}^*C_{2j} - C_{4i}^*C_{3i}C_{2j}^*C_{3j})$ $\Delta_{ij} = -t \sum_\sigma ((-1)^\sigma C_{2i}^*C_{1i}C_{3j}^*C_{4j} - C_{3i}^*C_{4i}C_{2j}^*C_{1j})$

  3. 单体费米子哈密顿量矩阵构建：

    • 将$H_{eff}[\Omega, d]$表示为$2N_{site} \times 2N_{site}$的矩阵形式。例如，如果用$(d_1, \dots, d_{N_{site}}, d_1^\dagger, \dots, d_{N_{site}}^\dagger)^T$作为基态向量，哈密顿量矩阵的块结构将包括$T_{ij}$和$\Delta_{ij}$。

  4. 对角化和费米子配分函数：

    • 使用标准数值线性代数库对$2N_{site} \times 2N_{site}$的复Hermitian矩阵进行对角化。
    • 根据费米子本征值$\epsilon_k$计算$Z_d(\Omega')$。这个步骤至关重要，因为$Z_d$是玻色子场采样的权重因子。

  5. 物理量计算：

    • 粒子数$\langle n_{i\sigma} \rangle$、双占据$\langle D_i \rangle$、自旋关联$\langle S_i \cdot S_j \rangle$和跳跃振幅等物理量，需根据论文公式（32）-（37）进行计算。这些量通常是$C_a(\Omega_i)$系数和费米子格林函数（$\langle d_k^\dagger d_l \rangle$）的组合。

论文中未提及任何特定的软件包或开源代码库。 这一点是需要特别注意的。因此，任何复现工作都需要从头开始实现上述算法，或利用通用的数值计算库。可能的通用数值计算工具包括：

  • Python生态系统：
    • NumPy 和 SciPy：用于矩阵操作、对角化和基本数值计算。特别是SciPy.linalg模块提供了高效的线性代数函数。
    • Matplotlib：用于数据可视化和绘图。
    • Numba 或 Cython：为了提高蒙特卡罗模拟中计算密集的部分的性能（例如矩阵对角化），可能需要对Python代码进行JIT编译。
  • C++生态系统：
    • Eigen 或 Armadillo：高性能C++模板库，用于线性代数运算和矩阵对角化。
    • Boost.Random：用于生成高质量的随机数。
  • Fortran生态系统：
    • LAPACK 或 BLAS：底层的线性代数库，提供矩阵对角化等核心功能。

开源仓库链接：

由于论文未提供，因此目前没有可以直接引用的开源仓库。若要复现，研究者需自行构建。考虑到蒙特卡罗模拟和矩阵对角化是计算瓶颈，选择合适的语言和库对于高效实现至关重要。

3.4 性能考量

  • 计算成本： 二格点系统的蒙特卡罗模拟已经需要一定的计算资源。对于更大尺度的系统（例如二维哈伯德模型），玻色子场的采样空间呈指数增长。每一步蒙特卡罗都需要对$2N_{site} \times 2N_{site}$矩阵进行对角化，其计算复杂度为$O((2N_{site})^3)$。因此，对于数百甚至数千个格点的系统，这将是巨大的计算挑战。
  • 收敛性： 蒙特卡罗模拟的平衡和收敛速度对于结果的可靠性至关重要。特别是在低温度或存在相变时，系统可能会遭遇临界慢化，导致需要大量的蒙特卡罗步才能获得统计上可靠的样本。
  • 高级采样算法： 论文提及“更复杂的更新算法可以提高模拟效率”。这可能包括使用Hybrid Monte Carlo (HMC)、并行回火（Parallel Tempering）或其他高级蒙特卡罗技术来克服局部最小值和加速收敛。

4. 关键引用文献，以及你对这项工作局限性的评论

4.1 关键引用文献

这篇论文的引用文献列表反映了其在强关联电子系统研究中的深厚背景，并将其工作置于现有理论框架中。以下是一些关键引用的重要性：

  1. [1] D. P. Arovas et al., The Hubbard Model, Annual Review of Condensed Matter Physics 13, 239 (2022).
    • 重要性： 这是一篇关于哈伯德模型的年度综述，奠定了模型的基本背景和当前研究进展。它强调了哈伯德模型的核心地位和复杂性，为本文提出新方法提供了动机。
  2. [2] D. J. Griffiths and D. F. Schroeter, Introduction to Quantum Mechanics, 3rd ed. (Cambridge University Press, Cambridge & New York, 2018).
    • 重要性： 量子力学教科书，介绍了摄动论、变分法和半经典近似这三种基本的近似方案。本文的工作正是对半经典近似框架的创新性探索，与这些经典方法形成了对比和补充。
  3. [4-11] 涉及哈伯德-斯特拉托诺维奇变换和半经典近似的工作。
    • 重要性： 这些文献代表了传统半经典近似方法，尤其是在HS变换后忽略辅助玻色子场虚时依赖性的工作。论文通过指出这些方法的局限性（例如HS场引入的任意性和无法捕捉非局部关联），凸显了其非常规相干态方法的必要性和优势。
  4. [14, 15] P. B. Wiegmann and A. Ferraz et al., t-J model路径积分表述。
    • 重要性： 论文将自身的方法视为对t-J模型和其他具有投影效应的模型（例如排除双占据态或费米子自由度）路径积分形式的扩展。这些工作为如何处理强关联模型中的投影或约束提供了参考，但论文的方法直接应用于保留了双占据态和费米子自由度的哈伯德模型。
  5. [16-19] H. Zhang, C. D. Batista et al., 多轨道莫特绝缘体和SU(N)相干态的工作。
    • 重要性： 这些工作探讨了多轨道系统和SU(N)相干态，为本文将非常规相干态推广到多轨道系统和利用SU(2)η × SU(2)spin对称性提供了背景和灵感。特别是大M极限的引入，与这些工作中利用SU(N)相干态处理大N极限的方法有异曲同工之妙。
  6. [20-27] S. Östlund et al., 非线性费米子变换。
    • 重要性： 这些文献讨论了哈伯德模型中电子算符的非线性变换，包括自旋-电荷分离和复合费米子算符。论文最终证明了其路径积分表述与Ref. [22]中引入的非线性费米子变换的等价性，从而为这种非常规相干态构造提供了坚实的理论基础和深刻的物理洞察。
  7. [33-36] A. M. Tsvelik, P. Coleman, J. Nilsson et al., 马约拉纳表示和近藤格点模型。
    • 重要性： 这些文献涉及马约拉纳费米子的概念以及它们在近藤格点模型中的应用。论文将哈伯德模型映射到马约拉纳-近藤格点模型，并利用马约拉纳表示来处理局部量子自旋，这与这些工作的思想紧密相连，揭示了哈伯德模型在马约拉纳基底下的新结构。

4.2 对这项工作局限性的评论

尽管本文提出的半经典近似框架具有显著的创新性和潜力，但其作为一种近似方法，也存在一些固有的局限性。

  1. 半经典近似的固有缺陷：

    • 忽略量子涨落： 核心的半经典近似在于忽略了玻色子场$\Omega_i$的虚时依赖性（即将其视为静态场）。这意味着玻色子场的量子涨落被完全忽略。这些涨落对于描述某些量子临界现象和相变至关重要。例如，在低维系统或临界点附近，量子涨落可能占据主导地位，此时该近似的准确性会降低。
    • 连续态密度 vs. 离散谱： 忽略量子涨落导致系统能谱从精确解的离散谱变为连续态密度。这种根本差异在一格点模型中表现为粒子数和双占据等物理量在渐近行为上的代数衰减与指数衰减的区别。这可能限制其在精确预测能隙、激发谱等方面的能力。
    • 能量尺度的偏差： 论文在二格点模型中明确指出，自旋关联的特征能量尺度（例如海森堡耦合常数$J$）在半经典近似中会系统性地降低一个常数因子（例如$J_{cl} = J/4$）。这表明该近似可能需要对结果进行适当的缩放才能与实验或更精确的量子计算进行定量比较。
    • 人为产物和权重损失： 在一格点模型中，高温下双占据的非单调行为被认定为半经典近似中连续态密度采样导致的人为产物。此外，单粒子格林函数的权重损失（$\langle |c_{i\sigma}|^2 \rangle < 1$）也暗示了近似的缺陷。这些细节表明，虽然方法在定性上表现良好，但在定量细节上可能引入伪影。

  2. 蒙特卡罗模拟的计算挑战：

    • 对于大于一或二格点的系统，对玻色子场$\Omega_i$的积分需要通过蒙特卡罗模拟完成。即使每个格点只有四个玻色子变量（两个SU(2)群的参数），对于大尺度的格点，状态空间依然巨大，蒙特卡罗采样效率可能成为瓶颈。论文也承认需要“更复杂的更新算法”来提高模拟效率，这暗示了当前实现的计算成本较高。
    • 模拟的收敛速度和平衡态的获取，尤其是在低温度下，可能是一个实际的挑战。

  3. 近似的系统性改进：

    • 论文提出了大M构造，使其成为一个受控的半经典极限，并暗示可以通过1/M展开来系统地引入量子涨落。然而，在实际应用中，超越静态近似来计算高阶1/M修正项的复杂性可能非常高，这本身就是另一个研究挑战。

  4. 非常规相干态的特异性：

    • 尽管非常规相干态的设计旨在最小化格拉斯曼变量并更好地处理自旋-电荷耦合，但这种特定选择可能对某些特定的量子态或物理相变存在偏向性或局限性。例如，它是否能够完美描述所有类型的超导配对或拓扑相，仍需在未来的工作中进一步验证。

  5. 缺乏具体代码和性能基准：

    • 论文没有提供任何代码实现、所用软件包或计算性能数据。这使得其他研究人员难以直接评估该方法在实际复杂系统中的适用性、可扩展性和效率，也增加了复现的难度。对于一个旨在提供实用工具的理论框架来说，缺乏这些信息是一个重要的实际局限。

总的来说，该工作在理论概念上取得了显著进展，为理解哈伯德模型及其推广模型提供了一个新的视角和工具。然而，作为一种半经典近似，它在预测细节和计算效率方面仍有改进空间，需要进一步的理论和数值工作来克服上述局限性。

5. 其他你认为必要的补充

5.1 理论意义与创新点

本文提出的哈伯德模型半经典近似框架具有显著的理论意义和多方面的创新点，为强关联电子系统的研究开辟了新途径。

  1. 创新性的相干态构造： 论文的核心突破在于提出了一种非常规的分级相干态表示。与传统的费米子相干态引入多个格拉斯曼变量不同，该方法通过单一的格拉斯曼变量来编码费米子宇称，而自旋和电荷自由度则由玻色子变量描述。这种构造巧妙地将局部希尔伯特空间分级为偶宇称和奇宇称扇区，并通过格拉斯曼变量作为两者之间的媒介。这种处理方式使得路径积分中的相互作用项自然地转化为格拉斯曼场的双线性形式，极大地简化了其处理。

  2. 避免哈伯德-斯特拉托诺维奇（HS）变换的任意性： 传统的半经典近似通常依赖于HS变换引入辅助玻色子场。HS场的选择具有一定的任意性，并且静态HS场难以捕捉非局部关联。本文的方法通过直接构造相干态，避免了这种任意性，并为处理非局部关联提供了一种更自然的框架，因为它在有效哈密顿量中保留了格点间的跳跃项，并通过玻色子场影响这些项的系数。

  3. 统一处理自旋和电荷自由度： 玻色子变量自然地编码了局部希尔伯特空间的自旋和电荷自由度，并清晰地展现了半满哈伯德模型的SU(2)η × SU(2)spin对称性。这使得该框架特别适用于描述涉及自旋和电荷纠缠序（如磁序与超导或电荷密度波共存）的系统，这是许多复杂材料的关键特征。

  4. 受控的半经典极限： 引入大M构造将哈伯德模型嵌入到一个M个副本的对称表示中，提供了一个受控的半经典极限。在M→∞的极限下，量子涨落被抑制，半经典近似变得精确。这为系统地改进近似和未来引入量子涨落修正项奠定了基础，使其成为一个可扩展的理论工具，类似于量子场论中的1/N展开。

  5. 自然推广到多轨道系统： 该相干态框架可以非常自然地推广到多轨道系统。通过为每个自旋轨道块分配一个格拉斯曼变量，可以处理具有更高维局部希尔伯特空间的复杂材料，而无需引入大量格拉斯曼变量，从而保持了计算的可处理性。

  6. 与非线性费米子变换和马约拉纳-近藤格点模型的联系： 论文不仅构建了一个新的近似框架，还深入探究了其与现有理论的深层联系。它证明了其有效作用量与非线性费米子变换（特别是Östlund等人提出的变换）是等价的，并且可以将哈伯德模型映射到马约拉纳-近藤格点模型。这种多重联系不仅验证了方法的正确性，也为从不同角度理解哈伯德模型的物理性质提供了新的视角。

5.2 未来研究方向与应用前景

鉴于该方法的创新性和所展现的潜力，未来有多个富有前景的研究方向：

  1. 超越静态近似： 当前工作忽略了玻色子场的虚时依赖性。未来的工作可以尝试系统地纳入这些量子涨落。例如，可以在大M展开中计算1/M阶修正，或者使用更高级的路径积分技术（如有效作用量的方法）来处理这些涨落。这将有助于捕捉更精细的量子效应，特别是低维系统和量子临界现象。

  2. 应用于更复杂的系统： 将该方法应用于更大尺度和更高维度的系统，如二维哈伯德模型，甚至是具有多轨道结构的实际材料模型。这将需要高效的蒙特卡罗采样算法和并行计算技术来克服计算挑战。特别是对于那些被认为具有纠缠序的材料（如铜氧化物超导体、重费米子系统等），该方法可能提供新的见解。

  3. 探索特定物理相： 利用该框架研究哈伯德模型中的各种物理相，如反铁磁序、电荷密度波、超导（特别是d波超导）和莫特绝缘体等。由于该方法能够自然处理自旋和电荷自由度，因此在研究这些纠缠序的共存和竞争机制方面具有独特优势。

  4. 与DMFT等先进方法的结合： 可以探索将该半经典框架与动力学平均场理论（DMFT）或其他先进的数值方法相结合。例如，可以将该半经典方法作为DMFT的杂质求解器，或者将DMFT的自洽性概念引入到玻色子场的平均中。这有望在保留局域量子涨落的同时，有效处理非局域关联。

  5. 数值算法优化与软件开发： 针对玻色子场的蒙特卡罗采样效率，开发更先进的更新算法，例如Hybrid Monte Carlo (HMC) 或并行回火。同时，将该方法封装成开源软件包，并提供清晰的文档和示例，将极大地促进其在凝聚态物理社区的广泛应用和进一步开发。

  6. 进一步探索与马约拉纳-近藤模型的联系： 深入研究哈伯德模型与马约拉纳-近藤格点模型之间的映射，这可能为理解强关联系统中的涌现现象（如自旋-电荷分离）提供新的理论工具。马约拉纳费米子在拓扑量子计算等领域具有重要意义，这种联系可能在交叉学科研究中产生意想不到的发现。

5.3 比较与展望

与现有方法相比，本文提出的半经典方法具有以下显著特点：

  • 优于传统静态HS方法： 相比于在HS变换后直接将辅助玻色子场视为静态的传统方法，该方法通过精心设计的相干态，使得相互作用项在格拉斯曼变量中变为双线性，从而在玻色子变量被冻结后，费米子问题仍然可以精确求解。这避免了HS变换的任意性，并能够更好地处理自旋和电荷的纠缠。
  • 与变分法和摄动法互补： 该方法提供了一个非微扰的框架，可以在摄动理论失效的强耦合区域提供见解。同时，它也比完全量子蒙特卡罗等方法计算成本更低，在物理直观性上更强。

总而言之，这项工作为研究哈伯德模型这类强关联电子系统提供了一个新的、系统性的理论工具。尽管半经典近似带来了定量上的偏差，但其在定性上捕获关键物理行为的能力、处理纠缠序的潜力以及与现有理论框架的深层联系，使其在凝聚态物理领域具有重要的研究价值和广阔的应用前景。