高级搜索

结合高压热弹性模拟数据的Hugoniot弹性极限实验解读

,  ,  ,  ,  ,  , 

downloadPDF

上一篇

下一篇

尹俊磊, 高兴誉, 王毅, 刘海风, 宋海峰, 李金山, 张平祥. 结合高压热弹性模拟数据的Hugoniot弹性极限实验解读[J]. 物理学报, 2026, 75(12): 120801. doi: 10.7498/aps.75.20260047
引用本文: 尹俊磊, 高兴誉, 王毅, 刘海风, 宋海峰, 李金山, 张平祥. 结合高压热弹性模拟数据的Hugoniot弹性极限实验解读[J]. 物理学报, 2026, 75(12): 120801. doi: 10.7498/aps.75.20260047
shu
Junlei YIN, Xingyu GAO, Yi WANG, Haifeng LIU, Haifeng SONG, Jinshan LI, Pingxiang ZHANG. Interpretation of Hugoniot elastic limit experiments based on high-pressure thermoelastic simulations[J]. Acta Physica Sinica, 2026, 75(12): 120801. doi: 10.7498/aps.75.20260047
Citation: Junlei YIN, Xingyu GAO, Yi WANG, Haifeng LIU, Haifeng SONG, Jinshan LI, Pingxiang ZHANG. Interpretation of Hugoniot elastic limit experiments based on high-pressure thermoelastic simulations[J]. Acta Physica Sinica, 2026, 75(12): 120801. doi: 10.7498/aps.75.20260047
shu

结合高压热弹性模拟数据的Hugoniot弹性极限实验解读

  • 1.

    西北工业大学, 凝固技术全国重点实验室, 西安　710072

  • 2.

    北京应用物理与计算数学研究所, 计算物理全国重点实验室, 北京, 100088

    • 作者简介: 尹俊磊: junlei_yin@mail.nwpu.edu.cn .
    通讯作者: E-mail: gao_xingyu@iapcm.ac.cn.;  E-mail: song_haifeng@iapcm.ac.cn.;  E-mail: pxzhang@nwpu.edu.cn.

