(Na2O)0.1-(V2O5)0.6-(TeO2)0.3玻璃结构的定量评估：基于第一性原理和机器学习的混合价氧化物体系原子间势研究

《Artificial Intelligence Chemistry》：Quantitative assessment of the structure of (Na2O)0.1-(V2O5)0.6-(TeO2)0.3 glass: An ab initio and machine learning interatomic potential study of oxide systems with mixed-valence

　　费拉斯·肖艾布（Firas Shuaib）、雷米·皮奥特罗夫斯基（Remi Piotrowski）、盖尔·德拉伊齐尔（Gaelle Delaizir）、皮埃尔-马里·杰弗罗伊（Pierre-Marie Geffroy）、史蒂夫·戴夫·万西·温吉（Steve Dave Wansi Wendji）、卡洛·马索布里奥（Carlo Massobrio）、毛罗·博埃罗（Mauro Boero）、吉多·奥里（Guido Ori）、菲利普·托马斯（Philippe Thomas）、奥利维尔·马松（Olivier Masson）、阿西尔·布兹德（Assil Bouzid）
法国利摩日大学（Université de Limoges）陶瓷研究机构（IRCER），隶属CNRS UMR 7315项目；欧洲陶瓷中心（Centre Européen de la Céramique），地址：利摩日Atlantis街12号，邮编87068

**摘要**
本研究采用第一性原理分子动力学（FPMD）和机器学习原子间势（MLIP）方法，探究了玻璃态系统的原子尺度结构。基于MLIP生成的结构模型与FPMD结果高度吻合，同时具备显著的计算效率，能够模拟远大于FPMD所能处理的系统规模，且保持相当的精确度。通过与实验测得的X射线总结构因子和配对分布函数的比对，验证了FPMD和MLIP生成结构的可靠性。研究发现玻璃网络由VO4、VO5、VO6多面体以及TeO3、TeO4、TeO5单元组成。通过最大局域化Wannier函数（MLWF）基的键合分析、配位数分析及键角分布研究了局部结构特征。此外，MLIP模型还被用于研究不同V5+和V4+含量对玻璃结构的影响，以探讨氧化态依赖性的变化。实验结果表明，利用MLIP基分子动力学创建并通过单次DFT计算放松的熔融急冷玻璃模型能够高精度地捕捉局部拓扑结构、配位关系及V–O键合模式的变化，从而证明了MLIP模型在描述混合价态氧化物玻璃系统中的有效性。这一现象可归因于过渡金属氧化物中的电荷自调节机制——即氧化态的变化通过局部金属–O电子重分布及相关结构重组实现。研究表明，MLIP能够准确描述具有多种成分和氧化态的复杂氧化物玻璃，进一步扩展了其在化学多样性无序体系中的应用范围。

**1. 引言**
含有过渡金属的氧化物玻璃因其多功能性和广泛的应用潜力而受到广泛关注，应用领域涵盖固态器件、能源储存、智能窗户及光学组件等。其中，基于二氧化碲（TeO2）的玻璃因具有高折射率[1]、高介电常数[4]、优异的红外透明性能[12, 13, 14, 15, 16]等特性而独树一帜。TeO2玻璃还具有出色的化学稳定性和机械稳定性[6, 17, 18, 19]，熔点相对较低，并且不像磷酸盐或硼酸盐玻璃那样具有吸湿性[18, 20, 21, 22, 23, 24, 25, 26]。这些特性使得TeO2玻璃在基础科学研究和技术应用中都极具吸引力。虽然纯TeO2本身难以形成玻璃态，但与其他氧化物（尤其是V2O5、MoO3、P2O5）结合后可形成稳定的玻璃[23, 27, 28, 29, 30]。在纯TeO2玻璃中，Te原子通常以TeO4双锥构型存在[5, 10, 11, 16, 27, 31, 32]，但加入Na2O等改性剂后会改变这种局部结构。例如，Na2O会促使TeO4双锥结构转化为TeO3三角锥结构[10, 21]，从而导致玻璃网络解聚并产生非桥接氧原子（NBOs）[10, 21]，这些NBOs显著影响材料的电子结构和导电性[33]。
钒氧化物（V2O5）作为典型的过渡金属氧化物，兼具玻璃形成能力和导电性能。在二元系统中，V2O5引入多价性（V5+和V4+），促进极化子跃迁并使其具备半导体特性[27, 34, 35, 36, 37]。当V2O5含量低于20%时，玻璃主要形成连续的碲酸盐网络[27]；而超过该浓度后，结构转变为钒酸盐为主的网络[27, 34, 38, 39, 40]，此时VO4四面体与VO5或VO6多面体共存，影响材料的弹性、光学吸收和导电性[20, 25, 27, 34, 35, 38, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54]。
掺杂碱金属（M）氧化物（如Na2O、LiO2）的三元钒碲酸盐玻璃展现出混合的电子-离子导电性[49, 50, 55]。在V/V比率较低时，玻璃靠极化子跃迁实现电子导电；而在碱金属含量较高时，离子导电性占主导[35, 56]。此前研究指出，这种导电性最小值出现在V/Na比值为1附近[50]，这一现象源于电子局域化（由NBO形成）与离子扩散（受碱金属改性剂增强）之间的相互作用。TeO2的存在为这些系统提供了额外的结构可调性。

相比之下，钠钒磷酸盐（NVP）玻璃在储能应用方面的结构研究更为深入[57, 58, 59, 60, 61, 62]，其结构基于通过PO4单元连接的VOn多面体，形成开放的网络结构，有利于Na离子的迁移[63, 64, 65, 66, 67, 68]。NVP玻璃也具有很强的玻璃形成能力和混合价态导电性。然而，直到最近，通过BOMD和MLIP模拟才得以详细解析其复杂的局部环境、键合特性及离子迁移规律[69, 70, 71]。尽管如此，由于磷酸盐基网络的吸湿性，其结构稳定性相对较弱，不如钒碲酸盐玻璃。

NVP玻璃的双重导电特性对化学成分引起的结构变化敏感[35, 72]，但这种转变的微观机制仍不明确。结构变化（尤其是氧桥的形成或破坏）被认为起关键作用[19, 35, 73, 74]，然而目前针对这类三元体系的先进计算方法（如红外/拉曼光谱或XRD）研究仍然有限[49, 50, 55]。尽管已有大量使用经典分子动力学（CMD）和Born–Oppenheimer分子动力学（BOMD）对NVP玻璃进行的模拟[75, 76, 77, 78, 79, 80, 81]，但针对其复杂结构特性的系统研究仍十分匮乏。尤其是将第一性原理分子动力学（FPMD）与机器学习原子间势（MLIP）及X射线衍射（XRD）结合的研究尚属首次。近期FPMD模型的进展突破了传统力场的局限[82, 83, 84]，实现了对无序体系的量子力学精确模拟[82, 85]，尽管计算成本较高（尤其是对于大系统及长模拟时间[86]。MLIP模型通过学习从头算数据，能够在保持接近DFT精度的同时模拟数千个原子的大规模系统[70, 94, 95]，例如高分辨近似势（GAP）和基于神经网络的方法已成功模拟了非晶材料的高压相变和成核过程，并预测了非晶态的局部电荷密度和键合环境[88, 96, 97, 98, 99, 100, 101, 102]。

尽管MLIP技术日渐成熟，但在包含三种不同阳离子和氧的三元非晶氧化物体系中的应用仍有限。目前仅有少数基于GAP的MLIP模型用于这类材料，且仅适用于特定成分范围。主要挑战在于这类体系具有较高的化学多样性和结构变异性，导致核函数模型的训练和验证难度较大。本研究首次结合FPMD、大规模GAP基分子动力学模拟及实验XRD分析，揭示了玻璃系统的详细结构特征，并探讨了V2O4/V2O5比例变化对玻璃结构的影响[0.0至0.15]。研究发现，仅基于V5+环境训练的GAP模型能够准确捕捉整个成分范围内的混合价态（V5+/V4+）结构特征，这一发现通过DFT自旋密度分析和实验验证得到证实，展现了过渡金属氧化物中局部电荷自调节的基本机制。

**2. 计算细节和模型生成流程**
2.1 第一性原理分子动力学模拟与模型构建
模拟工作基于密度泛函理论（DFT），使用CP2K模拟软件包[103]进行。电子结构计算采用高斯平面波（GPW）方法，其中Kohn–Sham轨道采用双ζ原子中心的高斯基组表示[104]，电子密度投影到动能截止为800 Ry的辅助平面波格子上。核-价层相互作用通过Goedecker–Teter–Hutter（GTH）赝势描述[105]，广义梯度近似（GGA）通过Perdew–Burke-Ernzerhof（PBE）交换相关泛函实现[106]。由于系统无晶态且体积庞大，布里渊区采样仅限于Γ点。

