凝聚相系统与外部电场的相互作用在自然和技术的无数过程中至关重要,包括细胞的场定向运动(趋电性),到地球化学和行星上冰相的形成,到场定向化学催化和能量存储和转换系统,包括超级电容器、电池和太阳能电池。电场存在下的分子模拟将为这些过程提供重要的原子洞察力,但最准确的方法(例如从头算分子动力学(AIMD))的应用范围因其计算费用而受到限制。在这里,我们引入扰动神经网络势分子动力学(PNNP MD)来推迟此类模拟的可访问时间和长度尺度。我们证明了液态水的重要介电特性,包括场引起的弛豫动力学、介电常数和场相关红外光谱,可以通过机器学习达到约 0.2 V 1 的令人惊讶的高场强,与 ab 相比,精度不会损失。-从头分子动力学。这是值得注意的,因为与大多数以前的方法相比,PNNP MD 所基于的两个神经网络专门针对从零场 MD 模拟中采样的分子配置进行训练,这表明这些网络不仅可以插值,而且可以可靠地推断场响应。PNNP MD 基于严格的理论,但它简单、通用、模块化且可系统改进,使我们能够获得对各种凝聚相系统与外部电场相互作用的原子洞察。
电场在自然和技术中无所不在。它们的存在引导大黄蜂寻找花蜜1,并且据信它们在金星2上引起了超离子冰VII相。它们在众多电子和能量转换设备中发挥着核心作用,包括场效应晶体管、(超级)电容器、电池和太阳能电池。在化学中,电场可用于控制催化中的选择性3和控制反应性4,而在物理学中,电场可用于将粒子加速到接近光速。这些例子中的场强跨越了非常大的范围。晴朗天气下的大气场约为 106 V 15,充当花卉提示的电场可高达 105 V 11。在电气元件中,场强可能会根据实际设计和设想的应用而发生显着变化。约 1010 至 103 V 1. 粒子加速器通常在高达 5 103 V 的场强下运行 1. 在带电电极处,可以发现场强约为 0.1 V 1 6,这也标志着化学键激活的开始7和水的电冷冻8,9。化学键的形成和断裂可能需要超过 1 V 110,11。
分子动力学 (MD) 或蒙特卡罗模拟是理解和预测液体、电解质和液体相互作用的最合适方法/固体与外部电场的界面,因为它们对分子构型的正确平衡分布进行采样。许多有限电场模拟12都是用经典力场13、14、15进行的,这些模拟通常可以很好地描述纯溶剂,但在电荷转移和极化效应变得重要的更复杂系统(例如电极/电解质)上往往表现不佳接口。从头开始 MD (AIMD) 模拟16从每个 MD 时间步长的第一原理(通常在密度泛函理论 (DFT) 的理论水平)求解系统的电子结构,从而为此类情况提供最准确的描述。事实上,有限电场的 AIMD 模拟已经在许多系统17上进行,范围从晶体12和纯液态水18,19到电极/电解质界面13,20,21,22。AIMD 模拟的一个严重缺点是它们对计算的要求仍然很高,模拟时间通常限制在几十皮秒,系统大小限制在几百个原子,具体取决于所选的密度泛函。因此,对重要的场诱发现象(例如介电弛豫、离子电导率15、双电层形成23或电容充电)进行 AIMD 模拟的成本高昂,因为需要大量的原子来忠实地模拟此类过程。
<机器学习 (ML) 的出现改变了 MD 模拟领域。现在,通过仅从几百个显式电子结构计算中训练 ML 势,就可以进行纳秒 MLMD 模拟,与 AIMD 相比,精度几乎没有损失24,25,26。最初,这些机器学习模型旨在仅计算势能和力,包括将系统总能量分解为内部和环境贡献的方案27,28。在过去的几年里,已经发布了一系列机器学习模型来预测偶极矩和极化率29,30,31,32,33,34,35,36,以及其他响应属性,包括原子极张量 (APT)37。随着这些发展,最近引入了一些机器学习模型,这些模型明确地描述了分子系统与外部电场的相互作用,并成功应用于此类系统和液体38,39,40,41,42。在大多数这些方法中,势能表面的场依赖性明确或隐含地是 ML 势的一部分38,39,40,42。因此,电场是 ML 势的输入参数,因此需要不同场强的训练数据,这在 AIMD 级别的计算成本很高。在这里,我们介绍了一种简单而鲁棒的方法电场 MLMD 的替代方案,专门从零场分子配置中学习场响应,即,它不需要通过使用电场运行 AIMD 模拟来生成训练配置。我们从标准势能面开始,通过一阶截断的级数展开以微扰方式解释与电场的相互作用,注意到一阶力项可以用 APT37,43 来写。训练两种 ML 模型:一种用于未受扰动势能面的标准 ML 势(此处为委员会神经网络势 (c-NNP)24,26,44),一种用于 APT 的 ML 模型(此处为等变图神经网络)表示为 APTNN37)。它们结合起来形成了扰动的 ML 势(此处为扰动的神经网络势 (PNNP))。这两个网络不了解电场效应(无论是明确的还是隐含的),因为它们都是专门针对零场采样的配置进行训练的,并且电场本身不是机器学习模型的输入。因此,MD 模拟过程中与场的相互作用完全是由图神经网络表示的 APT 与外部场的乘积给出的级数展开中的一阶项造成的。据我们所知,使用 APT 及其 ML 表示来计算原子上的场致力以前从未被探索过,这代表了本文的主要概念进展。
PNNP 基于基于严格的物理原理,并精确遵循量子化学计算中电场的处理,其中总能量在零场处扩展45。因此,可以通过在 ML 表示中添加高阶贡献(例如极化率和超极化率)来系统地提高准确性。由于我们没有引入或修改 ML 模型,因此继承了两种采用的 ML 模型的优点和准确性,但也继承了其局限性。重要的是,PNNP 遵循现代极化理论的精神46,因为 APT 与极化的变化相关,因此不受周期性边界条件下极化的多值性的影响。这使得 PNNP 适用于通常在周期性边界条件下建模的广泛凝聚相系统,包括固体、液体和界面。
在下一节中,我们将详细描述 PNNP 方法。经过验证,我们应用该方法来模拟纯液态水的介电响应。我们展示了随着场强的升高和降低,液态水的可逆极化和去极化。然后详细分析介电弛豫动力学,然后根据极化相对于场强的响应计算介电常数,结果与实验值非常吻合。此外,我们表明 PNNP 正确预测了 O-H 拉伸模式中场引起的红移和液态水解放模式中场引起的蓝移,与 AIMD 的结果非常一致。在讨论部分,我们将 PNNP 与之前介绍的用于计算电场分子系统的 ML 方法进行比较,并讨论当前的局限性以及克服这些局限性的方法。
在我们的方法中,原子系统与均匀外部电场 E 的相互作用通过场中一阶截断的级数展开进行微扰处理13,15,20,45
其中 \({{{{\mathcal{H}}}}}_{0}({{{{\bf{r}}}}}^{N},\,{{{{\bf{p}}}}}^{N})\) 是总的未扰动哈密顿量,由动量为 pN 的 N 个原子核的动能和取决于所有核位置 rN 的电子势能组成,\({{{{\mathcal{H}}}}}_{0}({{{{\bf{r}}}}}^{N},\,{{{{\bf{p}}}}}^{N})\,=\,{E}_{{{{\rm{kin}}}}}({{{{\bf{p}}}}}^{N})+{E}_{{{{\rm{pot}}}}}({{{{\bf{r}}}}}^{N})\),E M(rN) 是电场 E 作用引起的扰动零场时系统的总偶极矩 M(rN)。对于本工作中研究的弱和中强电场,场中的一阶截断预计是准确的,这将在下面进一步证明。偶极矩的场依赖性或等效的高阶项(极化率、超极化率)可以在更强的场处添加,以解释电子结构的场依赖性扰动47。应用哈密顿运动方程,我们得到作用在原子 i 上的力
其中 = x, y, z 和 = x, y, z 代表三个笛卡尔坐标,E 是 E 的 -分量。第一项是在没有电场的情况下作用在原子核上的力,第二项是场引起的贡献,可以用原子 i 的 APT 转置来写,\({{{{\mathcal{P}}}}}_{i}\)37,43(不要与极化 P,方程(6)混淆),带有元素
在我们的方法中,我们训练两个 ML模型,一个用于势能 (Epot(rN)),一个用于 APT (\({{{{\mathcal{P}}}}}_{i}\)) 并使用相应的力,方程 1(2)、在存在外部电场的情况下进行MD模拟。我们强调,这两个量都是在没有外部场的情况下进行训练的,这与之前建议的方案38、39、40、42不同。我们使用第二代高维神经网络势 24,26 (c-NNP) 的委员会 44 来对 Epot(rN) 进行建模,并使用 E(3) 等变图神经网络来对 APT (APTNN) 进行建模,正如最近由一位我们37。有关强制实现的详细信息,请参阅方法部分。下面,电子势能的 c-NNP 和 APTNN 组合模型,包括场项 Epot(rN) E M(rN),简称为扰动神经网络势 (PNNP)。
PNNP MD 模拟可以了解固体、液体和离子溶液的许多重要介电特性。总偶极矩的时间导数可以通过将所有 APT 乘以原子核各自的速度 vi37 求和来获得,因此,可以从 \(\dot{{{{\bf{M}}}}}\) 沿 PNNP 轨迹采样(参见下面的等式(9))。时间积分给出了电池的总偶极矩,
和极化 P,
其中 V 是模拟盒的体积。顺便提及,我们参考了现代固体极化理论47,其中 \(\dot{{{{\bf{M}}}}}(t)/V\) 是瞬态电流密度,当随时间积分时,给出流动偶极矩。方程式中的极化。(6) 非常重要,可以计算相关介电特性,包括介电常数(或相对介电常数,参见下面的方程(7))、电容和离子电导率。
我们通过在室温下模拟纯液态水(约 0.002 至 0.2 V 1 的不同电场强度)来演示我们的方法。有关液态水的 c-NNP 和 APTNN 训练的详细信息,请参见方法部分。从质心动量为零的平衡水构型初始化轨迹,并在 0.0129 V 1 的中间场强下运行 PNNP MD,我们获得了 2.2 109 Hartreeatom1 ps1 守恒总能量的漂移以及守恒总能量的平均幅度质心动量为 8.4×109 au,分别参见补充图 S1 A 和 B。在所有其他施加的电场强度下进行模拟时都会获得类似的值,并且它们也是无场模拟的典型值,证明了正确的力实施。
接下来,我们通过验证来评估 PNNP 的质量针对未见过的测试集预测的沿轨迹的力。在每个场强下,从 100 ps PNNP 轨迹等距提取 101 个配置。然后,使用相同场强下的 DFT 参考计算来计算每个原子上的总力,并与 PNNP 的总力进行比较(图 1A 中显示的数据为 0.0129 V 1 的示例性场强)。然后,我们通过在相同配置但没有外部电场的情况下执行额外的 DFT 计算,将 DFT 参考力分解为未受扰动和电场引起的贡献。有场和无场的力之间的差异隔离了场引起的力。这使我们能够单独验证未受扰动的 c-NNP 力贡献,即方程右侧的第一项。(2),具有未受扰动的 DFT 贡献(图 1B)和电场引起的 APTNN 贡献,方程右侧的第二项。(2),具有场诱导的 DFT 贡献(图 1C)。显然,对于总力以及未扰动和场诱导的贡献,我们获得了 PNNP 和 DFT 之间的强相关性。
原子力中的均方根误差 (RMSE) 作为函数场强如图1D所示。总力 RMSE 约为 90 meV 1,因此与先前发布的各种系统上的高维 NNP 和 c-NNP 一致44,48,49,50,但比最近报告的值稍大(参见讨论)。值得注意的是,我们发现即使在场强高达 0.2057 V 1 的情况下,未扰动贡献的 RMSE 也几乎保持恒定在约 90 meV 1 ,尽管训练集仅由未扰动的液态水(即零场平衡配置)组成。场感应力的贡献比未受扰动的贡献小大约两个数量级,它们的 RMSE 也是如此。值得注意的是,场感生力 RMSE 与场强成比例地增加,从 0.0026 V 1 时的 8.6×102 meV 1 增加到最大模拟场 0.2057 V 1 时的 9.05 meV 1。将 RMSE 除以场感生力的范围一在 0.0129 V 1 处获得 1.1% 的误差,在 0.2057 V 1 处获得 1.2% 的误差。或者,将 RMSE 除以参考 DFT 力的均方根,如文献中常见的做法50(另请参阅补充说明 4),我们获得了在场强约为 0.02 V 1 时场感贡献为 8.2% 的相对 RMSE,在 0.2057 V 1 场强时场感贡献为 9.8%(图 1E)。因此,两个误差指标上的相对力误差在所研究的场强范围内表现出显着的恒定性。
最后,我们验证了通过积分获得的沿 PNNP MD 轨迹的总偶极矩的计算。根据方程计算偶极矩的时间导数。(5)。我们发现偶极矩可以很好地跟踪从显式 DFT 计算中获得的参考偶极矩,持续时间约为 10100 ps,具体取决于场强(参见补充图 S1C)。然而,在较长的时间内,即使准确地再现了 \(\dot{{{{\bf{M}}}}}(t)\),由于在有限时间步长上积分时误差的累积,偏差也会变得更大。这里通过以周期性间隔计算沿 PNNP MD 轨迹的参考 DFT 偶极矩并将从 APTNN 获得的偶极矩的时间导数仅从一个 DFT 参考值积分到下一个来解决此问题。使用此积分过程,我们计算了 100 ps PNNP MD 轨迹上的平均偶极矩 M(t),作为两个 DFT 参考值之间的时间间隔的函数。图 1F 显示了在相同的 100 ps PNNP MD 轨迹(使用 1 ps1 的采样频率计算)上平均的参考 DFT 平均偶极矩的相对误差。我们发现,对于所有场强,对于 10 ps 的两个 DFT 参考偶极矩之间的间距,误差迅速减小到低于 5%。该积分程序的间距为 10 ps,用于本工作中报告的总偶极矩。因此,在实践中,只需要极少量的额外 DFT 计算即可以仅受 MD 时间步长限制的分辨率精确计算沿 PNNP MD 轨迹的偶极矩。
根据 DFT 参考数据验证了力和极化后,我们现在应用 PNNP 来模拟水取向因与电场相互作用而松弛。我们从使用无场 c-NNP MD 平衡 1 ns 的液态水样本中提取了 20 个独立配置,并将它们用作 40 ps PNNP MD 模拟的起始配置,每个配置均在 0.0257 V 1 的场强下运行。我们采用相同的初始配置,并在相同的场强和积分时间步长 (1 fs) 下执行显式 AIMD 模拟。给定时间 t 下水分子的方向由所施加电场的方向与所有水分子平均的两个分子内 OH 键向量 (t) 之间的平分线之间的角度来描述。结果如图 2 所示。在 t = 0 时,初始平均取向非常接近 90,对应于平衡时随机取向水分子的期望值。当场在 t = 0 时打开时,从 PNNP MD 获得的弛豫动力学与 AIMD 的结果非常一致。指数拟合得出 PNNP 的时间常数 = 5.9 ps1 (R2 = 0.99),而 AIMD 的时间常数 = 6.6 ps1 (R2 = 0.99)。此外,PNNP 和 AIMD 模拟分别收敛到几乎相同的最终平均角度:60 度和 61 度。由于两种方法的标准偏差重叠(图 2 中的阴影区域),我们将剩余的小差异归因于统计不确定性。需要大量总计几纳秒的轨迹来澄清这一点,但由于 AIMD 的计算成本,这是不可行的。
在下面我们及时改变场强并监测从 PNNP 模拟获得的水样的结构和介电响应。结果如图 3 所示。图 A 中显示了施加电场强度的时间分布,图 B 中显示了平均水取向 (t) 以及施加电场方向上的总偶极矩 Mz(t),在面板 C 中。当场从 0 增加到 0.0154 V 1(以 0.0026 V 1 的增量)时,我们在平均方向上获得了近似线性的响应,并且偶极矩随之线性增加。这与先前发布的线性范围上限估计值 0.03 至 0.07 V 1 17,51 一致,但略低。在较大的场强下,介电响应变弱,表明达到了非线性状态。在 0.2057 V 1 处,观察到沿着场方向的强定向排列,= 20。请注意,在这些高场强下,PNNP 仍然预测相当准确的力(见图 1D、E),尽管没有电子极化项和任何与现场相关的训练数据。在 0.4114 V 1 处,我们观察到水分裂成质子和氢氧根离子,这与之前关于在这些场强下发生化学键激活的报道一致7。在这里,人们可以进一步训练 c-NNP 以包含质子和氢氧根物质,从而测试更大的场强。然而,对这些物种的准确描述需要明确包含核量子效应52,53,这超出了本工作的范围。相反,运行了一组额外的模拟,其中场强从线性状态结束时的值 0.0154 V 1 反转到 0,增量为 0.0026 V 1。补充图 S2 中详细介绍了该扫描。我们获得了与前向扫描非常相似的平均取向和偶极矩,证明样品可以可逆地极化和去极化。
偶极矩的时间序列如图所示图 3C 对每个施加的电场强度进行时间平均,以获得沿场方向 Pz 的极化平均值(方程(6)),作为场强的函数。前向和后向扫描的结果如图 4 所示(数据分别为蓝色和红色)。线性区域中的数据点在图 4 的插图中放大显示。它们适合直线,并且根据方程 1 和用于获得静态介电常数 r 的斜率。(7),
其中0是真空介电常数。我们获得前向和后向扫描的值 r = 77.9 2.7 (R2 = 0.99)、80.6 2.7 (R2 = 0.96),通过与前向加权平均数据的线性拟合得到 r = 79.3 2.2 (R2 = 0.99,绿色数据)和向后扫描,与r = 78.454的实验值一致。由于稳健的线性相关性,仅使用一个数据点即可获得类似的值。值得注意的是,每个施加电场仅经过约 175ps 的模拟时间后,介电常数就收敛了(参见补充图 S3),类似于之前报道的采用有限场哈密顿量的经典 MD 模拟13。这比根据极化波动计算介电常数所需的仿真时间少了一个数量级13,55,56,57,
其中总偶极矩 M 在零处采样电场 和 = 1.7258 是液态水的光学介电常数。事实上,沿着零场 c-NNP 轨迹对总偶极矩 M 进行采样,我们仅在 3 ns 后就获得了收敛值 r = 90.6 1.7(参见补充图 S4),与之前的模拟类似59。场扫描和零场模拟介电常数存在差异的原因尚不清楚,但可能有多种原因。极化的波动可能对模拟细节比平均值更敏感,例如使用的恒温器、有限的系统尺寸和 PNNP 的剩余误差。此外,根据方程 2 的极化波动推导出介电常数的近似值。(7)60 也可能造成这种差异。
沿着 PNNP 轨迹的 APT 和核速度可以直接获取偶极矩的时间导数, \(\dot{{{{\bf{M}}}}}\),根据方程:(4)、IR 光谱中与频率 () 相关的比尔-朗伯吸收系数 (),
其中 n() 是与频率相关的折射率,c 是真空中的光速, kB 是玻尔兹曼常数,T 是温度。在图 5 中,我们展示了在零场强下从 c-NNP MD 和在有限场强下从 PNNP MD 获得的计算频谱。正如之前报道的那样,使用 RPBE-D3 泛函可以很好地再现零场的实验红外光谱 37。光谱对线性区域及更高范围内的电场的存在非常不敏感,最高可达约 0.0514 V 1。对于较大的场,我们观察到分子内 OH 伸缩振动在 3500 cm1 左右发生系统红移,大约为0.2057 V 时为 100 cm1 1. OH 伸缩的红移通常归因于形成更强的分子间氢键 62,63,64,此处由外部电场引起。这一观点得到了 0.2057 V 1 处约 600 cm1 处约 200 cm1 的解放带系统蓝移的进一步支持,这表明与场的相互作用导致水分子旋转运动潜力的加强65。这与之前的一项研究一致,该研究描述了在这些场强下水更像冰8,9。红外光谱的趋势与之前发布的从其他方法获得的结果非常吻合,这些方法在许多方面与我们的方法不同,例如使用外部电场的 AIMD 模拟8 和使用 FIREANN-wF 模型的 MLMD 模拟42。因此,计算出的能带偏移是局部电场的敏感指纹。反过来,此类计算可以提供一种根据测量的红外光谱来估计样品中局部场强的方法。
在这项工作中,我们扩展了 MLMD 模拟以包括通过将场引起的扰动添加到未扰动的哈密顿量的一阶来与外部电场相互作用(方程(1))。由此产生的力方程(方程(2))允许我们分别计算来自标准机器学习势(此处:c-NNP)的未扰动力和来自 APTNN 的场诱导力。因此,力的计算严格遵循静电定律,无需任何额外的近似。
由此产生的称为 PNNP MD 的方案具有多个优点:首先,该方法是模块化的,因为电场贡献与所采用的电场无关。未受扰动的势能面。因此,APTNN 可以与任何 ML 势(不仅是 NNP)耦合,甚至可以与力场耦合以包括与外部电场的相互作用。根据未扰动势和扰动势的精度要求,可以不同地选择用于APTNN和ML势训练的参考电子结构计算的理论水平。此外,模块化方法允许人们单独评估每个组件的准确性。其次,APT 对于任何原子系统都是明确定义的,因此,APTNN 可以开箱即用地进行训练,不需要任何概念定制或调整或使用任意定义的代理(例如原子电荷、分子偶极矩或类似的代理)。)37。重要的是,APT 在使用周期性边界条件时是唯一定义的,因为它们量化了总偶极矩或极化的变化。这与总偶极矩形成对比,总偶极矩在周期性边界条件下是多值的46。第三,该方法不需要根据所施加的电场对训练配置进行采样。然而,在接近 DFT 参考精度、高达 0.2 V 1 量级的高场强下,可以预测与电场的相互作用。第四,可以通过在扰动中包含更高阶项(例如极化率和超极化率)来系统地改进该方法45。这将允许在超过 0.2 V 1 的极高场强下进行更准确的模拟。
这些功能使我们的 PNNP 方法不同于在 ML 模拟中对与外部电场的相互作用进行建模的其他方法38,39,40,42。FIREANN-wF 模型42 在存在电场的情况下对原子力进行训练,以执行液态水的有限场模拟,因此需要与场相关的训练数据。该模型向每个真实原子添加一个伪原子,该伪原子明确依赖于电场,从而捕获有效原子偶极矩中对施加电场的响应。它已被非常成功地用于计算响应函数,例如偶极矩和极化率,随后,根据这些量计算液态水的场相关振动谱。原则上,这种方法允许准确地模拟任意大的电场和场梯度,因为相关的极化配置明确包含在训练集中。与 PNNP 相比,在 FIREANN-wF 中,场引起的对总能量和力的贡献不是明确学习的,而是隐含在仅力训练中42。隐式学习是否是一种可以应用于广泛系统的稳健策略还有待观察。
FieldSchNet39 通过引入原子环境的向量值表示(例如利用虚构的原子)来合并外部场效应偶极矩。然后,它通过添加偶极场和偶极-偶极相互作用项,使用矢量场对分子与任意外部环境之间的相互作用进行建模。与 PNNP 不同,它接受的数据包括各种介电环境和场强中的分子结构。场的静电响应特性可作为网络的解析导数使用,并已用于计算零场下进行的 MD 模拟的振动谱。最近有人指出了这种方法的某些局限性,例如,当原子偶极子与所施加的场方向正交时,场系统相互作用的描述不完整。
作为上述方法的补充,SCFNN40 依赖于训练万尼尔的位置以水分子为中心。引入该方法的主要目标是结合远程静电学,但也允许人们在存在外部电场的情况下评估配置。该方法确实可以被推广,以便可以处理水以外的其他分子,但是必须在Advance 40中知道一组可能的分子。此外,需要定义围绕每种分子的参考框架。在水分子的情况下,这是通过将四个Wannier中心与该分子相关联并检查其局部分子框架内的坐标来完成的。因此,当需要处理更多分子时,这种结构可能引入开销,无法构建分子框架(例如,简单的离子),或者在MD模拟过程中发生反应(例如,质子在水中的质子转移)。>返回PNNP MD,我们已经表明,与AIMD相比,我们的方法可以预测准确的偶极响应动力学(图2)和电场依赖性IR光谱(图5),以及与实验性相比的准确介电常数值(图4)。对于介电常数获得的几乎定量一致性是显着的,因为通常不能期望DFT和实验水平的电子结构计算之间的完美一致。此外,在我们的模拟中忽略了核量子效应52,53。已知RPBE-D3功能与经典MD模拟结合使用,可以描述水,水溶液和界面在广泛的温度和压力上非常好37,59,65,65,66,67,67,68,69,70。在PNNP模拟的帮助下,我们可以证明RPBE-D3的强劲性能在描述径向分布功能,自扩散,定向放松过程以及液态水相对于实验数据的振动光谱扩展到其介电响应特性。这部分是因为该功能的其余缺陷往往会因缺失的核量子效应有效地补偿。实际上,我们在零视野处计算的IR谱与实验数据几乎是定量的一致性。此外,在存在的电场存在下,我们在IR光谱中计算出的峰值变化几乎与dileann-WF Model42报告的结果几乎是一致的。其中采用了更准确且计算上更昂贵的RevPBE0-D3混合功能,并结合路径积分MD模拟来明确考虑核量子效应。
PNNP的性能可能令人惊讶在0.01 V 1的中等场强度下,对总力的贡献(图1C)的贡献甚至小于未扰动的力贡献和总力的RMS,约为90 MeV 1(图。1d)。我们通过指出不受干扰的力贡献中的误差大约是高斯分布式的,并且没有方向性优先偏好,而野外诱导的力贡献当然是沿场上的净方向。因此,当对许多构型和原子进行平均时,未受干扰的力贡献中的误差取消了,而野外诱导的力贡献则不会。因此,不受干扰的力贡献中的RMSE不会损害电场响应的准确性。
我们想指出,我们的PNNP的总力量(约90 MeV 1)具有先前已报告的ML模型的典型值,例如参考文献。26,48,49,50。然而,一些最近的研究报告说,RMSE大约较小23倍,例如,据报道,火星WF Model42的RMSE为39.4 MeV 1。此处用于C-NNP的训练集并非非常详尽,仅由260个128米分子盒的配置组成。因此,它比参考文献中使用的训练集要小得多。42.同样,APTNN仅接受了27种随机选择的液体水的配置37。另一个重要方面是基础参考DFT计算的质量。出于一致原因,我们采用了与先前关于液体水的AIMD研究中使用的DFT设置完全相同的67,68。特别是,我们采用了600 ry的动能截止,而在dileann-wf Model42的情况下,使用了非常紧密的截止时间为1200 ry。众所周知,电子结构计算的更严格的收敛性导致ML Model26中的噪声和较小的RMS。因此,我们预计,通过增加训练集的结合以及使用更紧密的收敛标准进行DFT计算,可以进一步降低我们的ML模型的RMS。与AIMD和实验数据相比,我们的PNNP模型的强大性能表明,目前的RMSE就本研究而言。我们想在下面讨论。首先,通过根据等式集成电流密度来获得巡回的细胞偶极矩。(5)。因此,与先前发表的技术相反,细胞偶极矩仅通过数值整合间接获得。39,40,42。为了补偿积分误差的积累,通过通过基础电子结构方法明确计算偶极矩来定期重新启动集成(图1F)。尽管这大大减少了所需的电子结构计算的量,但对于非常大的系统大小的额外DFT计算可能会限制性能,而该系统尺寸不能再在DFT级别(例如1000s水分子)进行常规处理。但是,DFT计算不限制PNNP MD可以访问的可访问时间尺度。如图1F所示,足以在PNNP MD模拟(数十NS)上访问的时间尺度上以周期性间隔(10 ps)进行DFT计算。值得注意的是,这些DFT计算仅是PNNP MD的后处理,因此可以并行计算它们。此外,许多相关特性不需要了解巡回细胞偶极矩的知识,例如振动光谱或传输特性,其中电流密度(等式(4))就足够了。后者在PNNP方案中很容易获得,并且不需要数值集成。
,而恰当地说明了所有远程静电和非局部电荷传递效应,而我们的神经网络表示APT(我们的神经网络表示)(APTNN)仅使用本地描述符,即,它不包含明确的远程静电信息。这在精神上与这里用于代表原子间潜力的第二代C-NNP相似。带有本地描述符的ML模型必定会因不再有效筛选长时间静电的蒸气阶段或Apolar培养基而失败,或者在电场引起远距离电荷转移的情况下。但是,这是我们PNNP方法的模块化提供明显优势的地方,因为它使我们可以根据需要替换基础的ML模型。例如,当前对原子间电位进行建模的第二代C-NNP可以被第四代C-NNP26,71,72取代,该C-NNP26,71,72明确地考虑了远程静电和非局部电荷转移。同样,带有本地描述符的APTNN最终可以用包括非本地描述符在内的版本代替。此处引入的想法通常是使用APT使用电场的原子间电位,这通常是有效的。
我们还注意到电场等式中的扰动系列扩展。(4)仅包含偶极期。场诱导的力贡献中的低相对RMSE表明,不需要高阶项(偶极偏振性和超极化能力),而高磁场强度约为0.2 V 1(图1E)。然而,此时,RMSE缓慢而稳定地增加,表明偶极子极化在这种电场状态下变得越来越重要。原则上,可以以与APT相似的方式学习两极分化的张量,并且我们目前正在为此目的探索有效的方案。最后,当前方案的另一个局限性是我们仅考虑均匀的电场。为了在诱导的不均匀磁场下(例如,通过STM尖端或电极界面)进行模拟,需要补充多极扩展等式。(4)具有四极相互作用项(因为现场梯度与四极矩相互作用)。这可以使用从APT获得的针对四极力矩的计算的出生指控来完成。
总而言之,我们已经实施了一种ML方法,表示PNNP,在存在外部的情况下运行MD模拟同质电场。这种发展的关键是APT,它是由ML模型训练的,用于计算场引起的力。后者与描述原子间电位的标准ML电位获得的力相结合。该方法是模块化的,它利用了在周期性边界条件下定义得很好的APT,它不需要训练作为应用电场的函数,并且可以系统地改进。
我们进行了PNNP模拟在几种不同的场强度上,我们的结果与参考AIMD计算和实验文献数据进行了比较。通常,我们发现对所研究的所有属性都非常达成共识。我们证明了随着电场的上下走动,水分子的可逆定向极化和去极化,并在AIMD数据的统计准确性中具有定向的松弛时间(图2)。这允许在模拟时间尺度上成功地计算相对介电常数R,该计算比在零场处极化波动的非常缓慢收敛的更慢的标准方法中短得多(图4)。获得的值与实验验证了我们在线性响应(低电场)制度中的方法非常吻合。我们还计算了红外光谱,并将它们与以前的显式AIMD模拟进行了比较。高电场上的库和分子内OH弹力振动的系统峰值移动很好地匹配了AIMD参考数据,从而验证了我们的非线性响应(高电场)状态的方法。
我们期望PNNPMD成为与外部电场相互作用的广泛凝结相系统的AB级模拟的有用工具,以及在这项工作中未研究的电场依赖性属性的计算,包括离子电导率和电容。相对于电场的能量扩展中高阶项的实施将在更高的场强度下进行准确的模拟,并在APT框架内访问浓缩相系统的拉曼和总和 - 频率生成光谱。此外,由于我们方法的模块化,人们可以很容易地利用非本地ML电位的未来发展,这对于在介电筛选不如水中的媒体中模拟可能是必不可少的。
实现了不受干扰的势能,EPOT和相应的力,这是等式右侧的第一项。(2),由最近在CP2K软件包44中实施的第二代高维神经网络电位的委员会建模。与APT有关的磁场依赖性贡献,第二项在等式的右侧。(2),由APTNN建模,该APTN基于基于Pytorch Library 74的E(3) - 等级图神经网络E3NN73。在CP2K中添加了一个新的力量评估环境,该环境将基于Fortran的CP2K和客户端/服务器方法中基于Python/C ++的Pytorch链接。受I-PI实现75CP2K的启发,并连接到等待接收配置并寄回APTNN模型预测的APTS的Python服务器。然后,APT在CP2K中使用,以评估相应的力贡献。使用内置混合力环境将它们添加到C-NNP力中,以获取使用速度 - verlet算法传播原子的总力。
使用N2P276对C-NNP的委员会成员进行了培训。神经网络的参数来自先前关于液体水的ML潜在研究59,并与通用对称函数结合使用44。使用RPBE-D3级别从DFT计算获得的能量和力对C-NNP进行了训练,如下所述,使用参考文献中报道的主动学习程序。44.关于APT,我们使用了最近发表的APTNN由US37之一开发的APTN,其中包含APT,该APT在其训练集中仅在RPBE-D337级别的训练集中仅27个平衡的128水分子盒的快照。通过单点有限差计算获得了最终的10,368个APT,这需要大量的计算工作。请注意,本工作的目的是介绍PNNP方法,并提供原则证明,而不是优化培训协议。在这里,我们只是使用了以前的工作中已经可以使用的最大培训套件37。这使我们能够基准基准PNNP的准确性,同时确保使用太小的训练集可能会遇到的错误。
我们非常有信心生成aptnn的培训数据的计算工作能力可以很大程度上显着在不久的将来减少。首先,我们先前表明,IR光谱仅需要9个训练配置,即128个水分子盒(即3456个APTS)才能准确地重现参考光谱,从而将所需的参考计算减少了三个。我们预计这也适用于其他感兴趣的可观察物。其次,我们目前针对APTNN的培训策略尚未得到优化。在给定配置中计算所有原子的APT可能包括冗余数据,这是由于相邻原子的APT之间的高相关性。以原子为基础的训练仅针对原子的一部分计算,预计将更有效。该方法也可以与主动学习结合使用,在该学习中,网络迭代为训练集选择原子。例如,通过使用C-NNP44中的查询方法使用查询方法。我们预计,这两种策略使我们能够减少参考APT计算的数量,甚至更多。最后,例如,通过使用分析密度函数扰动理论77而不是有限的差异计算,也可以使计算APT的方法变得更有效。探索这些策略超出了本手稿的范围,但将在以后的工作中解决。
使用参考培训和测试数据的电子结构计算计算CP2K78版本2023.1和QuickStep Module79。RPBE功能通过LIBXC软件包80评估,并补充D3分散校正81。我们采用了混合的高斯轨道/平面波(GPW)基set82,平面波截止为600 ry,相对临界值为20 ry,高斯轨道是使用三重质量的TZV2P基set83构建的,包括偏振功能,以前65,66,67,68。核心电子通过稳定的Goedecker-Hutter(GTH)伪电位84,85来描述。通过Umari和Pasquarello86在CP2K中实施的方法处理均匀的电场。通过最大局部的Wannier函数从DFT计算获得了模拟电池的偶极矩16。
电荷保守性通过pNNP严格执行,以通过pNNP实施。原子极张量77,87,88,
这是通过计算每个组件的总和除以原子数并在系统中所有原子中分布任何小的多余物品来完成的。请注意,此校正类似于纠正原子电荷,就像先前在文献中所做的那样,在整个System 89,90,91中均匀分配任何超额费用。
所有计算,包括训练,测试和生产运行,用于使用周期性边界条件(密度0.996 kg L1)的128水分子盒15.6627进行。使用CSVR恒温器在NVT集合中进行有限场PNNP MD模拟,其时间常数为1 ps92,MD积分时间步长为1 fs,零场C-NNP MD模拟在NVT中使用NONES在NONENMEMEN中进行 - 恒温器93,94,时间步长为0.5 fs,除非另有说明。在NVE集合中进行了定向放松动力学和红外光谱的模拟。像以前的工作一样,计算了啤酒 - 叶片吸收系数。65。对于每种电场强度(0.0514、0.1028和0.2057 V 1),从场扫描PNNP模拟中选择了20个独立的配置。对于NVE集合中的每个配置中的每一个简短(20 ps)PNNP模拟。计算了20个模拟中的每个模拟,然后进行平均。
有关研究设计的更多信息,请提供与本文有关的自然投资组合报告摘要。
源图数据可在FIGSHARE上获得以下https://doi.org/10.6084/m9.figshare.2671663195。补充数据也可以在Figshare上获得,可提供以下https://doi.org/10.6084/m9.figshare.2666935096.26666935096.
使用所有现场模拟进行了所有现场模拟CP2K-2023.1软件包,已在内部进行自定义。这是一个开发版本,可在以下公开场所公开可用:https://github.com/kjaj98/cp2k-apt-pnnp-paper,https://doi.org/10.5281/zenodo.136278497。原子极性张量神经网络和用于训练和预测原子极性张量的所有脚本均可公开提供:https://github.com/pschienbein/atomicpolartensor,https://doii.org/10.5281/zenodo.13323558798。N2P2的2.1.0版用于培训委员会神经网络。
Clarke,D.,Whitney,H.,Sutton,G。