高级搜索

玻璃形成液体弛豫结构特征的机器学习研究

,  ,  ,  , 

downloadPDF

上一篇

下一篇

潘强强, 杨阳, 罗健, 王琦, 彭海龙. 玻璃形成液体弛豫结构特征的机器学习研究[J]. 物理学报, 2026, 75(11): 110701. doi: 10.7498/aps.75.20260233
引用本文: 潘强强, 杨阳, 罗健, 王琦, 彭海龙. 玻璃形成液体弛豫结构特征的机器学习研究[J]. 物理学报, 2026, 75(11): 110701. doi: 10.7498/aps.75.20260233
shu
Qiangqiang PAN, Yang YANG, Jian LUO, Qi WANG, Hailong PENG. Machine learning study of structural features for the relaxations of glass-forming liquids[J]. Acta Physica Sinica, 2026, 75(11): 110701. doi: 10.7498/aps.75.20260233
Citation: Qiangqiang PAN, Yang YANG, Jian LUO, Qi WANG, Hailong PENG. Machine learning study of structural features for the relaxations of glass-forming liquids[J]. Acta Physica Sinica, 2026, 75(11): 110701. doi: 10.7498/aps.75.20260233
shu

玻璃形成液体弛豫结构特征的机器学习研究

  • 1.

    中南大学材料科学与工程学院, 粉末冶金全国重点实验室, 长沙　410083

  • 2.

    中国工程物理研究院材料研究所, 绵阳　621908

    • 通讯作者: E-mail: qwang_mse@caep.cn;  E-mail: hailong.peng@csu.edu.cn

  • 中图分类号: 64.70.P-, 02.70.-c, 05.10.-a, 61.20.-p

