DP还能干这个？DP结合Quantum Thermal Bath描述SrTiO3中奇妙的晶格量子效应

作者 | 吴宏宇，何日，钟志诚，卢毅

近日，宁波材料所钟志诚研究员课题组通过深度势能模型结合Quantum Thermal Bath，提出一种作用于大尺度分子动力学模拟研究晶格量子效应的策略。相关工作发表于Phys. Rev. B 106, 224102，该工作以SrTiO₃（STO）为例证明了DP+QTB方法的有效性，该例子展示了由量子效应引起的几种现象，如被抑制的结构相变温度、低温下的量子顺电基态和介电常数的量子临界行为。

研究背景

模拟计算中结合晶格量子效应的重要性以及困难

尽管在很多问题中原子核可以近似地被看作具有确定位置和速度的经典粒子，但是考虑到其量子本质，在某些问题中需要考虑他们的量子效应。比如说，基于能量量子化假设和玻色-爱因斯坦统计，原子的集体激发模式可以被看作声子。基于声子的理论，可以算一系列热力学性质，其中最著名的例子是爱因斯坦固体比热理论。从图1中可以看到即使是室温300K下，STO的比热与基于能量均分定律的经典极限（虚线）也差之甚远。这表明不仅在极低温、质量轻的元素这两种情况下要考虑原子核的量子效应，而且在室温乃至更高的温度以及更多凝聚态物质中也需要考虑。

图1. STO的比热，点线为实际情况（量子效应），虚线为经典能量均分定律结果。https://doi.org/10.18280/mmc_b.870403

目前可以在计算中考虑原子核量子效应的主流模拟方法有：基于路径积分的一系列方法比如路径积分分子动力学（PIMD）【Z. Phys. B Condens. Matter 95, 143 (1994)】，还有随机自洽简谐近似方法（SSCHA）【Nature 578, 66 (2020)】。但是这两种方法都存在一个问题就是计算成本非常高，这会导致计算非常耗时且模拟原胞的空间尺度会非常小，模拟的时间尺度会非常短。对于做统计相关的问题，例如热力学性质、相图、复杂结构（缺陷、位错、表面）、非平衡态等等，由于尺度的限制这些问题的计算很难达到热力学平衡。因此，有必要发展一种高效的能够描述核量子效应的大尺度原子模拟方法。基于上述问题，本课题组提出一种Deep-Potential（DP）+Quantum Thermal Bath（QTB）【Phys. Rev. Lett. 103, 190601 (2009)】的全新策略，既可以实现晶格量子效应的修正，又能实现具有密度泛函理论（DFT）精度的大尺度原子模拟。

SrTiO₃中的晶格量子效应

需要一个典型的例子来验证DP+QTB这种策略是否有效，这里我们选取了SrTiO₃（STO）作为例子。在STO中有一个与晶格量子效应相关的现象，在低温下STO会进入一个具有极高的恒定介电常数的相而不会进入到铁电相，这种现象叫做量子顺电现象（如图2所示），且其介电常数在有限温度内和温度的关系是  （图2），这些现象都是经典理论无法解释的。【Phys. Rev. B 19, 3593 (1979)，Nat. Phys. 10, 367 (2014)】

图2. 实验测得的STO介电常数和温度的关系 左右分别引用自 Phys. Rev. B 19, 3593 (1979)，Nat. Phys. 10, 367 (2014)

可以从一个简单的图像定性理解这一现象。在钙钛矿材料中电偶极矩的产生源自于中心原子的位置偏移。当施加电场后中心Ti原子虽然会产生一定偏移，但是由于零点涨落的存在，Ti原子的位置存在不确定性。因此，Ti原子的实际时间平均位置会更靠近原本中心，也就是量子涨落会压制电场带来的原子偏移（图3a）。

STO高精度DP势函数

STO中的的四方相到立方相的结构相变（图2bc）的能量差非常微小，平均到每原子约1meV。如此小的能量差距超出了以往MD原子间相互作用势的精度可以描述的范围，因此对于STO中的MD模拟是急需一种新的方法获取高精度的势函数的。

在之前的工作中【Phys. Rev. B 105, 064104 (2022)】，我们使用同步学习代码DP-GEN 【Comput. Phys. Commun. 253, 107206 (2020)】来有效地获得相关的构型空间。同步学习以四方相和立方相作为初始结构，对温度在50-400K、压强在0.001-50kbar范围内的构型进行探索，最终获得约2600训练集构型。利用该训练集训练得到的DP模型，计算得到的能量与密度泛函理论（DFT）结果的误差可以达到meV/atom以内，表明其具有超高精度足以准确地描述STO中的相转变。

图3. (a)STO中的量子涨落现象，以及由(b)四方相到(c)立方相的结构相变

量子热浴QTB

为了实现对量子效应的描述，很自然地会联想到从薛定谔方程出发，严格地用量子力学去求解每个原子受力。但这显然是不太可能实现的，严格求解量子多体问题直到今天仍是一个艰巨的挑战。换一个角度思考，分子动力学（MD）的基础是经典的牛顿力学方程，有没有可能通过直接修改牛顿力学方程来实现对量子效应的描述呢？原子核的量子效应是原子核的集体激发行为，而牛顿力学方程是以质点为研究对象的，怎么描述集体行为呢？这些问题爱因斯坦和朗之万已经给出了一个解决思路，他们在研究布朗运动理论的时候发现颗粒们的无规则运动是有时间关联的，而这种关联和涨落耗散定理有关，颗粒的运动方程可以表达为：

这个方程叫做朗之万方程，其中是代表空间分量， 是代表第  个原子，是随机力，是阻尼系数，这后两者分别代表了介质对颗粒的涨落力和耗散力，通过它们来描述了颗粒集体行为对单个颗粒的影响。这两部分可以通过涨落耗散定理联系起来：

