量子多体纠缠计算新范式：抽样约化密度矩阵

如果有一种数值方法，可以跳过求解基态而直接得到约化密度矩阵，将大大拓展纠缠计算的范围。经过多番尝试，我们通过量子蒙卡算法达到了这个目的。

参加过几次学校的博士生面试，某吴老师常问学生：“你觉得量子力学最有趣在哪里？”于我自己而言，我觉得量子纠缠肯定是其中之一了。量子纠缠如此神奇，以至于在民间被神化和佛化了。所谓的“量子佛学”，“几世纠缠”……在各种神棍口中被发扬光大，混吃混喝。然而，思而不学则殆，抛开量子力学和线性代数而空想量子纠缠实则是非常不科学乃至迷信的。建议大家有钱有闲，不妨买本量子物理的科普书，而不是买量子袜子、量子水等智商税产品。

物理人比较实在，我们更关注量子纠缠的计算、探测以及背后的机制等。近期我们在量子多体系统中如何计算纠缠信息有了一番新的认识和探索，在此与大家分享。

要谈量子纠缠，首先需要量化它。数学上有许多反映纠缠信息的量，比如，Von Neumann纠缠熵、Renyi纠缠熵、纠缠谱等等。而它们的定义都往往离不开约化密度矩阵。所谓约化密度矩阵，就是将系统分为两部分：关心的子系统A和不关心的环境B。如果我们对总系统的密度矩阵做部分求迹操作，把B的自由度都积掉，那么我们就得到了关于A的约化密度矩阵，即 ρA =tr_B(ρ)。这个矩阵中隐藏了A和B之间的量子纠缠。

换句话说，如果有一个数值方法，能够求得所需要的ρ_A，那么绝大部分纠缠的物理量就都可以得到了。然而，对于一个量子多体系统，其自由度本身随着尺寸而指数发散。要求解基态本身就是一件极其困难的事情，更何况现在还要将密度矩阵部分求迹，这大大限制了纠缠信息的计算和提取。

目前常用的计算纠缠方法就是，先用严格对角化或者密度矩阵重整化群一类的方法，得到体系的基态，然后得到对应的密度矩阵，接着对其部分求迹，得到ρ_A。受限于现代计算机的内存，这些方法能处理的总体系本身尺寸很小，或者局限于一维。要想再得到其中的A部分，尺寸就得变得更小了。

如果有一种数值方法，可以跳过求解基态而直接得到约化密度矩阵，那计算机的全部内存就只需要处理ρ_A的自由度，这将大大拓展纠缠计算的范围。经过多番尝试，我和学生丁一茗以及康复大学 (筹) 的毛斌斌老师，通过量子蒙卡算法达到了这个目的^[1]！关心细节的朋友们且听我慢慢道来。

传统的量子蒙卡计算着重于得到一个数字，比如配分函数。而约化密度矩阵则是一个矩阵： ρA =tr_B(ρ)。若将此矩阵用一组基展开，则对应的矩阵元满足：，这个公式用路径积分或者投影的方式来理解十分容易，此处不多展开。为了方便比较，我们把配分函数和约化密度矩阵都写成便于蒙卡计算的路径积分形式，则其示意图可如图1所示。

简单来说，如果你要trace掉某部分自由度，则在虚时间上将该部分作循环边界条件。如果你想要对某部分自由度保留矩阵形式，则在虚时间上作开放边界处理。通过这种形式，我们就在路径积分的形式上完成了部分求迹。不难得到，每个矩阵元的大小应该满足：，这里表示抽样到上/下时间边界构型是的数量，N_total表示总抽样数。就这样我们得到了约化密度矩阵ρ_A，记得要把它归一化：tr(ρ_A)=1。

图1. 配分函数(上)和约化密度矩阵(下)的路径积分形式。

上下边界相同颜色的格点表示此处虚时间为循环边界条件，反之则为开放边界条件。

这样一处对于传统量子蒙卡计算方式小小的改动，赋予了计算对象新的含义。从配分函数变成了约化密度矩阵。回过来看很简单，事实上为了追求约化密度矩阵的信息，蒙卡先贤们付出了巨大的努力和尝试。比如，通过构造多腿裤子结构加上解析延拓来尝试得到纠缠谱函数的信息^{[2， 3]}（约化密度矩阵的本征值取ln后就是纠缠谱）；通过计算不同n阶的Renyi纠缠熵来反推纠缠谱的能级^[4]。其中，前者受复本流形和解析延拓所限，只能得到非常低能范围内的谱函数，无法得到能级和简并度等重要信息^[5]。

而我们知道，在共形场论、拓扑序等信息提取时，纠缠谱的简并和精细能级是极其重要的。后者方法虽然能给出精细的能级，但是计算代价是极为巨大的。我们需要计算尽可能多的不同阶数的Renyi纠缠熵来反解纠缠谱。而高精度的Renyi熵计算本身就是极其困难的问题^{[6， 7]}，更何况要穷尽去算不同阶数的熵值。我们的新方法显然是以极低代价，给出了新的解决方案。

图2. (a) 在L=8的海森堡梯子中计算两腿之间的纠缠谱。红色的是严格对角化结果，黑色点为蒙卡抽样结果。可以看到随着抽样步数增加，能级从低到高收敛。(b) 相同系统中，更大尺寸下，蒙卡得到的纠缠谱数据。

