量子处理器上大型多体哈密顿量的Krylov对角化

多体系统的低能量估计是计算量子科学的基石。变分量子算法可用于在容错前量子处理器上准备基态，但它们缺乏收敛保证和不切实际的成本函数估计数量阻止了将实验系统地扩展到大型系统。在预容错设备上进行大规模实验需要变分方法的替代方案。在这里，我们使用超导量子处理器在多达 56 个位点的二维晶格上计算量子多体系统的本能，使用 Krylov 量子对角化算法，该算法类似于著名的经典对角化技术。我们使用在量子处理器上执行的 Trotterized 酉进化来构建多体 Hilbert 空间的子空间，并在这些子空间内经典地对角化多体交互哈密顿量。这些实验证明了对基态能量估计的指数收敛，并表明量子对角化算法已准备好在量子系统计算方法的基础上补充其经典算法。

介绍

求解量子多体系统的薛定谔方程是凝聚态物理学、量子化学和高能物理学等领域许多计算算法的核心。这项任务的量子优势将对自然科学产生深远的影响。在使用量子计算机进行特征态计算的方法中，迄今为止有两个主要讨论对象：量子相位估计 （QPE）1,2，包括它最近的进展（例如，refs.3、4、5） 和变分量子特征求解器 （VQE）6.前容错器件的实验实现主要集中在 VQE 上，VQE 已在各种实验平台上针对各种问题（例如 refs.6、7、8、9).然而，到目前为止，参数优化的瓶颈阻碍了它扩展到小型实例之外。另一方面，QPE 具有理论精度保证，但量子纠错对于达到有价值问题所需的电路深度是必要的，尽管已经实现了一些小例子10、11、12.

这些结果在迄今为止已经执行的小型演示与在容错量子计算机上使用 QPE 或相关方法的大规模、高精度模拟之间留下了特征状态估计方法的差距。在这项工作中，我们证明了 Krylov 量子对角化 （KQD）13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28，一种量子子空间对角化13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41，42，可以填补更普遍问题的空白。

KQD 的主要思想是使用量子计算机来近似地将哈密顿量投影到由初始参考状态的各种时间演变跨越的 Krylov 空间中。然后将得到的低维矩阵进行经典对角化，以获得近似的低位能量特征态13.这种方法与 VQE 共享变分特性（最高可达噪声的影响），但不需要迭代参数优化，而是依赖于单轮电路执行，然后进行经典的后处理。此外，只要噪声可以量化，该方法的准确性就可以在理论上受到限制27,28 元，这意味着 KQD 在过渡到容错时代时可以继续保持价值。在短期内，对于现有的量子计算机来说，精度要求不那么严格的模拟的时间演变并不是令人望而却步的。虽然在这种情况下量化所有误差源的影响可能具有挑战性，但 KQD 仍然作为一个启发式方法，具有潜在的指数收敛，朝向基态能量的噪声估计。

在这项工作中，我们使用 KQD 来估计海森堡模型在重六边形晶格上的基态能量。我们表明，尽管噪声对高精度构成了重大障碍，但即使使用高级错误缓解功能也是如此43,44 元，我们可以在多达 56 个量子比特上获得到基态能量的收敛。

结果克雷洛夫量子对角化理论

KQD 包括两个主要步骤。第一个是构造矩阵的量子子例程

$${tilde{H}}_{jk}=langle {psi }_{j}|H|{psi }_{k}rangle，qquad {tilde{S}}_{jk}=langle {psi }_{j}|{psi }_{k}rangle，$$(1)

它对应于哈密顿量投影到子空间 （{{mathcal{K}}}={{rm{Span}}}{leftvert {psi }_{0}rightrangle，ldots，leftvert {psi }_{D-1}rightrangle }）的重叠（Gram）矩阵。第二步是经典求解投影到子空间中的与时间无关的薛定谔方程，该方程由下式给出

$$tilde{H}c=Etilde{S}c，$$(2)