Interpretation of Hugoniot elastic limit experiments based on high-pressure thermoelastic simulations

  • 1.

    State Key Laboratory of Solidification Processing, Northwestern Polytechnical University, Xi’an 710072, China

  • 2.

    National Key Laboratory of Computational Physics, Institute of Applied Physics and Computational Mathematics, Beijing 100088, China

    • Corresponding authors: E-mail: gao_xingyu@iapcm.ac.cn.;  E-mail: song_haifeng@iapcm.ac.cn.;  E-mail: pxzhang@nwpu.edu.cn.

  • 摘要: 针对冲击加载实验中波前弥散及非稳态效应制约Hugoniot 弹性极限(HEL)精确解读的问题, 本文提出了一种结合理论模拟与实验数据的自洽解读方法. 该方法基于宏观粒子动量守恒定律与微观材料本征纵波声速响应的物理关联, 构建了以纵波声速和轴向应力为核心变量的自洽方程. 在求解策略上, 以实验测得的弹性前驱波粒子速度为约束变量, 并利用平均场势(MFP)耦合准静态近似(QSA)生成的宽压域热弹性数据集, 构建了连续的纵波声速-压力响应关系. 本文以金刚石为研究对象验证该方法的有效性, 首先构建了室温下0—1000 GPa压区的金刚石热弹性数据集, 其与15 GPa压力范围内静加载实验数据的相对偏差不超过0.4%. 在此基础上, 整合了峰值应力高达1 TPa的金刚石冲击实验数据, 应用上述方法重新解读了HEL. 结果表明, 重新解读的HEL与标准气炮加载实验值的相对偏差在2%以内, 并定量表征了其各向异性. 此外, 针对激光加载实验中因非稳态衰减导致的系统性偏差, 本方法将其HEL值向上修正了8%—13%. 本文数据集可在https://doi.org/10.57760/sciencedb.j00213.00274中访问获取.
    Abstract: The Hugoniot elastic limit (HEL) serves as the critical stress threshold demarcating the transition from purely elastic to elastoplastic response under dynamic loading. Accurate determination of the HEL is essential for understanding dynamic mechanical behavior. Traditionally, the HEL is derived from the “double-wave” structure of free-surface velocity profiles. However, wavefront dispersion of the elastic precursor and non-steady-state effects in short-pulse laser-driven experiments introduce systematic deviations in longitudinal sound velocity measurements, thereby limiting precise HEL determination. To address these challenges, a self-consistent interpretation method integrating macroscopic conservation laws with microscopic thermoelastic simulations is established.Based on momentum conservation, this work establishes a self-consistent equation that relates axial stress to longitudinal sound velocity and particle velocity. This equation corresponds to the process in which high-stress perturbations catch up with low-stress wavefronts. Using the experimentally measured elastic-precursor particle velocity as the constraint, the HEL and the corresponding longitudinal sound velocity can be determined self-consistently by iteratively solving this equation, without relying on ambiguous wavefront arrival-time measurements. To support this procedure, a continuous sound velocity-pressure constitutive relationship is derived from a high-pressure thermoelastic dataset spanning 0–1000 GPa, constructed using density functional theory combined with the mean-field potential (MFP) method and the quasi-static approximation (QSA).The validity of this approach is verified using diamond as a benchmark material. The predicted elastic moduli agree with static compression data within 0.4% relative deviation up to 15 GPa. For gas-gun experiments with peak stresses reaching 1 TPa, the reinterpreted HEL values show less than 2% relative deviation from standard experimental benchmarks. These results quantitatively characterize the crystalline anisotropy, where the longitudinal sound velocity follows the sequence [111] > [110] > [100]. For laser-driven experiments where non-steady-state attenuation previously caused systematic underestimation, the HEL values are corrected upward by 8%–13%, reconciling discrepancies between different loading platforms.This hybrid approach enables the retrospective correction of non-ideal experimental effects through solely through data reanalysis, thereby providing high-precision dynamic constitutive parameters under extreme conditions without requiring specialized experimental configurations. The datasets presented in this paper are openly available at https://doi.org/10.57760/sciencedb.j00213.00274.
  • 加载中

  • 图 1  (a)弹塑性响应理想双波结构的速度剖面示意图; (b)铁的实验速度剖面[17]; 在理想速度剖面(a)中, 自由面速度依序呈现弹性波(橙色)和塑性波(蓝色)的陡峭波前与平台特征. HEL可由粒子速度$ {u}_{\text{E}} $和纵波声速$ {C}_{\text{L}} $确定. $ {u}_{\text{E}} $近似为弹性波平台处自由面速度$ {u}_{\text{f}} $的一半. 纵波声速$ {C}_{\text{L}} $通常被定义为d/$ {t}_{\text{a}} $, 其中$ d $为靶板厚度; $ {t}_{\text{a}} $为弹性前驱波穿过样本所需时间, 即弹性波走时. 然而, 在铁的实验测量剖面(b)中, 由于复杂的动态响应, 平台和上升沿均显著呈现弥散现象. 此时, 自由面速度$ {u}_{\text{f}} $通常取自平台起始处. 弹性前驱波上升沿的弥散导致其到达时间的判读存在较大的不确定度$ ∆t $, 这直接导致了纵波声速取值的不确定性

    Figure 1.  (a) Schematic of the ideal double-wave velocity profile for elastic-plastic response and (b) experimental velocity profiles of iron[17]. In the ideal profile (a), the free surface velocity successively exhibits the steep wavefronts and plateau features of the elastic wave (orange) and plastic wave (blue). The HEL can be determined by the particle velocity $ {u}_{\text{E}} $ and the longitudinal sound velocity $ {C}_{\text{L}} $. Specifically, $ {C}_{\text{L}} $ is defined as the ratio of the target thickness $ d $ to the transit time of the elastic precursor $ {t}_{\text{a}} $. $ {u}_{\text{E}} $ is approximately half the free surface velocity $ {u}_{\text{f}} $ at the elastic plateau. However, in the experimental profiles for iron (b), both the plateaus and rising fronts exhibit significant smearing due to complex dynamic responses. In such cases, the free surface velocity is typically sampled at the onset of the plateau. The smearing of the elastic precursor’s rising front leads to considerable uncertainty in determining its arrival time, which directly translates into uncertainty in the determined longitudinal sound velocity.

    下载: 全尺寸图片 幻灯片

    图 2  融合实验与理论模拟的 HEL 混合解读方法计算流程图. (左侧) 实验数据处理模块: 从速度剖面中提取弹性前驱波粒子速度, 作为自洽方程求解的边界条件. (右侧) 理论模拟模块: 构建有限温度-压力空间下的纵波声速数据集. (中间) 自洽迭代求解流程: 结合上述实验边界条件与理论声速数据, 通过应力-声速自洽迭代, 确定临界屈服条件下的轴向应力(即 HEL)

    Figure 2.  Computational flowchart of the hybrid experimental-theoretical framework for HEL determination. (Left Panel) Experimental data processing module: Extraction of the elastic precursor particle velocity from velocity profiles to serve as the boundary condition for the self-consistent interpretation equation. (Right Panel) Theoretical simulation module: Construction of the longitudinal sound velocity dataset within the finite temperature-pressure domain. (Middle Panel) Self-consistent interpretation iterative solution process: Integration of the aforementioned experimental boundary conditions and theoretical sound velocity data to determine the axial stress at the critical yield condition (i.e., HEL) via stress-sound velocity self-consistent interpretation iteration.

    下载: 全尺寸图片 幻灯片

    图 3  模拟计算与实验测量的金刚石状态方程对比　(a)零压下体积随温度的变化(等压膨胀线), 散点数据分别引自文献[46](红色, Reeber)和文献[47](蓝色, Thewlis); (b) 900 K(红色)与常温(蓝色)下的体积随压力的变化(等温压缩线), 散点数据引自文献[48]利用电阻加热金刚石压砧(DAC)技术获得的X射线衍射测量结果; 图中实线为本研究基于QHA+MFP方法的模拟结果, 需要说明的是, 经典MFP模型由于没有考虑量子效应, 在极低温下的模拟会出现非物理的预测, 因此本研究只保留了298 K以上的数据

    Figure 3.  Comparison of the calculated equation of state (EOS) for diamond with experimental measurements: (a) Volume versus temperature at zero pressure (isobaric expansion), where the scatter points correspond to experimental data from Ref. [46] (red, Reeber) and Ref. [47] (blue, Thewlis). (b) Volume versus pressure at 900 K (red) and ambient temperature 298 K (blue) (isothermal compression). The solid lines represent the simulation results from this study based on the QHA+MFP method, while the scatter points are X-ray diffraction measurements by Ref. [48] using resistive-heating diamond anvil cell (DAC) techniques.

    下载: 全尺寸图片 幻灯片

    图 4  金刚石热弹性性质的理论计算与实验对比　(a)零压下弹性常数 C11(蓝色)和C44 (红色) 随温度的变化(上图)及相对偏差分析(下图), 散点为 Zouboulis等[49]的布里渊散射数据; (b)室温下弹性模量(上图)与纵波声速(下图)随压力的变化, 上图中实线与虚线分别代表由弹性常数导出和由 EOS 二阶导数导出的体积模量; 实验数据引自文献[50] (Hu)

    Figure 4.  Comparison of theoretical calculations and experimental measurements of the thermoelastic properties of diamond: (a) Temperature dependence of elastic constants C11 (blue) and C44 (red) at zero pressure (upper panel) and the corresponding relative deviation analysis (lower panel). The scatter points represent Brillouin scattering data from Ref. [49]. (b) Pressure dependence of elastic moduli (upper panel) and longitudinal sound velocity (lower panel) at room temperature. In the upper panel, the solid and dashed lines represent the bulk modulus derived from elastic constants and from the second derivative of the EOS, respectively. Experimental data are taken from Ref. [50] (Hu).

    下载: 全尺寸图片 幻灯片

    图 5  金刚石冲击波速与HEL的理论预测及实验验证　(a)纵波声速随粒子速度的非线性演化关系, 实线代表本研究针对多晶(黑色)及单晶[100](蓝色)、[110](绿色)、[111](红色)晶向的理论计算结果, 散点数据分别引自McSkimin和Andreatch[52]的准静态测量、Lang等[22]、Winey等[54]、Kondo和Ahrens[53]及Mashimo等[55] 的气炮冲击实验, 以及 McWilliams等[23]的激光冲击实验; (b)基于本文混合方法重评估的HEL与实验原始解读HEL 的对比, 图例符号及颜色定义与图(a)保持一致

    Figure 5.  Theoretical prediction and experimental validation of shock wave velocity and HEL for diamond: (a) The non-linear evolution of longitudinal sound velocity with particle velocity. Solid lines represent theoretical calculations from this study for polycrystalline (black) and single-crystal [100] (blue), [110] (green), and [111] (red) orientations. Scatter points are derived from experimental measurements, including quasi-static measurements by McSkimin and Andreatch[52], gas gun shock experiments by Lang et al.[22], Winey et al.[54], Kondo and Ahrens[53] and Mashimo et al.[55], and laser shock experiments by McWilliams et al.[23]. (b) Comparison between HELs re-evaluated by the proposed hybrid method and those originally interpreted from experiments. Symbols and colors are consistent with panel (a).

    下载: 全尺寸图片 幻灯片

    表 1  金刚石冲击实验数据汇总

    Table 1.  Summary of diamond impact experiment data.

    文
    献    		 加载
    方式    		 发表年 晶向与
    组织    		 粒子速度/
    (km·s–1)    		 纵波声速/
    (km·s–1)    		 文
    献    		 加载
    方式    		 发表年 晶向与
    组织    		 粒子速度/
    (km·s–1)    		 纵波声速/
    (km·s–1)
    [53] 气炮 1983 多晶 0.9±0.34 20.1±0.8 [55] 气炮 2025 $ \left\langle 100\right\rangle $ 0.703±0.011 18.1±0.3
    气炮 1983 多晶 0.899±0.063 19.61±0.23 气炮 2025 $ \left\langle 100\right\rangle $ 0.755±0.012 18.2±0.3
    [22] 气炮 2018 $ \left\langle 100\right\rangle $ 0.849±0.011 18.24±0.15 气炮 2025 $ \left\langle 100\right\rangle $ 0.552±0.009 18.1±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.829±0.011 18.19±0.14 气炮 2025 $ \left\langle 100\right\rangle $ 0.606±0.009 18±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.941±0.011 18.46±0.11 气炮 2025 $ \left\langle 100\right\rangle $ 0.833±0.013 18.2±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.917±0.011 18.29±0.11 气炮 2025 $ \left\langle 100\right\rangle $ 0.554±0.009 18±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.847±0.011 18.41±0.19 气炮 2025 $ \left\langle 100\right\rangle $ 0.562±0.009 18.1±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.939±0.012 18.49±0.17 气炮 2025 $ \left\langle 100\right\rangle $ 0.647±0.01 18±0.3
    气炮 2018 $ \left\langle 100\right\rangle $ 0.917±0.012 18.28±0.17 气炮 2025 $ \left\langle 100\right\rangle $ 0.576±0.009 18±0.3
    气炮 2018 $ \left\langle 110\right\rangle $ 0.867±0.011 20.4±0.19 气炮 2025 $ \left\langle 100\right\rangle $ 0.572±0.009 18±0.3
    气炮 2018 $ \left\langle 110\right\rangle $ 1.128±0.011 20.54±0.11 气炮 2025 $ \left\langle 100\right\rangle $ 0.304±0.005 17.8±0.3
    气炮 2018 $ \left\langle 111\right\rangle $ 0.951 ± 0.011 20.67±0.17 气炮 2025 $ \left\langle 100\right\rangle $ 0.328±0.005 17.8±0.3
    [54] 气炮 2020 $ \left\langle 100\right\rangle $ 1.25858 18.4±0.2 气炮 2025 $ \left\langle 100\right\rangle $ 0.366±0.006 17.9±0.3
    气炮 2020 $ \left\langle 100\right\rangle $ 0.84563 18.1±0.4 气炮 2025 $ \left\langle 111\right\rangle $ 0.819±0.012 20±0.3
    气炮 2020 $ \left\langle 110\right\rangle $ 1.28161 20.6±0.2 气炮 2025 $ \left\langle 111\right\rangle $ 0.44±0.007 19.3±0.3
    气炮 2020 $ \left\langle 110\right\rangle $ 1.78142 21.4±0.4 [23] 激光 2010 多晶 2.6±0.44 21.9±1
    气炮 2020 $ \left\langle 110\right\rangle $ 1.93615 21.6±0.5 激光 2010 $ \left\langle 100\right\rangle $ 1.33±0.13 19.29±0.8
    气炮 2020 $ \left\langle 110\right\rangle $ 2.0785 21.9±0.3 激光 2010 $ \left\langle 100\right\rangle $ 1.4±0.05 19.65±0.68
    气炮 2020 $ \left\langle 110\right\rangle $ 2.5806 22.6±0.5 激光 2010 $ \left\langle 110\right\rangle $ 2.58±0.32 21.5±1.6
    气炮 2020 $ \left\langle 110\right\rangle $ 1.48373 20.9±0.4 激光 2010 $ \left\langle 110\right\rangle $ 2.75±0.14 21.6±1.7
    气炮 2020 $ \left\langle 110\right\rangle $ 2.37292 22.3±0.3 激光 2010 $ \left\langle 110\right\rangle $ 1.107±0.071 20.72±0.68
    气炮 2020 $ \left\langle 111\right\rangle $ 0.97336 20.2±0.2 激光 2010 $ \left\langle 110\right\rangle $ 2.44±0.24 21.5±1.4
    气炮 2020 $ \left\langle 111\right\rangle $ 1.39708 20.8±0.4 激光 2010 $ \left\langle 111\right\rangle $ 0.88±0.11 20.74±0.3
    气炮 2020 $ \left\langle 111\right\rangle $ 0.91596 20.1±0.2 激光 2010 $ \left\langle 111\right\rangle $ 0.828±0.043 20.75±0.36
    气炮 2020 $ \left\langle 111\right\rangle $ 2.31926 22.1±0.5 激光 2010 $ \left\langle 111\right\rangle $ 0.932±0.093 20.81±0.4
    气炮 2020 $ \left\langle 111\right\rangle $ 0.95815 20.1±0.4 激光 2010 $ \left\langle 111\right\rangle $ 1.83±0.11 22.43±0.62
    气炮 2020 $ \left\langle 111\right\rangle $ 1.27341 20.6±0.2
    注: 1) Lang等[22]的单波结构实验数据未被收录, 因其无法明确界定是否达到弹性极限状态;
    2) McWilliams等[23]的激光加载实验中的双层堆垛金刚石的数据未被收录, 因为界面处的应力松弛使得两块样品中应力载荷不一致.
    下载: 导出CSV
  • [1] Tang Y, Wang R X, Xiao B, Zhang Z R, Li S, Qiao J W, Bai S X, Zhang Y, Liaw P K 2023 Prog. Mater. Sci. 135 101090 doi: 10.1016/j.pmatsci.2023.101090
    [2] Flanagan R M, Fensin S J, Meyers M A 2022 Appl. Phys. Rev. 9 011305 doi: 10.1063/5.0053693
    [3] Jiang Z J, He J Y, Wang H Y, Zhang H S, Lu Z P, Dai L H 2016 Mater. Res. Lett. 4 226 doi: 10.1080/21663831.2016.1191554
    [4] Chen M W, McCauley J W, Dandekar D P, Bourne N K 2006 Nat. Mater. 5 614 doi: 10.1038/nmat1689
    [5] Li Y, Cheng X W, Ma Z L, Li X H, Wang M 2021 J. Mater. Sci. Technol. 93 119 doi: 10.1016/j.jmst.2021.03.052
    [6] Euser V K, Jones D R, Martinez D T, Valdez J A, Trujillo C P, Cady C M, Fensin S J 2023 Acta Mater. 250 118874 doi: 10.1016/j.actamat.2023.118874
    [7] Zhao S P, Feng Z D, Li L X, Zhao X J, Lu L, Chen S, Zhang N B, Cai Y, Luo S N 2023 J. Mater. Sci. Technol. 134 178 doi: 10.1016/j.jmst.2022.05.060
    [8] Dai Z H, Lu L, Chai H W, Xiao X H, Gong X L, Luo S N 2020 Mater. Sci. Eng. : A 789 139690 doi: 10.1016/j.msea.2020.139690
    [9] 孙占峰, 贺红亮, 李平, 李庆忠 2012 物理学报 61 096201 doi: 10.7498/aps.61.096201 Sun Z F, He H L, Li P, Li Q Z 2012 Acta Phys. Sin. 61 096201 doi: 10.7498/aps.61.096201
    [10] Weng J D, Tan H, Wang X, Ma Y, Hu S L, Wang X S 2006 Appl. Phys. Lett. 89 111101 doi: 10.1063/1.2335948
    [11] Cheng Y H, Song C M, Tan Y Z, Yue S L, Li G, Qiu Y Y, Xu T H, Liu Y, Gao F, Yu Q 2021 Shock Vib. 2021 6677077 doi: 10.1155/2021/6677077
    [12] Meyers M A 1994 Dynamic Behavior of Materials (New York: Wiley) p668
    [13] 李俊, 吴强, 于继东, 谭叶, 姚松林, 薛桃, 金柯 2017 物理学报 66 146201 doi: 10.7498/aps.66.146201 Li J, Wu Q, Yu J D, Tan Y, Yao S L, Xue T, Jin K 2017 Acta Phys. Sin. 66 146201 doi: 10.7498/aps.66.146201
    [14] Graham R A 1993 Solids under High-Pressure Shock Compression (New York: Springer-Verlag) pp 32
    [15] Liu X, Mashimo T, Ogata K, Kinoshita T, Sekine T, Zhou X M, Nellis W J 2013 J. Appl. Phys. 114 243511 doi: 10.1063/1.4854435
    [16] Taylor J W, Rice M H 1963 J. Appl. Phys. 34 364 doi: 10.1063/1.1702614
    [17] Barker L M 1975 J. Appl. Phys. 46 2544 doi: 10.1063/1.321931
    [18] 冯晓伟, 李俊承, 王洪波, 常敬臻 2016 物理学报 65 166201 doi: 10.7498/aps.65.166201 Feng X W, Li J C, Wang H B, Chang J Z 2016 Acta Phys. Sin. 65 166201 doi: 10.7498/aps.65.166201
    [19] Fratanduono D E, Munro D H, Celliers P M, Collins G W 2014 J. Appl. Phys. 116 033517 doi: 10.1063/1.4890014
    [20] Zhao R, Xu R Q, Yang B, Shen Z H, Lu J, Ni X W 2006 Appl. Opt. 45 265 doi: 10.1364/AO.45.000265
    [21] Yu H, Yoh J J 2016 AIP Adv. 6 055314 doi: 10.1063/1.4950855
    [22] Lang J M, Winey J M, Gupta Y M 2018 Phys. Rev. B 97 104106 doi: 10.1103/PhysRevB.97.104106
    [23] McWilliams R S, Eggert J H, Hicks D G, Bradley D K, Celliers P M, Spaulding D K, Boehly T R, Collins G W, Jeanloz R 2010 Phys. Rev. B 81 014111 doi: 10.1103/PhysRevB.81.014111
    [24] Smith R F, Volz T J, Celliers P M, Braun D G, Swift D C, Gorman M G, Briggs R, Fernandez-Pañella A, Kirsch L E, Marshall M C, McNaney J M, Eggert J H, Fratanduono D E, Ali S J 2023 Rev. Sci. Instrum. 94 083003 doi: 10.1063/5.0150031
    [25] Shui M, Chu G B 2015 Chin. Phys. B 24 94701 doi: 10.1088/1674-1056/24/9/094701
    [26] Oleynik I I, Landerville A C, Zybin S V, Elert M L, White C T 2008 Phys. Rev. B 78 180101 doi: 10.1103/PhysRevB.78.180101
    [27] Zhou X M, Li J, Nellis W J, Wang X, Li J B, He H L, Wu Q 2011 J. Appl. Phys. 109 083536 doi: 10.1063/1.3575330
    [28] Yang Z, Xian J W, Gao X Y, Tian F Y, Song H F 2024 J. Chem. Phys. 161 194101 doi: 10.1063/5.0231337
    [29] Liu Z J, Qi J H 2007 Chin. Phys. B 16 499 doi: 10.1088/1009-1963/16/2/035
    [30] Yang Z, Xian J W, Xu Y, Lin D Y, Wang Q C, Gao X Y, Tian F Y, Song H F 2025 Comput. Phys. Commun. 316 109814 doi: 10.1016/j.cpc.2025.109814
    [31] Yang Z, Xian J W, Gao X Y, Tian F Y, Song H F 2024 Chin. Phys. B 33 116401 doi: 10.1088/1674-1056/ad73b0
    [32] Wang Y, Wang J J, Zhang H, Manga V R, Shang S L, Chen L Q, Liu Z K 2010 J. Phys. Condens. Matter 22 225404 doi: 10.1088/0953-8984/22/22/225404
    [33] Robert G, Legrand P, Bernard S 2010 Phys. Rev. B 82 104118 doi: 10.1103/PhysRevB.82.104118
    [34] Wang Y 2000 Phys. Rev. B 61 R11863 doi: 10.1103/PhysRevB.61.R11863
    [35] Song H F, Liu H F 2007 Phys. Rev. B 75 245126 doi: 10.1103/PhysRevB.75.245126
    [36] Gao X Y, Wang W Y, Chen X, Chong X Y, Xian J W, Tian F Y, Wang L F, Chen H J, Liu Y, Huang H B, Song H F 2025 Mater. Genome Eng. Adv. 3 e70029 doi: 10.1002/mgea.70029
    [37] Kresse G, Furthmüller J 1996 Comput. Mater. Sci. 6 15 doi: 10.1016/0927-0256(96)00008-0
    [38] Kresse G, Furthmüller J 1996 Phys. Rev. B 54 11169 doi: 10.1103/PhysRevB.54.11169
    [39] Perdew J P, Burke K, Ernzerhof M 1996 Phys. Rev. Lett. 77 3865 doi: 10.1103/PhysRevLett.77.3865
    [40] Blöchl P E 1994 Phys. Rev. B 50 17953 doi: 10.1103/PhysRevB.50.17953
    [41] Perdew J P, Ruzsinszky A, Csonka G I, Vydrov O A, Scuseria G E, Constantin L A, Zhou X, Burke K 2008 Phys. Rev. Lett. 100 136406 doi: 10.1103/PhysRevLett.100.136406
    [42] Blöchl P E, Jepsen O, Andersen O K 1994 Phys. Rev. B 49 16223 doi: 10.1103/PhysRevB.49.16223
    [43] Methfessel M, Paxton A T 1989 Phys. Rev. B 40 3616 doi: 10.1103/PhysRevB.40.3616
    [44] Hu Y K, Gao X Y, Zhao Y F, Liu X, Song H F 2022 Phys. Rev. B 106 104101 doi: 10.1103/PhysRevB.106.104101
    [45] Wang P, Zhang C, Jiang S E, Duan X X, Zhang H, Li L L, Yang W M, Liu Y G, Li Y G, Sun L, Liu H, Wang Z B 2021 Matter Radiat. Extrem. 6 035902 doi: 10.1063/5.0039062
    [46] Reeber R R, Wang K 1996 J. Electron. Mater. 25 63 doi: 10.1007/BF02666175
    [47] Thewlis J, Davey A R 1956 Philos. Mag. 1 409 doi: 10.1080/14786435608238119
    [48] Dewaele A, Datchi F, Loubeyre P, Mezouar M 2008 Phys. Rev. B 77 094106 doi: 10.1103/PhysRevB.77.094106
    [49] Zouboulis E S, Grimsditch M, Ramdas A K, Rodriguez S 1998 Phys. Rev. B 57 2889 doi: 10.1103/PhysRevB.57.2889
    [50] Hu Q Y, Li B S, Gao X, Bi Y, Su L, Mao H K 2021 Proc. Natl. Acad. Sci. 118 e2118490118 doi: 10.1073/pnas.2118490118
    [51] Feng J, Xiao B, Zhou R, Pan W 2013 Acta Mater. 61 7364 doi: 10.1016/j.actamat.2013.08.043
    [52] McSkimin H J, Andreatch P 1972 J. Appl. Phys. 43 2944 doi: 10.1063/1.1661636
    [53] Kondo K, Ahrens T J 1983 Geophys. Res. Lett. 10 281 doi: 10.1029/GL010i004p00281
    [54] Winey J M, Knudson M D, Gupta Y M 2020 Phys. Rev. B 101 184105 doi: 10.1103/PhysRevB.101.184105
    [55] Mashimo T, Liu X, Tokuda M, Kawai N, Isobe H, Tsushida M, Yamamuro T, Kurio F, Sumiya H 2025 Phys. Rev. B 112 224111 doi: 10.1103/gs3m-hbb4
  • 加载中
