张晓清, 李 强,2*

(1. 攀枝花学院 钒钛学院, 四川 攀枝花 617000; 2. 四川省材料腐蚀与防护重点实验室, 四川 自贡 643002)

低热导材料具有优良的隔热作用,因此在航空航天、能源、冶炼等领域得到了广泛应用[1-3].开发和设计具有优良隔热性能的新型低热导材料是当前材料研究的重要方向之一.传统的耐火隔热材料(如陶瓷)通常具有一定的高温强度,但隔热性能有时并不优异.由于其孔隙率过高,存在高温腐蚀渗透的风险.为此,许多研究人员希望提高陶瓷的致密性[4-7].一般隔热材料具有优良的保温性能,但缺乏高温强度,因此提高隔热材料的高温强度和保持较低的热导系数是非常重要的工作.然而,通过实验开发低热导材料的成本高昂,需要经过反复调整实验方案,而理论预测低热导材料将大大降低实验探索低热导材料的成本.一些研究人员采用理论方法已经预测了一些新型低导热材料,如Seko等[8]通过第一性原理方法发现了221种低导热材料,其中2种具有较小的带隙,成为热电材料的潜在候选对象;Liu等[9]使用第一性原理计算预测的二维PbSe材料表现出较低的声速值和超低的导热系数(约0.26 W·m-1·K-1);Mortazavi等[10]使用机器学习原子间势方法结合第一性原理研究了graphene/borophene异质结构,发现该结构具有低晶格热导率.由于第一性原理在研究材料的力学和热学等方面具有很多优势[11-13],因此本文通过机器学习和第一性原理相结合的方法,预测2种新型低热导材料Ba2O2S的结构,并研究其电子特性和声子特性.

1 计算方法

为了寻找新型Ba2O2S低导热材料,本文使用CrySPY软件[14]进行机器学习搜索,能量计算由Quantum-Espresso软件[15]完成.采用GBRV赝势[16]处理电子和原子核之间的相互作用力.选择不同的价电子组态O(2s22p4)、S(3s23p4)和Ba(5s25p66s2)来处理系统能量,采用广义梯度近似和PBE泛函处理电子间的相互作用.选择612 eV的截止能量来优化和驰豫电子.倒易空间采样k点间距取为0.4 nm-1,使系统能量能够收敛到4.8×10-7eV/atom,原子间的力收敛到1.4×10-3eV·nm-1,应力收敛到0.02 GPa.使用Phono3py软件[17]计算声子谱和热导率.采用VESTA软件显示晶体结构,选择应力应变法计算弹性常数.

2 结果与讨论

2.1 Ba2O2S的结构和电子性质为了获得Ba2O2S的低能结构,将原子类型、原子符号和收敛标准作为输入参数,通过机器学习算法实现全局结构搜索.搜索过程结束后,找到能量最小的2种结构.通过对称性分析,发现这2种结构均为六方晶体,空间群分别为P-3m1和P63/mmc,具体结构如图1所示.为方便后续讨论,将这2种结构分别命名为B1和B2结构.表1列出了Ba2O2S2种结构的相应原子坐标、晶格参数和wyckoff位置.

图1 预测的2种Ba2O2S结构,空间群分别为 (a) P-3m1和(b) P63/mmc

从图1可以看出,2种类型的Ba2O2S具有六方结构,并且按照O—Ba—S—Ba—O原子层堆叠形成.这2种结构的主要区别在于它们具有不同的Ba原子堆垛顺序,这导致B1-Ba2O2S的能量比B2-Ba2O2S约低0.5 meV.

为了获得Ba2O2S的成键信息,运用Bader量子理论[18]计算键临界点电荷密度(ρBCP)和拉普拉斯电荷密度(∇2ρBCP),结果如表2所示.可以看出,Ba—O和Ba—S键的临界点电荷密度非常小,为0.001 5～0.004 0,这意味着Ba—O与Ba—S键为离子键.获得的Ba—O和Ba—S键的拉普拉斯电荷密度大于零,进一步说明这2种键具有明显的离子键特征,这与MgO的情况相同[19].

表 1 计算的Ba2O2S的结构参数

表 2 计算的Ba2O2S的键临界点电荷密度和拉普拉斯密度

由于B1-Ba2O2S和B2-Ba2O2S具有几乎相同的电子结构,因此仅在图2中呈现了B1-Ba2O2S的能带.可以看出,B1-Ba2O2S属于直接半导体,最小带隙位于Γ点位置,带隙约为4.6 eV.Ba 6s态位于导带顶部,O 2p和S 3p态位于价带边缘.与B1-Ba2O2S相比,B2-Ba2O2S的带隙略小,约为4.5 eV.

图2 计算的B1-Ba2O2S的能带

2.2 Ba2O2S的弹性性质六角形晶体[20-21]的弹性稳定性通过以下公式进行评估:

C44>0,

(1)

M1=C11-|C12|>0,

(2)

(3)

通过应力应变方法获得的B1-Ba2O2S的C11、M1和M2计算值分别为17.9、26.5和46.2 GPa,B2-Ba2O2S的相应值为13.4、19.8和52.5 GPa.这些值都是大于零的,符合方程(1)～(3)给出的弹性稳定性标准.因此,Ba2O2S的2种结构在零压条件下是力学稳定的.利用获得的弹性常数Cij,可以通过Voigt-Reuss-Hill(VRH)近似[22]计算Ba2O2S的体积模量(B)、剪切模量(G)和杨氏模量(Y).结果表明B1-Ba2O2S的B、G和Y的计算值分别为38.6、14.7和39.2 GPa,大于B2-Ba2O2S的对应值37.5、12.6和33.9 GPa.这意味着B1-Ba2O2S应该是低温相,B2-Ba2O2S是高温相.

