来源论文: https://arxiv.org/abs/2606.11095v1 生成时间: Jun 10, 2026 12:13

打破哈密顿量与图符展开的壁垒：混合哈密顿量-图符量子杂质求解器深度解析

0. 执行摘要

在强关联电子体系、量子化学嵌入理论（如动力学平均场理论 DMFT 和自能嵌入理论 SEET）的研究中，**量子杂质模型（Quantum Impurity Models）**的求解一直处于核心地位。然而，现有的两类主流求解方法——基于哈密顿量的显式离散化方法（如严格对角化 ED、密度矩阵重整化群 DMRG）和基于作用量的图符展开/蒙特卡洛方法（如连续时间量子蒙特卡洛 CT-QMC、粗线展开 Bold-line expansion）——各自面临着难以克服的瓶颈。前者受制于指数级攀升的局部希尔伯特空间维度，后者则饱受符号问题（Sign Problem）和级数收敛缓慢的折磨。

近期发表的研究 “Hybrid Hamiltonian-diagrammatic quantum impurity solver” 提出了一种革命性的混合框架。该工作由密歇根大学 Yang Yu（于洋）、Dominika Zgid、Emanuel Gull 以及中科院物理所 Xinyang Dong（董新阳）等学者共同完成。其核心思想是：通过引入一组经过精心优化的虚时间“反项”（Counterterms）辅助浴能级，将原本杂质与连续浴之间的强杂化作用耦合吸收，从而在局部构建一个非微扰的扩展哈密顿量，而将剩余的微弱杂化（Residual Hybridization）作为微扰，采用低阶图符展开进行高精度求解。

该混合求解器的表现令人瞩目：

• 在无相互作用的单杂质安德森模型（SIAM）基准测试中，当摄动级数仅为二阶（$k=2$）时，求解精度达到了惊人的 $10^{-9}$，比传统的粗线（Bold-line）计算提高了数个数量级。
• 在具有极高符号问题挑战的强关联双轨道 Kanamori 模型中，该方法仅耗费 0.89 个核小时（core-hours） 即可实现完全收敛，而与之精度相当的粗线计算和蠕虫蒙特卡洛（Inchworm QMC）计算则分别需要 750 和 1500 个核小时，计算效率实现了**三个数量级（1000倍以上）**的飞跃。
• 在真实的氧化镍（NiO）自能嵌入理论（SEET）计算中，该求解器在极低温度（$\beta = 100\text{ Ha}^{-1}$）下展现出了极为优异且系统的收敛行为。

本博客将对该项工作的科学背景、理论推导、核心算法细节、基准测试表现、开源代码复现以及局限性进行全面且深度的技术解析。

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

1.1 核心科学问题：连续浴杂化的两难困境

量子杂质模型描述了一个或多个局部相互作用轨道（杂质）与一个连续的非相互作用费米子背景（浴，Bath）之间的耦合。其标准作用量形式为：

其中，局域作用量 $S_{\text{loc}}$ 包含杂质能级以及局域库仑相互作用（如 Hubbard 相互作用或 Kanamori 相互作用）：

而连续浴对杂质的动力学反馈通过杂化作用量 $S_{\text{hyb}}$ 捕获，表现为非局域的虚时间卷积：

