在纠缠中窥见自然的奥秘
1. 引子
“身无彩凤双飞翼，心有灵犀一点通”。
唐代诗人李商隐的这首诗，道尽了人世间痴情男女虽天各一方、但彼此心意相通、似有一条丝线相连的缠绵。
其实，电子的世界亦如此。配对的电子，一旦形成一个纠缠态，无论相距多远，它们都能感受到彼此的“心意” (自旋)。它们一个向上，一个向下，仿佛约好了似的。只要有一方改变状态，那么另一方就会立刻感受到，并且调整成与对方相反的状态。
以上的爱情故事，用数学的语言来写就是这样：
其中，|ψ> 这个符号指代两个电子所处的状态。这是一种约定俗成的简便写法，物理学家习惯称它为“态”或者“波函数”。上面公式 (1) 中的这个波函数的意思是：两个电子有 1/2 的可能性处于“一上一下”的状态，还有 1/2 的可能性处于“一下一上”的状态。|ψ> 和距离以及时间无关，也就是说两个电子能随时随地知道到对方的“心意”。这是一种紧紧相连的纠缠，如此大概是所有异地爱侣们心向往之的美满爱情吧。
可是，自然界中电子的遭遇和人世间饮食男女的故事其实也并无两样。那种美满的、只有两两配对的状态是极少见的。实际上，两个电子这样配对的状态并不稳定。几乎只要很小的外界干扰，就能摧毁这种状态，甚至连观测本身也会导致这种状态发生可怕的坍缩。难以想象，自然界怎么会有两个电子从产生到湮灭，只和对方有纠缠，而不会改变心意呢？！
如此说来，电子的世界和人的世界其实没什么两样。在笔者研究的多体物理范畴里，多电子之间的纠缠和相互作用，构成了物质世界千奇百怪的物态。比如导体、半导体、绝缘体，仅仅因为固体中电子所处的环境以及电子间的相互作用形式不同，就在宏观上造成不同的能带结构。而不同的能带结构又决定了导电能力的强弱。如图 1 所示，即表达了不同能带结构所对应的完全不同的导电能力及物质状态。
图 1. 不同导体的能带示意图。
这正如人与人不同的交往相处方式构成了不同的文化，不同的文化又组成了不同的国家、种族与文明。关乎到个人生活的种种小事的差异，其实在宏观上最终导致了截然不同的文化与文明。它们共同组成了这个纷繁复杂而又充满魅力的世界。
物质的态，当然不止图1所示的那些。现已口耳相传的超导体、石墨烯等材料，不能用上述那么简单的能带结构图来解释。量子材料的物态，要更加复杂和丰富，所对应的多电子之间的纠缠与关联方式也千奇百怪。如何去衡量纠缠强弱、如何理解电子之间纠缠强弱的程度与不同物质状态之间的关系，现在已经成为物理人极其关注的一个问题。
2. 纠缠熵
在谈论多电子纠缠前，我们暂时先回到配对电子的问题。公式 (1) 只是纠缠态的一种，两个电子还可以有不同的纠缠态，比如公式 (2) 所表达的也是一种纠缠方式。正如感情的浓烈程度会有差异，不同纠缠态的纠缠强度也不同。于是，物理人用纠缠熵来刻画纠缠的强度，其中一种常用的纠缠熵是冯 • 诺伊曼熵，其定义如公式 (3) 所示：
冯 • 诺伊曼熵刻画了体系包含的子系统 A 与系统其它部分纠缠的强弱程度。其中 ρ 是系统的密度矩阵，而约化密度矩阵 ρA 表示把密度矩阵中有关 B 的自由度积掉后剩余的部分。
图 2. 一个有子系统 A 和 B 的多体系统示意图。
根据这一定义，|ψ>1 中两个电子之间的纠缠熵为 0.69, 而 |ψ>2 中两个电子之间的纠缠熵为 0.56。所以，我们可以简单理解为：|ψ>1 这小两口相比于 |ψ>2 要更“恩爱”一些。
而我们关心的量子材料问题，需要考虑的往往是很多电子 (实际上，这个数目大到惊人)的相互作用和纠缠问题。此时，公式 (3) 代表的其实是子系统 A 作为一个整体和系统中其它部分的纠缠熵，如图 1 所示。如果前面考虑的是个人间纠缠的问题，那么对于量子材料的研究来说，要考虑是类似于河南省与湖北省的关联与纠葛的故事了：难度和复杂性可想而知、可想而不知！
对于公式 (3) 所定义的纠缠熵，其实已经有了很多有趣的研究结果。最有名的一个结果是 A、B 之间纠缠熵是与 A – B 之间的边界长度呈正比的！当然，这个结果仅限于有限关联长度的系统 [1, 2, 3, 4]。其实，这也很容易理解：想象一下在几百年甚至上千年前那个交通十分落后、信息闭塞的年代，两个地区之间的关联，其实仅限于边境人民的交流。边境线越长，交流的机会也就越多。所以，两个地区之间纠缠的强度，就和边境的长度成正比。不过，这个结果看起来很显然，但证明起来一点也不轻松。从公式 (3) 出发，要在场论中进行一通操作。在普通人的眼里，这个证明就好像变魔术一样证出来了。不过，作为普通人的一员，我们现在只满足于知道这个结论就好了。
事实上，从纠缠熵中可得到的信息不止于此。实际的纠缠熵与体系、时空维度、纠缠边界的形状都有关。在纠缠熵和边界长度简单成正比的规律之外，还有一些修正项。自然可以想到，这个修正项刻画的是体系里长程的纠缠行为 (因为短程的纠缠已经包含在和长度成正比的项里了)。物理人发现，这些长程的行为，有时候更加本质一些。它像是量子材料的“护照”，可以区分它们的“国籍”。比如，我们考虑纠缠熵 SA 和边界长度 l 的关系，如公式 (4) 所示。这里，ln(l) 项可以来自于连续自发性对称破缺相的 Goldenstone mode 数目，也可以来自于不光滑纠缠边界 (也就是说有 A 和 B 的边界是折线、而不是光滑曲线) 的贡献 [4, 5, 6]。
接下来，我们会介绍几种相和相变点。我们会发现，公式 (4) 中参数 s 和 c 的取值，确实可以作为不同类别相的“护照”。
图 3. SU(2) 对称性破缺的反铁磁海森堡模型的纠缠熵。
3. 对称性自发破缺的奈尔相、O(3) 量子临界点、拓扑序量子自旋液体
首先，我们不妨考虑一些具体的物理模型。最简单的莫过于二维反铁磁海森堡模型了。这个模型，在零温时 SU(2) 对称性破缺，它的 Goldstone mode 数目是 2。可以证明此时公式 (4) 中 s = -Ng/2，Ng 是 Goldstone mode 数目。并且，这个关系对于所有自发连续性对称破缺相都成立。所以，通过 s 的取值就可以判断不同相的 Goldstone mode 数目 (也就是起到了区分这些相的目的)。
用改进的非平衡累加算法(后文会详细描述) 来计算图 3中插图里的 A 区域和其它部分的纠缠熵。从图中可以看到，A 和它的补集的边界是光滑的，所以 ln(l) 项的系数唯一的贡献来自 Goldstone mode。用计算的结果来拟合公式 (4)，即得到SA(2)= 0.092(1)l + 1.00(9)ln(l) - 1.63(3), 即 s = -1.00(9)。这与理论值符合得很好，说明了本算法对于探测 Goldstone mode 数目 (连续自发性对称破缺相的特征) 是十分有效的。
对于量子相变点，比如从 Néel state 到 singlet state 的 O(3) 相变点，在边界是折线的情况下，用 s 的取值可以区分不同类别的相变。比如 Ising 相变点的 s 值就与 O(2)、O(3) 相变点的 s 值不同，它们大致有个 1:2:3 的关系 [5]。通过 s 取值判断相变类型的这一性质对研究多体问题来说是十分有用的工具。如果说回 O(3) 相变点，因为此时没有 Goldstone mode (没有连续自发对称性破缺)，ln(l) 项的系数 s 由折线上的 corner 来贡献。如图 4 所示，选择 A 区域为一个正方形，有四个 90 度角，理论计算得知它们会对系统贡献大概 - 0.08 ln(l) 左右的纠缠熵。通过非平衡累加算法，我们测量了不同取值下的纠缠熵，并且根据公式 (4) 来拟合，得到的结果为 SA(2) = 0.167(2)l - 0.081(4) ln(l) - 0.124(7)，即 s = 0.081(4)。s 的值再一次和理论值很好吻合，再一次说明我们的算法对于区分不同的相变点是十分有效且精确的。
图 4. O(3) 相变点的纠缠熵。
上述两类相和相变点的区分都是以 s 的数值来判定的，而 c 在上述的例子里都是很平凡的常数。而对拓扑序相，就该轮到 c 来大显身手了。我们挑选了比较受关注的自旋液体相来测试我们的算法。此时常数 c 等于拓扑纠缠熵的取值，它是普适的，可以用来区分不同类别的拓扑相。如果 A 和 B 之间的边界是光滑的，那么在自旋液体相，s 的取值为0。而理论预言在周期性边界条件下 c = 2lnD，其中 D 是体系总的量子维度。在 Z2 自旋液体相中，D= 2，所以在公式 (4) 右边的 al 项以外，系统的纠缠熵还有一个修正项 -2ln2。
我们这里考虑一个典型的 Z2 自旋液体——二维笼目晶格上的 Balents – Fisher – Girvin (BFG) 模型 [7，8]。因为这个模型十分复杂，模拟难度较高，修正项 2ln2 这个数值以前还从未被真正观测到过。基于我们发展的非平衡累加算法，因为计算的并行度得到极大提高，纠缠熵测量的效率和精度也大幅度提高了，所以精确测定拓扑纠缠熵再也不是一个遥不可及的目标。
图 5. 两种不同长宽比下的笼目晶格以及对应的 BFG 模型的纠缠熵的测量值。
图 6. 配分函数 ZA(2) 的时空构型。
如图 5 所示，我们在长宽比分别为 1:1 和 2:1 的两种晶格上分别测量了不同边界长度下纠缠熵的取值，并且按照 SA = al - c 拟合，得到的结果为c = 1.4(2)，和理论值 c = 2ln2 ~ 1.386 吻合得很好。我们的结果，不仅从数值计算的角度验证了自旋液体理论，而且说明了非平衡累加算法的有效性。我们相信，我们发展的非平衡累加算法及结果，对其它更复杂模型纠缠熵的数值计算 (即计算测量)，如费米子体系纠缠熵的测量，将有很大的示意作用。
4. 算法
最后，我们来具体描述我们发展的非平衡累加算法。
公式 (3) 给出了冯 • 诺伊曼熵的定义。但是，实际算法用的是另一种广义的 n 阶 Rényi 熵，其定义如公式 (5) 所示。
在 n → 1 的时候，Rényi 熵会收敛到冯 • 诺伊曼熵。Rényi 熵有一个优点，就是它可以写成两个配分函数的商的形式，便于我们用蒙特卡洛方法进行计算，如公式 (6) 所示，其中 Z(n) 是 n 个独立的配分函数 Z 的集合，而 ZA(n) 的构型是一个 n 片状的黎曼面 (n - sheeted Riemann surface)。当 n = 2 的时候，ZA(2) 在蒙卡中的时空构型如图 6 所示，这是一个类似于裤子的几何结构 (实际上我们考虑的是有空间周期性的情况，所以这个构型应该是把两个裤腿塞进裤裆里的感觉)。当 A 是一个空集的时候，ZA(2) 会回到 Z(2)。
我们知道，- lnZ 在统计物理里等于自由能 F 的统计平均。所以，如果可以设计一个从 Z(2) 到 ZA(2) 的演化过程，那么纠缠熵的值就可以转化为在这个过程中总的做功了，即公式 (7) 所定义，其中 WA(n) 就是我们要测量的总的做功的值 [9]。注意，当这个演化过程是非平衡时，公式 (7) 的证明要用到 Jarzynski 等式 [10]，这里就不详细讲证明过程了。
图 7. 非平衡累加算法示意图。
那么，如何设计这样一个从 Z(2) 到 ZA(2) 的演化过程呢？
注意到，前文说到 Z(2) 其实就是 ZA(2) 中 A 是空集的一个特殊情形。所以，只要缓缓地改变子集 A 中格点连接的方式 (也就是说像一个裁缝一样缓慢的调节裤裆和裤腿的长度比)，就能够把 Z(2) (只有两条裤腿的几何构型) 缓慢的变成图 6 所示的一个尺寸合适穿着舒服的裤子了！好吧，还是回归数学描述，这样的定义可见公式 (8)，其中，B 是 A 的子集，而 NA 和 NB 分别代表 A 和 B 中格点的数目。不难发现，ZA(2)(λ = 0) = Z(2)、ZA(2)(λ = 0) = ZA(2)。利用这个以 λ 为参数的配分函数，缓慢地把 λ 的取值从 0 调节到 1，并且测量这个过程中所做的功，就可以根据公式 (7) 测量出纠缠熵了。
我们发展的非平衡累加算法，主要的新意是算法将这个缓慢的演化过程并行化了。如图 7 所示，该算法将 0 到 1 的过程分成许许多多的小段，同时分别计算每一小段演化过程中所作的功，最后将这些功加起来，就可以得到整个过程的总功了。这个并行化的过程，显著提高了算法的效率。如图 8 所示，并行度 N = 24 时，我们的算法随演化时间收敛到正确值的速度远远大于不并行 (也就是 N = 1) 的情形。
图 8. 并行和不并行的算法效率对比。
5. 小结与后记
行文至此，可以对这一工作进行几句小结。本工作中，我们发展了一种新的测量纠缠熵的蒙特卡洛算法，提升了纠缠熵测量的效率和精度。这种效率的提升，使得蒙特卡洛模拟中纠缠熵的测量变得更加简单。对一些很复杂的模型，比如前文提到的 BFG 模型，其纠缠熵精确测量也变成了现实。纠缠，真的变成了多体世界一个可以量化的物理量和工具，可以用来分辨不同相和相变。由此，我们也加深了对多体世界关联和纠缠的认识！
这一工作刚刚发表在《npj QM》上，信息如下：
Measuring Rényi entanglement entropy with high efficiency and precision in quantum Monte Carlo simulations
Jiarui Zhao, Bin-Bin Chen, Yan-Cheng Wang, Zheng Yan, Meng Cheng & Zi Yang Meng
npj Quantum Materials volume 7, Article number: 69 (2022)
https://www.nature.com/articles/s41535-022-00476-0
小结之外，笔者也不免有些感慨。我们感悟到，科研之难，难于上青天。在这个要求研究必有新进展、发现必有突破性的时代，物理人凭借兴趣本身去探索自然的难度很大。在这样的苦闷里，有一些人放下了执念，将自己的全部身心投身于奉献里，终于求得了内心的安宁。笔者有幸，在与身边诸君的交流与相处中，也获得了一丝奋发向上之气。可以说，全世界科学家的心都是紧紧纠缠在一起的，只要对科学有兴趣，总会在某个时刻突然感受到那种“身无彩凤双飞翼，心有灵犀一点通”的感动。
最后，以笔者导师卡洛君答谢 《npj QM》 主编 Ising 君的一篇短文结尾，想必能读到这里的朋友里定会有和他心有灵犀之人吧！
********
Ising 先生足下。N 日不见，如隔 M 秋。确数不详，聊寄敬意。港大赵学士家瑞君，严博士正君，陈博士斌斌君，承贵刊厚爱，有秋裤纠缠文章发表，属予做一 teaser 以释之。余嘉其能不拘于时，以钻研为乐，寄人生于格物，戏为一篇以贻之。
辞曰：同心倩女，至离枕上之魂；千里良朋，犹识梦中之路。而况量子纠缠，执多体系统之精髓；基态简并，通拓扑物态之本原也哉！指数高墙，纠缠难做。前期算法，多系花拳。确数落落，唯靠 hand - waving。叹普适拓扑纠缠之难得，来高维量子多体之艰深。不意三位以秋裤配分，携 Jarzynski 非平衡做功定理，极量子蒙卡重要性抽样之算力，一击中地，攻克高维多体纠缠熵之难题。信可乐也。
噫！花面逢迎，世情如鬼。内卷之风，南北一辙。古今痛哭之人，卞和唯尔；颠倒逸群之物，伯乐伊谁？余窃喜诸君子尚能寄情格物，跌宕自喜，不随世风逐流，不为外物所动。人生世上，只须合眼放步，余唯愿与诸君子齐听造物之低昂而已矣。
参考文献：
[1]. Srednicki M. Entropy and area [J]. Physical Review Letters, 1993, 71(5): 666.
[2]. Holzhey C, Larsen F, Wilczek F. Geometric and renormalized entropy in conformal field theory [J]. Nuclear Physics b, 1994, 424(3): 443-467.
[3]. Bombelli L, Koul R K, Lee J, et al. Quantum source of entropy for black holes [J]. Physical Review D, 1986, 34(2): 373.
[4]. Calabrese P, Cardy J. Entanglement entropy and quantum field theory [J]. Journal of Statistical Mechanics: Theory and Experiment, 2004, 2004(06): P06002.
[5]. Kallin A B, Stoudenmire E M, Fendley P, et al. Corner contribution to the entanglement entropy of an O (3) quantum critical point in 2 + 1 dimensions [J]. Journal of Statistical Mechanics: Theory and Experiment, 2014, 2014(6): P06009.
[6]. Laflorencie N. Quantum entanglement in condensed matter systems [J]. Physics Reports, 2016, 646: 1-59.
[7]. Balents L, Fisher M P A, Girvin S M. Fractionalization in an easy-axis Kagome antiferromagnet [J]. Physical Review B, 2002, 65(22): 224412.
[8]. Wang Y C, Cheng M, Witczak-Krempa W, et al. Fractionalized conductivity and emergent self-duality near topological phase transitions [J]. Nature Communications, 2021, 12(1): 1-7.
[9]. D’Emidio J. Entanglement entropy from nonequilibrium work [J]. Physical Review Letters, 2020, 124(11): 110602.
[10]. Jarzynski C. Nonequilibrium equality for free energy differences [J]. Physical Review Letters, 1997, 78(14): 2690.
备注：
(1) 本文作者赵家瑞，现系香港大学物理系博士生三年级。
(2) 封面图片即文中的图 7，系非平衡算法示意图。
精选文章
欢迎订阅 npj Quantum Materials 的 Email Alerts
点击 https://idp.nature.com/register/natureuser 完成免费订阅