冲击加载下氟化锂晶体热导率的第一性原理研究*

杨 陈，王 艳，陈南迪，余柏树，曾召益，刘 勋

(1. 重庆师范大学 物理与电子工程学院，重庆 401331； 2. 南方科技大学 理学院物理系，广东 深圳 518055；3. 武汉理工大学 理学院物理系，武汉 430070)

0 引 言

氟化锂(LiF)是一种化学性质稳定的透明绝缘体，熔点为1118 K[1]。在冲击和静态压缩实验中广泛应用.当冲击压力达到200 GPa时，LiF仍然具有良好的透光性[2]。这使其成为动高压速度测量和温度测量重要的窗口材料[3]。因此,研究其在动高压下的热导率具有重要的应用价值。

高压下材料的热导率研究,在理论方面并不完善.Leibfried和Schlomann在德拜模型基础上提出了Leibfried-Schlomann方程(L-S方程)[4],认为晶体中绝大部分能量是由声学声子传递，而晶体缺陷对声学声子散射效应影响很小,以及光学声子对热传导的贡献也可以忽略不计.Roufosse和Jeanloz等[5]先在常压条件下考虑温度变化对晶体热导率的作用,然后在等温条件下改变压力研究热导率的变化情况,推导出了晶体高温高压Roufosse -Jeanloz 热导率方程(R-J方程)。孙峪怀等[6]认为L-S方程中只考虑随温度升高声子平均自由程减小,是考虑不周的。由于动高压加载技术和测温技术限制,对于高温高压LiF热导率的实验研究,相关数据不多见。在静高压加载的热导率测量实验中,由于温度加载多属于对热导率压力效应和低温压力效应的研究。高压下热电偶的标定开始影响实验的测量结果。为提高实验测量精度,近年来发展了以光学测量为主导的热学光栅法[7]。由于该方法受样品尺寸的限制,加载的最高压力不能超过40 GPa, 而且高温高压条件下,传压介质的热力学性质可能对测量结果有影响。Goncharov小组[8]一直致力于光脉冲加热技术的高压热导率研究工作，由于该方法对激光能量的建模存在较大的不确定度，以及高温热辐射也对实验结果产生较大影响。胡金彪等选用液体CHBr3作为样品测量高温高压下的热导率,但是实验数据分析中热导率系数为理论值,没有得到实验验证。早期美国Ahrens等提出了冲击压缩下夹心法热导率测量技术[9-10]，获得LiF晶体高压的热导率数据,但实验压力范围很窄。赵万广等[11]对R-J方程进行参数修正后,在实验压力达250 GPa时,测得晶格热导率与原R-J方程数值相差一倍。由于推导过程不同和实验条件限制,上述理论变化趋势和热导率关系式互相存在明显差异[9]。

目前,计算晶格热导率的早期方法有经典分子动力学模拟[12]和弛豫时间近似[13]。其中,经典分子动力学方法主要是采用经验势模拟。势函数的准确性直接影响到计算结果,而势函数在高温高压条件能否适用仍未可知。弛豫时间近似是通过引入弛豫时间(τ)参量来描述散射过程，而忽略了真实的声子散射[14]。因此,这两种方法均无法准确的获得晶格热导率。本文基于第一性原理计算,结合原子间力常数,采用迭代法求解声子玻尔兹曼输运方程,获得LiF晶体冲击高压的热导率。

1 计算方法

本文中对LiF晶体(空间群Fm-3m)计算总能量的电子结构时,采用的是VASP软件包中的PAW方法,交换关联选择广义梯度近似的PBE形式。我们将Li的2S1和F的2S22P5作为价电子来处理。在计算中,为了获得准确的计算结果,我们对k点和q点以及其他相关参数都进行了严格的收敛性测试。平面波截断能均取为500 eV,在静态能量计算和声子计算中,分别采取了15×15×15和10×10×10的布里渊区K点取样。体系总能量的收敛值取了1×10-6eV,原子间力的收敛精度为1×10-4eV/(0.1 nm)。声子贡献的热导率就是晶格热导率。晶体的自由能F(T.V)可以写成如下形式：

(1)

其中，U(V)是T=0K时晶格的内能,只与晶体体积有关,而与温度(或晶格振动)无关。后项与晶格振动有关,即与温度有关。声子频率的几何平均定义为：

(2)

其中，ωqj是波矢q处的第j支声子的频率，Nqj是总的声子的振动模式数目。经测试,以上参数可使声子的几何频率收敛到1 cm-1。在结构优化和自洽计算中采取15×15×15的布里渊区网格点。自洽计算收敛和结构优化完成后，构建4×4×4的超晶胞。采取有限位移法进行声子谱计算,进而确定二阶原子间力常数[15]；在计算三阶力常数时,原子间相互作用的截断半径经过严格测试,为确保结果的准确性,选取原子间六阶最近邻.三声子散射过程必须满足能量和动量守恒:

ωj(q)±ωj′(q′)=ωj″(q″)

(3)

q±q′=q″+G

(4)

ωj(q)是(j，q)空间的频率，当G=0时，称为正规过程；当G≠0时，称为翻转过程。采用ShengBTE[14]程序,迭代求解声子玻尔兹曼输运方程[16],从而获得LiF的晶格热导率等相关物理性质。晶格热导率κ可以表示为：

(5)

其中，Ω是原胞的体积，N为原子个数。Cq为热容量，νq为群速度，τq为声子寿命。

2 计算结果与讨论

2.1 物态方程及声子色散关系

物态方程是固体能量随体积的变化关系,描述了固体在压缩或者拉伸作用下的行为。通过对LiF的静态能量的计算获得了物态方程,结果如图1所示。Liu等[17]在300 K温度实验测量了37 GPa范围内体积随压力的变化情况,与Smirnov等[2]的理论结果对比,我们的计算结果与实验数据更为吻合。物态方程的参数如表1所示,零压下LiF的平衡体积V0为0.01662 nm3,与实验值0.01642 nm3接近[17]，体积模量B为72.52 GPa,与实验值73.0 GPa接近[17]，体积模量对压力的偏导数B′为4.03，与实验值3.90吻合[17].弹性常数C11、C12和C44与Briscoe等[18]的实验数据吻合。

图1 LiF零压0 K 和300 K的物态方程,空心圆为Liu等[17]的实验结果,点划线和圆点分别为Smirnov等[2]在0 K和300 K的理论结果,实线和虚线分别为本文计算结果Fig.1 Equation of state of LiF at zero pressure 0 K and 300 K ,comparison with experimental data (open circles) of Liu et al [17], and the theoretical results (dot and dash) of Smirnov et al [2], the solid and dotted lines are the results of this paper

表1 零温零压LiF平衡晶格常数a(nm)、平衡体积V0(nm3)、体模量B(GPa)、体模量对压力的导数B′(GPa)和弹性常数C11，C12，C44。

本文采用有限位移方法计算二阶力常数,我们获得零温下LiF晶体的声子色散关系，如图2(a)所示。可以看出声子曲线没有出现虚频,表明晶体结构是稳定的。本文的计算结果与Dolling等[22]实验结果基本吻合。我们也分别计算了高压下的LiF声子色散关系,如图2(b)所示。从图中可以看出,随着压力的增加声子频率峰值上移,没有出现虚频现象,表明在本文计算的压力内LiF一直保持体心立方的稳定结构。通过声子色散曲线,我们也获得了高温高压下的自由能。

图2 LiF在零压(a)和高压(b)的声子色散关系,实心球为Dolling等[22]的实验结果Fig.2 Phonon dispersion relation of LiF at zero pressure (a) and high pressure (b), comparison with the experimental data (solid circles) of Dolling et al [22], the solid line is the result in this paper

2.2 Hugoniot物态方程及声速

Hugoniot物态方程是材料的重要性质之一,它可以同时反映出材料对压力和温度的响应特性。在冲击波物理中,冲击压力PH内能EH和体积VH之间满足 Rankine-Hugoniot关系:

(6)

式中，E0,P0,V0分别是冲击加载前初态的摩尔内能、压力和摩尔体积。EH,PH,VH,分别是冲击加载后终态的摩尔内能、压力和摩尔体积。

根据上式,我们计算了LiF的Hugoniot物态方程,如图3(a)所示。从图3(a)可以看出,PH-VH随Hugoniot与实验结果十分吻合。PH-TH曲线如图3(b)所示,我们的结果与实验测量及理论计算的结果基本吻合。图中空心圆为Kormer等[23]的实验结果,三角形表示赵万广等[11]在改进并完善液体夹心法测热导率实验技术后获得的实验数据,黑色实线为本文计算结果。我们的结果略低与实验测量的结果,主要原因在于,在冲击加载下LiF保持透明性，动高压实验无法直接测量冲击温度，实验结果中的冲击温度通过热力学计算获得，而热力学计算中相关参数的选择可能导致冲击温度估算过高。

采用直接的应力应变方法,我们计算了LiF随Hugoniot压力的弹性常数,如图4所示。计算结果与Smirnov等[2]的理论结果符合,随着压力增大，C11和C12呈线性增加趋势明显，C44变化趋势相对较小。表明在纵向应力的影响下晶胞的拉伸或压缩程度明显,体积膨胀或压缩明显。将弹性常数通过Voigt-Reuss-Hill 近似计算得到体积模量(B)和剪切模量(G)，压缩波声(VP)和体波(VB)声速。由Hugoniot弹性常数进一步计算得到的Hugoniot声速如图5所示,在0 GPa时，LiF的压缩波和体波声速分别为7.1 km/s和5.1 km/s,其中压缩波声速与刘前程等[26]的实验数据吻合.随着压力增大,压缩波声速和体波声速均增加,且与实验数据基本一致[26-27]。