Machine learning study of structural features for the relaxations of glass-forming liquids

  • 1.

    State Key Laboratory of Powder Metallurgy, School of Materials Science and Engineering, Central South University, Changsha 410083, China

  • 2.

    Institute of Materials, China Academy of Engineering Physics, Mianyang 621908, China

    • Corresponding authors: E-mail: qwang_mse@caep.cn;  E-mail: hailong.peng@csu.edu.cn

  • MSC: 64.70.P-, 02.70.-c, 05.10.-a, 61.20.-p

  • 摘要: 不均匀弛豫是过冷液体最典型的动力学特征, 其随时间的演化通常呈现非指数形式. 这种非指数行为使得液体的弛豫时间在体系慢动力学出现时以超Arrhenius的形式增加, 而液体结构仅有细微的变化. 在该过程中, 急剧变化的动力学性质与静态结构是否存在关联仍具有争议. 本文以典型二元玻璃形成体系$\mathrm{Cu}_{50}\mathrm{Zr}_{50}$为研究对象, 采用经典分子动力学模拟结合监督机器学习方法, 从原子尺度系统研究了静态结构与 α 弛豫动力学之间的内在关联, 定义了具有物理可解释性的原子尺度软度序参量$S_i$. 研究结果表明, 随着温度降低, 过冷液体中原子迁移率分布由单峰逐渐演化为双峰结构, 明确标志着低温区动力学不均匀性的显著增强. 由机器学习预测得到的软度场显示, 高$S_i$原子与在后续 α 弛豫过程中发生显著位移的原子高度相关. 并在空间上形成清晰“高软度-低软度原子群”双群体分布. 对于给定$S_i$的原子群体, 其动力学激活概率在过冷液态区间内服从Arrhenius行为, 对应的激活能与$S_i$呈近似线性关系; 相比之下, 其α弛豫时间在低温区表现出明显的非Arrhenius特征, 反映了协同运动的增强. 与动力学行为相对应, 软度结构序参也呈现双峰分布, 表明在过冷液态区间, 原子局域结构的激活倾向经历了显著的重分配. 本研究为揭示过冷液态过程中结构演化规律以及动力学不均匀性所对应的结构特征提供了新的视角.
    Abstract: Dynamic heterogeneity is one of the hallmarks for the dynamics of glass-forming liquids, strongly associated with non-exponential relaxation. The non-exponential relaxation can lead to the structural relaxation increasing as a super-Arrhenius formula with quenching, while with only minor change in structures. Whether there is a relationship between the statistic structure and the drastically changed dynamics is still in debate in glassy community. In this study, we combine large-scale molecular dynamics simulations with the supervised machine learning to investigate the dynamics of the binary glass-forming system Cu50Zr50. The local atomic environment is characterized using symmetry functions, and atomic mobility, i.e., the logarithm of the migration rate $\ln (m_i)$. We develop an interpretable atomic-scale structural parameter, termed “softness” ($S_i$), using an XGBoost model to elucidate the structure-dynamics relationships at the atomic level. Our results show that the probability distribution of the atomic mobility evolves from a unimodal to a bimodal distribution with quenching, consistent with the emergence of dynamical heterogeneity. The softness mapping reveals the coexistence of high- and low-softness atomic populations. For atoms with a given softness $S_i$, the activation of local structural rearrangements follows Arrhenius behavior, with the activation energy approximately linearly correlated with $S_i$. In contrast, the softness dependent α relaxation time exhibits a clear non-Arrhenius behavior, indicating the strong collectivity between different structural units. This distinct response behavior indicates that $S_i$ effectively separates the “intrinsic structural activation energy” from the “collective cooperative effects”, providing new insights into the structure-dynamics coupling during the glass transition. The high and low softness atomic clusters identified by $S_i$ can be regarded as distinct structural units with different local activation tendencies. With temperature decreasing, the distribution of softness parameter exhibits a transition from a high-softness state to a low-softness state, indicating a pronounced redistribution of the activation propensity of local atomic structures in the supercooled liquid regime. These results provide new insights on the structural evolution and the structural significance for the dynamic heterogeneity during glass transition.
  • 加载中

  • 图 1  原子迁移率的概率分布函数　(a)不同温度下原子迁移率$ \ln(m_i) $的概率密度分布函数(PDF); (b)温度从1800 K到780 K过程中原子迁移率$ \ln(m_i) $的分布情况

    Figure 1.  Probability distribution functions of atomic mobility: (a) Probability density functions (PDFs) of atomic mobility, $ \ln(m_i) $, at different temperatures; (b) total PDF of atomic mobility, $ \ln(m_i) $, in the temperature range from 1800 K to 780 K.

    下载: 全尺寸图片 幻灯片

    图 2  机器学习训练所得静态结构参量“软度”对过冷液体不均匀动力学性质的预测　(a) XGBoost 模型在不同温度下的预测精度; (b)—(f)原子软度$ S_i $的空间分布图, 其中图(b), (c)切片厚度为3 Å, 图(d)—(f) 切片厚度为6 Å, 温度依次为1600, 1200, 1000, 900和780 K. 白色圆圈为后续α弛豫时间τ内发生显著位移($ \ln(m_i) \geqslant 1.1 $)的原子

    Figure 2.  Prediction of the machine-learning-derived structural descriptors softness on the heterogeneous dynamics of supercooled liquids: (a) Prediction accurancy of the XGBoost model at different temperatures; (b)–(f) spatial distributions of atomic softness $ Si $, with the slice thickness is 3 Å for panels (b) and (c), and 6 Å for panels (d)–(f). The temperatures corresponding to panels (b)–(f) are 1600, 1200, 1000, 900, and 780 K, respectively. The white circles represent atoms that undergo significant displacement ($ \ln(m_i) \geqslant 1.1 $) during the subsequent α relaxation time τ.

    下载: 全尺寸图片 幻灯片

    图 3  结构序参软度与宏观动力学行为的关联　(a)不同软度$ S_i $结构的弛豫时间τ随温度的变化, 其中实线为整个体系的弛豫时间, 红色圆圈为偏离Arrhenius行为的转折点; (b)不同软度$ S_i $下原子被激活为快原子的概率$ P_{\mathrm{f}}(S) $与温度的关系, 虚线表示Arrhenius拟合曲线; (c)从Arrhenius拟合中提取出的依赖于原子软度$ S_i $的激活能$ \Delta E $和Σ演变

    Figure 3.  Relationship between the atomic softness and the macroscopic dynamical properties in the supercooled liquids. (a) Temperature dependence of τ for atoms with different $ S_i $. The solid line denotes the relaxation time of the system, and the red circle indicates the turning point where the behavior deviates from Arrhenius. (b) Probability of a particle being activated as a mobile one, $ P_{\mathrm{f}}(S) $, as a function of temperature for different $ S_i $ values. Dashed lines represent the Arrhenius fits. (c) Activation energy, $ \Delta E $, and entropy factor, Σ, extracted from the high-temperature Arrhenius fits, as a function of $ S_i $

    下载: 全尺寸图片 幻灯片

    图 4  结构序参软度$ S_i $的概率分布及其温度演化　(a)不同温度下$ S_i $的概率密度函数$ P (S) $; (b)—(c)对780 K与1800 K的$ P (S) $进行双峰高斯拟合. 红色表示“高软度原子群”, 蓝色表示“低软度原子群”

    Figure 4.  Probability distribution and temperature evolution of atomic softness $ S_i $: (a) Probability density functions $ P (S) $ of atomic softness $ S_i $ at different temperatures; (b)–(c) bimodal Gaussian fits of $ P (S) $ at 780 K and 1800 K. Red curve indicates the “high softness atomic group”, while the blue one indicates the “low softness atomic group”.

    下载: 全尺寸图片 幻灯片
  • [1] Du Y, Hua D, Zhou Q, Pei X, Wang H, Ren Y, Wang H, Liu W 2024 Composites Part B: Engineering 272 111222 doi: 10.1016/j.compositesb.2024.111222
    [2] Zhang Z, Ding Q, Gao Y, Wei X, Zhang Z, Bei H 2025 Mater. Futures 4 011001 doi: 10.1088/2752-5724/ad9bce
    [3] Wagner H, Bedorf D, Küchemann S, Schwabe M, Zhang B, Arnold W, Samwer K 2011 Nat. Mater. 10 439 doi: 10.1038/nmat3024
    [4] Bin S J B, Fong K S, Chua B W, Gupta M 2022 J. Magnes. Alloys 10 899 doi: 10.1016/j.jma.2021.10.010
    [5] Ma H, Zhao Y, Su Y, Luo J, Xiang J, Zheng T, Wu H, Liaw P K, Wu Y 2026 Mater. Futures 5 015003 doi: 10.1088/2752-5724/ae293f
    [6] Wu Y, Pan Q Q, Ning R, Peng H 2025 Chin. Phys. B 34 076402 doi: 10.1088/1674-1056/adcb25
    [7] 王峥, 汪卫华 2017 物理学报 66 176103 doi: 10.7498/aps.66.176103 Wang Z, Wang W H 2017 Acta Phys. Sin. 66 176103 doi: 10.7498/aps.66.176103
    [8] Angell C A, Ngai K L, McKenna G B, McMillan P F, Martin S W 2000 J. Appl. Phys. 88 3113 doi: 10.1063/1.1286035
    [9] Ediger M D, Harrowell P 2012 J. Chem. Phys. 137 080901 doi: 10.1063/1.4747326
    [10] Berthier L 2011 Physics 4 42 doi: 10.1103/Physics.4.42
    [11] Starr F W, Sastry S, Douglas J F, Glotzer S C 2002 Phys. Rev. Lett. 89 125501 doi: 10.1103/PhysRevLett.89.125501
    [12] Hu Y C, Li F X, Li M Z, Bai H Y, Wang W H 2015 Nat. Commun. 6 8310 doi: 10.1038/ncomms9310
    [13] Ding J, Cheng Y, Ma E 2013 Acta Mater. 61 3130 doi: 10.1016/j.actamat.2013.02.004
    [14] Ding J, Cheng Y Q, Ma E 2014 Acta Mater. 69 343 doi: 10.1016/j.actamat.2014.02.005
    [15] Wang F T, Liu X, Cheng J 2024 Mater. Futures 3 041001 doi: 10.1088/2752-5724/ad7619
    [16] Cai G, Cao Z, Xie F, Jia H, Liu W, Wang Y, Liu F, Ren X, Meng S, Liu M 2024 Mater. Futures 3 025601 doi: 10.1088/2752-5724/ad19e2
    [17] Behler J, Parrinello M 2007 Phys. Rev. Lett. 98 146401 doi: 10.1103/PhysRevLett.98.146401
    [18] Schoenholz S S, Cubuk E D, Sussman D M, Kaxiras E, Liu A J 2016 Nat. Phys. 12 469 doi: 10.1038/nphys3644
    [19] Wang Q, Jain A 2019 Nat. Commun. 10 5537 doi: 10.1038/s41467-019-13511-9
    [20] Fan Z, Ding J, Ma E 2020 Mater. Today 40 48 doi: 10.1016/j.mattod.2020.05.021
    [21] Fan Z, Ma E 2021 Nat. Commun. 12 1506 doi: 10.1038/s41467-021-21806-z
    [22] Bapst V, Keck T, Grabska-Barwińska A, Donner C, Cubuk E D, Schoenholz S S, Obika A, Nelson A W R, Back T, Hassabis D, Kohli P 2020 Nat. Phys. 16 448 doi: 10.1038/s41567-020-0842-8
    [23] Cubuk E D, Schoenholz S S, Kaxiras E, Liu A J 2016 J. Phys. Chem. B 120 6139 doi: 10.1021/acs.jpcb.6b02144
    [24] Smessaert A, Rottler J 2013 Phys. Rev. E 88 022314 doi: 10.1103/PhysRevE.88.022314
    [25] Cubuk E D, Schoenholz S S, Rieser J M, Malone B D, Rottler J, Durian D J, Kaxiras E, Liu A J 2015 Phys. Rev. Lett. 114 108001 doi: 10.1103/PhysRevLett.114.108001
    [26] Mendelev M I, Kramer M J, Ott R T, Sordelet D J, Yagodin D, Popel P 2009 Philos. Mag. 89 967 doi: 10.1080/14786430902832773
    [27] Plimpton S 1995 J. Comput. Phys. 117 1 doi: 10.1006/jcph.1995.1039
    [28] Yu H B, Wang Q 2024 Next Mater. 3 100168 doi: 10.1016/j.nxmate.2024.100168
    [29] Wang B, Gao X, Su R, Guan P 2024 Sci. China Phys. Mech. Astron. 67 236111 doi: 10.1007/s11433-023-2273-6
    [30] Lam C H 2018 J. Chem. Phys. 149 164909 doi: 10.1063/1.5052659
    [31] Peng H, Liu H, Voigtmann T 2022 Phys. Rev. Lett. 129 215501 doi: 10.1103/PhysRevLett.129.215501
    [32] Wang L M, Richert R 2007 Phys. Rev. B 76 064201 doi: 10.1103/PhysRevB.76.064201
    [33] Feng Y, Feng Y, Li Z, Peng H 2023 Comput. Mater. Sci. 226 112256 doi: 10.1016/j.commatsci.2023.112256
    [34] Jin C, Wu Y, Wang J, Han F, Tan M, Wang F, Xu J, Yi J, Li M, Zhang Y, Huo J, Wang J, Gao M 2024 J. Mater. Sci. Technol. 194 63 doi: 10.1016/j.jmst.2024.01.026
  • 加载中