为了研究玻璃态结构演变，采用了Born–Oppenheimer分子动力学（BOMD）[107, 108]，时间步长为1 fs，确保总能量的稳定守恒。动力学过程中采用正则系综（NVT）方案，温度控制通过Nosé–Hoover恒温器链实现[109, 110]。模拟体系为三元玻璃组合NVT60，立方晶胞各方向周期性，反映了参考文献[111]中报告的3.514 g/cm3密度，晶胞尺寸为18.34 ?，包含432个原子（16个Na、96个V、24个Te和296个O），随机分布在晶胞内，所有V原子均呈五价态（V5+）。初始配置在0 K下放松后，经过6000 K的高温平衡阶段（5 ps），随后在4000 K下进行短暂平衡（5 ps），再经历多阶段温度退火：T=2000 K（10 ps）、T=1800 K（15 ps）、T=1500 K（20 ps）、T=1200 K（30 ps）、T=900 K（20 ps）、T=600 K（25 ps），最后在T=300 K下进行35 ps的统计平均。

在T=1200 K的熔化阶段，所有原子均表现出显著扩散行为，扩散系数约为10?5 cm2 s?1。Te和V原子的位移远超典型原子间距，表明系统已完全失去初始结构记忆。

玻璃中的键合结构通过最大局域化Wannier函数（MLWFs）[112, 113]进行分析。MLWFs wn(r)及其中心rnWFC=(xn,yn,zn)通过对Kohn–Sham轨道ψi(r)的幺正变换得到[114]，变换公式为wn(r)=∑i=1NoccAi,n(p)ψi(r)，其中Ai,n(p)是依赖于迭代的幺正矩阵[114]。通过最小化扩散函数Ω=∑i=1Nocc〈i|r2|i〉?〈i|r|i〉2来优化变换（量化轨道的空间范围）。Wannier中心坐标例如沿x轴计算为xn=?Lx2π?ln〈wn|e?i2πxLx|wn〉，Lx表示x方向的晶胞长度；yn和zn的计算方法类似。实际计算中，Wannier中心是基于T=300 K时分子动力学轨迹最后25 ps内的100个构值得出的，每个中心代表一对电子，可能关联到化学键（WB）或原子α上的孤对电子WαLP。如果两个原子α和β共享一个Wannier中心，并且满足|dαβ?dαW?dβW|≤0.05?（其中dαβ是原子间距离，dαW和dβW分别是每个原子到共享WFC的距离），则认为它们是键合的。为了排除非键合贡献，还规定了O、Te和孤对中心WTeLP之间至少73°的角截止条件，如参考文献[115]、[116]中所展示的。利用这一框架，氧原子可以被分类为桥接（BO），形成如TeOTe、TeOV或VO这样的连接结构，或者被分类为非桥接（NBO），形成如TeO或VO这样的端Bond结构。这种分析方法能够在不依赖于任意距离截止条件的情况下，化学上一致地确定配位数和结构基元。

下载：下载高分辨率图片（388KB）
下载：下载全尺寸图片

图1. 从GAP MLIP和FPMD计算得到的能量（上图）和力分量（下图）在0.1、0.6、0.3玻璃系统测试数据库上的比较。

2.2. 机器学习建模
基于高斯近似势（GAP）框架开发的机器学习原子间势（MLIP）用于补充对0.1、0.6、0.3玻璃系统的研究。GAP采用基于核的方法，从DFT数据中高保真地插值能量和力，从而能够有效地对无序系统进行建模，同时保持接近从头算的准确性[88]、[117]。尽管GAP已被证明对纯玻璃和二元玻璃有效[70]、[94]、[95]、[118]，但其应用于三元非晶系统的能力仍然有限。在我们之前对不同钒氧化态的NVP玻璃的研究中，基于GAP方案和SOAP描述符的MLIP模型被证明特别可靠且稳健[69]、[70]、[71]。

在GAP框架内，总能量表示为短程两体贡献和基于原子位置平滑重叠（SOAP）描述符的多体项之和[85]、[88]、[117]、[119]、[120]。在这种形式主义中，原子i的局部环境由所有满足rij≤rcut的邻近原子j定义。相应的物种分辨邻居密度表示为（1）ρiα(r)=∑j∈Ni,Zj=αfcut(rij)exp?|r?rij|2/2σatom^2，其中α标记化学物种，fcut是一个平滑截止函数，σatom是高斯模糊宽度。该密度展开为径向函数gn(r)和球谐函数Y?m(r?)的基，公式为（2）ρiα(r)=∑n?mcn?mi,αgn(r)Y?m(r?)，由此构造出旋转不变的SOAP功率谱（3）pnn′?i,αβ=∑mcn?mi,αcn′?mi,β*。所有pnn′?i,αβ的集合定义了与原子i相关的SOAP描述符pi。然后总能量表示为（4）E=∑i
mlip的超参数（见补充表s1）是通过结合对navteo玻璃网络的物理直觉以及我们之前为化学相关复杂玻璃开发的gap-soap模型的经验来选择的[70]、[94]、[95]、[121]。从这个基础出发，通过检查插值误差和得到的势能的整体行为，迭代地改进了超参数。实际上，我们通过调整nmax和lmax以及稀疏表示和正则化参数来变化描述符的分辨率，并保留了仍然能再现dft参考数据并且在考虑的结构和温度范围内表现稳定的最简单模型。这一点很重要，因为对于基于核的模型（如gap），可转移性不仅由描述符的大小决定，还受到描述符定义、稀疏化、正则化以及训练数据库中代表的局部环境多样性的综合影响。

特别注意了soap的分辨率。采用的值nmax=8和lmax=9，被发现能够提供足够的径向和角度细节，以区分这种玻璃中存在的扭曲和化学多样化的局部基元，包括nao环境以及vo和teo多面体。在试验模型中，较低的值导致这些环境之间的区分能力明显降低，力预测也更差。尽管更大的稀疏集可以稍微提高固定描述符的回归质量，但它无法完全恢复从分辨率不足的soap基础中缺失的信息。因此，我们选择了一个具有足够内在分辨率的描述符，然后使用稀疏表示来控制回归的有效复杂性。将nmax和lmax增加到选定的值之外并没有带来任何有意义的改进，所以最终选择代表了准确性と效率之间的实际折中。

截止半径rcut=4.27 ?和高斯模糊宽度σatom=0.6 ?被选中，以提供平滑的描述符，同时保留玻璃网络的第一层配位壳和最相关的近邻角相关性。截止过渡宽度rδ=1.0 ?确保了描述符在局部环境边界处的平滑衰减。soap能量尺度参数δ=0.88ev设置了由soap贡献表示的局部能量变化的特征幅度。稀疏集大小nsparse=1200用于soap，对于两体项则为70，一旦测试误差不再随进一步增加而显著改善，同时仍允许以合理的成本进行长时间分子动力学模拟。对于两体贡献，高斯长度尺度固定为θ=1.0 ?，提供了一个平滑的短程排斥基线。最后，正则化参数σvenergy=0.001 ev/原子、σvforce=0.05 ev/?和σvvirial=0.05 ev/原子被选中，以保持与dft能量景观的紧密一致，同时避免过度拟合fpmd快照中存在的有限温度噪声。

gap mlip模型的训练数据是从t=2000 k到t=300 k的宽温度范围内进行的fpmd模拟中提取的，以确保结构多样性。具体来说，从每个模拟平台采样了70个快照，总共得到490个配置。为了提高数据集的保真度，使用更高的平面波截止值1000 ry，重新计算了相关的dft能量、力和维里应力。这个数据集在训练集和测试集之间平均分配。我们注意到，基于gap的mlip即使在小型训练数据集上也能够达到化学准确性，通常只需要几百个配置就可以处理无序材料[94]。这种数据效率来源于gap的基于核的性质，它能够从稀疏的、有代表性的原子环境中进行稳健的外推。得到的gap模型表现出强大的预测准确性，在测试集上力学的rmse值为0.38 ev/?（见图1），能量的rmse值为2.91 mev/原子，完全在常见的5 mev/原子阈值[98]、[99]范围内。

为了评估mlip强再现dft计算的能量和力的能力是否也扩展到准确的结构描述，我们在不同大小的系统上进行了模拟。使用一个432原子的模型（mlip_s）与fpmd结果进行了直接比较，同时还构建了一个包含多达3456个原子的更大系统（mlip_b），以研究fpmd模拟计算上无法获得的潜在大小效应。对于所有模拟，模型都是用随机初始化的结构和实验上一致的密度（3.514 g/cm^3）构建的。分子动力学模拟是在lammps代码[122]中实现的，采用nvt系综和0.5 fs的时间步长。所有系统都经历了多个温度平台的热退火：在t=300 k、t=500 k、t=750 k和t=1000 k时各5 ps，然后在t=2000 k、1800 k、1500 k、1200 k、900 k和300 k时各100 ps。结构性质是从t=300 k时轨迹的最后50 ps中分析的。