图( 5) 表( 1)
计量
  • 文章访问数:  288
  • HTML全文浏览数:  288
  • PDF下载数:  1
  • 施引文献:  0
出版历程
  • 收稿日期:  2026-01-11
  • 刊出日期:  2026-06-20

结合高压热弹性模拟数据的Hugoniot弹性极限实验解读

    通讯作者: E-mail: gao_xingyu@iapcm.ac.cn.; 
    通讯作者: E-mail: song_haifeng@iapcm.ac.cn.; 
    通讯作者: E-mail: pxzhang@nwpu.edu.cn.
    作者简介: 尹俊磊: junlei_yin@mail.nwpu.edu.cn
  • 1. 西北工业大学, 凝固技术全国重点实验室, 西安　710072
  • 2. 北京应用物理与计算数学研究所, 计算物理全国重点实验室, 北京, 100088
  • 收稿日期:  2026-01-11

摘要: 针对冲击加载实验中波前弥散及非稳态效应制约Hugoniot 弹性极限(HEL)精确解读的问题, 本文提出了一种结合理论模拟与实验数据的自洽解读方法. 该方法基于宏观粒子动量守恒定律与微观材料本征纵波声速响应的物理关联, 构建了以纵波声速和轴向应力为核心变量的自洽方程. 在求解策略上, 以实验测得的弹性前驱波粒子速度为约束变量, 并利用平均场势(MFP)耦合准静态近似(QSA)生成的宽压域热弹性数据集, 构建了连续的纵波声速-压力响应关系. 本文以金刚石为研究对象验证该方法的有效性, 首先构建了室温下0—1000 GPa压区的金刚石热弹性数据集, 其与15 GPa压力范围内静加载实验数据的相对偏差不超过0.4%. 在此基础上, 整合了峰值应力高达1 TPa的金刚石冲击实验数据, 应用上述方法重新解读了HEL. 结果表明, 重新解读的HEL与标准气炮加载实验值的相对偏差在2%以内, 并定量表征了其各向异性. 此外, 针对激光加载实验中因非稳态衰减导致的系统性偏差, 本方法将其HEL值向上修正了8%—13%. 本文数据集可在https://doi.org/10.57760/sciencedb.j00213.00274中访问获取.