图( 4)
计量
  • 文章访问数:  358
  • HTML全文浏览数:  358
  • PDF下载数:  0
  • 施引文献:  0
出版历程
  • 收稿日期:  2026-02-11
  • 刊出日期:  2026-06-05

玻璃形成液体弛豫结构特征的机器学习研究

    通讯作者: E-mail: qwang_mse@caep.cn; 
    通讯作者: E-mail: hailong.peng@csu.edu.cn
  • 1. 中南大学材料科学与工程学院, 粉末冶金全国重点实验室, 长沙　410083
  • 2. 中国工程物理研究院材料研究所, 绵阳　621908
  • 收稿日期:  2026-02-11

摘要: 不均匀弛豫是过冷液体最典型的动力学特征, 其随时间的演化通常呈现非指数形式. 这种非指数行为使得液体的弛豫时间在体系慢动力学出现时以超Arrhenius的形式增加, 而液体结构仅有细微的变化. 在该过程中, 急剧变化的动力学性质与静态结构是否存在关联仍具有争议. 本文以典型二元玻璃形成体系$\mathrm{Cu}_{50}\mathrm{Zr}_{50}$为研究对象, 采用经典分子动力学模拟结合监督机器学习方法, 从原子尺度系统研究了静态结构与 α 弛豫动力学之间的内在关联, 定义了具有物理可解释性的原子尺度软度序参量$S_i$. 研究结果表明, 随着温度降低, 过冷液体中原子迁移率分布由单峰逐渐演化为双峰结构, 明确标志着低温区动力学不均匀性的显著增强. 由机器学习预测得到的软度场显示, 高$S_i$原子与在后续 α 弛豫过程中发生显著位移的原子高度相关. 并在空间上形成清晰“高软度-低软度原子群”双群体分布. 对于给定$S_i$的原子群体, 其动力学激活概率在过冷液态区间内服从Arrhenius行为, 对应的激活能与$S_i$呈近似线性关系; 相比之下, 其α弛豫时间在低温区表现出明显的非Arrhenius特征, 反映了协同运动的增强. 与动力学行为相对应, 软度结构序参也呈现双峰分布, 表明在过冷液态区间, 原子局域结构的激活倾向经历了显著的重分配. 本研究为揭示过冷液态过程中结构演化规律以及动力学不均匀性所对应的结构特征提供了新的视角.