通过分析整个系统和每个组成元素（na、o、te、v）的振动密度态（vdos）来比较fpmd和mlip模拟的振动性质。vdos是从速度自相关函数（vacf）的傅里叶变换计算得到的。沿md轨迹的原子速度得到的vacf定义为c(t)=1/n∑i=1nvi(0)?vi(t)，其中n是原子数量，vi(t)是时间t时原子i的速度。然后vdos由g(ω)∝∫0∞c(t)cos(ωt)dt给出，或者等效地通过从轨迹评估的标准化vacf的离散傅里叶变换得到。如图2所示，mlip与fpmd结果非常一致，捕捉了所有元素和频率范围内的vdos的关键特征。特别是，主要峰和整体谱形都被很好地再现了，证实了mlip准确模拟了玻璃的振动响应。这验证了mlip在进一步的玻璃结构、能量和动力学研究中的可靠性。

作为对获得的mlip的额外验证，我们使用我们的mlip评估了体积模量，并将结果与fpmd和相关的碲酸盐-钒酸盐玻璃的实验数据进行比较。

下载：下载高分辨率图片（456kb）
下载：下载全尺寸图片

图2. 玻璃的振动密度态，比较了fpmd（黑色）和mlip（红色）模拟。总vdos及其元素分辨的贡献（na、o、te、v）是从速度自相关函数的傅里叶变换计算得到的。体积模量是通过在432原子模型上的10个独立非晶快照上执行的能量-体积计算获得的。对于每个快照，通过将晶格参数围绕参考值变化±3%的步骤为±0.5%来生成一组均匀缩放的立方晶胞，每个快照产生13个能量-体积点。在dft和mlip级别应用了相同的各向同性缩放协议，每个e(v)数据集都用birch-murnaghan状态方程[123]进行拟合，以提取平衡体积模量b0。报告的值对应于10个快照的平均值，而引用的不确定性表示相关的标准偏差。

使用这个程序，我们得到b0fpmd=24.85±1.47gpa和b0mlip=27.88±4.82gpa。这些结果表明，从fpmd和mlip计算得到的平均体积模量在估计的标准偏差范围内几乎相当。此外，考虑到弹性响应的固有快照间变异性，这种一致性对于目前的非晶氧化物玻璃来说是合理的。获得的值也与相关碲酸盐-钒酸盐玻璃的实验范围一致。在二元（1?x）teo2-xv2o5玻璃中，室温下的超声测量报告的体积模量根据成分的不同从大约27.2到37.3gpa不等，纯teo2的体积模量为31.7gpa[20]、[124]。在更复杂的teo2 bi2o3 v2o5 na2o tio2玻璃系统中，实验报告的体积模量范围为16.4–25.2gpa[125]。钠钒酸盐玻璃（成分xna2o (35?x)v2o5 65teo2）也显示出成分依赖的弹性模量，随着钠含量的增加而软化[126]。因此，这里得到的fpmd和mlip体积模量都在相关碲酸盐-钒酸盐氧化物玻璃的更广泛实验窗口内，这支持了计算值的物理合理性。

最后，我们想强调的是，当前的mlip是在对应于单一密度的配置上训练的，这本质上限制了其在更宽密度范围内定量再现体积依赖性质（如体积模量）的能力。从这个意义上说，体积模量是对模型特别严格的测试，因为它探测了能量景观相对于各向同形压缩和膨胀的曲率。因此，mlip和dft体积模量之间的一致性强调了我们潜在模型足够准确地描述潜在能量表面的能力。

如前所述，玻璃中的钒可以采取多种氧化态，从v2+到v5+[127]、[128]、[129]不等，这取决于成分和加工条件。v5+通常是主导的，即使在非有意还原条件下也常常存在不可忽视的v4+。这些状态与不同的局部原子环境相关[127]、[128]、[130]。 2θ^2，而多体项则以稀疏核形式表示为（6）?isoap=∑m=1MSOAPαmSOAPKSOAP(pi,pm)。这里，pm表示从训练数据库中选取的代表性环境的稀疏集合，SOAP核定义为归一化点积核的幂（7）KSOAP(pi,pm)=Δ^2pi?pmpi?pipm?pmζ，其中Δ设置SOAP贡献的特征能量尺度，ζ控制相似度测量的锐度。回归系数αm(2B)和αmSOAP通过正则化拟合DFT参考能量、力和维里应力获得。实际上，拟合是通过使用能量、力和维里应力的稀疏高斯过程回归来进行的，从而确保得到的势能表面在保持大规模分子动力学模拟所需的局部性、平滑性和计算效率的同时，再现了DFT势能表面[88]、[117]、[119]、[120]。 mlip的超参数（见补充表s1）是通过结合对navteo玻璃网络的物理直觉以及我们之前为化学相关复杂玻璃开发的gap-soap模型的经验来选择的[70]、[94]、[95]、[121]。从这个基础出发，通过检查插值误差和得到的势能的整体行为，迭代地改进了超参数。实际上，我们通过调整nmax和lmax以及稀疏表示和正则化参数来变化描述符的分辨率，并保留了仍然能再现dft参考数据并且在考虑的结构和温度范围内表现稳定的最简单模型。这一点很重要，因为对于基于核的模型（如gap），可转移性不仅由描述符的大小决定，还受到描述符定义、稀疏化、正则化以及训练数据库中代表的局部环境多样性的综合影响。 特别注意了soap的分辨率。采用的值nmax=8和lmax=9，被发现能够提供足够的径向和角度细节，以区分这种玻璃中存在的扭曲和化学多样化的局部基元，包括NaO环境以及VO和TeO多面体。在试验模型中，较低的值导致这些环境之间的区分能力明显降低，力预测也更差。尽管更大的稀疏集可以稍微提高固定描述符的回归质量，但它无法完全恢复从分辨率不足的SOAP基础中缺失的信息。因此，我们选择了一个具有足够内在分辨率的描述符，然后使用稀疏表示来控制回归的有效复杂性。将nmax和lmax增加到选定的值之外并没有带来任何有意义的改进，所以最终选择代表了准确性と效率之间的实际折中。 截止半径rcut=4.27 ?和高斯模糊宽度σatom=0.6 ?被选中，以提供平滑的描述符，同时保留玻璃网络的第一层配位壳和最相关的近邻角相关性。截止过渡宽度rδ=1.0 ?确保了描述符在局部环境边界处的平滑衰减。soap能量尺度参数δ=0.88eV设置了由SOAP贡献表示的局部能量变化的特征幅度。稀疏集大小nsparse=1200用于SOAP，对于两体项则为70，一旦测试误差不再随进一步增加而显著改善，同时仍允许以合理的成本进行长时间分子动力学模拟。对于两体贡献，高斯长度尺度固定为θ=1.0 ?，提供了一个平滑的短程排斥基线。最后，正则化参数σvenergy=0.001 ev 原子、σvforce=0.05 ev ?和σvvirial=0.05 ev 原子被选中，以保持与dft能量景观的紧密一致，同时避免过度拟合fpmd快照中存在的有限温度噪声。 gap mlip模型的训练数据是从t=2000 k到t=300 k的宽温度范围内进行的fpmd模拟中提取的，以确保结构多样性。具体来说，从每个模拟平台采样了70个快照，总共得到490个配置。为了提高数据集的保真度，使用更高的平面波截止值1000 ry，重新计算了相关的dft能量、力和维里应力。这个数据集在训练集和测试集之间平均分配。我们注意到，基于gap的mlip即使在小型训练数据集上也能够达到化学准确性，通常只需要几百个配置就可以处理无序材料[94]。这种数据效率来源于gap的基于核的性质，它能够从稀疏的、有代表性的原子环境中进行稳健的外推。得到的gap模型表现出强大的预测准确性，在测试集上力学的rmse值为0.38 ev ?（见图1），能量的rmse值为2.91 mev 原子，完全在常见的5 mev 原子阈值[98]、[99]范围内。 为了评估mlip强再现dft计算的能量和力的能力是否也扩展到准确的结构描述，我们在不同大小的系统上进行了模拟。使用一个432原子的模型（mlip_s）与fpmd结果进行了直接比较，同时还构建了一个包含多达3456个原子的更大系统（mlip_b），以研究fpmd模拟计算上无法获得的潜在大小效应。对于所有模拟，模型都是用随机初始化的结构和实验上一致的密度（3.514 g cm^3）构建的。分子动力学模拟是在lammps代码[122]中实现的，采用nvt系综和0.5 fs的时间步长。所有系统都经历了多个温度平台的热退火：在t=300 k、t=500 k、t=750 k和t=1000 k时各5 ps，然后在t=2000 k、1800 k、1500 k、1200 k、900 k和300 k时各100 ps。结构性质是从t=300 k时轨迹的最后50 ps中分析的。 通过分析整个系统和每个组成元素（na、o、te、v）的振动密度态（vdos）来比较fpmd和mlip模拟的振动性质。vdos是从速度自相关函数（vacf）的傅里叶变换计算得到的。沿md轨迹的原子速度得到的vacf定义为c(t)=1/N∑i=1Nvi(0)?vi(t)，其中N是原子数量，vi(t)是时间t时原子i的速度。然后VDOS由g(ω)∝∫0∞C(t)cos(ωt)dt给出，或者等效地通过从轨迹评估的标准化VACF的离散傅里叶变换得到。如图2所示，MLIP与FPMD结果非常一致，捕捉了所有元素和频率范围内的VDOS的关键特征。特别是，主要峰和整体谱形都被很好地再现了，证实了MLIP准确模拟了玻璃的振动响应。这验证了MLIP在进一步的玻璃结构、能量和动力学研究中的可靠性。 作为对获得的mlip的额外验证，我们使用我们的mlip评估了体积模量，并将结果与fpmd和相关的碲酸盐-钒酸盐玻璃的实验数据进行比较。 下载：下载高分辨率图片（456kb） 下载：下载全尺寸图片 图2. 玻璃的振动密度态，比较了fpmd（黑色）和mlip（红色）模拟。总vdos及其元素分辨的贡献（na、o、te、v）是从速度自相关函数的傅里叶变换计算得到的。体积模量是通过在432原子模型上的10个独立非晶快照上执行的能量-体积计算获得的。对于每个快照，通过将晶格参数围绕参考值变化±3%的步骤为±0.5%来生成一组均匀缩放的立方晶胞，每个快照产生13个能量-体积点。在dft和mlip级别应用了相同的各向同性缩放协议，每个e(v)数据集都用birch-murnaghan状态方程[123]进行拟合，以提取平衡体积模量b0。报告的值对应于10个快照的平均值，而引用的不确定性表示相关的标准偏差。 使用这个程序，我们得到b0fpmd=24.85±1.47GPa和B0MLIP=27.88±4.82GPa。这些结果表明，从FPMD和MLIP计算得到的平均体积模量在估计的标准偏差范围内几乎相当。此外，考虑到弹性响应的固有快照间变异性，这种一致性对于目前的非晶氧化物玻璃来说是合理的。获得的值也与相关碲酸盐-钒酸盐玻璃的实验范围一致。在二元（1?x）TeO2-xV2O5玻璃中，室温下的超声测量报告的体积模量根据成分的不同从大约27.2到37.3GPa不等，纯TeO2的体积模量为31.7GPa[20]、[124]。在更复杂的TeO2 bi2o3 v2o5 na2o tio2玻璃系统中，实验报告的体积模量范围为16.4–25.2gpa[125]。钠钒酸盐玻璃（成分xna2o (35?x)v2o5 65teo2）也显示出成分依赖的弹性模量，随着钠含量的增加而软化[126]。因此，这里得到的fpmd和mlip体积模量都在相关碲酸盐-钒酸盐氧化物玻璃的更广泛实验窗口内，这支持了计算值的物理合理性。 最后，我们想强调的是，当前的mlip是在对应于单一密度的配置上训练的，这本质上限制了其在更宽密度范围内定量再现体积依赖性质（如体积模量）的能力。从这个意义上说，体积模量是对模型特别严格的测试，因为它探测了能量景观相对于各向同形压缩和膨胀的曲率。因此，mlip和dft体积模量之间的一致性强调了我们潜在模型足够准确地描述潜在能量表面的能力。>
mlip的超参数（见补充表s1）是通过结合对navteo玻璃网络的物理直觉以及我们之前为化学相关复杂玻璃开发的gap-soap模型的经验来选择的[70]、[94]、[95]、[121]。从这个基础出发，通过检查插值误差和得到的势能的整体行为，迭代地改进了超参数。实际上，我们通过调整nmax和lmax以及稀疏表示和正则化参数来变化描述符的分辨率，并保留了仍然能再现dft参考数据并且在考虑的结构和温度范围内表现稳定的最简单模型。这一点很重要，因为对于基于核的模型（如gap），可转移性不仅由描述符的大小决定，还受到描述符定义、稀疏化、正则化以及训练数据库中代表的局部环境多样性的综合影响。