它表明涨落力乘积的平均值可以由耗散系数表示，以上是经典的布朗运动的特例也叫做白噪音。那么对应到凝聚态物质的量子情况又应该是如何呢？把原子核考虑成量子谐振子，对于谐振子而言最重要的性质是其振动的频率，因此应该在频率空间考虑涨落力和耗散力的关联。根据量子涨落耗散定理，此时随机力不再是白噪音而是一种“有色噪音”，它的功率频谱为：

其中

这个形式和玻色-爱因斯坦统计非常类似，但实际上他是服从量子统计的基本要求。特别需要强调的一点是第一项，它代表的是量子谐振子的零点能，这是经典理论中没有对应的一项。

所以简单来说，就是把描述质点运动的牛顿力学方程修改为可以描述集体行为的朗之万方程，而描述集体行为的是粒子受力的涨落与耗散之间的关联，这就是热浴。而恰好这种热浴能描述量子谐振子的行为，即量子热浴。【Phys. Rev. Lett. 103, 190601 (2009)】

研究结果

热力学性质的差异

本文先从热力学性质的结果出发，带大家直观地体验一下MD计算中经典的和考虑量子效应的区别。

图4. 经典和量子模拟的结果对比。（a）原子的平均能量 （b）Debye-Waller factor，也就是原子平均位移。

图4a中所指的经典理论即能量均分定律，也就是一条斜率为1的直线。而经典MD模拟出来的结果，恰好落在了这条线上。对能量做温度的一阶导，也就是比热，得出的结果是比热是不随温度变化的一个常数，这显然是不符合实际情况的。这表明在计算中忽略晶格的量子效应，那么连最基本的比热都无法正确求得，更不要说其他热力学性质了。

基于量子谐振子的理论，随温度变化的能量（图4a中的红线）为：

其中  是STO的声子态密度。可以从图4a看到QTB模拟的结果能够复现量子谐振子理论能量，而且有在高温时回归经典极限的趋势。

为了衡量不同材料中的量子效应强弱，我们定义了零点温度

它也代表着零点涨落的强度，只不过是以温度形式来表示的。对于STO，大概是230K，这都是由STO的零点能贡献的，它和经典下的0K有着如此大的差距，这表明量子效应的贡献是非常巨大的。同时我们也可以确定其他材料的来求得它们的量子效应的强度，具体的数值请看下图。这表明在很多材料中量子效应同样显著，是不应该被忽略的。

图5. 不同材料的

另外图4b中展示了STO中原子平均位移和温度的关系，可以看到经典模拟下原子平均位移几乎是随着温度线性增加，这是由热涨落所贡献的。而在量子效应被考虑的情况下，在零温下就有原子位置的零点涨落（如下面动图所展示），且其幅度相比于经典情形来说是非常大的，这也是非常巨大的差距。

图6. STO的原子位置零点涨落

既然量子效应带来的热力学性质差异这么大，那它一定也会影响到STO的结构相变，我们同时研究了STO中的四方相-立方相相变与量子效应的关系，如图7所示。图中灰色的虚线是用经典的MD模拟出来的，可以看到氧八面体旋转角在200K左右才消失，表明经典计算的相变温度也在这附近，但这与实验的105K是有很大差距的。引入量子效应后，模拟的相变温度来到了150K附近，相对经典的结果已有很大的改善。但是这和实验还是有很大的差距，下文会讨论这点。

图7. STO的晶格常数随温度的变化

介电常数的差异

STO的介电常数，在经典和考虑量子效应两种情况同样存在巨大的差距。

图8. 由经典和QTB两种方案模拟的STO介电常数随温度的变化

从图8中可以看出来，经典的介电常数和温度的关系是始终满足  的，而这种关系是经典理论所给出的结果，说明经典的MD模拟只能给出符合经典理论的结果。而结合QTB的结果在低温附近时不再满足  的关系，但是它的量级和实验还对不上。这和上文相变温度的计算有着相同的问题，那就是DP势函数在训练的时候所依赖的DFT势函数它不能正确的描述STO的晶格常数，它在描述STO晶格常数时明显时偏小的，这等效于给STO在模拟的时候加了一个向内的净水压强（如图9所示）。所以，在模拟的时候我们应该施加一个等效的负压来平衡这种误差。

图9. 位移铁电体的理论相图

原本STO应该是位于图9中位移铁电体的临界压强  上，我们不断尝试施加负压找到了真正的  ，其模拟结果如图10所示。我们能看到和实验结果（图2a）一致的数量级，并且也能清晰地看到靠近高温时经典的  现象，以及零温时介电常数向平台转化的现象。更重要的是插图中所展示的在有限温下介电常数的  现象，在这里通过QTB的模拟很好地复现了出来。

图10. 靠近临界压强附近的STO介电常数

总结和展望

我们提出了一种新的基于第一性原理的策略DP+QTB，它可以以DFT的精度有效处理大尺度模拟问题，并且解决经典MD中无法处理的晶格量子效应问题。DP+QTB的优势之处在于高度的泛用性和高效性，理论上来说任何材料只要能通过DP训练出对应的模型，都可以高效地研究其中的晶格量子效应问题。这对我们解决更多领域内的难题提供了一个崭新的思路，即对于微观机理尚不清楚，经典理论无法解决的问题，我们可以通过做更大规模的模拟并考虑晶格量子效应，从微观出发探索凝聚态的宏观性质，正所谓“见微知著”。

本工作由中科院宁波材料所吴宏宇博士生、何日助理研究员、钟志诚研究员、南京大学卢毅教授共同合作完成。

参考文献

https://doi.org/10.1103/PhysRevB.106.224102

推荐阅读