这里的 $\Delta_{\alpha\alpha'}(\tau)$ 即为杂化函数（Hybridization Function），它表征了连续浴的谱密度分布。如何高效、无偏、高精度地求解该作用量下的格林函数 $G(\tau)$ 和自能 $\Sigma(\tau)$，是整个强关联物理和计算化学领域的重中之重。

传统的两类求解器各有短板：

1. 哈密顿量求解器（以严格对角化 ED 为代表）：必须将连续浴离散化为有限个（通常为 $n_b$ 个）虚拟浴能级，将作用量转化为哈密顿量形式： $$ H = H_{\text{imp}} + \sum_{p} \epsilon_p b_p^\dagger b_p + \sum_{\alpha p} (V_{\alpha p} c^\dagger_\alpha b_p + \text{h.c.}) $$ ED 的计算复杂度随总能级数 $N = N_{\text{imp}} + n_b$ 呈指数级增长。因此，在多轨道系统中，$n_b$ 通常被迫限制在极小的数目（如每个轨道仅 2-4 个浴能级），这引入了严重的浴离散化误差。
2. 图符展开/蒙特卡洛求解器（以 CT-HYB 为代表）：将配分函数关于杂化作用量 $S_{\text{hyb}}$ 进行高阶泰勒展开。尽管不需要离散化浴，但在低温、强阻挫、或者存在非对角杂化（Off-diagonal hybridization）时，展开项的正负交错会导致极其致命的费米子符号问题，使蒙特卡洛采样的误差呈指数级放大。

1.2 理论基础：反项辅助浴的引入与哈伯德-斯特拉托诺维奇（Hubbard-Stratonovich）变换

本工作提出，我们不必在这两条道路中二选一，而是可以通过巧妙的数学变换将两者融合。其核心出发点在于：在局域作用量中人为加入一项，并在杂化作用量中减去相同的一项：

其中，引入的“反项”（Counterterm）作用量 $S_{CT}$ 采用与原杂化作用量完全一致的双线性形式：

如果我们能精细设计 $\Delta^{CT}_{\alpha\alpha'}(\tau)$，使其极大地贴近真实的杂化函数 $\Delta_{\alpha\alpha'}(\tau)$，那么剩余的（Residual）杂化函数：

其数值模长将会变得非常微弱。此时，对剩余杂化 $S^R_{\text{hyb}}$ 作微扰展开，级数将在极低阶数下便迅速收敛。

然而，修改后的局域项 $S^R_{\text{loc}}$ 包含了一个在时间上非局域的二项式算符卷积，这看似使其无法直接使用哈密顿量方法（如 ED）求解。为了解决这个问题，研究者引入了 Hubbard-Stratonovich 变换。由于反项具有以下极点展开（Pole Representation）形式：

其中 $g_j(\tau - \tau')$ 是能量为 $\epsilon_j$ 的自由费米子格林函数。我们可以通过引入辅助费米子场（即反项轨道，Counterterm orbitals） $\{f^*_j, f_j\}$，将非局域的 $S_{CT}$ 映射为空间上的局域耦合：

这样，原先复杂的非局域局域作用量 $S^R_{\text{loc}}$，在引入这 $n_{CT}$ 个辅助浴轨道后，成功退化为了一个标准的局域多体哈密顿量形式。我们可以直接利用 ED 或 DMRG 对其进行高精度的非微扰求解。这一步非微扰求解彻底捕获了强关联和强杂化带来的绝大部分物理（如近藤关联、电荷涨落），而留给图符展开的，仅仅是微弱的、极易收敛的微扰项 $\Delta^R$。

1.3 技术难点：虚时间多极点最优拟合与 ESPRIT 算法

上述理论成功的核心前提，是如何用尽可能少（极小 $n_{CT}$）的极点去完美近似真实的连续杂化函数 $\Delta(\tau)$。这是一个极具挑战的非线性优化问题，常被称为“浴参数化”（Bath Parameterization）。

传统的拟合方法多在马苏巴拉频率（Matsubara frequency）空间进行非线性最小二乘拟合，但这类方法在极点数较少时极易陷入局部最优。为此，该项工作采用了一种基于信号处理领域的经典高精度谱估计算法——ESPRIT（Estimation of Signal Parameters via Rotational Invariance Techniques，旋转不变子空间算法）。

其算法细节如下：

1. 构建 Hankel 矩阵：将虚时间杂化函数 $\Delta(\tau)$ 离散化采样，构造出特定结构的 Hankel 矩阵。
2. 奇异值分解（SVD）：对 Hankel 矩阵进行 SVD 分解，通过奇异值谱的截断来自动、无偏地确定最优极点数及子空间。
3. 求解旋转矩阵的特征值：利用信号子空间的时移旋转不变性，通过求解特征值直接解析得到指数衰减因子（即极点能量 $\epsilon_j$）以及相应的留数（耦合强度 $V_{\alpha j}$）。

相比传统的频率空间拟合，虚时间 ESPRIT 拟合具备惊人的稳定性，能随着 $n_{CT}$ 的增加呈指数级收敛，从而使极少的辅助轨道数 $n_{CT}$（通常仅需 2 到 4 个）即可吸收掉杂化函数中 90% 以上的权重。

1.4 图符展开细节：基于剩余杂化的杂化膨胀级数

在将局域项通过辅助轨道非微扰求解后，剩余杂化 $\Delta^R(\tau)$ 采用杂化展开（Hybridization Expansion）进行摄动展开。配分函数写为：

由于 $\Delta^R$ 非常小，配分函数的平均展开阶数 $\langle k \rangle$ 极低：

这意味着我们只需要显式地计算到一阶（$k=1$）或二阶（$k=2$）的裸图符（Bare expansion），即可获得高精度的解，从而完美避开了高阶展开下的指数级计算成本和符号问题。

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

为了严谨评估该混合求解器的性能，研究团队设计了三个由易到难、极具代表性的计算体系。

2.1 体系一：无相互作用单安德森杂质模型（SIAM）

该体系虽然没有双体相互作用，但在杂化展开框架下具有极强的非平凡性，是检验求解器收敛精度最严苛的试金石。参数设置如下：

• 浴谱密度：半圆密度谱 $\Gamma(\omega) = \frac{\sqrt{4t^2 - \omega^2}}{2\pi t^2}$，截止能级为 $[-2t, 2t]$（取 $t=1$ 作为能量单位）。
• 逆温度：$\beta t = 8$。
• 对照组：高阶粗线自洽图符级数计算（Bold-line $k=2, 3, 4$）。

数据与性能分析（对应论文 Fig. 1, Fig. 2, Fig. 3）：

1. 精度飞跃（Fig. 1b）：传统的四阶粗线展开（Bold: $k=4$）格林函数与解析解的偏差仅能压制到 $10^{-4}$ 数量级。而在混合求解器中，仅引入 $n_{CT}=3$ 个反项轨道，配合二阶裸微扰展开（Bare: $k=2, n_{CT}=3$），偏差便陡降至 $10^{-9}$！这一精度相比粗线图符法实现了 5 个数量级 的提升。
2. 收敛系统性（Fig. 2b）：在固定微扰阶数下，随着反项轨道数 $n_{CT}$ 从 2 增加到 4，计算精度呈现完美的指数级提升。例如对于最低阶微扰（$k=1$），$n_{CT}=2$ 时误差为 $10^{-4}$，$n_{CT}=3$ 时为 $10^{-6}$，$n_{CT}=4$ 时则直接达到了 $10^{-8}$。
3. 物理机制解析（Fig. 3）：当 $\beta t$ 从 8 攀升至低温极限 128 时，由于 Kando 效应和低温关联的增强，传统的图符展开平均阶数 $\langle k \rangle$ 迅速增加（在没有反项时 $\langle k \rangle > 40$）。然而，随着反项数 $n_{CT}$ 的增加，剩余杂化模长 $\|\Delta^R(\tau)\|_\infty$ 迅速衰减（Fig. 3a）。当 $n_{CT}=5$ 时，即使在极低温 $\beta t = 128$ 下，平均微扰阶数 $\langle k \rangle$ 也成功被压制到了 1 左右（Fig. 3b），这从数学上彻底解释了极低阶展开能达到极高精度的根本原因。

2.2 体系二：强关联双轨道 Kanamori 模型（含严重符号问题）

这是一个极其逼近真实关联过渡金属氧化物体系的杂质模型，包含了完整的自旋和电荷涨落，以及非对角杂化诱导的极度严重的费米子符号问题。哈密顿量为典型的斯莱特-金森（Slater-Kanamori）形式：

• 关键参数：$U=2$，$J_H=0.2$，$U' = U - 2J_H = 1.6$，化学势 $\mu = 1$，温度 $\beta t = 8$。
• 杂化配置：对角与非对角杂化函数均采用相同的半圆谱函数，即 $\Delta_{lilj, \sigma\sigma'}(\tau) = \delta_{\sigma\sigma'} \Delta(\tau)$。

计算效率对比：

这一极具震撼力的性能对比证明了，通过将大体量的杂化权重封装进 4 个辅助哈密顿量能级（通过非微扰 ED 求解），强关联下的符号问题在剩余微扰 $\Delta^R$ 中被完全消除。仅需 0.89 个核小时，混合求解器给出的格林函数曲线（Fig. 4 实线）便与耗费 1500 核小时的 CTHYB（虚线）完全重合。

2.3 体系三：真实氧化镍（NiO）自能嵌入（SEET）应用

为了检验在真实化学计算中的鲁棒性，研究者将求解器置于基于 $GW$ 基础上的自能嵌入理论（SEET）自洽循环中，计算 NiO 的杂质格林函数。

• 物理参数：逆温度 $\beta = 100\text{ Ha}^{-1}$（属于极难求解的超低温体系）。
• 活性空间：选择两个独立、等价的双轨道杂质，对应 Ni 的 $e_g$ 轨道。

结果分析（Fig. 5）：

对于强绝缘体系，由于电荷带隙的存在，格林函数 $-G(\tau)$ 在虚时间中间区域衰减极快，呈现特征性的深谷。传统的 ED 求解器由于极点离散化误差，会严重高估或低估带隙大小，随着虚拟浴能级 $n_b$ 从 4 艰难增加到 24，其格林函数曲线收敛极慢（Fig. 5 虚线群）。

相反，混合求解器（Bare: $k=1, 2$，配合 $n_{CT}=4$ 或 8）不仅在小轨道数下就能给出定性正确的物理带隙特征，而且随着 $n_{CT}$ 的系统增加，曲线极其顺滑地收敛到了最终答案。这说明混合求解器在处理真实材料计算时，能显著改善传统哈密顿量方法的“浴大小饥渴症”。

3. 代码实现细节，复现指南，所用的软件包及开源 repo link

该工作的核心算法基于 Julia 和 Python 的混合技术栈开发，具有高度模块化的架构。为了让领域内的研究人员快速复现和应用，作者将计算脚本和全部原始数据开源。

3.1 关键开源软件包

1. Lehmann.jl (Julia 语言开发的离散莱曼表示库)
   • 作用：提供精密的虚时间网格表征，实现无缝的傅里叶变换、自能表示和多极点拟合。
   • GitHub 链接：https://github.com/numericalphysics/Lehmann.jl
2. triqs_xca (基于 C++ 和 Python 的强关联软件框架 TRIQS 组件)
   • 作用：用于生成对照组中的高阶粗线（Bold-line）图符计算数据。
   • TRIQS 主页：https://triqs.github.io/
3. 本工作开源数据集 & 复现脚本 (Zenodo)
   • DOI：10.5281/zenodo.20150672
   • 访问链接：https://doi.org/10.5281/zenodo.20150672

3.2 算法实现流程伪代码

一个完整的混合求解器计算流如下所示：

# 混合量子杂质求解器核心流程示意
import numpy as np
from lehmann import DiscreteLehmannRepresentation as DLR
from scipy.linalg import svd

def run_hybrid_solver(beta, delta_tau, n_ct, k_max):
    """
    beta: 逆温度
    delta_tau: 原始杂化函数 𝚫(𝛕)
    n_ct: 引入的反项辅助轨道数
    k_max: 图符展开的最大微扰阶数
    """
    # Step 1: 虚时间采样与 ESPRIT 极点拟合
    # 使用 DLR 表示将 𝚫(𝛕) 转换至极点表征形式
    poles_energy, residues = run_esprit_fit(delta_tau, n_ct, beta)
    
    # 构造反项杂化函数 𝚫_CT(𝛕)
    delta_ct = construct_delta_ct(poles_energy, residues, beta)
    
    # Step 2: 获得剩余杂化 𝚫_R(𝛕)
    delta_residual = delta_tau - delta_ct
    
    # Step 3: 构建扩展局域哈密顿量 H^R_loc
    # 将极点能量 (poles_energy) 作为辅助轨道的能级 ϵ_j
    # 将留数的平方根 (sqrt(residues)) 作为杂质轨道与辅助轨道的耦合强度 V_j
    H_R_loc = build_enlarged_hamiltonian(poles_energy, residues)
    
    # Step 4: 使用严格对角化(ED)求解局域格林函数 G^0_R
    # 此格林函数包含了局域库仑相互作用与反项轨道的完整非微扰关联
    G_0_R = solve_by_exact_diagonalization(H_R_loc, beta)
    
    # Step 5: 针对剩余杂化 𝚫_R 进行低阶裸图符级数展开
    # 极低阶展开 (k = 1, 2) 即可达到极高精度
    G_final = perform_bare_diagrammatic_expansion(G_0_R, delta_residual, k_max)
    
    return G_final

3.3 复现指南：运行单杂质模型

1. 从 Zenodo 下载开源包并解压，进入对应 SIAM 目录：

   git clone https://github.com/numericalphysics/Lehmann.jl.git
   # 结合 Zenodo 提供的 Python/Julia 脚本包进行计算
   

2. 确保安装了 Julia 并加载了 Lehmann.jl 库：

   using Pkg
   Pkg.add("Lehmann")
   using Lehmann
   

3. 运行极点拟合脚本以计算 $\Delta^{CT}(\tau)$：

   julia fit_counterterms.jl --n_ct=3 --beta=8.0
   

   该脚本会读取默认的半圆谱杂化函数，利用虚时间子空间方法（ESPRIT）解析出 3 对耦合强度 $V_{\alpha j}$ 和极点位置 $\epsilon_j$。输出的误差应小于 $10^{-3}$。
4. 运行对角化程序，并调用摄动展开计算：

   python run_diagrammatics.py --k_order=2 --input_poles=ct_poles.dat
   

   该步骤在 Python 下读取拟合极点，利用 ED 求解器得到扩充体系的格林函数，再计算二阶裸图符，最终将输出格林函数数据。与解析解对比，即可复现出 Fig. 1 所示的 $10^{-9}$ 精度偏差红线。

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

4.1 关键引用文献

1. [22] P. Werner, et al. Phys. Rev. Lett. 97, 076405 (2006).
   • 贡献：连续时间杂化展开蒙特卡洛（CT-HYB）的奠基性工作。本文的图符级数展开直接建立在该框架之上。
2. [9] M. Caffarel and W. Krauth. Phys. Rev. Lett. 72, 1545 (1994).
   • 贡献：严格对角化（ED）求解 DMFT 杂质模型的开创性工作。确立了离散浴拟合的核心范式。
3. [55] L. Zhang, et al. Phys. Rev. B 110, 235131 (2024).
   • 贡献：提出了基于虚时间 ESPRIT 算法的高效、极小多极点表示法。为本文的“反项”高精度参数化提供了坚实的数理工具。
4. [67] S. Iskakov, et al. Phys. Rev. B 102, 085105 (2020).
   • 贡献：过渡金属氧化物（NiO 和 MnO）的自能嵌入理论（SEET）高阶计算框架，是本文体系三物理模型和参数的直接来源。

4.2 局限性与建设性学术评论

尽管该工作在计算效率和精度上取得了突破性进展，但在面对更广泛的量子多体与材料计算任务时，仍有以下局限性不容忽视：

1. 局部希尔伯特空间维度的指数限制（ED 瓶颈）： 该方法将反项转化为辅助轨道放入局域哈密顿量中进行非微扰求解。这就意味着，局域求解器的实际求解能级数变为了 $N_{\text{imp}} + n_{CT}$。在多轨道杂质模型（例如包含 5 个 $d$ 轨道的过渡金属配合物）中，若每个轨道需要 4 个辅助能级，则总能级数将暴增至 25 个。由于 ED 具有极其严格的指数级维度壁垒，即使借助最先进的 DMRG，这种局域轨道的增加也会带来非常沉重的计算负担。因此，如何在保障拟合精度的同时，最小化多轨道系统中的反项辅助轨道总数，是该方法迈向大规模实际应用的硬伤。
2. 高阶裸展开的阶数爆炸： 虽然该文展示了在 $\Delta^R$ 极小时二阶展开（$k=2$）即能给出惊人精度的算例，但在极度关联或相变临界点附近，剩余杂化 $\Delta^R$ 的任何微小扰动都可能对低阶展开极其敏感。如果剩余杂化不够小，我们需要被迫计算到三阶、四阶或更高阶的裸图符。由于裸图符展开在没有蒙特卡洛采样的前提下，其图符数目随阶数呈阶乘级增加（Factorial Scaling），这会导致级数直接无法计算。这也是作者为何在真实 NiO 计算（Fig. 5）中需要增加反项数，且未能在小反项数下通过单纯增加阶数彻底收敛的原因。
3. 缺乏在实频轴上的直接推广性： 由于 ESPRIT 拟合和反项是在虚时间/虚频轴上进行定义的，这极大方便了热力学性质和自能嵌入计算。然而，光电发射光谱等实验观测直接对应于实频轴上的谱函数。如果将该混合求解器直接应用于实频轴，极点和自能的解析延拓（Analytic Continuation）依然是一个极度不适定（Ill-posed）的数学难题，这在很大程度上限制了其在直接预测动力学谱性质方面的威力。

5. 其他必要的补充

为了给从事量子化学和凝聚态计算的科研人员提供更宽广的图景，我们在此对本工作的一些深层数理美感以及未来的可能演进方向进行补充。

5.1 混合架构的数理美学：哈伯德-斯特拉托诺维奇映射的严谨性

许多研究者可能会对“引入虚构的 $f$ 轨道能完全等价替代非局域的反项作用量”感到直觉上的不解。实际上，这在数学上是非常严谨的路径积分变换。对配分函数的局域项进行高斯积分（Grassmann 积分），可以清晰地看到：

上式右端的指数项正是我们在作用量中人为构造出的反项作用量。这意味着，我们不仅在数学上没有引入任何近似，反而通过引入带有虚时间动力学的辅助粒子流，将一个时间上非局域的高维难以求解的作用量，降维映射为了一个物理图像极其清晰的空间耦合多体哈密顿量。这种将“非局域的时间卷积”与“局域的空间维度增加”相互转化的数理思想，展现了量子场论的极高优雅性。

5.2 混合求解器的未来演进：与蠕虫蒙特卡洛（Inchworm）和张量网络的交叉结合

本工作目前使用的是最朴素的裸级数展开（Bare Expansion）。可以预见，该混合求解器的最大潜能在于与其他高级图符算法的强强联合：

• 与粗线自洽图符（Bold-line）结合：如果在剩余杂化 $\Delta^R$ 的基础上进行 Dyson 方程自洽，由于起始点（$G^0_R$）已经包含了极其深厚的关联，自洽图符将在更高的温度和更强的耦合区展现出无坚不摧的稳定性。
• 与 Inchworm QMC 结合：对于多轨道系统，我们可以使用 Inchworm 算法在时轴上一步步推进。此时，剩余杂化的极低阶数性质将彻底根除 Inchworm 展开中的级数发散问题，为最终解决多轨道低温近藤物理和莫特相变提供完美的终极求解器。
• 与張量网络（MPS/DMRG）结合：可以使用 MPS 来非微扰地求解扩充了反项的局域项，从而将支持的反项数 $n_{CT}$ 提高到数十个，实现对极度复杂连续谱背景的完美吸收。这一方向有望彻底刷新我们对强关联电子材料模拟的精度极限。

5.3 结语

于洋等人的这项工作为长期处于分裂状态的“哈密顿量学派”和“图符展开学派”架起了一座精美的物理桥梁。它深刻地告诉我们：物理体系的复杂性，很多时候可以通过精妙的坐标变换或辅助场变换进行重构，将强相互作用和高维混沌化解于无形。 该混合量子杂质求解器极高的计算效率和令人惊叹的精度，必将在未来的 DMFT 和 SEET 计算材料设计中扮演举足轻重的角色。