特别注意了soap的分辨率。采用的值nmax=8和lmax=9，被发现能够提供足够的径向和角度细节，以区分这种玻璃中存在的扭曲和化学多样化的局部基元，包括nao环境以及vo和teo多面体。在试验模型中，较低的值导致这些环境之间的区分能力明显降低，力预测也更差。尽管更大的稀疏集可以稍微提高固定描述符的回归质量，但它无法完全恢复从分辨率不足的soap基础中缺失的信息。因此，我们选择了一个具有足够内在分辨率的描述符，然后使用稀疏表示来控制回归的有效复杂性。将nmax和lmax增加到选定的值之外并没有带来任何有意义的改进，所以最终选择代表了准确性と效率之间的实际折中。

截止半径rcut=4.27 ?和高斯模糊宽度σatom=0.6 ?被选中，以提供平滑的描述符，同时保留玻璃网络的第一层配位壳和最相关的近邻角相关性。截止过渡宽度rδ=1.0 ?确保了描述符在局部环境边界处的平滑衰减。soap能量尺度参数δ=0.88ev设置了由soap贡献表示的局部能量变化的特征幅度。稀疏集大小nsparse=1200用于soap，对于两体项则为70，一旦测试误差不再随进一步增加而显著改善，同时仍允许以合理的成本进行长时间分子动力学模拟。对于两体贡献，高斯长度尺度固定为θ=1.0 ?，提供了一个平滑的短程排斥基线。最后，正则化参数σvenergy=0.001 ev/原子、σvforce=0.05 ev/?和σvvirial=0.05 ev/原子被选中，以保持与dft能量景观的紧密一致，同时避免过度拟合fpmd快照中存在的有限温度噪声。

gap mlip模型的训练数据是从t=2000 k到t=300 k的宽温度范围内进行的fpmd模拟中提取的，以确保结构多样性。具体来说，从每个模拟平台采样了70个快照，总共得到490个配置。为了提高数据集的保真度，使用更高的平面波截止值1000 ry，重新计算了相关的dft能量、力和维里应力。这个数据集在训练集和测试集之间平均分配。我们注意到，基于gap的mlip即使在小型训练数据集上也能够达到化学准确性，通常只需要几百个配置就可以处理无序材料[94]。这种数据效率来源于gap的基于核的性质，它能够从稀疏的、有代表性的原子环境中进行稳健的外推。得到的gap模型表现出强大的预测准确性，在测试集上力学的rmse值为0.38 ev/?（见图1），能量的rmse值为2.91 mev/原子，完全在常见的5 mev/原子阈值[98]、[99]范围内。

为了评估mlip强再现dft计算的能量和力的能力是否也扩展到准确的结构描述，我们在不同大小的系统上进行了模拟。使用一个432原子的模型（mlip_s）与fpmd结果进行了直接比较，同时还构建了一个包含多达3456个原子的更大系统（mlip_b），以研究fpmd模拟计算上无法获得的潜在大小效应。对于所有模拟，模型都是用随机初始化的结构和实验上一致的密度（3.514 g/cm^3）构建的。分子动力学模拟是在lammps代码[122]中实现的，采用nvt系综和0.5 fs的时间步长。所有系统都经历了多个温度平台的热退火：在t=300 k、t=500 k、t=750 k和t=1000 k时各5 ps，然后在t=2000 k、1800 k、1500 k、1200 k、900 k和300 k时各100 ps。结构性质是从t=300 k时轨迹的最后50 ps中分析的。

通过分析整个系统和每个组成元素（na、o、te、v）的振动密度态（vdos）来比较fpmd和mlip模拟的振动性质。vdos是从速度自相关函数（vacf）的傅里叶变换计算得到的。沿md轨迹的原子速度得到的vacf定义为c(t)=1/n∑i=1nvi(0)?vi(t)，其中n是原子数量，vi(t)是时间t时原子i的速度。然后vdos由g(ω)∝∫0∞c(t)cos(ωt)dt给出，或者等效地通过从轨迹评估的标准化vacf的离散傅里叶变换得到。如图2所示，mlip与fpmd结果非常一致，捕捉了所有元素和频率范围内的vdos的关键特征。特别是，主要峰和整体谱形都被很好地再现了，证实了mlip准确模拟了玻璃的振动响应。这验证了mlip在进一步的玻璃结构、能量和动力学研究中的可靠性。

作为对获得的mlip的额外验证，我们使用我们的mlip评估了体积模量，并将结果与fpmd和相关的碲酸盐-钒酸盐玻璃的实验数据进行比较。

下载：下载高分辨率图片（456kb）
下载：下载全尺寸图片