其中 c 是 Krylov 空间中的坐标向量。在整个希尔伯特空间或对称扇区内的近似基态能量为 （2） 的最低特征值。两个不同的分量会影响近似的精度27,28 元：将完整特征值问题向下投影到子空间的内在误差，这与足够低能量的状态与子空间的重叠有关，以及任何其他算法、统计和硬件错误。

子空间对角化方法的主要区别在于子空间的选择。在经典计算中，一种常见的方法是通过局部运算符（例如费米子的跳跃项）生成相关性来构造子空间，就像在多引用配置交互中一样45.或者，可以使用全局运算符。例如，经典的 Lanczos 方法采用哈密顿量的幂级数来构造子空间 （{{{mathcal{K}}}}_{P}={{rm{Span}}}{{H}^{，j}leftvert {psi }_{0}rightrangle }），也称为幂或多项式 Krylov 空间。这种结构的主要优点是，解的精度随着子空间大小 D 的出现呈指数级提高46,47,48 元.经典 Lanczos 和相关方法的限制因素是，由于需要表示纠缠量子态，它们不可避免地会受到内存消耗的影响，内存消耗会随着系统大小呈指数级增长。

虽然已经提出了这种方案对量子计算机的各种改编13,14,15,16,17,18,21,22,23,24,25,26,27,28,29,30,31,32,35,36,37,38,39,40,41,42，最适合近期量子计算机的是使用实时演化作为全局算子来生成 Krylov 空间：

$${{{mathcal{K}}}}_{U}={{rm{Span}}}{{U}^{， j}leftvert {psi }_{0}rightrangle }，quad j=0,1，ldots，D-1，$$(3)

其中 U = e−iH dt是某个时间步 dt 的时间演化算子13,14,15,16,17,18,21,22,23,24,25,26,27,28.这样做的优点是双重的：首先，时间演化可以通过足够短的深度电路来近似，以便在现有量子设备上实现。其次，可以证明，即使在存在噪声的情况下，投影到这个酉 Krylov 空间引起的误差也会与 Krylov 维度呈指数级快速收敛，就像在经典的 Krylov 算法中一样。噪声只是贡献一个额外的误差项，只要它不是太大以至于完全压倒信号27,28 元.这意味着有可能达到近似基态能量与有限维的 Krylov 空间的收敛。

虽然在量子计算机上评估 Krylov 矩阵解决了内存问题，这是经典方面扩展的主要障碍，但量子方面的主要障碍是噪声。两个主要贡献是有限散粒采样导致的统计噪声和器件中的硬件噪声。来自时间演变近似的算法误差也进入了，但下面我们用数字表明它的影响低于硬件误差的水平。另一方面，为了扩大模拟的规模，抑制和减轻这些硬件错误被证明是至关重要的：我们为此目的应用了实验技术（详见补充说明 4），同时保持量子电路尽可能浅，同时保持 Krylov 空间的全局耦合结构。

为了简化我们的电路，我们利用了许多凝聚态模型所具有的 U（1） 对称性，包括我们关注的海森堡模型。作为量子比特运算符，U（1） 对称性可以表示为汉明权重守恒;就 spin-1/2 算子而言，它对应于总自旋的 z 分量守恒。等效地，我们可以将对称子空间视为 k 粒子子空间，将 ↑（↓） 自旋视为粒子的不存在（存在）。

图 1 显示了原则上可用于计算矩阵元件的电路序列 （1）。面板 （a） 显示了标准 Hadamard 检验，这将是此类计算的默认工具。图 （b） 说明了我们如何使用自旋守恒来避免实现传统 Hadamard 测试中存在的受控时间演化：相反，我们实现了参考状态 （leftvert {psi }_{0}rightrangle） 的受控初始化，然后依赖于时间演化保持“真空态”（leftvert 00ldots 0rightrangle） 直到经典可计算阶段的事实。

图 1：Krylov 量子对角化示意图。