接下来我们先看看这个方法的有效性。以一个反铁磁海森堡模型梯子链为例，我们计算了两条链之间的纠缠谱。如图2(a)所示，可以看到，随着蒙卡步数的增加，越低能的纠缠谱能级越早收敛，其中红线是标准的严格对角化结果，而黑色的数据点则是蒙卡结果。此处的梯子长度为L=8。进一步，我们尝试了更大尺寸的纠缠谱数据，如图2(b)所示，这在过去的计算中是难以达到的^[8]。

从这里我们可以简单分析下此方法的优劣。由于我们对环境部分B用蒙卡做了trace计算，所以该部分自由度的计算代价是代数的，故而该方法并不受限于B的大小。但是因为我们要对所关心系统A的自由度做矩阵操作，事实上导致了A部分计算复杂度是指数增长的。

幸运的是，同时，由于蒙卡自带的重要性抽样特性，低能的纠缠谱会被优先得到。而在凝聚态物理中，大家普遍关注的都是低能纠缠谱的普适行为，所以高能谱可以被丢掉。所以如果我们限制蒙卡抽样数量的话，相当于可以对A的态空间做有效截断，只提取其低能信息。故而A部分的计算复杂度也可以根据关注的低能范围做优化。最极端的例子是，如果我们只关注纠缠谱是不是有能隙（所谓的Schmidt gap^[9]），事实上只需要极少数的抽样就可以捕捉该信息。

我们在文章中展示了几个例子^[1]，分别是反铁磁/铁磁海森堡梯子两条腿之间的纠缠谱，以及二维反铁磁海森堡模型的纠缠谱。此处不再一一赘述。我们着重展示下二维系统的例子。如图3(a)所示，我们在二维反铁磁海森堡模型的基态中，选择一条链(A)和一个方块(A’)为子区域。计算了他们各自的纠缠谱，如图3的(b)和(c)所示。其中我们都发现了反映自发对称性破缺特征的塔状（tower of states， TOS）纠缠结构^[10]。这类计算在之前是难以做到的。如果大家对于其他两个例子也感兴趣的话，欢迎阅读我们近期的文章^[1]。

图3. (a) 自旋1/2的二维正方晶格反铁磁海森堡模型中，选择一条链(A)和一个方块(A’)为子区域。计算他们各自的纠缠谱。(b) A区域的低能纠缠谱和(c) A’区域的低能纠缠谱。两个图中都出现了反映连续对称性破缺的塔状结构(tower of states， TOS)。

文末，在写这篇小文时，看到了我老东家卡洛组里贴出了类似精神的文章，用约化密度矩阵计算了纠缠负度（entanglement negativity）等物理量^[11]。赞英雄所见略同之时，深感卡圈太卷！望各位读者朋友继续净化职场！

参考文献：

[1] Bin-Bin Mao， Yi-Ming Ding， and Zheng Yan. "Sampling reduced density matrix to extract fine levels of entanglement spectrum." arXiv preprint arXiv:2310.16709 (2023).

[2] Fakher F. Assaad， Thomas C. Lang， and Francesco Parisen Toldin， “Entanglement spectra of interacting fermions in quantum monte carlo simulations.” Phys. Rev. B 89， 125121 (2014).

[3] Zheng Yan and Zi Yang Meng， “Unlocking the general relationship between energy and entanglement spectra via the wormhole effect.” Nature Communications 14， 2360 (2023).

[4] Chia-Min Chung， Lars Bonnes， Pochung Chen， and Andreas M. Läuchli， “Entanglement spectroscopy using quantum monte carlo.” Phys. Rev. B 89， 195147 (2014).

[5] 严正，好蛋大爷的纠缠谱猜想——路径积分虫洞效应揭示纠缠谱与能谱的迷离关系 (qq.com)

[6] Jiarui Zhao， et al. "Measuring Rényi entanglement entropy with high efficiency and precision in quantum Monte Carlo simulations." npj Quantum Materials 7， 69 (2022).

[7] Jiarui Zhao， et al. "Scaling of entanglement entropy at deconfined quantum criticality." Phys. Rev. Lett. 128， 010601 (2022).

[8] Didier Poilblanc， “Entanglement spectra of quantum heisenberg ladders，” Phys. Rev. Lett. 105， 077202 (2010).

[9] Abolfazl Bayat， Henrik Johannesson， Sougato Bose， and Pasquale Sodano， “An order parameter for impurity systems at quantum criticality，” Nature communications 5， 3784 (2014).

[10] F. Kolley， S. Depenbrock， I. P. McCulloch， U. Schollwöck， and V. Alba， “Entanglement spectroscopy of su(2)-broken phases in two dimensions，” Phys. Rev. B 88， 144426 (2013).

[11] Ting-Tung Wang， et al. "Entanglement Microscopy: Tomography and Entanglement Measures via Quantum Monte Carlo." arXiv preprint arXiv:2402.14916 (2024).

本文受科普中国·星空计划项目扶持

出品：中国科协科普部

监制：中国科学技术出版社有限公司、北京中科星河文化传媒有限公司

本文来自微信公众号：返朴 （ID：fanpu2019），作者：严正（西湖大学物理系）