图2. 玻璃的振动密度态，比较了fpmd（黑色）和mlip（红色）模拟。总vdos及其元素分辨的贡献（na、o、te、v）是从速度自相关函数的傅里叶变换计算得到的。体积模量是通过在432原子模型上的10个独立非晶快照上执行的能量-体积计算获得的。对于每个快照，通过将晶格参数围绕参考值变化±3%的步骤为±0.5%来生成一组均匀缩放的立方晶胞，每个快照产生13个能量-体积点。在dft和mlip级别应用了相同的各向同性缩放协议，每个e(v)数据集都用birch-murnaghan状态方程[123]进行拟合，以提取平衡体积模量b0。报告的值对应于10个快照的平均值，而引用的不确定性表示相关的标准偏差。

使用这个程序，我们得到b0fpmd=24.85±1.47gpa和b0mlip=27.88±4.82gpa。这些结果表明，从fpmd和mlip计算得到的平均体积模量在估计的标准偏差范围内几乎相当。此外，考虑到弹性响应的固有快照间变异性，这种一致性对于目前的非晶氧化物玻璃来说是合理的。获得的值也与相关碲酸盐-钒酸盐玻璃的实验范围一致。在二元（1?x）teo2-xv2o5玻璃中，室温下的超声测量报告的体积模量根据成分的不同从大约27.2到37.3gpa不等，纯teo2的体积模量为31.7gpa[20]、[124]。在更复杂的teo2 bi2o3 v2o5 na2o tio2玻璃系统中，实验报告的体积模量范围为16.4–25.2gpa[125]。钠钒酸盐玻璃（成分xna2o (35?x)v2o5 65teo2）也显示出成分依赖的弹性模量，随着钠含量的增加而软化[126]。因此，这里得到的fpmd和mlip体积模量都在相关碲酸盐-钒酸盐氧化物玻璃的更广泛实验窗口内，这支持了计算值的物理合理性。

最后，我们想强调的是，当前的mlip是在对应于单一密度的配置上训练的，这本质上限制了其在更宽密度范围内定量再现体积依赖性质（如体积模量）的能力。从这个意义上说，体积模量是对模型特别严格的测试，因为它探测了能量景观相对于各向同形压缩和膨胀的曲率。因此，mlip和dft体积模量之间的一致性强调了我们潜在模型足够准确地描述潜在能量表面的能力。

如前所述，玻璃中的钒可以采取多种氧化态，从v2+到v5+[127]、[128]、[129]不等，这取决于成分和加工条件。v5+通常是主导的，即使在非有意还原条件下也常常存在不可忽视的v4+。这些状态与不同的局部原子环境相关[127]、[128]、[130]。>从计算角度来看，传统的分子动力学（MD）必须使用预定义的指标来表示钒的物种，并且力场参数用于区分不同的氧化态。相比之下，在基于BOMD和DFT的模拟中，钒的物种分布是根据电子结构自动生成的，没有任何人为的限制，这一点在之前对NVP玻璃的研究中已经得到证实[69][70][71]，其中氧化态由组成和局部环境决定。另外，基于MLIP的MD可以在同一组成下提供具有FPMD级别准确性的0.1 0.6 0.3玻璃结构模型，但这种情况下只考虑了V2O5单元。然而，由于V5+和V4+会促进(n = 4–6)多面体的重叠分布，尽管它们的相对数量不同，并且VO和VO键长也非常接近[71]，因此MLIP模型可以有效地探索V5+/V4+比例变化对玻璃结构和性质的影响，从而以较低的成本实现结构和性质的优化。为了评估这一可能性，所开发的MLIP模型也被用来制备不同V2O4浓度的玻璃系统（z = 0.0、0.05、0.10和0.15）。这一范围允许我们进一步分析V2O4的影响，超出了之前FPMD模拟的范围。这些玻璃系统是通过保持与原始NVT60玻璃相同的密度，并采用类似于小系统和大型NVT60模型的熔融-淬火热处理策略生成的。轨迹的最后50皮秒在T=300 K时被用于统计分析，从而能够有效地研究这些玻璃内部的结构和组成变化。

图3. 在T=300 K时，0.1 0.6 0.3玻璃系统的总X射线结构因子（SX(q)。实验得到的X射线结构因子（黑色虚线）与通过实空间中成对相关函数的傅里叶变换计算得到的FPMD SX(q)（红线）以及分别含有432和3456个原子的MLIP玻璃系统（ green线和blue线）进行了比较。

3. 0.1 0.6 0.3玻璃的结构

3.1 模型验证
我们模型的质量是通过将计算得到的结构函数SX(q)和成对分布函数（PDF）GX(r)与参考文献[111]中合成的、具有相同组成的玻璃的实验数据进行比较来评估的。实验PDF是在专门的实验室设备上获得的，该设备能够测量高达Qmax值21.796 ?^-1的SX(q)[111]。
计算得到的PDF通常如下获得：(8) GX(r) ≈ ∑α=1... ∑β=1... cαcβfα(q0)fβ(q0)|〈f(q0)〉2|[gαβ(r)?1]。其中，α和β表示原子种类（Te、Na、V或O），cα是原子浓度，fα(q)是在q0=0时评估的原子X射线散射因子，gαβ(r)是部分成对相关函数。尽管这种近似方法能够捕捉到G(r)的一般形式，但它对第一个峰的强度描述并不完美。为了解决这个问题，我们采用了参考文献[131]中更准确的公式，其中PDF是通过修改后的部分PDF的加权和来计算的。
图3显示了NVT60玻璃系统的计算X射线总结构函数SX(q)与实验数据之间的比较。实验得到的SX(q)在1.90 ?^-1处显示出一个主要峰，随后在3.95 ?^-1、5.5 ?^-1和7.75 ?^-1处出现较宽且强度较低的峰，然后在10 ?^-1附近有衰减振荡。FPMD和MLIP计算得到的SX(q)都在整个q范围内再现了实验峰的位置和强度，并且具有典型的统计波动。图4显示了来自FPMD和MLIP模型的计算总PDF与实验数据的比较，表1总结了从FPMD和MLIP模型得到的PDF峰位置，并与实验结果进行了对比。实验PDF显示出大约在1.9 ?、3.5 ?和4.4 ?处的三个主要峰。第一个主要峰实际上由两个重叠的峰组成，分别对应于典型的短双键/钒酰VO和长单键VO（1.58–1.83 ?）以及TeO（1.89–2.16 ?）的键长[132][133]。我们还注意到在2.63 ?附近有一个小峰，这可能与NaO和OO之间的距离有关[133][134]。第二个主要峰的范围是3.04–3.78 ?，第三个峰的范围是3.95–4.90 ?，这些峰主要代表了第二和第三邻居之间的相互作用，这些相互作用基本上归因于阳离子之间的相关性。

表1. NVT60玻璃系统的总X射线成对分布函数G(r)中中心峰位置的比较：FPMD和MLIP模型计算值与实验数据对比。对于这些模型，还添加了作为部分PDF gαβ(r)中第一个最大值位置的键长rαβ。表中还包含了基于MLWF方案的FPMD和MLIP模型计算得到的配位数及其相应的误差范围。我们注意到(rcut)是相对于系统中原子β的数量给出的。括号中的值是基于将PDF积分到第一个最小值所得到的配位数。

下载：高分辨率图像（173KB）
下载：全尺寸图像

图4. 在T=300 K时，0.1 0.6 0.3玻璃系统的总X射线成对分布函数。FPMD计算得到的PDF（红线）与实验结果（黑色虚线）以及含有432和3456个原子的MLIP玻璃系统（绿色和蓝色线）进行了比较。

表1. NVT60玻璃系统的总X射线成对分布函数G(r)中中心峰位置的比较：FPMD和MLIP模型计算值与实验数据对比。键长rαβ也被列为部分PDF gαβ(r)中第一个最大值的位置。表中还包含了基于MLWF方案的FPMD和MLIP模型计算得到的配位数及其相应的误差范围。我们注意到(rcut)是相对于系统中原子β的数量给出的。括号中的值是基于将PDF积分到第一个最小值所得到的配位数。此外，在FPMD模型中TeO键长略有延长，在MLIP模型中则被进一步高估。表2显示了基于MLWF形式主义计算得到的l-fold百分比，其中Te的rcut = 2.46 ?，V的rcut = 2.50 ?，以及NVT60玻璃中Na的rcut = 2.93 ?，这些数据来源于FPMD和MLIP模型。l-fold值小于约2%的未予报告。括号内的数值是基于积分PDF至第一个最小值所得到的配位数。αl（%）分别为：FPMD为19.9，MLIP_S为13.8，MLIP_B为10.9。

玻璃态系统的局部环境可以通过配位数nαβ来描述，实际计算时是指在截止距离rcut内围绕物种α的物种β的数量。特别是，平均配位数可以通过积分gαβ(r)至其第一个最小值来估算，这标记了第一配位壳层的边界。然而，对于TeO对来说，这种方法并不直接，因为(r)没有明确的最低点。为了解决这个问题，我们采用了第2.1节中解释的MLWF形式主义。应用到小系统时，通过将TeWB和OWB的平均距离相加，得到了TeO的最佳键长截止值为2.46 ?，这与我们之前的研究[115]、[116]一致。对于VO距离，截止值设定为2.50 ?，对应于(r)中观察到的第一个最小值（见图5）。同时，对于NaO键，选择了2.93 ?的截止距离，这与Na2Te4O9和NaVTeO5晶态系统报告的最长NaO键长相一致[133]、[134]。