用于计算 〈 形式的矩阵元素的 Hadamard 电路ψ我∣P∣ψj〉，它依赖于 Krylov 基态的受控单一实现。b 通过利用对称性（如粒子数守恒）来简化电路。c 这项工作中采用的结构。只需要一个时间演进电路，第二个受控制备电路被吸收到测量的基础上。d 经典后处理构造矩阵 （tilde{H}） 和 （tilde{S}），这会产生一个广义特征值问题。矩阵对于 a、b 中所示的电路是 Hermitian，对于 c，它们是 Toeplitz Hermitian。请注意，由黑线括起来的对角线元素可以经典地计算。

作为第二个简化，我们注意到对于确切的时间演变，（langle {psi }_{0}|{U}_{j}^{{dagger} }H{U}_{k}|{psi }_{0}rangle=langle {psi }_{0}|H{U}_{k-j}|{psi }_{0}rangle） 中，它为我们提供了两种形式上等效的方法来测量相同的矩阵元素，第二种方法产生了更简单的电路，因为它只涉及一次时间演化。然而，一旦时间演变通过 Trotterization 近似化，这两个表达式就不再相等。在图 .1c，我们显示了与后一个版本相对应的电路。

目前尚不清楚人们是否应该更喜欢图 b 或 c 中所示的电路。1，纯粹从 Trotter 错误的角度来看。图 .1b 是它仍然对应于子空间中的变分优化，因为每个矩阵元素仍然具有 （1） 形式。然而，在存在有限样本和器件噪声的情况下，即使如此28.图 1c 是显式强制执行 Toeplitz 结构的版本，从电路深度的角度来看更可取，原因有两个：它只需要一次演化，因此，第二次受控初始化可以作为 Clifford 变换应用于哈密顿量中的 Pauli 可观察对象，而不是在电路中显式实现。在实践中，由于使用了正则化技术来避免特征值问题 （2） 的不良调节，我们没有看到这种方法对变分的显著破坏（详见补充注释 5）。例如，我们比较了图 1 中电路的精确经典仿真结果。1b、c，它们将在实验中编译，如图 1 所示。2 （有关详细信息，请参阅下一节）。图 20 量子比特系统的能量曲线。3 表示变分随着 Trotter 步数的增加而迅速恢复。这些发现使用图 1 中所示的电路版本进行激励。1c.

图 2：Krylov 量子对角化的量子电路。

a 每个回路在目标粒子扇区内执行初始状态的受控准备，然后是 Trotterized 时间演化。b 受控制备准备一个计算基态，其中汉明权重对应于给定实验的粒子数，在辅助量子比特上进行控制。由于重六边形晶格可以是边缘三色的（图中颜色由红色、绿色和蓝色表示），因此可以使用三个独特的双量子比特门层的序列来实现受控准备和 Trotterized 时间演化，这些层与单量子比特旋转交错。有关详细信息，请参阅正文。c 每层双量子比特门都经过 Pauli 旋转，以便根据稀疏 Pauli-Lindblad 噪声模型 Λ 定制噪声43,44 元，前面是其扩增 ΛG对于 PEA。请注意，源自源电路、旋转层或噪声放大层的单量子比特门的相邻层始终组合在单个层中;为了清楚起见，它们未合并到图中。CZ 层的 d （12 + 1） 量子比特示例。

图 3： 算法误差的数值研究。

用于数值模拟的海森堡模型的 （20 + 1） 量子比特布局，绿色和红色圆圈表示控制和激发量子比特。b 能量与克雷洛夫空间维度。虚线和实线表示图 1 中电路的结果。1b、c 分别。c 使用 4 个二阶 Trotter 步骤，k = 5 粒子扇区具有不同 dt 和 D 的基态能量误差 ΔE 的热图。白色箭头表示 π/||H||.

大型实验演示

在我们的实验中，我们研究了自旋 1/2 反铁磁海森堡模型，该模型将一组边 E 定义为