图4 LiF的Hugoniot弹性常数,实心球为Smirnov等[2]理论结果，实线为本文的结果Fig.2 Hugoniot elastic constants of LiF, comparison with the theoretical result (solid circles) of Smirnov et al[2], the solid line is the calculation result in this paper

图5 LiF的Hugoniot声速实心球和实心三角形分别为Yu等[27]和Liu等[27]的实验数据,实线为本文的计算结果Fig.5 Hugoniot acoustic velocities of LiF, comparison with the experimental data (solid circles of Yu et al[27], and solid triangles of Liu et al[26]), the solid line is the calculation result in this paper

2.3 晶格热导率

LiF晶体作为典型的窗口材料,其高温高压下的热导率研究一直是研究热点。零压下300～1 000 K范围内LiF的晶格热导率如图6所示。结果显示，300 K时，结果为12.8 W/(m·K)。2010年Elhadj等[28]实验得到低吸收系数时晶格热导率为13.2 W/(m·K), Singh等[29]采用Callaway-Holland模型获得的结果为13.6 W/(m· K), Liang等[30]采用迭代解求解玻耳兹曼输运方程的方法得计算值13.89 W/(m·K)。本文计算结果与以上数据符合较好。图6显示，随着温度升高晶格热导率降低,当温度升高到1 000 K时,热导率降为3.86 W/(m·K).本文计算结果实验值和其他理论值趋势一致。

图6 LiF常压下300～1 000 K的晶格热导率，实心正方形为Liang等[30]的理论结果，三角形和空心正方形分别为Men团队[31]和Elhadj等[28]的实验数据，实心球为本文结果Fig.6 Lattice thermal conductivity of LiF at 300-1 000 K at atmospheric pressure, comparison with theexperimental data (open square of Elhadj et al[28]and open triangles of Men et al [31]), solid square is the theoretical result of Liang et al [30], solid sphere is the result calculated in this paper

图7 (a)零压300 K下声子散射率,(b)不同Hugoniot压力下声子散射率Fig.7 Phonon scattering rates at 300 K at zero pressure and different Hugoniot pressures

赵万广等[11]采用液体夹心法热导率测量技术，结合理论推导给出了修正后的R-J方程。通过数值拟合获得～39、～70和～100 GPa 3个压力值的LiF单晶冲击高压热导率数据，分别为13.2±1.0、7.1±1.4和4.8±0.6 W/(m·K) (这3个数据点是实测数据，先有的实测数据，再根据实验结果修正的R-J方程)。与本文的计算结果的变化趋势一致,如图8所示。结果表明，LiF晶体在动高压冲击100 GPa范围内,热导率先缓慢增加后快速降低,40 GPa附近达峰值。出现这种现象主要是受声子非谐散射效应的影响,散射率越低，声子寿命越长,从而具有较大的晶格热导率。 压力较低时，我们的计算结果远高于实测结果，压力较高时，二者则趋于一致，造成这种差异的具体原因尚不清楚。需要指出的是，虽然赵万广等[11]采用的液体夹心法热导率测量技术比早期Ahrens等[9]的实验方法有所改进，但在得出LiF热导率的过程中还是使用了一些经验参数，如高压下三溴甲烷、LiF的定容比热及LiF的格林艾森参数等，这些热力学参数的选取会影响最终的实验结果。

图8 Hugoniot热导率,实心球为赵万广等[11]的实验数据，方形为本文的计算结果，实线为线性拟合结果Fig.8 Obtain lattice thermal conductivity along the Hugoniot, comparison with the experimental data (solid sphere) of Zhao et al[11], the square is the result of calculation and the solid line is the result of linear fitting

3 结 论

基于第一性原理晶格动力学方法，通过求解声子的玻尔兹曼输运方程，获得了高压下LiF晶体的Hugoniot物态方程和Hugoniot晶格热导率。高压下声子色散关系结果表明LiF晶体在计算的压力温度范围内一直保持体心立方稳定结构。采用准谐近似的方法,计算获得的Hugoniot物态方程与早期的理论结果及实验数据符合较好。采用迭代解方法求解声子的玻耳兹曼输运方程，准确计算了100 GPa,2000 K范围内LiF晶体的晶格热导率,在3个压力值～39、～70和～100 GPa的实验数据分别为13.2±1.0、7.1±1.4和4.8±0.6 W/(m·K)，与计算结果的变化趋势一致.最后预测的LiF沿Hugoniot线的晶格热导率将为动高压温度测量提供重要的参考数据。