关注Te的平均配位数，通过积分部分(r)至接近2.46 ?，我们发现平均TeO配位数从FPMD模型中的约4.26增加到MLIP_S模型中的4.82，然后在大型MLIP_B模型中略微下降至4.77。这一趋势直接反映在表2中呈现的l-fold分布上。在FPMD模拟中，网络主要由TeO4单元组成（49.6%），还有相当比例的TeO5（33.6%）和TeO3（13.8%）。相比之下，MLIP_S模型显示出向更高配位的显著转变，其中TeO5单元增加到62.6%，TeO6单元出现为11.0%，同时TeO4减少到22.4%，TeO3的存在最少，仅为3.6%。这一趋势在大型MLIP_B模型中持续存在，但过配位的程度稍低：TeO5占主导地位，为58%，TeO6占10%，而TeO4增加到28.3%，TeO3几乎可以忽略不计（3%）。Te配位数的增加与MLIP_S和MLIP_B中观察到的更强烈且略微偏移的(r)峰值一致，这些峰值出现在约2.04 ?，而FPMD中为2.00 ?。

在FPMD和MLIP_S系统中，MLWF分析通过明确计算成键电子对提供了更准确的化学描述。我们观察到FPMD的绝对值略低，约为3.97，而MLIP_S为4.38。这种相对于基于PDF的积分的减少是因为MLWF自然排除了在积分截止范围内仍然有贡献的较长或相互作用较弱的TeO对。尽管存在这种差异，趋势仍然保持一致，MLIP_S再次显示出比FPMD更高的Te配位数。表1显示，MLIP_S中增强的Te配位数主要涉及桥接氧原子（TeBO），使得?4.06，而FPMD为3.67，这表明网络更加互连，而与非桥接氧原子（TeNBO）键合的Te比例在两种模型中都保持较低（?0.30）。我们注意到，纯TeO2玻璃的平均TeO配位数通常在约3.6至3.9之间，主要由TeO4单元构成，其中包含一些三角形的TeO3[115]、[116]、[138]、[140]、[141]。在我们的系统中添加高浓度的V2O5已知会略微增加Te的配位数，这与二元和三元碲酸盐-钒酸盐玻璃报告的值（约4.0–4.12）一致[142]、[143]。因此，我们的FPMD值完全在这个范围内，而MLIP模型略微高估了Te的配位数，这强调了捕捉Te4+上立体化学活性孤对结构效应的持续困难。

NVT60玻璃中的钒局部环境显示出主要由四面体和金字塔单元构成的结构特征。根据表2中的l-fold统计，FPMD模型显示约66%的V原子处于四面体VO4结构中，约30%处于VO5结构中，只有4%的VO6单元。MLIP模型紧密再现了这一分布，显示出约62%–63%的VO4和约28%–31%的VO5，反映了ML势能通常略微偏好更高配位的倾向。V的配位数从FPMD中的4.38增加到MLIP模型中的4.44，与VO5比例的小幅增加一致。关键的是，将(r)分解为桥接和非桥接氧原子的贡献后，发现VNBO的贡献在FPMD中几乎保持不变，约为0.84，而在两个MLIP模型中为0.80。这些发现与类似玻璃的文献结果一致，这些文献一致报告VO的配位数在4.1–5.4之间，主要由VO4单元构成，同时含有相当比例的VO5单元[138]、[142]、[144]、[145]。

从MLIP_B和FPMD模型获得的配位数相同，均为6.05，而MLIP_S显示的值略高，约为6.25，这是由于形成了较少比例的NaO9（2.6%）和更多比例的NaO8单元（约8%）。FPMD和MLIP_B显示更紧凑的Na分布，以NaO6和NaO7为中心。进一步分析Na与BO和NBO的键合情况，如表1所示，FPMD给出约3.3的Na配位数，而MLIP_S显示约2.8的Na配位数，表明Na与NBO和BO的相互作用大致相等。MLIP_S高估了（约为3.4），与其相对于FPMD更大的平均配位数一致。在类似玻璃中 also 发现了类似的广泛Na环境，其中Na通常与约6个氧原子配位，有利于离子传输。在Na-磷酸盐-钒酸盐玻璃中，NMR和建模研究报道Na的配位数在5.5–6.9之间，反映了适应局部电荷需求[81]、[146]、[147]。总体而言，我们的结果证实了Na修饰剂的作用，保持了在BO和NBO之间平衡的约6倍配位数。

综合来看，我们的结果表明，在测试范围内（432–3456个原子），结构属性（包括总和对分布函数及配位数）对系统大小的敏感性非常小。不同大小的FPMD和MLIP模型之间几乎完全一致的结果证实了熔化-淬火协议和局部结构特征（例如VO4、TeO4单元）主要受化学因素而非系统大小控制，这在其典型的统计波动范围内是如此。

下载：下载高分辨率图像（403KB）
下载：下载全尺寸图像

图6. NVT60玻璃态系统中Te、V和O原子的键角分布。我们使用与计算配位数时相同的截止距离来计算键角分布。

3.4. 键角分布
图6展示了三种模型的OVO、VOV、TeOTe、OTeO和TeOV键角分布（BADs）。其他BADs在补充图S2中提供。OVO BAD显示出复杂的分布，具有三个显著特征：第一个峰值大约在86.8°，对应于有缺陷的八面体和金字塔几何结构；一个尖锐而强烈的峰值位于109°，代表近乎完美的四面体配位；以及一个较宽且强度较低的峰值位于约165°，表明存在有缺陷的正方形金字塔或高度扭曲的三棱柱形结构。重要的是，这种分布在三种计算模型（FPMD、MLIP_S和MLIP_B）中都得到了一致的再现，表明所有模型都能稳健地捕捉到V周围的局部角环境。关注FPMD的分布，将OVO角分布分解为VO4、VO5和VO6多面体的贡献，如图S3所示，我们观察到109°的尖锐峰值主要由VO4四面体引起，而85–90°的肩部主要来自VO5和VO6单元。165°的宽峰也与较少见的扭曲VO5/VO6连接相关。在纯V2O5玻璃[136]、[142]中已经显示，VO4四面体占主导，导致角度约为109°的典型四面体角度。在二元玻璃中，添加碲氧化物并不会显著改变四面体峰值，但会引入更多扭曲的几何结构。在我们之前关于含钠钒酸盐玻璃的研究中，我们报告了类似的OVO角分布，VO4的尖锐峰值在108–110°之间，而VO5和VO6的较低角度肩部低于90°，这受到钠作为电荷补偿剂和网络修饰剂的作用[71]、[81]。

VOV键角分布进一步描述了钒酸盐单元的中程有序性和互连性。FPMD得到的NVT60中的VOV分布显示出一个以约135°为中心的宽主峰，以及一个约为90°的小的肩部。约135°的峰值表明VO4和VO5单元之间主要的角共享[55]、[136]、[142]，与开放的聚合玻璃网络一致。同时，约90°的肩部暗示了较少的边共享VOV连接[55]、[136]。在MLIP_S和MLIP_B模型中，VOV分布显示出更宽的特征，峰值向较低角度移动，强度增加，表明更扭曲和紧凑的多面体连接。

至于OTeO BADs（图6），我们观察到三种结构模型中的类似特征：一个主要在87°附近的高强度峰值，以及一个更宽、强度较低的峰值在165°附近。第一个峰值很好地对应于在扭曲的TeO4双锥体或有缺陷的八面体单元内的OTeO equatorial 连接。这些角度通常在约70°至110°之间，与纯无定形TeO2报告的角度相符[116]、[144]、[148]。第二个较宽的特征在150°至180°之间，以165°为中心，归因于轴向的OTeO键，反映了TeO4或TeO5多面体沿一个轴的伸长[116]。