English Abstract

    全文HTML

    1.   引　言

    • 金属玻璃(metallic glasses, MGs)是指金属熔体在快速冷却过程中成功规避晶化, 从而将液态原子排布“冻结”于非晶态的一类合金材料. 由于具有“短程有序、长程无序”的复杂拓扑结构特征, 金属玻璃在力学强度[1,2]、弹性极限[3]、耐腐蚀性[4,5]以及软磁性能[6]等方面表现出优异的性能. 然而, 其形成机理及动力学本质至今仍未得到统一理解. 普遍认为, 非晶态的本质问题深植于其形成过程中的过冷液体(supercooled liquid)阶段. 当金属熔体被冷却至熔点$ T_{\mathrm{m}} $以下而尚未发生结晶时, 体系进入一种热力学亚稳、动力学高度迟滞的状态[7]. 在该温度区间内, 体系黏度随温度降低呈指数级增长, 变化幅度可跨越十余个数量级, 而局域结构却仅发生微弱变化. 这种“结构近似保持而动力学剧烈演化”的反常行为, 构成了玻璃转变(glass transition)的核心特征[8]. 需要强调的是, 玻璃转变并非传统意义上的一级或二级热力学相变, 其既不伴随热力学量的跃变, 也不存在全局序参量的突变, 而是一种本质上非平衡、非遍历的动力学冻结过程.

      玻璃转变过程中最显著的动力学特征之一, 是随温度降低而逐渐增强的动力学不均匀性(dynamical heterogeneity). 在相同外部条件下, 不同空间区域或不同原子可表现出截然不同的动力学行为—部分原子仍能快速重排, 而邻近区域则近乎被“固化”[9,10]. 这一现象表明, 在过冷液体中可能存在某种隐匿的、具有结构起源的动力学前驱体, 它们主导着体系由液态向玻璃态的渐进式演化. 然而, 长期以来, 如何在缺乏长程序的非晶体系中识别这类能够被称为“结构缺陷”或“动力学易激活单元”的动力学前驱体, 始终是玻璃物理中的核心难题.

      早期研究主要从自由体积[11]、五重对称性[12]、几何优先团簇[13,14]等结构指标出发, 试图建立结构-动力学关联, 但这些指标往往具有较强的体系依赖性, 且预测能力有限. 近十余年来, 随着计算能力的提升与数据科学的发展, 监督式机器学习(supervised machine learning, ML)为解决这一问题提供了新的研究范式[15,16]. 由Behler和Parrinello[17]提出的原子中心对称函数(symmetry functions, SFs)最初被用于高维神经网络势函数的构建, 该函数天然满足平移、旋转不变性, 并能够表征局域径向结构信息. 这使其成为表征非晶局域结构的理想工具. 基于此, Schoenholz等[18]率先将SFs引入无序体系研究, 以原子跳跃重排度$ p_{\mathrm{hop}} $作为监督信号, 定义了原子尺度的“软度”(softness), 并证明该参量能够有效预测玻璃及过冷液体中局域结构重排的倾向性. 随后, 相关研究不断拓展这一思想, 发展出结构柔性[19]、密度涨落驱动序参量[20]以及基于深度学习[21,22]的预测模型, 持续加深了人们对非晶结构-动力学关联的理解. 然而, 现有多数基于softness的研究仍主要依赖离散且稀有的动力学事件作为监督信号[18,23-25], 实际上此类参量对时间窗口选择、体系笼效应及统计采样高度敏感, 亦在一定程度上模糊了局域结构本征动力学与体系整体协同运动之间的物理区分. 尤其是在过冷液体向玻璃态演化的关键温区, 如何构建一种既连续、又能直接反映α弛豫尺度动力学的结构序参量, 仍是一个亟待解决的问题.

      基于上述背景, 本文以典型二元金属玻璃形成体系$ \mathrm{Cu}_{50}\mathrm{Zr}_{50} $为研究对象, 系统研究其从高温熔体到过冷液体阶段的结构与动力学演化过程. 这里引入以原子α弛豫时间τ为基础定义的原子迁移率$ {\mathrm{ln}} (m_{i})$作为动力学判据, 构建了一种新目标的监督式机器学习框架, 从而在原子尺度上建立静态局域结构与α弛豫动力学之间的直接联系. 该方法不仅避免了对稀有事件的依赖, 也使结构参量能够自然区分本征结构激活与集体协同效应. 通过这一框架, 进一步揭示了过冷液体中高软度/低软度原子群分布由温度诱导的重分配过程, 以及其在玻璃转变过程中的关键作用.

    2.   模拟方法

      2.1.   分子动力学模拟

    • 采用经典分子动力学(molecular dynamics, MD)模拟方法制备$ \mathrm{Cu}_{50}\mathrm{Zr}_{50} $金属玻璃过冷液体样品. 原子间相互作用通过Mendelev等[26]开发的嵌入原子方法(embedded atom method, EAM)势函数描述. 模拟体系为立方盒子结构, 包含5000个原子, 其中Cu原子和Zr原子各2500个. 模拟盒子的初始尺寸在三个空间方向上均约为46 Å. 初始构型中, 原子随机分布在模拟盒内, 并在2000 K高温下充分弛豫, 以获得均匀的液态结构. 随后, 采用1 fs的时间步长, 以1012 K/s的降温速率将熔体依次冷却至目标温度(1800, 1600, 1400, 1200, 1100, 1000, 900, 850, 810和780 K), 从而获得不同温度下的过冷液体结构. 冷却过程在等温等压(NPT)系综下进行, 通过动态调整模拟盒子尺寸, 使体系在三个方向的压力始终保持在零附近. 在完成快速冷却后, 每个目标温度下的样品盒子尺寸由10个独立构型的平均值计算得出. 随后, 在等体积等温(NVT)系综下进行充分弛豫, 以获得热力学稳定的过冷液体结构. 最终, 轨迹采样采用等构型(iso-configurational)采样方法: 在选定的参考构型下, 保持原子空间坐标不变, 仅对原子速度按照玻尔兹曼分布进行随机重新赋值, 并对每个构型重复进行200次独立动力学演化, 以获得单个原子平均运动轨迹. 时间采样采用指数递增的时间间隔进行数据记录. 所有分子动力学模拟均使用大规模原子/分子并行模拟器(LAMMPS)软件包完成[27].

      • 2.2.   机器学习框架

      2.2.1.   分类参量

    • 针对机器学习任务, 利用单原子的动力学性质来构建分类参量. 具体而言, 采用原子迁移率$ \ln(m_{i})=-\ln(\tau_{i})$ 作为衡量指标, 其中$ \tau_{i} $为原子i的α弛豫时间, 定义为该原子的自中间散射函数(self-intermediate scattering function, SISF)衰减至$ 1/{\mathrm{e}} $所需的时间. 每个原子的SISF表达式为

      其中$ \boldsymbol{r}_i(t)$为原子i在时刻t的位置, 波矢模长$ \boldsymbol{q} $取为静态结构因子第一峰对应的位置, $ \langle \cdots \rangle_{{\mathrm{iso}}} $表示等构型(iso-configurational)系综平均. 为消除原子热运动的随机性, 在相同的初始静态构型下生成200条独立动力学轨迹, 每条轨迹的初始速度均从Maxwell-Boltzmann分布中随机抽取, 从而获得统计意义上稳定且可靠的动力学倾向(dynamic propensity), 这是后续机器学习模型训练的关键基础. 在此基础上, 进一步选取$ \ln(m_i)=1.1 $作为区分快慢原子的物理阈值. 当$ \ln(m_i)\geqslant 1.1 $时, 原子被判定为“快”原子, 反之则为“慢”原子. 该阈值确定的依据见后续分析.

      • 2.2.2.   局域结构描述符

    • 采用一组原子对称函数(symmetry functions, SFs)对原子局域结构环境进行描述. 该方法最早由 Behler 和 Parrinello 提出[17]. 该类描述对体系的平移和旋转保持不变, 能够同时捕捉邻近原子的径向结构信息, 广泛用于非晶和复杂结构体系的机器学习建模. 在$ \mathrm{Cu}_{50}\mathrm{Zr}_{50} $体系中, 为每个原子i构建径向对称函数, 其定义为

      其中, $ \displaystyle\sum\nolimits_{j \in X} $表示对类型为X(Cu/Zr)的邻近原子求和, $ r_{ij} $表示原子i与j之间的距离, r和σ分别为高斯函数的中心位置和宽度参数. 截断半径取Rc = 10 Å, $ \sigma$ = 0.1 Å, 参数r在2—10 Å范围内以0.1 Å的步长取值. 由此, 每个原子可获得一个包含160维分量的结构特征向量, 其中来自Cu邻居和Zr邻居的径向对称函数各80个. 该特征向量作为机器学习模型的输入, 用以表征原子的局域结构环境.

      • 2.2.3.   模型训练

    • 每个粒子的软度($ S_i $)由监督式机器学习分类模型的预测结果定义. 这里首先根据迁移率阈值$\ln({m_i})= 1.1$对原子进行分类标注. 对于每个独立的过冷液体样品, 随机选取400个“快”粒子和400个“慢”粒子作为训练样本, 数据涵盖10个温度, 每个温度包含10个独立样品. 为了确保模型在不同温度区间上的预测能力, 将所有温度的数据合并, 构建统一数据集, 每次采样得到的数据集规模为80000个原子. 为保证结果具有统计上的稳定性, 重复随机采样200次, 每一次训练都采用5折交叉验证(K-fold cross validation). 具体地, 在构建训练集和测试集时, 先在单一温度的数据内部进行5折划分, 再将各温度对应的折次数据合并, 从而保证训练集与测试集中来自各温度的数据量保持一致. 机器学习模型采用XGBoost, 并通过网格搜索(grid search)对超参数进行优化. 模型性能以受试者工作特征曲线下面积(area under the receiver operating characteristic curve, AUC-ROC)作为评估指标. 最终, 对 200 次独立训练模型所得的“快”原子预测概率取平均, 并将该平均值作为每个原子的软度$ S_i $.

    3.   结果与讨论

    • 为直观对比过冷液态中原子动力学的差异, 将$ \ln(m_i)$定义为原子迁移率(详见方法部分), 以此作为动力学指标. 该指标能够自适应反映原子完成局部结构弛豫所需的时间, 和基于固定时间窗口的位移指标相比[28,29], 该参数能够更加精确地区分高温与低温具有不同特征弛豫时间时原子可动性的差异. 图1(a)展示了不同温度下$ \ln(m_i)$的概率分布函数(probability distribution function, PDF). 在高温区($ T\geqslant 1400\; \mathrm{K} $), 分布呈单峰且近似高斯形态, 反映体系动力学相对均匀; 当温度降至1200 K以下, PDF开始显著展宽, $ \ln(m_i)$分布逐渐演化为明显的双峰结构. 主峰对应的是占多数的低迁移率原子, 即仍然受到局域笼困或刚性环境限制的慢原子背景, 其位置随温度降低逐渐向更低迁移率方向移动; 与此同时, 在高迁移率一侧逐渐出现肩峰并发展为可分辨的次峰, 对应在后续α弛豫时间尺度内更容易发生显著位移的高迁移率原子群. 需要指出的是, 这里的双峰结构反映的是在给定观测时间窗口内不同原子局域动力学活性的统计分群, 更是体系慢动力学出现的证明, 而不是时间相关函数或频谱中的弛豫过程分解. 因此, 该分布更适合理解为快、慢动力学原子群的出现及其分离, 体现了动力学不均匀性的增强, 其主峰与次峰可能与弛豫谱上α 弛豫和β弛豫相关, 但不是严格等同. 类似的单原子迁移率分布双群体特征在 Lam[30]的工作中也曾观察到. 为进一步量化动力学差异, 统计了所有温度$ \ln(m_i)$的分布(图1(b)). 图1(b)中由于低温下动力学不均匀性的增加, 导致高温下的主峰占据主导作用, 基于此, 将$ \ln(m_i)= 1.1 $设为“快”原子判据(该值对应PDF主峰拐点, 且与图1(a)中快弛豫亚峰位置一致, 该阈值会随温度最大值进行微调, 但不会影响最终的结论), 并以此作为后续$ P_{\mathrm{f}}(S)$计算的激活阈值.

      为探究过冷液体动力学不均匀性的结构起源, 采用SFs编码每个原子的局部环境, 采用XGBoost模型将原子局域环境与较大的迁移率原子进行关联. 通过分类器预测所有原子属于“快”类的概率, 该概率即定义为原子i的软度$ S_i $, 其表征在给定静态结构下原子发生α弛豫的倾向性. 图2(a)展示了模型在各温度下的交叉验证性能: AUC值均大于93%, 该值验证了模型分类的高度可靠性. 模型预测精度随温度呈现非单调变化, 并在1200 K附近出现极小值, 该温度与图1(a)中原子迁移率分布由单峰向双峰演化、体系慢动力学开始出现的温区相对应. 在该温度点, 不同类型原子之间的局域结构环境相似. 这可能导致模型对快慢原子的分类能力相较于其他温度较弱, 模型精度出现最小值. 图2(b)—(f)直观展示了$ S_i $场的空间分布及其与真实动力学事件的关联. 在1600—780 K的温度范围内, $ S_i $的空间分布图中, 高$ S_i $区域(暖色)与后续α弛豫时间τ内发生显著位移的原子(白圈)高度重合. 同时, 一个有趣的现象出现, 高温下$ S_i $场普遍较高, 且在不同的温度下, $ S_i $场都呈现出典型的红蓝相间图样或“高软度-低软度原子群”双群体分布(图2(d)—(f)), 这与动力学不均匀性的温度依赖性相一致. 这表明, $ S_i $不仅是一个稳健的预测指标, 更是一种可解释的结构序参量—它将微观局域结构映射到宏观动力学.

      在软度序参量能够精确描述原子局部弛豫行为的基础上, 进一步考察该参量与宏观动力学行为的关联. 通过软度, 可以将不同局域结构进行分类, 进而分析其宏观弛豫行为. 图3(a)给出了不同$ S_i $的弛豫时间τ的温度依赖性. 可以看到, 对于任意给定$ S_i $子集, 在高温区域($ T \gt 1200\; \mathrm{K} $), $ \tau(S)$呈典型的Arrhenius行为; 而在低温区域($ T \lt 1200\; \mathrm{K} $), $ \tau(S)$相对于高温Arrhenius关系出现了明显的上偏离, 这种偏离与动力学不均匀性的增强相伴出现, 该偏离温区与金属玻璃体系中扩散系数脱离Arrhenius行为的温度区间相似[31]; 其与分子体系具有相类似的偏离关系, 但在具体弛豫时间值上存在差异[32]. 类似现象也出现在近期 Feng等[33]对AlTi液体协同运动的第一性研究中. 这种非Arrhenius转变可以从协同运动的角度加以理解. 作为一种热激活过程, 结构弛豫时间可表示为$ \tau=\tau_0\exp[E_{\mathrm{A}}/(k_{\rm B} T)] $, 其中$ E_{\mathrm{A}} $表示激活能, 对应于弛豫过程中被激活的一定体积内的原子群. 在高温区间($ T \gt 1200\; \mathrm{K} $), 体系弛豫主要由局域事件主导, 因此Arrhenius行为成立. 随着温度降低($ T \lt 1200\; \mathrm{K} $), 局域弛豫开始受到周围环境的限制. 激活能随之增大, 从而导致Arrhenius关系失效, 并反映出激活体积或动态关联长度的增大. 因此, $ \tau(S)$的非Arrhenius行为可视为玻璃形成液体中局域协同运动增强的标志.

      如图1和图2所示, 玻璃形成液体中存在两种不同的动力学状态. 若将PDF中次峰右侧、具有较高迁移率的原子视为“激活态”原子(见图1(b)), 则可以在给定结构条件下计算原子的激活概率. 例如, 对于具有特定软度$ S_i $的原子, 其激活概率$ P_{\mathrm{f}}(S)$定义为在观测时间窗口$ \Delta t=\tau $内满足$ \ln(m_i)\geqslant 1.1 $的概率. 图3(b)给出了$ P_{\mathrm{f}}(S)$在不同温度下的演变规律. 结果表明, 在整个温度范围(780—1800 K)内, $ P_{\mathrm{f}}(S)$随温度变化始终服从Arrhenius关系: $ P_{\mathrm{f}}(S)=P_0 \exp [\varSigma -\Delta E/T] $. 该行为与整体弛豫时间τ随温度变化的超Arrhenius关系(图3(a))形成鲜明对比. 该结果表明, 如果将原子结构重排简化为“是否被激活”的事件, 并在固定局域结构条件下进行统计, 则不同结构环境下的局域激活过程可以用近似单一的有效能垒来描述. 换言之, 通过软度这一结构序参量, 可以将整体体系复杂的非Arrhenius动力学分解为多个由局域结构表征的、近似Arrhenius的激活通道. 从这一意义上讲, 该结果与动力学不均匀性理论中的一个观点具有一定相似性: 当整体体系中复杂的异质性被分解为若干更均一的局域子群后, 每个子群的动力学行为会更加接近单尺度描述. 但需要强调的是, 本文图3(b)所讨论的是温度域上的局域激活统计, 而动力学不均匀性理论中所谓“纯指数行为”通常是指时间相关函数在时间域上的衰减形式, 因此两者虽然相关, 但是不完全一致. 它们的严格一致性需要进一步验证. 类似利用原子尺度结构描述动力学激活行为的研究, 在以原子跳跃重排度作为监督信号的机器学习工作中也曾得到类似结果[23,18].

      对不同软度S的数据进行拟合可得到激活能$ \Delta E(S)$与$ \varSigma(S)$. 二者均随S呈近似线性变化(如图3(c)). 低温区域少数低S原子的偏离主要源于样本数量不足, 导致其拟合的激活能$ \Delta E(S)$偏离线性趋势. 总体而言, 这一线性关系与Liu等[23,18]在Lennard-Jones经典体系中的发现一致. 这表明, 软度$ S_i $可以作为有效结构参量, 其值越高, 原子越易被热激活. 针对$ \tau(S)$与$ P_{\mathrm{f}}(S)$截然不同的温度行为, 我们提出以下物理图景: $ \tau(S)$的非Arrhenius行为来源于体系协同效应的增强. 随着温度降低, 黏度急剧上升, 使得有效弛豫过程需要更大体积的原子协同参与. 这体现在图3(a)中$ \log(\tau)$随$ 1/T $斜率的上升. 因此弛豫对应的激活能$ E_{\mathrm{A}} $随温度降低而增大. $ P_{\mathrm{f}}(S)$的Arrhenius行为则体现了局域结构的本征激活特征. 对于固定$ S_i $的原子群, 其激活能$ \Delta E(S)$完全由局域结构环境决定, 与体系宏观黏度无关. 换言之, 软度$ S_i $将“局域结构的本征激活倾向”与“宏观协同行为”这两类动力学贡献有效解耦, 为理解玻璃形成液体的多尺度动力学提供了清晰的结构化框架.

      为进一步解释这种激活概率和弛豫时间随温度演变的差异, 计算了软度的概率密度分布$ P(S)$随温度的变化(如图4(a)). 可以很明显地看出, 高温时$ P(S)$在高软度位置出现一个主要的峰. 随着温度降低, 该峰逐渐向较小$ S_i $方向移动, 最终在低$ S_i $处形成主峰. 这一演化表明, 体系中不同局域结构环境的统计权重随温度发生持续变化. 为了更加定量地描述这一过程, 将$ P(S)$拟合为两个高斯分量之和, 分别对应“高软度原子群”(软区, 红色)和“低软度原子群”(硬区, 蓝色), 具体结果如图4(b)和图4(c)所示. 需要指出的是, 此处定义的“软区”与“硬区”本质上是基于模拟与机器学习识别出的结构-动力学统计分群, 其物理含义是两类具有不同局域激活倾向的结构基元, 而并非严格意义上的热力学两态. 在物理图景上, 它们可以与实验中常讨论的类液区、类固区以及局域弹性异质性和动态等概念[34]进行定性对应, 但当前工作并未建立软区与硬区与具体实验结构单元之间的一一映射关系, 因此这里仅将其作为与已有实验图景相一致的一种物理解释. 高温下(1800 K), 分布集中于高$ S_i $处, 体系以“软区”为主; 随着温度降低, 该分布逐渐向低$ S_i $处演变; 低温下, 具有较低$ S_i $的“硬区”逐渐暂居主导地位. 当双峰的相对权重与峰位变化随降温持续演化并逐渐趋于稳定时, 体系进入更强的玻璃态前驱区间. 值得注意的是, 在$ T\sim1200\; \mathrm{K} $附近, $ P(S)$关于$ S_i\approx0.5 $附近呈现近似对称的双峰分布, 表明软区与硬区的统计权重在该温区发生了显著的重分配过程. 该温区同时对应迁移率分布开始出现明显双峰化, 以及整体弛豫时间偏离高温Arrhenius关系的起始区域. 上述结果表明, 在这一温区内局域结构的统计重分配与动力学不均匀性的增强同步发生, 体现出明显的结构-动力学协同演化特征, 即一个显著的结构-动力学交叉过渡温区(structural-dynamical crossover). 随着温度进一步降低, 体系由动力学相对均匀的液体状态逐渐进入具有显著动力学不均匀性的玻璃态前驱区间. 类似的结构分布转变也曾在胶体玻璃体系中被观察到[24], 表明局域结构环境的统计重分配可能是玻璃化过程中普遍存在的重要结构特征.

    4.   结　论

    • 本文基于对$ \mathrm{Cu}_{50}\mathrm{Zr}_{50} $玻璃形成液体从1800 K到780 K热力学路径的分子动力学模拟, 结合以原子迁移率$ \ln(m_i)$为标签的监督机器学习方法, 构建了可解释的原子尺度软度序参量$ S_i $, 并系统揭示了其与动力学行为之间的关联. 研究表明, 随着温度降低至约1200 K, 原子迁移率分布由单峰演变为双峰; 同时, 软度$ S_i $概率分布$ P(S)$也在该温区呈现近似对称的双高斯结构, 表明高软度与低软度原子群的统计权重在这一温区发生了显著的由温度诱导的重分配过程. 该温区同时对应协同运动的出现、弛豫时间关于温度偏离Arrhenius关系, 以及局域结构环境双峰分布近似的临界现象. 这表明玻璃化过程并非简单的原子个体冻结, 而是伴随着局域结构激活倾向统计分布的持续演化, 并表现出显著的结构-动力学交叉过渡行为(structural-dynamical crossover). 此时, 体系由动力学相对均匀的液体逐渐进入具有强动力学不均匀性的玻璃态前驱区间. 进一步分析表明, $ S_i $具有解耦多尺度动力学的关键能力, 对于固定$ S_i $值的原子, 其热激活概率$ P_{\mathrm{f}}(S)$遵循Arrhenius关系, 且激活能$ \Delta E(S)$与$ S_i $呈近似线性关系; 对应的弛豫时间$ \tau(S)$的非Arrhenius行为则源于低温下协同运动增强导致的团簇激活能升高. 这种差异化的响应机制说明, $ S_i $有效分离了“本征结构激活能”与“集体协同效应”两个层次的动力学贡献, 从而为理解玻璃转变中的结构-动力学耦合关系提供了新的视角. 本文识别出的高软度与低软度原子群可视为两类具有不同局域激活倾向的结构基元, 其统计分布与动力学不均匀性的增强密切相关, 也为从微观结构角度理解非晶材料的动力学不均匀性提供了新的分析框架, 并为未来通过调控非晶合金的局域结构分布来优化材料性能提供了可能的理论启示.

    参考文献 (34)

目录

/

返回文章
返回