$$H={sum}_{（i，j）in E}{J}_{ij}（{X}_{i}{X}_{j}+{Y}_{i}{Y}_{j}+{Z}_{i}{Z}_{j}）$$(4)

带均匀联轴器 J我J= 1，其中 X我、 Y我、 Z我表示第 i个站点上的 Pauli 矩阵。相互作用集 E 是重六边形晶格的子集（见图 D）。4）. 请注意，虽然重六边形晶格是二分的，因此可以使用路径积分蒙特卡洛方法有效地模拟整个希尔伯特空间中的基态49，符号问题通常存在于激发态中。在激发态中，我们关注几个 k 粒子子空间中能量最低的特征态。k 粒子子空间的维度为 O（Nk).请注意，电路构造依赖于 U（1） 对称性，而不是 SU（2） 对称性，因此我们的方法也直接适用于 XXZ 模型。

图 4：多体哈密顿量的实验对角化。

a 在系统大小分别为 N = 56、44 和 42 的情况下，粒子数 k = 1、3 和 5 的海森堡模型的每位点能量。误差线表示通过 bootstrap 估计的标准差。虚线表示无噪声经典模拟的能量，黑色实线表示给定 k 粒子子空间中的确切最低能量。b–d 量子比特布局图。绿色和红色圆圈分别表示粒子的控制和初始位置。例如，单个粒子数 k = 1、3 和 5 的能量曲线。h–j 误差矩阵 （Delta tilde{H}/N ：=|{tilde{H}}_{exp }-{tilde{H}}_{{{rm{num}}}}|/N） 和 （Delta tilde{S} ：=|{tilde{S}}_{exp }-{tilde{S}}_{{{rm{num}}}}|），其中下标 “exp” 和 “num” 分别表示来自实验和数值计算的数据

我们在三个不同的 k 粒子扇区中运行了实验：k = 1、3、5。所有三种情况下的初始状态都是计算基态，其数字为 （leftvert 1rightrangle） s，由 k 给出：例如，在单粒子情况下，（leftvert {psi }_{0}rightrangle=leftvert 00ldots 1ldots 0rightrangle）。因此，k的不同值的电路实现在受控准备中有所不同（见图）。1 和 2）。k = 1 情况对应于在初始状态下仅生成一个粒子，这可以通过控制量子比特（Hadamard 测试中的辅助）和相邻量子比特之间的 CX 门轻松实现。对于 k > 1，我们为粒子选择了大致均匀分布在量子比特图上的位置。

重六边形晶格允许对其边缘进行三种着色，其中每种颜色对应于一层可以同时实现的双量子比特门（见图 D）。2）. 由于每个不同的双量子比特层都需要自己的噪声学习以进行概率误差放大 （PEA—见下文），因此最大限度地减少电路中不同层的数量是有利的。受控准备电路可以使用一组双量子比特层来实现，这些层对应于重六边形中边缘的三种颜色，与任意层相比，开销是恒定的（有关详细信息，请参见补充说明 2）。对于我们的 Trotterized 时间演化，我们将哈密顿项划分为同一组层。因此，对于每个实验，我们只需要学习总共三个独特层的噪声模型。

电路的受控初始化部分的深度与量子比特图中两个相距最远的初始粒子之间的距离成正比。我们使用两个二阶 Trotter 步骤来近似我们所有实验中的时间演变。具有三个哈密顿项交换组的 r 二阶 Trotter 步骤需要 4 个r +1 个双量子比特层（参见图 b 中的面板 b）。2），在我们的例子中，电路的时间演化部分产生了 9 层。

为了测量对应于 （tilde{H}） 和 （tilde{S}） 中矩阵元素的实部或虚部的可观察量，我们将可观察量划分为尽可能少的局部交换集（测量基数），因为这些集合是可共同测量的7.缩短的电路，如图 3 的第三行所示。1，需要第二个受控初始化电路共轭哈密顿项 （4），因为它没有物理实现。这产生了相同数量的 Pauli 可观察对象，因为受控初始化是一个 Clifford 电路，并且可以证明这些可观察对象可以划分为 2（k + 2） 个测量基;见补充说明 2。