转向TeOTe角度，图6显示出三种模型之间的差异。在FPMD模型中，TeOTe BAD显示出两个主要峰值，分别位于约110.8°和约128.5°，表明多面体之间的连接程度不同。相比之下，MLIP_S模型显示出三个峰值，分别位于约102°、119°和133°，而在MLIP_B模型中它们变成了两个较宽的峰值，分别位于约102°和118°。MLIP_S和MLIP_B的整体略微更宽的分布表明角无序程度略有增加。这与MLIP_S和MLIP_B中五倍Te配位数比例较高的情况相符（见表2），这内在地改变了桥接角度以适应增加的局部堆积。将TeOTe角度分解为不同碲配位单元的贡献可以更清楚地看出这一点（见图S5）。在FPMD模型中，第一个尖锐峰值主要由和贡献引起，两者在约110–112°显示出强烈的峰值，并在约125–130°处有显著的强度最大值，特别加强了总TeOTe分布中的第二个峰值。我们还看到在这里和也起着明显的作用，达到约111°的最大值，并有助于BAD的更尖锐的第一个峰值。值得注意的是，FPMD几乎没有连接，反映了五倍Te单元较少（见表2）。在MLIP_S模型中，的贡献显示出显著減弱的强度，并且向较小的角度约102°移动。它还显示出一个从约119°延伸到约140°的宽肩部，有助于总Te-O-Te BAD中的第二个和第三个最大值，并表明在MLIP_S中TeO4单元采用了更多样的桥接几何结构。此外，与FPMD不同，在MLIP_S中，和都显示出两个峰值，分别位于约119°和约133°，直接导致了总TeOTe分布中的第二个和第三个峰值。总体效果是MLIP_S显示出较小的、强度较低的第一个峰值，约为102°，与FPMD的尖锐峰值约111°相比。

最后，观察TeOV键角分布（图6），我们看到FPMD显示出一个以约130°为中心的强烈峰值，分布范围从90°到165°。MLIP_S和MLIP_B具有相似的形状但更宽的分布，同时在130°附近的峰值强度较低，而在约95°附近的强度略高。这表明FPMD展示了一个更加有序的网络，具有明确定义的TeOV桥接，并暗示了TeO4（或略微畸变的TeO5）与VO4或VO5单元之间的稳定角共享连接，而MLIP模型则捕捉到了更多种类的角度环境，反映了更高的无序性。MLIP_S和MLIP_B中更广泛的角度分布主要是由于它们含有更高比例的过协调TeO5和VO6单元，正如之前在表2中所看到的。这种过协调影响了TeOV连接，使得角度向后移动到95°，并减少了排列整齐的桥接数量。

下载：下载高分辨率图像（583KB）
下载：下载全尺寸图像

图7. 玻璃模型片段中的电子局域化函数表示。紫色、棕色、绿色、红色和白色分别代表Na、V、Te、O和H。我们注意到，这里的H原子仅用于钝化提取片段的末端O原子。

4. 电子局域化函数分析
从离子和共价贡献的角度分析化学键合，可以进一步洞察玻璃网络的连通性。为了表征这一方面，我们在玻璃基体的一个代表性片段上计算了电子局域化函数（ELF）。ELF的值范围是从0到1。接近1的值表示电子对的高度局域化，例如孤对或共价键合区域，而0.5的值对应于均匀电子气的参考状态。在不均匀系统中，大约0.5的值通常表示中等程度的局域化，而接近0的值则对应于电子密度非常低的区域。因此，ELF提供了电子结构的实空间描述，特别有助于区分共价键、离子键和非局域化键的特征。

在实际操作中，我们提取了一个包含四个VO4单元和一个TeO5单元的片段，这些单元围绕一个Na原子。为了防止悬挂键，片段边界的末端氧原子被H原子覆盖。然后，在保持Te、V和O原子固定的情况下，对H原子的位置进行了PBE理论水平的放松。随后，在PBE0理论水平上计算了ELF。图7显示了这个结果片段以及ELF在片段平面上的投影。

ELF分析揭示了无定形片段中Te、V和Na的局部键合特征存在明显差异。TeO键合区域的ELF值约为0.5，表明电子密度在该区域有中等程度的积累和局域化。这与极性共价键的特性一致，其中部分电子密度在Te和O之间共享，这与通常报道的碲酸盐环境中的情况相符[115]、[149]。至于VO键合区域，其ELF值较低，约为0.3，表明电子密度在键合区域的局域化程度较低。这同样表明了极性共价键的特性，但具有更大的离子贡献，电子密度更倾向于集中在氧原子上，而不是作为共享的键合对，这与氧化物系统中VO键的一般描述一致[150]。因此，我们发现在这种情况下TeO键的共价特性比VO键更强。

Na的表现则明显不同。Na原子周围的ELF值非常低，通常为0.2或更低，看起来基本上与邻近原子隔离。在NaO区域没有观察到显著的ELF积累，表明Na与周围氧化物框架之间没有局域化和方向性的键合。这一预期结果表明NaO键主要是离子性质的，电荷主要转移到了氧化物网络中。这种行为与Na作为网络修饰剂的作用一致，与直接参与玻璃键合框架的Te和V不同。

下载：下载高分辨率图像（347KB）
下载：下载全尺寸图像

图8. 不同z值的玻璃系统的总X射线对相关函数G(r)：0.0、0.05、0.10和0.15，使用MLIP计算，并与z = 0.0时FPMD计算的G(r)（红线）以及实验结果（黑色虚线）进行比较。

5. ML-GAP处理氧化态能力的扩展
通过使用MLIP模型，我们可以超越FPMD模拟在计算上的可行性，探索一系列玻璃组成（z = 0.05、0.10和0.15），以研究混合钒氧化态如何影响玻璃结构。图8给出了这三个系统的总对分布函数与实验数据和FPMD结果（z = 0.0）的比较。z = 0.05的MLIP模型与实验数据的一致性稍好，而更高的V2O4比例（例如z = 0.15）则导致显著差异。这一结果与合成样品中存在少量V2O4的实验证据相符[111]。

5.1. 基于MLIP的钒物种鉴定
为了确定模型中的钒氧化态，我们遵循了参考文献[71]中概述的方法。使用PBE0理论水平的自旋极化单点计算来确定钒的物种[151]。每个V中心的局部自旋种群定义为σspin=Nα?Nβ，用于量化未配对电子。这里，α和β表示与给定V相关的电子种群。对于钒，V5+ ([Ar]3d0) 预期显示σspin≈0，V4+ ([Ar]3d1) 显示σspin≈1，V3+ ([Ar]3d2) 显示σspin≈2，尽管共价相互作用和自旋离域可能会导致小的偏差。计算是针对三个系统（z=0.05、0.10和0.15）进行的，使用的是从T=300 K轨迹中提取的五个快照。

图9(a)中绘制的自旋密度σspin的分布显示了两个峰值，其中一个峰值约为0.5，表明沿键合区域的电子密度有一定程度的积累和局域化。这与极性共价键的特性一致，其中部分电子密度在Te和O之间共享，这与通常报道的碲酸盐环境中的情况一致[115]、[149]。至于VO键合区域，其ELF值较低，约为0.3，表明电子密度在键合区域的局域化程度较低。这也表明了一种极性共价键，但具有更大的离子贡献，电子密度更倾向于集中在氧原子上，而不是作为共享的键合对，这与氧化物系统中VO键的一般描述一致[150]。因此，我们发现在这种无定形环境中，TeO键的共价特性比VO键更强。

Na的表现则明显不同。Na原子周围的ELF值非常低，通常为0.2或更低，基本上与邻近原子隔离。在NaO区域没有观察到显著的ELF积累，表明Na与周围氧化物框架之间没有局域化和方向性的键合。这一预期结果表明NaO键主要是离子性质的，电荷主要转移到了氧化物网络中。这种行为与Na作为网络修饰剂的作用一致，与直接参与玻璃键合框架的Te和V不同。

下载：下载高分辨率图像（417KB）
下载：下载全尺寸图像

图9. 在PBE0理论水平上计算的每个V位点的自旋密度σspin（以e为单位），对于z = 0.05、0.10和0.15的玻璃系统。结果以(a)直方图分布，(b)原始数据（每个模拟模型在不同配置下的单个数据点分布）的形式呈现，对于从MLIP模型获得的结果，以及(c)PBE放松后的MLIP模型的原始数据。(d)在z = 0.15模型上使用不同放松方案计算的自旋密度σspin。

为了克服这一限制，我们首先在PBE0理论水平上重新计算了五个PBE放松后的MLIP结构的自旋种群。图9(c)中显示的结果揭示了两个定义清晰的V5+和V4+ σspin位点，其自旋值分别接近0和1，这与之前关于钒磷酸盐玻璃的研究一致[71]。使用σspin=0.5的截止值确定了V4+和V5+的比例，并在表3中报告。我们观察到随着z值的增加，V5+的比例减少，而V4+的比例增加。值得注意的是，根据上述分析，在与NVT60样品的实验PDF最一致的情况下获得的模型中，V4+的比例约为8.54%，与实验测量结果（约11%）非常接近[111]。

通过对结构在PBE0理论水平上进行放松，可以实现钒物种鉴定的进一步精细化。这可以通过两种不同的方式实现：一种是直接且成本较高的从MLIP结构开始的PBE0放松；另一种是成本较低的PBE放松后跟随PBE0放松（PBE+PBE0）。这两种放松方案都应用于z = 0.15系统的最终配置，自旋密度始终在PBE0理论水平上计算。图9(d)展示了结果。在两种情况下，我们都观察到由两种自旋密度a组成的钒位点之间的清晰分离，这两种自旋密度围绕0和1紧密分布，证明了PBE0在达到高精度方面的优越性，并且尽管计算成本有很大差异，这两种放松方案是等效的。

