来源论文: https://arxiv.org/abs/2606.31117v2 生成时间: Jul 10, 2026 05:57

基于自由基碎片多体展开(MBE2)的线性烷烃量子化学计算:突破 NISQ 时代硬件限制的变革性技术

0. 执行摘要

在嘈杂中等规模量子(NISQ)时代,如何将高精度的量子化学计算扩展到具有实际化学意义的大分子体系,是量子计算领域最核心的挑战之一。传统的全分子量子化学模拟(如全构型相互作用 FCI)所需的量子比特数和线路深度随体系规模呈指数或高阶多项式增长,使得在现有硬件上直接计算中大型分子变得完全不可行。

本文深度解析了一项突破性的研究工作:基于自由基碎片的多体展开方法(Radical-Fragment Many-Body Expansion, MBE2)。该方法专为线性烷烃体系设计,通过对 C-C 键进行均裂(Homolytic Cleavage),将长链烷烃拆分为极小规模的开壳层(Open-shell)自由基单体(甲基自由基 $\text{CH}_3^\bullet$ 和亚甲基双自由基 $^\bullet\text{CH}_2^ullet$)及对应的双体。与传统的碎片分子轨道(FMO)方法不同,该方法完全摒弃了复杂的静电嵌入潜能和自洽电荷(SCC)迭代,而是直接在孤立状态下使用受限开壳层哈特里-福克(ROHF)方法处理碎片。

由于线性烷烃的高度平移对称性,无论分子链有多长,MBE2 组装公式均可将其简化为**仅需 4 个独特碎片(2 个单体,2 个双体)**的独立计算。在基组 STO-3G 下,最大碎片($\text{CH}_3-\text{CH}_2$ 双体)的比特需求被永久锁死在 30 个量子比特(通过活动空间 CASSCF 进一步压缩至 8 个量子比特)。在对己烷($\text{C}_{26}\text{H}_{54}$)的基准测试中,该方法实现了高达 **12.3 倍的量子比特数消减(从 368 个降至 30 个)**以及 12.8 倍的计算量消减。本文不仅在经典算法(RHF, CCSD)和理想量子模拟器(VQE, ADAPT-VQE)上验证了其极高的数值稳定性,更在 **IBM 的真实超导量子处理器(ibm_pittsburgh Heron)**上利用样本量子对角化(SQD)算法成功复现了关联能,展示了该框架在当前硬件下的非凡实用性与高精度。


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

1.1 核心科学问题与技术瓶颈

量子化学是量子计算最被看好的杀手级应用领域之一。然而,直接在量子计算机上求解大分子的电子薛定谔方程面临着严峻的资源瓶颈。以线性烷烃 $\text{C}_n\text{H}_{2n+2}$ 为例,在最基础的最小基组 STO-3G 下,每个碳原子贡献 5 个原子轨道,每个氢原子贡献 1 个原子轨道。在 Jordan-Wigner 映射下,每个空间轨道对应 2 个自旋轨道,进而需要 2 个自旋轨道量子比特。因此,全分子的比特需求量公式为:

$$q = 2 \times (5n + 2n + 2) = 14n + 4$$

对于链长仅为 26 个碳原子的己烷($\text{C}_{26}\text{H}_{54}$),直接计算需要高达 368 个物理量子比特。更致命的是,分子体系的希尔伯特空间维度随着电子和轨道数的增加而呈组合数爆炸增长(斯莱特行列式数量急剧膨胀)。在 NISQ 时代的真实量子芯片上,高噪声、低相干时间和有限的门保真度使得直接模拟超过 20 个比特的分子体系变得极其困难。因此,如何将一个庞大的全分子体系在物理和数学上不失真地拆解为当前量子硬件可承受的微小单元,是解决量子化学可扩展性问题的核心科学问题。

1.2 理论基础:多体展开(MBE)与传统 FMO 的对比

多体展开(Many-Body Expansion, MBE)是一种将复合体系的总能量表达为碎片及其相互作用贡献之和的通用数学框架。在二体截断(MBE2)级别,体系的总能量可表示为:

$$E_{\text{MBE2}} = \sum_{I} E_I + \sum_{I其中 $E_I$ 为独立单体 $I$ 的能量,$\Delta E_{IJ} = E_{IJ} - E_I - E_J$ 为相邻键合双体 $IJ$ 的非加和性校正能。为了控制计算复杂度,此处实施了“仅限键合(bonded-pair-only)”的物理限制,即只计算在化学键上相邻的碎片双体,忽略了空间中不相邻碎片之间的长程弱相互作用。

传统 FMO 方法的局限性

传统的碎片分子轨道(FMO)方法为了处理共价键断裂,通常采用**异裂(Heterolytic Cleavage)**的方式。为了维持碎片闭壳层(Closed-shell)的稳定状态,FMO 引入了“键脱离原子(Bond-Detached Atoms, BDA)”或氢原子饱和修饰(Hydrogen-capping),并且必须在一个极其复杂的外部静电嵌入势(Electrostatic Embedding Potential)中对各碎片进行自洽电荷(SCC)迭代循环。这种处理带来了三大技术痛点:

  1. 电荷人工引入风险:异裂会导致碎片带上非自然的局部电荷,扭曲原本非极性共价键的电子云分布。
  2. SCC 迭代收敛困难:静电嵌入潜能的自洽迭代在量子计算机上的实现极为繁琐,且容易遇到不收敛问题。
  3. 量子硬件兼容性差:外加电场的静电势算符会显著增加分子哈密顿量中 Pauli 项的数量,大幅提升量子测量开销。

1.3 技术创新:基于均裂的自由基碎片 MBE2 框架

针对上述痛点,本研究提出了一个颠覆性的方案:对 C-C 键进行均裂(Homolytic Cleavage)。其物理和化学合理性在于:对于非极性的 C-C 键,两个参与成键的碳原子电负性完全相同($\Delta \chi = 0$)。因此,均裂是最符合物理真实的断键描述。

具体断键和碎片化策略如下:

  • 单体生成:将每一个成键电子对均等地一分为二,每个碎片在断键处保留一个未成对的自由电子。这产生了两类电中性、开壳层的自由基单体:
    • 甲基自由基 ($\text{CH}_3^\bullet$):含 9 个电子,单态双重态(Doublet, $S = 1/2$),具有 8 个空间轨道。
    • 亚甲基双自由基 ($^\bullet\text{CH}_2^\bullet$):含 8 个电子,开壳层单重态或三重态(本文采用开壳层处理,具有 7 个空间轨道)。
  • 双体生成:相邻的单体重新结合为双体,依然保留其固有的自旋特性,避免了人工闭壳层化带来的自旋污染:
    • $\text{CH}_3-\text{CH}_2$ 双体:含 17 个电子,双重态($S = 1/2$),15 个空间轨道。
    • $\text{CH}_2-\text{CH}_2$ 双体:含 16 个电子,三重态($S = 1$),14 个空间轨道。

这些碎片完全处于**孤立状态(In isolation)**进行计算,无需静电嵌入、无需加氢饱和、无需 SCC 迭代,极大地净化了量子哈密顿量。

对于包含 $n$ 个碳原子的线性烷烃,利用平移对称性,MBE2 组装公式可被精简压缩为:

$$E_{\text{MBE2}}(n) = 2 E_{\text{CH}_3^\bullet} + (n-2) E_{^\bullet\text{CH}_2^\bullet} + 2 \Delta E_{\text{CH}_3-\text{CH}_2} + (n-3) \Delta E_{\text{CH}_2-\text{CH}_2}$$

该公式的伟大之处在于:无论烷烃链多长(即使 $n=1000$),都只需要进行这 4 个独特碎片的能量计算。量子计算资源的最大需求直接被限制在最大碎片($\text{CH}_3-\text{CH}_2$)的固有尺寸上。在 STO-3G 基组下,其最大空间轨道数为 15,映射后仅需 30 个量子比特,彻底打破了 $14n + 4$ 的线性增长红线!

1.4 量子求解器管道(Quantum Solver Pipelines)方法细节

为了在不同的计算约束下求解这 4 个碎片的基态能量,本工作设计了三种不同的量子求解器管道(如图 1 所示):

                      +---------------------------------+
                      |    线性烷烃 C_n H_{2n+2} (STO-3G) |
                      +---------------------------------+
                                       |
                                       | 均裂 C-C 键
                                       v
                      +---------------------------------+
                      | 4个独特自由基碎片 (CH3*, *CH2*...)|
                      +---------------------------------+
                                       |
                                       +-----------------------+
                                       |                       |
                                       v (VQE / ADAPT-VQE)     v (SQD 管道)
                        +------------------------------+ +--------------------+
                        |  ROHF 求解自洽分子轨道空间     | | ROHF 分子轨道生成  |
                        +------------------------------+ +--------------------+
                                       |                           |
                                       v                           | (无活性空间压缩)
                        +------------------------------+           |
                        | CASSCF(n,4) 活性空间压缩到8比特 |           |
                        +------------------------------+           |
                                       |                           |
                     +-----------------+-----------------+         | 14-30比特
                     |                                   |         v
                     v (VQE 管道)                        v (ADAPT) v (SQD 硬件执行)
       +--------------------------+          +-----------------------+ +---------------------+
       | EfficientSU2 硬件高效电路 |          | GSD算符池动态迭代构建 | | LUCJ电路 + 经典S-CORE|
       | SPSA  noisy 模拟优化测量  |          | 状态矢量精确梯度更新  | | 真实超导芯片采样对角化 |
       +--------------------------+          +-----------------------+ +---------------------+

1.4.1 变分量子本征求解器(VQE)管道

  • 活性空间压缩:为使 VQE 和 ADAPT-VQE 线路在极浅的深度内收敛,本工作引入了经典完全活性空间自洽场(CASSCF)进行预处理。通过将碎片压缩至统一的 $\text{CAS}(n, 4)$ 活动空间(即 4 个空间轨道,对应 8 个自旋轨道比特),将复杂的碎片哈密顿量压缩为标准的 8 比特算符。
  • 变分电路(Ansatz):采用硬件高效的 EfficientSU2 启发式线路,包含交替的单比特旋转门($R_Y, R_Z$)和线性拓扑的双比特 $CZ$ 纠缠门。层数控制在 1-3 层。
  • 经典优化:鉴于测量的噪声,使用联立摄动随机逼近算法(SPSA)进行参数更新。每次能量评估使用 10,000 次测量采样(shots),迭代 1500 至 2500 步。

1.4.2 自适应导数装配伪 Trotter 变分量子本征求解器(ADAPT-VQE)

  • 算符池构建:为了规避固定硬件高效电路的表达力局限和贫瘠高原(Barren Plateaus)问题,ADAPT-VQE 采用动态生长机制。算符池使用广义单双激发(Generalized Singles and Doubles, GSD)算符:

    $$\hat{G}_{pq}^{(1)} = i \left( \hat{a}_p^\dagger \hat{a}_q - \hat{a}_q^\dagger \hat{a}_p \right)$$$$\hat{G}_{pqrs}^{(2)} = i \left( \hat{a}_p^\dagger \hat{a}_q^\dagger \hat{a}_r \hat{a}_s - \hat{a}_s^\dagger \hat{a}_r^\dagger \hat{a}_q \hat{a}_p \right)$$
  • 动态生长:在当前状态 $|\Psi\rangle$ 下,计算算符池中每个算符 $\hat{A}_k$ 的能量梯度:

    $$g_k = \left\| \frac{\partial E}{\partial \theta_k} \right\|_{\theta_k=0} = \left| \langle\Psi| [\hat{H}, \hat{A}_k] |\Psi\rangle \right|$$

    选择梯度绝对值最大的算符 $k^*$,将其对应的单参数幺正算符 $e^{i\theta_{k^*}\hat{A}_{k^*}}$ 插入电路中,随后对所有电路参数进行联合重新优化。当最大梯度低于 $10^{-3}$ 时宣告收敛。

1.4.3 样本量子对角化(SQD)

  • 大尺寸硬件运行:SQD 无需活性空间压缩,直接操作完整的 STO-3G 碎片哈密顿量(14 至 30 个物理量子比特)。

  • 量子状态准备:利用经典得到的精确保连簇单双激发(CCSD)振幅 $t_1, t_2$,通过局部幺正耦合簇 Jastrow(LUCJ) Ansatz 构建量子态准备电路:

    $$|\Psi\rangle = \prod_{\mu=1}^L e^{\hat{K}_\mu} e^{i \hat{J}_\mu} e^{-\hat{K}_\mu} |\Phi_0\rangle$$
  • 量子测量与后选择:在真实的量子计算机上执行电路,获取计算基下的比特串样本。对收集到的原始比特串进行粒子数后选择,剔除不守恒 $(n_\alpha, n_\beta)$ 的嘈杂自旋扇区状态。

  • S-CORE 优化与 Davidson 对角化:将保留的优势配置集合 $S$ 输入经典自洽轨道恢复与扩展(S-CORE)程序,精炼投影哈密顿量:

    $$\hat{H}_S = \hat{P}_S \hat{H} \hat{P}_S, \quad \hat{P}_S = \sum_{x \in S} |x\rangle\langle x|$$

    利用 Davidson 迭代算法求解 $\hat{H}_S$ 的最低特征值,得到最终的关联能估算。该方法彻底避免了量子处理器上的变分循环,极其适合 NISQ 硬件。


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

2.1 基准测试体系设置

研究团队构建了涵盖 11 种线性烷烃的完备测试集:从丁烷($\text{C}_4\text{H}_{10}$)到己烷($\text{C}_{26}\text{H}_{54}$),具体包含 $n = 4, 5, 6, 7, 8, 9, 10, 15, 20, 25, 26$。出于规避边界伪影的考虑,排除了丙烷($n=3$),因为在其 MBE2 分解中,两个终端 $\text{CH}_3$ 的边界效应严重,无法代表长链烷烃的真实物理化学特性。

所有烷烃分子均被设置为理想的全反式(all-anti)构象,碳骨架完全平铺在 $xz$ 平面上,其几何参数固定为:碳碳键长 $d_{\text{CC}} = 1.54 \text{ Å}$,碳氢键长 $d_{\text{CH}} = 1.09 \text{ Å}$,键角 $\angle\text{CCC} = 112^\circ$。基组统一采用 STO-3G。

2.2 核心性能数据与分析

2.2.1 固有碎片化误差(Intrinsic Fragmentation Error)分析

首先,为了度量二体键合多体展开(MBE2)本身的近似精度,研究团队在经典层次对比了 MBE2 与全分子直接计算(Full-molecule)的能量差异。数据如表 1 所示(能量误差单位为 $\text{kcal/mol}$):

体系分子碳原子数 $n$MBE2-RHF 相对 Full-RHF 误差MBE2-CCSD 相对 Full-CCSD 误差
丁烷 ($\text{C}_4\text{H}_{10}$)4+6.6+8.6
己烷 ($\text{C}_6\text{H}_{14}$)6+33.8+38.1
癸烷 ($\text{C}_{10}\text{H}_{22}$)10+95.6+105.1
二十烷 ($\text{C}_{20}\text{H}_{42}$)20+263.9+288.4
己烷 ($\text{C}_{26}\text{H}_{54}$)26+364.6+396.5

理论解读

  1. 误差的物理来源:MBE2 方法高估了总能量(误差为正值)。这是由于“仅限键合键”近似忽略了非相邻碳原子之间的长程弱排斥/吸引作用,以及三体及以上的多体极化效应。
  2. CCSD 误差偏大:MBE2-CCSD 的碎片化误差略高于 MBE2-RHF。这符合物理直觉——电子关联能具有比均值场能更长程的物理特征,截断二体相互作用会损失稍多一些的关联能校正。
  3. 单碳原子平均误差的收敛性:虽然总绝对误差随链长增长而线性累积,但每个碳原子的平均误差表现出完美的收敛特性
    • MBE2-RHF 每碳误差从 $n=4$ 时的 $+1.7 \text{ kcal/mol/C}$ 快速渐近收敛至 $+14.0 \text{ kcal/mol/C}$。
    • MBE2-CCSD 每碳误差从 $+2.2 \text{ kcal/mol/C}$ 收敛至 $+15.3 \text{ kcal/mol/C}$。 这证实了该均裂自由基碎片框架具有极佳的大小一致性(Size-extensivity)。只要该每碳误差是常数,即可通过简单的线性修正在热力学极限下予以消除。

2.2.2 量子求解器精度评估

1. MBE2-VQE 表现

MBE2-VQE 在 8 比特 CAS 空间中表现极为稳健。其总能量误差曲线与经典的 MBE2-RHF 基准线几乎完美重合(如图 2 所示)。对于丁烷,总误差为 $+4.9 \text{ kcal/mol}$;对于最长的己烷($n=26$),总误差为 $+365.6 \text{ kcal/mol}$。相对于全分子计算,MBE2-VQE 的相对能量偏差在丁烷处仅为 0.005%,在己烷处也仅为 0.058%。这证明 VQE 在硬件高效电路上通过 SPSA 成功且稳定地找到了各碎片的活动空间基态。

2. ADAPT-VQE 的不平衡收敛与正向漂移

令人意外的是,ADAPT-VQE 在长链烷烃中表现出了显著的正向误差积累(如图 2 紫色曲线所示)。在丁烷处,其总能量显著低于 MBE2-RHF 约 $24.8 \text{ kcal/mol}$(表明它在活动空间内找到了更低、更精确的能级)。然而,到了己烷($n=26$),其积累误差居然飙升至 $+457.5 \text{ kcal/mol}$,大幅偏离了 MBE2-RHF 基准。

不平衡收敛定理推导: 这一异常现象源于多体展开公式对碎片能量收敛不一致的敏感性。令 $\epsilon_{1}, \epsilon_{2}, \epsilon_{12}, \epsilon_{22}$ 分别为单体 $\text{CH}_3^\bullet$、单体 $^\bullet\text{CH}_2^\bullet$、双体 $\text{CH}_3-\text{CH}_2$、双体 $\text{CH}_2-\text{CH}_2$ 在 ADAPT-VQE 迭代中的微小残差。将带误差的碎片能量代入 MBE2 组装公式:

$$E_{\text{MBE2}}^{\text{calc}}(n) = 2(E_1 + \epsilon_1) + (n-2)(E_2 + \epsilon_2) + 2(E_{12} + \epsilon_{12} - E_1 - \epsilon_1 - E_2 - \epsilon_2) + (n-3)(E_{22} + \epsilon_{22} - 2E_2 - 2\epsilon_2)$$

对上式进行同类项消去与化简,得到最终的总装配误差公式:

$$\epsilon_{\text{total}} = 2\epsilon_{12} + (n-3)\epsilon_{22} - (n-4)\epsilon_2$$

当 ADAPT-VQE 运行于不同碎片时,由于单体 $\epsilon_2$ 尺寸小(14 比特)、对称性高,其电路极易被优化至极高精度($\epsilon_2 \approx 0$)。而双体 $\epsilon_{22}$(28 比特)由于哈密顿量复杂,优化残差较大($\epsilon_{22} > 0$)。此时,上式中的误差项化简为:

$$\epsilon_{\text{total}} \approx (n-3)\epsilon_{22}$$

这意味着,双体收敛不足的残差 $\epsilon_{22}$ 将会随着链长 $n$ 呈线性放大! 本文深刻地指出了这一点:量子求解器在多体展开应用中,不能仅追求单个碎片的“绝对极低能量”,而必须保证所有碎片之间的精度均衡收敛,否则多体组装公式中的减法项会成为灾难性的误差放大器。

3. 真实硬件上的 MBE2-SQD 运行数据

SQD 算法直接部署在 IBM 拥有 156 个物理量子比特的 “ibm_pittsburgh” Heron 超导处理器上。每个碎片使用 25,000 个测量采样。所得数据堪称惊艳:

  • 无活性空间压缩:直接计算 14 到 30 比特的哈密顿量。
  • 数据追踪:MBE2-SQD 能量曲线几乎与经典的 MBE2-CCSD 理论上限平行。对于丁烷,SQD 误差为 $+10.0 \text{ kcal/mol}$(仅比经典 CCSD 碎片误差 $+8.6$ 高出 $+1.4$);对于己烷($n=26$),SQD 测得的总误差为 $+401.7 \text{ kcal/mol}$(仅比经典高出 $+5.2$)。
  • 每碳额外硬件惩罚:随着碳链增长,在真实的超导硬件上,SQD 每多处理一个 $\text{CH}_2$ 单元,仅产生 $+0.17 \text{ kcal/mol}$ 的极其微小的测量噪声偏差。这直接证实了在真实硬件上,采用量子采样对角化策略能够完全免疫由于变分线路参数发散带来的灾难性硬件噪声。

2.3 关键资源节省指标

本研究最显著的技术成就体现在对量子硬件资源的巨大节省:

  • 量子比特数极限削减(图 3 & 图 4)
    • 丁烷($n=4$):从全分子直接计算所需的 60 比特 降至最大碎片(双体)所需的 30 比特,实现 2.0 倍 缩减。
    • 十烷($n=10$):从 144 比特 降至 30 比特,实现 4.8 倍 缩减。
    • 己烷($n=26$):从 368 比特 降至 30 比特,实现 12.3 倍 的惊人缩减(比特数缩减率高达 91.8%)。
    • 若使用 CAS 压缩,所有烷烃的求解器物理比特需求更是被统一降至固定的 8 比特
  • 计算唯一性与对称性红利
    • 对于己烷($n=26$),全分子包含 26 个单体和 25 个相邻双体(共计 51 个碎片)。
    • 基于平移对称性,只需要计算 4 个独特碎片。其余 47 个碎片的计算被完全免除,计算效率提升了 12.8 倍。这一红利会随着链长的增加呈线性无限放大,使得计算超大聚合物分子成为可能。

3. 代码实现细节、复现指南与开源工具链

3.1 软件栈与开源仓库链接

要完整复现该论文的研究成果,需要构建以下经典的经典-量子混合软件栈:

  1. 经典量子化学计算PySCF (Python-based Simulations of Chemistry Framework)。用于分子的几何构建、ROHF 自洽场计算、活性空间 CASSCF 选择、经典 CCSD 振幅提取及双电子积分生成。
  2. 量子算法与模拟器Qiskit。用于执行 Jordan-Wigner 变换、构建 VQE/ADAPT-VQE 变分电路、状态矢量模拟以及调用 SPSA 优化器。
  3. 硬件执行平台Qiskit IBM Runtime。用于在 IBM Quantum 真实超导量子处理器上提交和执行 SQD 采样任务。

3.2 核心算法复现:自由基碎片生成与哈密顿量构建管道

以下是一个完整的、可执行的 Python 复现脚本示例。该脚本展示了如何使用 PySCF 构建甲基自由基($\text{CH}_3^\bullet$)的开壳层 ROHF 分子轨道,进行 $\text{CAS}(5, 4)$ 活性空间压缩,并利用 Qiskit 将其映射为 8 比特的费米子哈密顿量算符:

import numpy as np
from pyscf import gto, scf, mcscf
from pyscf.symm import param
from qiskit_nature.second_q.drivers import PySCFDriver
from qiskit_nature.second_q.mappers import JordanWignerMapper
from qiskit_nature.second_q.operators import FermionOperator

# ======= 第一步:经典分子几何定义与 ROHF 自由基计算 =======
# 甲基自由基 CH3* 具有 9 个自旋电子 (5个alpha, 4个beta),自旋多重度为双重态 (2S+1 = 2)
mol = gto.M(
    atom='''
    C    0.000000    0.000000    0.000000
    H    0.000000    1.080000    0.000000
    H    0.935307   -0.540000    0.000000
    H   -0.935307   -0.540000    0.000000
    ''',
    basis='sto-3g',
    charge=0,
    spin=1, # 2S = 1 (未成对电子数为1)
    verbose=0
)

print("Executing Classical Restricted Open-Shell Hartree-Fock (ROHF)... ")
# 必须使用开壳层 ROHF 避免自旋污染
my_rohf = scf.ROHF(mol)
my_rohf.kernel()
print(f"ROHF Energy = {my_rohf.e_tot:.8f} Hartree")

# ======= 第二步:完全活性空间自洽场 (CASSCF) 压缩 =======
# 选择活跃空间:5个活性电子,4个活性空间轨道 -> 对应 8 个自旋轨道 (Qubits)
n_active_electrons = 5
n_active_orbitals = 4

print(f"Performing CASSCF({n_active_electrons}, {n_active_orbitals}) optimization...")
my_cas = mcscf.CASSCF(my_rohf, n_active_orbitals, n_active_electrons)
my_cas.kernel()
print(f"CASSCF Energy = {my_cas.e_tot:.8f} Hartree")

# ======= 第三步:提取分子积分并进行费米子哈密顿量构建 =======
# 从 CASSCF 提取单体核心有效势 (Frozen Core) 以及活性空间内的单、双电子积分
h1_active = my_cas.get_hcore()
# 进行轨道系数变换得到活性空间下的MO表象积分
mo_coeff = my_cas.mo_coeff[:, my_cas.ncore:my_cas.ncore+n_active_orbitals]
# 计算活性空间一电子积分
h1_mo = np.einsum('pi,qj,pq->ij', mo_coeff, mo_coeff, my_rohf.get_hcore())
# 计算活性空间二电子积分 (使用 PySCF 内部高效转换)
iri = mol.ao2mo(mo_coeff)
eri_mo = iri.reshape((n_active_orbitals,)*4)
# 冻结核能量校正
e_core = my_cas.energy_nuc() + my_cas.e_casdepth # 包含闭壳层核轨道能量贡献

print(f"Active space integrals generated. Core energy = {e_core:.6f} Ha")

# ======= 第四步:Qiskit Nature 2nd-Quantization 与 Jordan-Wigner 映射 =======
# 构建量子费米子哈密顿量
fermionic_hamiltonian = FermionOperator()

# 1. 填入单电子算符项 t_ij a^i_dagger a_j
for p in range(2 * n_active_orbitals):
    for q in range(2 * n_active_orbitals):
        # 保持自旋守恒 (p 和 q 必须同为 alpha 或 beta)
        if p % 2 == q % 2:
            p_spatial = p // 2
            q_spatial = q // 2
            coeff = h1_mo[p_spatial, q_spatial]
            if abs(coeff) > 1e-12:
                fermionic_hamiltonian += FermionOperator(((p, 1), (q, 0)), coeff)

# 2. 填入双电子算符项 0.5 * V_pqrs a^p_dagger a_q_dagger a_s a_r
for p in range(2 * n_active_orbitals):
    for q in range(2 * n_active_orbitals):
        for r in range(2 * n_active_orbitals):
            for s in range(2 * n_active_orbitals):
                if (p % 2 == r % 2) and (q % 2 == s % 2):
                    p_sp, q_sp, r_sp, s_sp = p//2, q//2, r//2, s//2
                    val = 0.5 * eri_mo[p_sp, q_sp, r_sp, s_sp]
                    if abs(val) > 1e-12:
                        fermionic_hamiltonian += FermionOperator(((p, 1), (q, 1), (s, 0), (r, 0)), val)

# 加上恒等算符对应的核排斥与冻结核能量
fermionic_hamiltonian += FermionOperator((), e_core)

# 3. 执行 Jordan-Wigner 映射生成 Qubit 算符
mapper = JordanWignerMapper()
qubit_op = mapper.map(fermionic_hamiltonian)

print(f"\nSuccessfully generated Qubit Hamiltonian!")
print(f"Number of Qubits = {qubit_op.num_qubits}")
print(f"Total number of Pauli terms = {len(qubit_op)}")

3.3 经典-量子多体展开重构工作流指南

为了完整复现整条能量重构曲线,研究者应遵循以下步骤:

  1. 坐标生成:编写一个脚本,按照全反式构象公式,生成 $n=4$ 到 $n=26$ 的线性烷烃的 .xyz 笛卡尔坐标文件。
  2. 碎片基础计算
    • 构建丁烷的几何结构。在 C2-C3 键中间截断,使用上文中的方法,分别运行 CH3**CH2*CH3-CH2**CH2-CH2* 碎片的自洽场计算。
    • 提取各碎片的电子总能量($E_I$ 和 $E_{IJ}$)。
  3. 量子模拟/硬件采样
    • 对 4 个碎片分别执行 VQE、ADAPT-VQE 模拟,或将哈密顿量提交至 IBM Quantum Experience 超导后端(如 ibm_pittsburgh)执行 SQD 线路运行与后处理。
  4. 多体组装
    • 得到 4 个碎片在相应量子求解器下的能量值:$E_{\text{CH}_3^\bullet}$、$E_{^\bullet\text{CH}_2^\bullet}$、$E_{\text{CH}_3-\text{CH}_2}$、$E_{\text{CH}_2-\text{CH}_2}$。
    • 按照 MBE2 二体解析重构公式(见第 1.3 节公式 (2)),计算任意碳原子数 $n$ 下的烷烃估算总能量 $E_{\text{MBE2}}(n)$。

4. 关键引用文献与局限性批判评述

4.1 关键引用文献

  1. Kitaura et al. (1999) [Chem. Phys. Lett. 313, 701]:首次提出了碎片分子轨道(FMO)方法,奠定了大分子拆分计算的基石。
  2. Fedorov & Kitaura (2007) [J. Phys. Chem. A 111, 6904]:系统阐述了多体展开在共价键断裂中的 FMO 应用,是本文二体截断(MBE2)公式的重要理论先导。
  3. Peruzzo et al. (2014) [Nat. Commun. 5, 4213]:变分量子本征求解器(VQE)的奠基性工作。
  4. Grimsley et al. (2019) [Nat. Commun. 10, 3007]:首次提出 ADAPT-VQE 算法,展示了动态生长量子变分算符池的优越性。
  5. Robledo-Moreno et al. (2025) [Sci. Adv. in press]:样本量子对角化(SQD)算法的理论源头,为 NISQ 硬件的高效无变分对角化提供了技术路径。

4.2 对本项工作的深度局限性批判

尽管本工作在打破量子比特线性增长红线、实现真实量子硬件计算长链分子方面取得了非凡的进展,但作为一名理性的量子化学研究人员,我们必须客观地指出其存在的数个底层技术局限性:

1. 物理层面的“截断误差”不可忽略

该方法将多体展开强行截断在二体键合级别(MBE2-bonded-only),虽然换取了超高的计算压缩比,但其带来的系统性物理误差(如己烷处的 $+396.5 \text{ kcal/mol}$)远远超出了化学精度($1 \text{ kcal/mol}$)的要求。尽管作者强调“每碳原子平均误差已收敛”,可以通过经典的经验线性公式进行系统性外推和校正,但这在本质上依然是一种半经验的折中。若要追求绝对的化学精度,必须将非键合双体相互作用(Non-bonded long-range Coulomb/Dispersion interactions)纳入考量,或者引入三体展开项(MBE3)。而这两者都会立即使碎片的唯一性数量(计算量)以及最大量子比特数显著攀升,破坏目前的 $O(1)$ 恒定资源开销优势。

2. 基组选择的局限(STO-3G 玩具模型)

所有计算均在最小基组 STO-3G 下完成。这在现代电子结构理论中被视为“玩具模型(Toy model)”,其物理描述和轨道极化能力极其有限。一旦升级到中大型计算常用的高精度双分裂价键基组(如 cc-pVDZ6-311+G(d,p)),单个碳原子的空间轨道数将呈数倍增长。届时,即使是二体碎片的比特需求也将迅速逼近甚至超越目前物理硬件的极限(例如,在 cc-pVDZ 下,$\text{CH}_3-\text{CH}_2$ 双体将需要接近 100 个量子比特),这将对活性空间压缩(CASSCF)技术提出更加严苛的稳定性挑战。

3. ADAPT-VQE 的残差不平衡漂移揭示了框架的脆弱性

论文中发现的 ADAPT-VQE 在长链烷烃下的能量正向漂移(+457.5 kcal/mol)极具警示意义(我们在第 2.2.2 节中给出了其不平衡收敛的数学证明)。多体展开公式中的减法抵消项(如 $\Delta E_{IJ} = E_{IJ} - E_I - E_J$)是一把双刃剑,它在消除多余能量贡献的同时,也极大地放大了量子求解器在不同大小碎片上表现出的精度不均一性。这表明,该框架极度依赖量子求解器精度的“绝对均衡”。如果未来的研究将此方法推广到更复杂的分子拓扑中,如何设计自适应的收敛标准以确保所有碎片的误差协同一致,将是一个极难调和的技术泥潭。

4. 体系特异性强,泛化难度极高

该均裂自由基碎片方案几乎是为非极性、高度平移对称的线性烷烃“量身定制”的。面对以下体系,其优势将荡然无存或面临底层重构:

  • 极性共价键体系(如聚醚、聚酯):均裂会导致电荷不均匀分配的物理失真,必须重新引入复杂的静电嵌入。而静电嵌入哈密顿量的引入会使 Qubit 算符中的 Pauli 项数呈几何级数爆炸,再次推高测量成本。
  • 非对称分子/支链烷烃:哪怕只是引入一个甲基支链,都会导致“独特碎片”的数量激增,对称性红利大打折扣。

5. 补充探讨:均裂自由基碎片的物理化学机制与硬件适配前瞻

为了给读者提供更加立体的物理化学图景,我们在此补充两个维度的高阶探讨。

5.1 均裂 vs 异裂:深层次的自旋与自旋污染(Spin Contamination)问题

在开放壳层体系的经典量子化学计算中,处理未成对电子一直是个难点。通常有两种处理方式:

  1. 无限制哈特里-福克(UHF):允许 $\alpha$ 自旋轨道和 $\beta$ 自旋轨道具有不同的空间部分。虽然能降低均值场能量,但其波函数不是自旋平方算符 $\hat{S}^2$ 的本征态,会引入严重的自旋污染(高自旋态掺杂),导致物理性质预测失真。
  2. 受限开壳层哈特里-福克(ROHF):强制所有配对电子共享相同的空间轨道,只有未成对的活性电子独占单独的轨道。ROHF 的波函数是严格的自旋本征态,彻底消除了自旋污染。

本研究之所以能够完全摆脱传统 FMO 的复杂外加氢饱和,核心就在于其选择在 ROHF 层面处理均裂产生的开壳层自由基。这保证了每一个被切断的共价键电子在物理上处于高度局域化的自由基自旋状态,保留了最自然的局部电子结构。实验证明,这种处理方式在与随后的关联能方法(如 CCSD 或量子 VQE、SQD)结合时,表现出了极佳的相容性,避免了异裂碎片因人为闭壳层化而不得不引入人工电荷与饱和氢的尴尬。

5.2 硬件适配前瞻:从 NISQ 碎片化到早期 FTQC 的黄金桥梁

随着量子计算技术从 NISQ(嘈杂中等规模量子)阶段向早期的 FTQC(容错量子计算)阶段迈进,物理比特与逻辑比特的转化率(Overhead)成为了最昂贵的系统开销。在容错量子计算的极早期,受限于物理纠错码的配比(如 Surface Code 中一个逻辑比特需要数千个物理比特),研究人员很可能在很长一段时间内只能操作数十个(例如 50 个)高保真度的逻辑量子比特

在这种背景下,本文提出的自由基碎片 MBE2 方案展现出了极具前瞻性的实用价值:

  • 逻辑比特资源完美适配:它将一万个碳原子的超长链聚乙烯分子的计算需求,直接压缩到了 30 个逻辑量子比特的极小物理上限内。对于早期只能提供 30-50 个逻辑比特的 FTQC 芯片,这无异于一场及时雨——我们可以在完全无需等待千比特级 FTQC 硬件成熟的前提下,直接在数十个逻辑比特上实现对无限长链聚合物体系具有完全关联能精度的量子化学模拟。
  • 量子并行化的终极形态:由于 4 个独特碎片的计算是完全物理孤立的,没有任何经典或量子数据的交叉依赖。这意味着它们可以被完美分发到不同的量子芯片上并行执行(Quantum Parallelism),进一步缩短整体计算壁垒,为云端量子化学计算服务(Quantum-Chemistry-as-a-Service)的商业落地铺平了道路。