我们在 Heron R1 处理器IBM_montecarlo上进行了实验。这是一个 133 量子比特设备，具有固定频率的 transmon 量子比特，通过可调谐耦合器相互连接。与前几代固定耦合器件相比，Heron 处理器具有更快的双量子比特门（持续时间与单量子比特门相似）和更低的串扰。为了进一步改进测量的可观察对象（参见图 .1），我们使用了概率误差放大 （PEA）44和旋转读出误差消除 （TREX）50，以减少 SPAM 错误，以近似无噪声期望值。我们还采用了错误抑制，特别是 Pauli 旋转和动态解耦。增补说明 4 中提供了错误缓解和抑制的详细信息。

在我们的实验中，对于每个测量基础，生成了一定数量的旋转实例，然后针对不同的噪声放大因子值重复测量每个实例。对于单粒子 （k = 1） 实验，我们使用了 300 个旋转实例，每个实例 500 次拍摄，噪声放大系数为 1、1.5、3。对于 k = 3、5，我们使用了 100 个旋转实例，每个实例有 500 次拍摄，噪声系数为 1、1.3、1.6。为较大的实验引入减少旋转实例是为了减少总运行时间，因为测量碱基的数量以及电路大小随着 k 的增加而增加。噪声放大因子的调整是由于较深电路中的噪声率增加。电路的受控初始化部分涉及创建控制量子比特的最大纠缠状态和初始粒子位置。随着粒子数量的增加，这转化为在电路开始时准备的更大的最大纠缠状态，这反过来又使结果更容易受到噪声的影响。

在所有实验中，Krylov 空间的大小固定为 D = 10，以便在设备重新校准过程的时间尺度 （24 h） 内实现算法的总运行时间。对于固定值 k ，实验在特定的量子比特子集上运行，该子集根据设备的当前状态使用启发式例程进行选择，以实现最佳量子比特映射51.k = 1 实验在 57 量子比特子集上执行，k = 3 实验在 45 量子比特子集上执行，k = 5 实验在 43 量子比特子集上执行（布局如图 4 所示）。4）. 后两个部分是手工选择的，以便在每种情况下有五个完整的重六角形。

尽管时间步长 dt 理论上具有 π/||H||27,28 元，则对低粒子数子空间的限制会改变这一点。因此，我们启发式地选择了时间步长，k = 1、3 和 5 的值分别为 0.5、0.022 和 0.1。

结果如图 1 所示。4. 图 a 在归一化能量尺度上总结了结果，而 e、f 和 g 显示了每个单独实验的收敛曲线。相应的量子比特图分别显示在面板 b、c 和 d 中。这些收敛曲线是评估噪声 KQD 实验结果的有用诊断工具。我们从理论分析中了解到，如果误差率足够低，可以解析信号，即区分 Krylov 空间中的最低能量状态和纯噪声，那么我们应该看到能量向与真实基态能量偏移一个常数的值呈指数衰减，具体取决于误差率27,28 元.我们的结果表明，这种行为在一定程度上会有所波动，这是意料之中的，因为理论结果只提供了指数衰减的上限。然而，如果噪声完全主导了信号，那么相对于系统大小，与子空间维度的收敛速率将呈指数级缓慢，因此我们不会在图 1 中看到最初的快速收敛。4. 有关更多详细信息，请参阅补充说明 4 和 5。

在我们的实验结果中，噪声和算法误差（由于 Trotter 近似以及有限的 Krylov 维数）仍然是重要的限制因素，最准确的估计能量（D = 10 时）与真实值之间的差异证明了这一点。我们使用自举法估计了实验能量的标准差，因为求解正则化、广义特征值问题 （2） 的后处理使直接误差传播变得困难。这产生了图 1 中的误差线。4;有关更多详细信息，请参阅补充说明 5。图 4 还显示了我们电路的理想经典仿真的能量收敛曲线，这些曲线可以通过仅在受限粒子数子空间中表示向量和运算符来实现。虽然由于实验结果的噪声导致误差线很大，但我们对两个较大 k 值的估计能量与几乎所有点的这些标准差的理想模拟曲线一致。