English Abstract

    全文HTML

    1.   引　言

    • 理解材料在爆炸或激光烧蚀等强载荷驱动下的力学响应, 是先进防护设计、航空航天结构安全、高能物理及行星科学等领域的关键科学问题[1-3]. Hugoniot弹性极限(HEL)作为材料从纯弹性响应向弹塑性响应转变的临界轴向应力阈值[4-6], 是描述材料动态力学行为的核心参数[7-9]. 通常, HEL的测定主要基于位移干涉仪系统(DISAR)或速度干涉仪系统(VISAR)所记录的自由面速度剖面特征[10]. 当冲击载荷超过材料的HEL时, 入射波在材料内部发生分离, 形成传播较快的弹性前驱波与传播较慢的塑性波. 理想的自由面速度剖面呈现出典型的“双波”结构(见图1(a)): 弹性前驱波率先到达, 驱动自由面速度跃升至第一平台, 该平台对应材料发生宏观屈服(如原子重排或缺陷成核)前的临界状态; 随后到达的塑性波则将速度进一步推升至更高平台. 该双波结构揭示了材料内禀力学响应与波速间的内在关联[11], 并为HEL的计算提供了物理依据, 即利用弹性前驱波对应的粒子速度与纵波声速进行反演[12,13](如无特殊说明, 本文中粒子速度与纵波声速均为HEL临界状态下物理量).

      然而, 实际测量中弹性前驱波上升沿的波前弥散影响了纵波声速的判定精度, 是HEL测量误差的来源之一. 与理想的阶跃函数响应不同, 实际的弹性前驱波表现为具有有限上升时间的连续变化. 这种波前弥散可能存在多种潜在的诱因, 例如材料内部的残余应力或微观结构(如成分)的非均匀性等[14,15]. 由于物理机制的复杂性, 弹性波抵达自由面的确切时刻存在定义上的模糊性, 在低冲击载荷下该问题更加突出. 对于低声速材料, 为避免塑性过驱, HEL的研究往往受限于较低的冲击载荷区间, 从而不可避免地放大了上述弥散带来的判读难题. 以图1(b)所示的铁的弹-塑性响应为例, 不同判读标准引入了声速计算偏差. Taylor和Rice[16]以自由面起跳时刻为基准, 计算得出纵波声速为6.18 km/s; Barker[17]则主张选取上升沿的中点时刻, 所得结果则为6.00 km/s. 此类由波形处理方法引起的纵波声速差异, 将直接影响HEL测量的准确性.

      除波前弥散外, 加载波形的非稳态特性也会影响声速判读[18], 这个问题在短时加载的激光驱动实验中尤为关键[19,20]. 除加载源自身扰动可能加剧波前弥散外, 激光加载产生的冲击波在传播中会发生衰减[21], 导致“波在介质中以恒定相速度传播”的假设失效, 通过走时测量的平均声速随之产生偏差. 例如, Lang等[22]指出, 受冲击波衰减的影响, McWilliams等[23]通过激光加载测得的金刚石HEL相较于气炮加载实验数据偏低7%—11%. 尽管存在上述局限, 激光加载凭借其极高的压力覆盖范围和效费比优势, 仍是冲击波物理研究的重要手段. 因此, 建立一套能够修正非稳态效应的声速解读方法, 对于消除实验误差、准确评估材料动态响应具有重要意义.

      针对此问题, 研究者们发展了多种实验修正技术. 例如, 通过多步采样方法[24]在时间演化维度上捕捉前驱波随传播距离的衰减规律, 从而消除应力松弛带来的判读偏差; 利用线成像VISAR[25]等技术在空间横向维度上解读冲击波阵面的平整度与倾斜度, 剥离了因几何非理想性导致的表观波前弥散. 但是上述纯光学的修正手段难以完全消除由材料内禀属性或动力学演化引起的信号失真. 针对非稳态冲击波引入的复杂性, Fratanduono等[19] 提出一种修正策略, 通过在样品旁设置透明窗口实时监测扰动的传播行为, 建立了可观测的冲击波速度变化与压力扰动之间的定量关联. 然而, 上述方法普遍依赖于复杂的实验构型或特殊的靶材制备工艺, 推高了实验成本, 限制了其在常规测试中的广泛应用.

      鉴于此, 本研究提出了一种融合理论模拟与实验数据的自洽解读方法. 该方法将材料的本征纵波声速-压力依赖关系与宏观动量守恒定律相结合, 构建了描述轴向应力的自洽方程. 在此基础上, 利用第一性原理计算, 结合平均场势(MFP)和准静态近似(QSA)理论, 构建了材料的高压热弹性基础数据集, 并通过插值处理为自洽方程的求解提供了连续的纵波声速-压力响应曲线. 最后, 选取测量精度较高的弹性前驱波粒子速度作为关键约束条件, 通过自洽迭代算法反演求解HEL值及其对应的纵波声速. 该方法建立了HEL与弹性前驱波粒子速度之间的直接映射机制, 有效避免了因非理想波形(如弥散)导致的声速判读误差, 为极端条件下材料动态力学行为的精确表征提供了新途径. 该方法的独特优势在于仅需后处理即可修正受非理想效应影响的数据. 这使得我们可以重新分析早期文献中的数据(如本文对McWilliams 数据的修正), 或在受限的实验条件下获得高精度的本构参数.

      本文选取金刚石作为模型材料验证了该方法的有效性. 我们梳理了现有的静态高压弹性数据(压强至15 GPa), 以及涵盖宽载荷区间的气炮(峰值应力至 1 TPa)和激光加载(峰值应力至1 TPa)实验速度剖面数据. 基于MFP与 QSA方法, 本文构建了室温下0—500 GPa压力区间的金刚石热弹性数据集, 其弹性模量预测值与静态实验数据的相对偏差在0.4%以内. 应用该数据集, 重新评估了金刚石在气炮加载下的HEL, 定量解读了其各向异性特征, 结果与实验基准值的偏差小于2%. 进一步利用该方法修正了激光加载实验中因非稳态冲击效应引入的系统性偏差, 为高精度解读激光冲击数据提供了理论依据.

    2.   模型与方法

      2.1.   应力-声速自洽模型

    • 在单轴应变加载条件下, 波阵面扫过前后的宏观状态遵循动量守恒定律, 其关系式可表示为

      其中, $ {\sigma }_{0} $和$ {u}_{0} $分别是波阵面抵达前的轴向应力和粒子速度, $ \sigma $和$ u $分别是波阵面扫过后的轴向应力和粒子速度, $ {\rho }_{0} $是材料的初始密度, $ {C}_{\text{L}} $是波阵面的传播速度, 即纵波声速[13]. 该守恒定律描述了平面冲击下轴向应力、粒子速度与波速之间的普适关系, 因此同样适用于弹性极限状态. 在理想情况下, 弹性前驱波被视为无限薄的间断界面, 波阵面前的粒子处于静止未扰动状态, 即轴向应力和粒子速度均为零, 波阵面后的物理量发生瞬时跳变. 由此可以推导出弹性极限处的轴向应力 ($ {\sigma }_{\text{HEL}} $) 与纵波声速和粒子速度($ {u}_{\text{E}} $)的乘积成正比: $ {\sigma }_{\text{HEL}}={\rho }_{0}{u}_{\text{E}}{C}_{\text{L}} $. 这也是传统实验中基于双波结构反演HEL的理论基础[12,14].

      然而, 实际实验观测到的弹性前驱波通常表现为有限厚度的应力扰动, 其内部物理量(应力、密度等)经历剧烈但连续的变化. 这个应力扰动持续向前传播, 其宏观轴向迁移速率就是冲击纵波声速. 当波前应力扰动扫过材料时, 晶格发生压缩致使纵向模量硬化[26]. 这意味着, 位于波阵面后方高应力区的扰动将以更高的声速向前传播, 并呈现出一种动态“追赶”效应, 驱动波形向陡峭化方向演化. 这一非线性硬化机制解释了为何冲击载荷下的宏观波速通常高于常压声速, 且表现出显著的轴向应力依赖性[27]. 受此机制驱动, 波阵面内部的动力学演化最终将达到一个宏观上波形稳定、速度恒定的自洽传播状态. 在此状态下, 波后的轴向应力与维持其高纵向模量及高纵波声速所需的应力水平相一致. 弹性极限状态下的自洽条件可表述为

      自洽方程(2)的建立是基于HEL临界状态下的纵波声速. 在达到HEL之前, 纵波声速主要依赖于冲击加载轴向应力; 在达到HEL之后, 基于连续介质假设, 纵波声速主要依赖于冲击波后的热力学状态. 考虑到冲击波前沿的传播速度远高于热扩散速度, 弹性前驱波的传播可视为绝热过程, 其热力学路径近似为从室温常压状态出发的等熵压缩[11]. 因此, 二者的等熵声速是相等的, 可以使用室温压缩线上的等熵声速作为方程(2)中的HEL状态下的纵波声速, 即: $ C_{\text{L}}^{\text{shock}}\left({\sigma }_{\text{HEL}}\right)=C_{\text{L}}^{\text{hydro}}\left(P\right) $. 通过从速度剖面读取粒子速度, 结合纵波声速关于压力的函数, 求解自洽方程(2)即可同时获得HEL及其对应的纵波声速.

      需要指出的是, 弹性前驱波粒子速度本质上是材料动态屈服行为的宏观体现, 其幅值与加载应变率密切相关. 在更高的应变率下, 材料可以维持更高的粒子速度而不发生原子重排或缺陷成核. 自洽方程(2)虽然形式上未显式包含应变率项, 但通过引入实验测得的弹性波粒子速度可反映应变率效应. 因此, 通过求解自洽方程解读的HEL结果已经包含了对应加载条件下的应变率效应.

      值得注意的是, 本研究提出的理论-实验自洽解读方法仅适用于密实的金属材料和无机非金属材料, 前提是实验能够获得具有双波结构的速度剖面. 对于一些声速较低的材料, 即便在较低的冲击载荷下也会发生塑性过驱动, 速度剖面会跳过HEL阶段呈现塑性单波; 对于一些冲击载荷较低的实验条件, 冲击波不足以达到HEL临界状态, 速度剖面会呈现具有弹性波平台的弹性单波. 而上述自洽方程(2)仅对双波剖面的HEL临界状态有效.

      • 2.2.   高压热弹性第一性原理计算

    • 纵波声速作为材料的本征物理属性, 可通过第一性原理模拟获得. 为了在模拟中平衡有限温度下高压热弹性的计算效率与精度, 本研究采用准静态近似(QSA)框架下的能量-应变方法[28,29]. 其核心假设在于: 材料在有限温度下的热弹性行为主要由热膨胀引起的体积变化主导, 而离子振动和电子热激发等对应变自由能的直接贡献则相对较小, 可以忽略. 基于此近似, 可以将有限温度和压力$ \left(P, T\right) $下由应变引起的亥姆霍兹自由能变化$ \text{Δ}F $近似为在对应体积 $ V\left(P, T\right) $下由相同应变引起的零温静态能量变化$ \text{Δ}E $:

      进而, 有限温下的弹性常数$ {C}_{ij}\left(P, T\right) $可近似为在体积$ V\left(P, T\right) $下的零温弹性常数:

      其中$ \epsilon $是Voigt标记法下的应变分量.

      根据2.1节的分析, 自洽方程(2)求解需要等熵声速, 为了从等温弹性常数$ C_{ij}^{T} $ 得到等熵弹性常数$ C_{ij}^{S} $, 本研究引入热膨胀对弹性常数的绝热修正项[30]:

      其中$ \lambda $和$ {C}_{\eta } $满足以下关系:

      $ {\alpha }_{j} $是沿$ j $方向的线热膨胀系数, $ {C}_{P} $为定压热容.

      在QSA框架下, 将不同静水压力映射到对应体积$ V\left(P, T\right) $的函数关系, 需要从精确的状态方程获得. 准谐近似(QHA)是基于第一性原理构建状态方程的广泛应用方法[31]. 例如Wang等[32]利用QHA方法构建了Al, Cu, Ni, Mo, Ta, NiAl和Ni3Al等多种金属及合金在近常压下的状态方程, 并通过QSA计算得到了与实验高度一致的热弹性常数. 然而, 完全的QHA需要在多个离散的体积点上进行完整的声子谱计算, 这使得QHA在构建宽压域状态方程时可能面临巨大的计算开销[33].

      为克服这一效率限制, 本研究采用平均场势(MFP)方法来计算离子的振动自由能贡献[34,35]. MFP方法的核心思想是, 将离子与周围所有其他原子的复杂多体相互作用, 近似为一个等效平均场势$ g $. 离子偏离其平衡位置$ V $达到距离$ r $时所受到的势能可通过第一性原理计算的基态能量$ E_{\rm c} $进行构建:

      将平均场势$ g\left(r,V\right) $代入玻尔兹曼因子, 并在空间中积分, 即可得到离子在该势阱中可活动的自由体积$ {v}_{\text{f}} $:

      最后根据粒子在势场中的配分函数分布, 离子的振动自由能$ {F}_{\text{ion}} $可以表示为

      其中$ M $为离子热运动的粒子质量. 该方法避免了声子计算, 提升了宽压域状态方程的构建效率.

      • 2.3.   HEL解读流程

    • 图2展示了该HEL解读方法的完整计算流程, 该流程已集成到ProME极端材料性能集成计算平台[36]. 通过三次样条插值离散的第一性原理热弹性数据, 构建在压力空间内连续且一阶导数光滑的纵波声速响应函数$ {C}_{\text{L}}(P) $, 将其代入方程 (2), 便可以迭代求解对应于特定粒子速度的HEL和纵波声速. 该迭代过程以0 GPa作为轴向应力初始值, 物理上模拟了后方高应力扰动追赶前方低压波前的动力学过程. 方程的自洽解在物理意义上对应了高速纵波的应力条件与波后轴向应力相匹配的状态. 自洽迭代收敛判据设定为前后两次迭代的轴向应力差$ \mid \text{Δ}P\mid \lt {10}^{-3}\;\text{GPa} $. 该判据的选择依据在于$ {10}^{-3}\;\text{GPa} $的精度高于目前实验技术对HEL的分辨能力[22].

      对于HEL解读流程涉及的第一性原理计算, 采用Vienna从头算模拟程序包 (VASP)[37,38], 并使用针对固体优化的Perdew-Burke-Ernzerhof (PBEsol)广义梯度近似(GGA)[39-41]. 波函数在Γ中心Monkhorst-Pack网格上采样, KSPACING设为0.12. 平面波截断能和电子自洽能量收敛判据分别设为600 eV和10–6 eV/atom, 以保证计算精度. 所有结构在固定体积下通过基于力的投影梯度下降方法进行完全弛豫[42,43], 而最终总能计算则采用包含Blöchl修正的四面体方法[44].

    3.   结果与讨论

    • 金刚石是高压物性研究实验中的基准材料[45], 屈服阈值高, 是目前公认的具有最高抗塑性流变阻力的材料之一. 在气炮加载高压实验中, 金刚石弹性前驱波与后续塑性波分离清晰, 弹性波上升沿弥散对波速解读的影响较小. 因此, 我们针对金刚石材料峰值应力达1 TPa的气炮冲击加载实验, 以基于速度剖面中弹性波走时计算得到的HEL作为参考, 对本文的HEL解读方法进行验证. 相对于金刚石声速较低的材料, 波前弥散可能受微观非均匀性、化学成分起伏等材料内禀属性的影响. 因此, Taylor和Rice[16]以及Barker[17]对于铁的HEL实验解读存在系统性的分歧, 不能提供直接的验证. 但是理论上本研究提出的方法对这种低声速材料是适用的.

      • 3.1.   物态方程数据集

    • 本研究中HEL的解读依赖于精确的状态方程和热弹性数据集. 本工作首先采用QHA+MFP的模拟框架, 构建了金刚石覆盖宽广温区(298—3000 K)与压区(0—1000 GPa)的P-V-T 数据集, 以充分验证模型在极端条件下的适用性. 针对后续的HEL分析, 截取了298—1000 K数据子集. 图3(a)描绘了零压下体积随温度的演化(等压膨胀线). 结果显示, 与文献[46,47]的实验测量结果相比, 理论模拟结果与实验数据最大偏差不超过0.3%. 对于冲击加载环境, 极端压缩条件下状态方程的适用性检验尤为关键. 图3(b)展示了常温(298 K)及高温(900 K)下的等温压缩曲线, 并对比了Dewaele等[48]利用电阻加热金刚石压砧(DAC)获得的X射线衍射数据. 在0—60 GPa的实验压力范围内, 本研究模拟得到的体积与实验基准的相对偏差在0.75%以内. 综上所述, 该数据集在宽广的热力学空间内展现了良好的一致性, 足以为后续热弹性模拟和HEL解读提供高精度输入.

      • 3.2.   高压热弹性数据集

    • 在精确构建P-V-T 状态方程数据集的基础上, 本研究进一步在QSA框架下采用能量-应变法, 计算了金刚石在宽温压范围内的弹性常数矩阵. 为全面评估数据集的可靠性, 分别从热软化效应与压致硬化效应两个维度进行严格的基准测试.

      针对弹性常数的温度依赖性, 图4(a)展示了零压下弹性常数C11和C44随温度的演化规律. 模拟结果表明, 随着温度升高至1500 K, 金刚石的晶格刚度呈现出单调的线性衰减趋势. 与Zouboulis等[49]的布里渊散射实验数据对比可见, QSA 方法捕捉了这一热软化行为. 在整个测试温区内, 理论预测值与实验测量值的相对偏差始终控制在4%以内. 考虑到高温实验中热辐射与样品非均匀性带来的测量不确定度, 这一吻合度验证了QSA方法在有限温度下的预测精度

      材料的宏观热弹性通过 Voigt-Reuss-Hill (VRH)近似由弹性常数导出. 图4(b)展示了室温下剪切模量、体积模量随压力的变化. 与Hu等[50]的高压超声测量数据对比, 剪切模量与体积模量的模拟结果与实验的相对偏差小于0.75%. 此外, 为验证计算体系的内部自洽性, 本研究还通过状态方程(EOS)对体积求二阶导数计算了体积模量. 结果显示, 这两种独立的方法所得体积模量高度一致, 偏差小于0.8%, 直接证实了本研究中热力学状态方程与力学本构关系之间的高度自洽性.

      更为关键的是, 纵波声速作为反演HEL的直接输入参数, 其精度决定了最终结果的置信度. 图4(b)同时展示了环境条件下等温压缩与等熵压缩的纵波声速. 由于金刚石极高的热导率与德拜温度, 其热膨胀系数极小, 其等温与等熵声速的差异仅约为0.5%. 对比Hu等[50]在50 GPa范围内的实验结果, 理论模拟预测的等温纵波声速曲线几乎完全落在实验误差带内, 最大相对偏差仅为0.4%. 这一精度比气炮加载实验中速度剖面解读的不确定度低约一个数量级.

      • 3.3.   HEL解读

    • 自洽方程的求解不仅需要实验测得的粒子速度作为边界约束, 还需要纵波声速与静水压力(或体积)之间的本构映射关系. 本研究利用三次样条插值高压热弹性数据集, 构建了高精度的纵波声速响应函数, 从而确保了声速及其对压力的一阶导数在整个热力学状态空间内的连续性与光滑性. 图5(a)展示了由此计算的纵波声速随粒子速度的演化关系. 值得注意的是, 单晶声速的各向异性由Christoffel方程的解导出[51]:

      其中, $ v $是声速, $ {u}_{j} $是位移矢量, $ {C}_{ikjm} $是刚度矩阵, $ {l}_{k} $和$ {l}_{m} $是方向余弦. 对于立方晶体, 存在3个独立的声波传播主方向: [100], [110]和[111]. 其对应的纵波声速分别为

      当粒子速度为0 km/s时(即准静态条件), 模拟的声速与实验值偏差不超过0.4%[52].

      为全面验证理论预测的准确性, 本研究在表1中收集了截至目前已发表的金刚石冲击实验数据[22,23,53–55], 并绘制在图5中. 可以看出, 纵波声速随冲击载荷的增强而非线性地增大. 在相同粒子速度下, 3个晶向的纵波声速呈现出清晰的各向异性: [111] > [110] > [100]. 这解释了为何由Winey等[54]以及Kondo和Ahrens[53]沿[100]方向进行的冲击实验更容易发生塑性过驱, 最终呈现为单波结构.

      将模拟结果与气炮实验数据对比, 代表多晶的计算曲线完全落在实验解读值的不确定度区间内, 与Kondo和Ahrens[53]的实验结果相比, 最大偏差不超过2%. 对比单晶金刚石的冲击实验, 模拟的单晶纵波声速基本落在Lang等[22]的气炮实验数据置信区间内, 与实验测量值的偏差不超过0.8%. 值得注意的是, 最近Mashimo等[55]采用大厚度(1.52—1.92 mm)高纯型金刚石单晶(氮含量<0.1 ppm)进行了系统的冲击压缩实验, 获得了稳态HEL值((20.5±2.0) GPa), 显著低于此前薄样品(0.08—0.51 mm)报道的高HEL值(53—144 GPa). 本研究理论预测的纵波声速与[100]和[111]晶向上的实验数据偏差分别不超过0.3%和1.2%. 然而, 此前HEL数据的分散性不仅源于样品厚度差异导致的非稳态效应, 更本质地反映了HEL对冲击载荷(或应变率)的固有依赖性. Mashimo等[55]实验结果相对集中, 部分得益于其冲击载荷范围相对较窄(44—108 GPa), 使得应变率效应的影响趋于稳定. 本研究通过建立HEL与粒子速度的直接映射, 将Mashimo等[55]揭示的厚度效应与HEL的本征载荷依赖性统一于同一理论框架. 二者并非互斥, 厚度效应解释了短传播距离内的非稳态衰减, 而载荷依赖性则揭示了材料在宽广动态范围内的本构响应.

      与气炮实验相比, McWilliams等[23]的激光加载数据显示出较高的不确定度和离散性[23]. 针对沿[111]和[100]晶向的低冲击载荷实验, 模拟结果基本位于实验置信区间内, 与实验点估计值的最大偏差分别为0.8%和4%. 对于沿[110]晶向且粒子速度超过1.3 km/s的冲击, Lang等[22]的气炮加载数据比McWilliams等[23]的激光实验结果高7%—11%. Lang等[22]将此偏差归因于激光加载中非稳态冲击波. 实验HEL测定依赖于弹性前驱波的走时. 在激光加载下, 波的非稳态衰减延长了测量走时, 从而低估了平均声速. 我们的方法基于动量守恒, 利用波阵面抵达时刻的局部粒子速度与材料纵波声速的映射关系, 独立于冲击波的传播路径, 解耦了波幅衰减效应与材料热力学响应. 经修正, 激光加载条件下的纵波声速提升了8%—13%.

      随后, 将得到的自洽纵波声速与实验粒子速度相结合, 我们重估了相应的 HEL 值, 如图5(b)所示. 由于重估的HEL与纵波声速是自洽求解的, 其与实验直接解读值的偏差也表现出一致性. 与气炮实验相比, 多晶的相对偏差不超过2%, 单晶不超过1.2%. 对于激光加载实验, 本方法将因冲击波不稳定而被低估的HEL 值向上修正了8%—13%. 总体而言, 重解读的纵波声速和HEL与气炮实验数据展现出高度的一致性.

    4.   结　论

    • 本文提出了一种结合第一性原理计算与实验数据自洽求解HEL的方法, 旨在减小传统冲击波实验解读中波前弥散或非稳态效应带来的数据解读歧义. 以峰值应力达1 TPa且具有清晰速度剖面的金刚石气炮加载实验为参照, 本方法解读的HEL与基于弹性波走时的解读结果相对偏差在2%以内, 并可定量分辨其晶向各向异性特征. 本方法初步应用于激光加载HEL实验的数据解读, 将此前文献中被低估的HEL向上修正了8%—13%, 降低因冲击波加载平台差异带来的数据分歧.

      另一方面, 本方法的热弹性计算流程已应用于多种材料体系的宽温压域弹性及声速模拟研究[30], 预测了金属铍冲击相变前后的声速跃变[31], 可为通过Hugoniot声速实验研究冲击相变提供理论参考.

    数据可用性声明

    • 支撑本研究成果的数据集可在科学数据银行https:// doi.org/10.57760/sciencedb.j00213.00274访问获取.

      感谢西北工业大学的赵婷婷、巴黎理工学院的胡雨宽和北京应用物理与计算数学研究所的王越超对本工作的有益建议与支持.

    参考文献 (55)

目录

/

返回文章
返回