表3. 每个系统中V2O5和V2O4的数量。在z = 0.05、0.10和0.15的情况下，PBE放松后的玻璃系统中V5+和V4+的比例，以及σspin = 0.5的截止值。

这些发现表明，尽管MLIP仅在V5+环境下进行训练，但在适当的氧含量驱动下，它能够生成结构一致的玻璃模型，这些模型与不同的钒氧化态兼容。虽然MLIP本身不显式区分不同的钒氧化态，但在DFT水平上对MLIP生成的配置进行后处理可以识别它们并进行准确的结构表征。然后可以通过一次性DFT放松MLIP衍生的结构来获得电子属性。这种方法允许高效地调整系统组成，以非常低的计算成本实现与实验观察结果的更紧密匹配，特别是在事先不知道精确钒物种的情况下。

5.2. 混合价态钒环境的结构特征
我们在这里关注V2O4含量最高的PBE+PBE0放松MLIP模型（z = 0.15），以突出玻璃网络对引入混合价态钒物种的结构响应。图10绘制了获得的VO部分PDF的首个峰值，这些峰值对钒的氧化态最为敏感，并与未经放松的MLIP模型的VO部分进行了比较。未经放松的MLIP模型的gVO(r)在1.62和1.76 ?处显示两个峰值，并且在较大的r值处有一个中等程度的尾部，如前所述。对于放松后的PBE+PBE0模型的gVO(r)，观察到两个峰值向较小距离方向小幅移动了约0.03 ?，并且在1.94 ?处出现了明显的肩部。我们注意到，未经放松模型的主要峰值分裂成两个贡献，分别位于1.69和1.77 ?，这显然是由于使用了单个放松快照导致的统计有限。为了理解这些特征的起源，我们将PBE0下放松模型的总gVO分解为V5+和V4+物种的贡献（见图S6）。可以观察到，两个主要贡献来自V5+单元，第一个峰值的小贡献来自V4+单元，而1.94 ?处的肩部则来自V4+单元。

下载：下载高分辨率图像（173KB）
下载：下载全尺寸图像

图10. z = 0.15模型的VO部分对分布函数：通过MLIP结构（未经放松）、用PBE放松（五个放松快照的平均值）以及用PBE后跟PBE0放松（PBE+PBE0，一个快照）获得的结果。

这表明，在PBE+PBE0理论水平上的放松可以清楚地区分V4+环境和V5+环境。这一结果与之前观察到的良好分辨的自旋种群一致。通过查看仅使用PBE放松获得的中间模型的gVO，可以进一步评估局部环境描述的质量与相应自旋值之间的联系。我们清楚地观察到1.9 ?处出现的肩部来自V4+单元。尽管如此，主要峰值的位置与未经放松的MLIP结构相同，导致肩部的分辨率较低。虽然PBE放松后的结构比MLIP作为获得的模型有所改进，但自旋种群仍然略偏离理想的0和1值（见图9c）。图11展示了V5+和V4+的l-fold（l = 4、5和6）单元的分布。我们观察到V5+基本上是四重配位的（83.3%），而V4+主要是五重配位的（80.3%）。这一观察与我们之前关于NVP玻璃的结果完全一致，其中较高的钒氧化态与较低的平均配位数相关。与在NVT60（未使用V4+）上获得的FPMD和MLIP_S结果相比，我们注意到四重单元的比例显著增加（在NVT60中为62%），而五重和六重单元的比例减少，导致六重单元的消失。下载：下载高分辨率图像（94KB）下载：下载全尺寸图像图11. MLIP模型中z = 0.15时V5+和V4+的l重V单元比例，按每种氧化态的原子数进行了归一化。这种环境在晶体相中观察到，其中要么存在Na[152]的混合V价态，要么在Te和Na氧化物[133]的存在下只有V5+环境。观察四重环境中的V5+，有两个短VO键长为1.64 ?和两个较长键长为1.81 ? [133]。至于五重环境中的V4+，有一个短VO键长约1.65 ?和四个较长键，平均长度为1.95 ? [152]。这些观察到的距离与部分gVO中观察到的峰值一致，特别是由于五重V4+环境造成的肩部。表4. 全部V原子群的Mulliken电荷（平均值±标准差），以及被分类为V5+-like和V4+-like的位点的Mulliken电荷，同时还计算了考虑半径2.5 ?范围内所有O原子时的平均V-O距离。z全部VV5+-likeV4+-like〈dV5+?O〉（?）〈dV4+?O〉（?）0.000.96±0.070.96±0.07–1.78–0.050.95±0.080.95±0.080.95±0.071.771.880.100.95±0.080.96±0.070.90±0.071.751.870.150.96±0.080.97±0.070.94±0.091.751.865.3. 关于过渡金属氧化物中的电荷自我调节为了进一步研究V形式氧化态的作用，我们分析了被分类为V5+-like和V4+-like的代表性局部环境中V位点的Mulliken种群及其自旋种群。尽管Mulliken电荷不应被解释为形式氧化态，但它们为在同一计算框架内比较相对电荷趋势提供了统一的度量标准。我们发现，局部在V上的Mulliken电荷在这两类环境之间只有轻微的变化。在所有组成中，表4显示全部V原子群的平均电荷几乎保持在0.95–0.96之间，而V5+-like和V4+-like位点之间的差异仅几百分之一电子。此外，检查电荷分布本身（见图12）显示这些分布有很强的重叠性，并且随着V4+-like和V5+-like位点比例的变化，它们的形状变化很小。这表明这两种环境之间的区别不是由明显分离的原子电荷来捕捉的，而是由适度的电荷变化和更明显的结构重排来实现的。相比之下，局部几何结构显示了更清晰的区别，其中V5+-like位点的平均VO距离较短，为1.75–1.78 ?，而V4+-like位点的平均VO距离较长，为1.86–1.88 ?。因此，这两种环境在结构上是不同的，但它们的电荷有很强的重叠性，远比纯粹的离子图像所暗示的更为接近。这种行为与过渡金属氧化物中的电荷自我调节[153]一致，在这种环境下，形式氧化态的变化只部分反映在金属中心上的电荷上，因为周围的氧网络适应了部分电子的重新分布。在这个背景下，MLIP的可转移性可能是由于它是由局部原子-环境描述符构建的，因此与不同形式氧化态分配相关的环境不需要对应于明显分离的局部电荷。下载：下载高分辨率图像（235KB）下载：下载全尺寸图像图12. 模型中所有V位点（总计）以及被分类为V5+-like和V4+-like的位点的Mulliken电荷分布。6. 结论本研究使用第一性原理分子动力学和机器学习原子间势研究了0.1 0.6 0.3玻璃态系统的原子尺度结构。这两种方法都通过实验数据进行了验证，包括X射线总结构因子和对分布函数，显示出良好的一致性，并表明它们在捕捉所研究玻璃态系统的基本结构特征方面的可靠性。MLIP使得在更大系统和更长时间尺度上进行高效模拟成为可能。结构分析显示，这两种方法都能再现关键特征，如总和对分布函数、键角分布以及局部配位环境，包括桥接和非桥接氧原子的分布。特别强调了TeO和VO的配位，MLIP的配位数和角度分布与FPMD和先前报道的实验结果非常吻合，捕捉到了基于碲酸盐的玻璃的混合网络结构和灵活性。MLIP通过引入受控的V5+/V4+比例来研究氧化态效应。模拟显示，V4+浓度的增加导致VO配位环境的变化，有利于形成五重单元，并扩展了VO键长分布。VO部分PDF的分解和配位数分析表明，V4+主要形成多面体，而V5+则有利于形成配位单元，这与实验观察到的混合价态钒玻璃的结果相符。总之，这项工作表明，仅基于V5+参考数据训练的MLIP可以外推到含有V4+的环境中，且准确度损失最小。这种能力可以通过过渡金属氧化物中的局部电荷自我调节机制来解释。这一发现为高效的大规模氧化还原活性玻璃态系统模拟打开了大门，否则这些系统将无法使用第一性原理方法进行计算。CRediT作者贡献声明Firas Shuaib：撰写 – 原始草稿，软件，研究，形式分析，数据整理。Remi Piotrowski：研究。Gaelle Delaizir：研究。Pierre-Marie Geffroy：研究。Steve Dave Wansi Wendji：研究。Carlo Massobrio：撰写 – 审查与编辑，研究。Mauro Boero：撰写 – 审查与编辑，研究。Guido Ori：撰写 – 审查与编辑，研究。Philippe Thomas：撰写 – 审查与编辑，研究。Olivier Masson：撰写 – 审查与编辑，监督，方法学，研究，资金获取，概念化。Assil Bouzid：撰写 – 审查与编辑，监督，方法学，研究，资金获取，概念化。