在 k = 1 实验中，结果偏离于真正的最低能量，表明噪声已经从 k = 1 子空间产生了有效的泄漏。这说明了依赖对称守恒来保留在特定子空间中的风险，尽管研究全局基态不会受到这种担忧。

精确对角化也可以在本实验中研究的 Hilbert 空间扇区中进行，但不能在整个 Hilbert 空间中进行。然而，除了受控初始化的电路深度减小外，实验并未以任何方式依赖于那些特定的粒子数扇区，因此缩放没有定性或结构障碍，只有噪声的影响。在我们关注的特定情况下，海森堡模型在二维重六边形晶格上的基态，仍然可以使用张量网络等经典技术计算精确的近似值。

有人可能会问，为什么使用 KQD 而不是最近为近期或早期容错设置中的基态能量估计而开发的其他各种算法之一，例如，3、4、5、52.一个主要原因，除了 KQD 相对广为人知的噪声容忍度21,25,27,28 元，这些替代方法都从时间演化中提取特征能，而不是直接从哈密顿量本身的投影中提取特征能。在我们这样的环境中，这是一个问题，其中 Trotter 电路随着时间步长数量的增加而保持固定（这对于最小化电路深度是必要的），因为 Trotter 电路的频谱与理想时间演变的频谱不同，实际上变得周期性，周期取决于时间步长和固定的 Trotter 电路。因此，在这种约束下，随着时间步数的增加，仅依赖于演化的算法将在某个时候停止收敛。

讨论

这里介绍的 KQD 方法丰富了在前容错量子处理器上进行基态估计的量子算法的前景，填补了 VQE 和 QPE 之间的空白。一种称为基于样本的量子对角化 （SQD） 的互补子空间算法，基于采样和复杂的经典后处理，使用以量子为中心的超级计算 （QCSC） 架构53最近被用于演示超越蛮力解决方案的化学量子模拟。这种 QCSC 方法产生经典可验证的能量，并且不需要近似时间演化，这使得它在短期内对于包含大量项（例如分子哈密顿量）的哈密顿量来说很容易处理。对于凝聚态应用，KQD 具有可证明的收敛保证，给定具有逆多项式重叠的初始参考状态，并且其电路在预容错处理器上是可行的，如本工作所示。

为了将这项工作的规模与先前的基态能量实验量子模拟进行比较，可以基于量子比特计数或基于使用的希尔伯特空间维度进行比较。两者都是相关的，因为，例如，即使我们的单粒子实验的子空间只有 56 维，它也会产生 57 量子比特量子电路的误差。我们最大的子空间维度是 5 粒子实验，其子空间维度为 850668，介于 19 和 20 个量子比特的完整希尔伯特空间维度之间。我们在表 1 中展示了用于基态能量模拟的端到端量子算法的最大（据我们所知）先前的实验演示，并通过上述两个指标进行评估。请注意，我们只包括实现整个算法的实验，而不是例如，经典地优化 VQE 的参数，然后只在量子计算机上估计单个能量。

1. QC-QMC 是指混合量子-经典量子蒙特卡洛算法。

如表 1 所示，我们的实验在量子比特方面比以前的工作高出两倍多，在希尔伯特空间维度上比以前高出两个数量级以上，但 ref 除外。53，它与上一段讨论的适用范围完全不同。人们可能还会注意到，这些实验所达到的精度差异很大，但正如我们的重点是如上所述，在规模上实现收敛，而规模的增加使我们的工作与小规模实验（包括一些未在表 1 中所示的实验）处于不同的状态，后者实现了更高的准确性。换句话说，这项工作代表了基态能量量子模拟技术的最新进展。