2.3 Ba2O2S的声子特性计算的B1-Ba2O2S和B2-Ba2O2S的声子谱和声子分态密度如图3所示.可以看出,对于2种结构,O原子的声子态密度主要位于5.8 THz以上的高频区,而S原子的主要贡献位于2～5 THz区域,这是由于O原子质量比S轻造成的.B1-Ba2O2S的最高声子频率出现在L点,约为9.27 THz,而B2-Ba2O2S的最高声子频率位于M点,约为9.58 THz.此外,B1-Ba2O2S和B2-Ba2O2S的声光耦合出现在0～4 THz,这将导致Ba2O2S的晶格热导率包含一些光学振动模的贡献,从而使声光耦合对晶格热导率产生一些积极的影响.

基于弛豫时间近似,B1-Ba2O2S和B2-Ba2O2S的晶格热导率可以通过以下公式[17]获得:

(4)

式中,V0是系统体积,Cλ是热容,vλ是群速度,τλ是声子寿命.为了避免声子寿命发散,平均自由程选择为1 m.结果如图4所示.

图4 计算的B1-和B2-Ba2O2S的热导率和群速度外积值

从图4可看出,由于B1-Ba2O2S和B2-Ba2O2S的xy平面各向同性,因此计算的[100]方向的热导率kxx等于[010]方向的热导率kyy.在室温下,B1-Ba2O2S的kxx为2.1 W·m-1·K-1,其约为[001]方向热导率kzz(～1.3 W·m-1·K-1)的1.6倍,这表明B1-Ba2O2S的热传输是各向异性的.对于B2-Ba2O2S,在300 K时,kxx和kzz分别为1.6和1.0 W·m-1·K-1,其低于B1-Ba2O2S的值,这是因为B2-Ba2O2S的弹性模量小于B1-Ba2O2S,以及声速较小造成的.如果Ba2O2S是多晶形态,则热导系数显示各向同性.计算表明,在300 K时B1-Ba2O2S和B2-Ba2O2S的各向同性热导率kiso分别为1.86和1.40 W·m-1·K-1,比MgO等氧化物小,因此可以用于隔热材料.随着温度的升高,热导率逐渐降低.最后,在1 000 K时2种结构的kiso分别变为0.55和0.42 W·m-1·K-1.通过图4还可以看到,8 THz时B1-Ba2O2S的群速度外积值达到62 THz2·nm-2,是B2-Ba2O2S的1.35倍,这是导致B1-Ba2O2S的热导率高于B2-Ba2O2S的主要原因.

B1-Ba2O2S和B2-Ba2O2S不同声子振动模对热导率kiso的贡献以及声学模和光学模的声子态密度如图5所示,可以看到B1-Ba2O2S和B2-Ba2O2S晶体热导率的最大贡献都源于3个声学模,分别达到73.4%和38.2%.光学模对B1-Ba2O2S热导率的贡献仅为26.6%.与B1-Ba2O2S不同,B2-Ba2O2S光学模对热导率的贡献达到61.8%,约为B1-Ba2O2S的2.3倍.同时,从图5(c)～(d)的声学模和光学模的声子态密度可以看到,B2-Ba2O2S的声光耦合强度要高于B1-Ba2O2S.因此声光耦合对Ba2O2S的热传输具有重要作用.

为了进一步分析声子寿命和热容对热导率的影响,计算B1-Ba2O2S模式热导率作为声子寿命和热容的函数,结果如图6所示.可以看出,B1-Ba2O2S的模式热导率主要集中在1.0×10-5W·m-1·K-1附近,当温度从300 K升高到1 000 K时,平均声子寿命从1.0 ps变为0.5 ps,说明对于B1-Ba2O2S,声子寿命对其热导率的影响不是太大.随着温度的升高,模式热容从8.0×10-5eV增加到8.5×10-5eV,而对应的热导率平均值从5×10-4变为1×10-5W·m-1·K-1,因此,热容对B1-Ba2O2S热导率有较大影响.

图6 B1-Ba2O2S模式热导率作为声子寿命和热容的函数

3 结论

基于机器学习算法和第一性原理方法预测了2种新型Ba2O2S结构.这2种结构都属于六方晶系,分别具有P-3m1和P63/mmc空间群,Ba原子的不同堆垛排布是造成2种结构不同的最大原因.声子谱和弹性系数的计算结果,显示这2种结构满足动力学和弹性稳定性的要求.通过Bader量子理论计算结果,发现Ba—S和Ba—O键都是离子键.室温下,2种结构的热导率存在各向异性,其中B1-Ba2O2S的kxx为2.1 W·m-1·K-1,约为kzz的1.6倍.

对于B2-Ba2O2S,300 K时kxx和kzz分别为1.6和1.0 W·m-1·K-1,其低于B1-Ba2O2S的值.同时,发现B1-Ba2O2S的热导率高于B2-Ba2O2S的主要原因是由于B1-Ba2O2S的群速度外积分值高于B2-Ba2O2S导致的.B1-Ba2O2S和B2-Ba2O2S晶体热导率的最大贡献都源于3个声学模,但是B2-Ba2O2S的声光耦合强度要高于B1-Ba2O2S,光学模对热导率的贡献要大于B1-Ba2O2S.

致谢攀枝花学院校级科研教研项目(035001731、JJ2261)和四川省大学生创新创业培训项目(S202211360049)对本文给予了资助,谨致谢意.