含有倾斜界面硅/锗超晶格的导热性能*

2021-06-18刘英光任国梁郝将帅张静文薛新强

物理学报 2021年11期

刘英光任国梁郝将帅张静文薛新强

(华北电力大学能源动力与机械工程学院，保定 071003)

采用非平衡分子动力学(NEMD)方法模拟含有倾斜界面的硅/锗(Si/Ge)超晶格在不同倾斜角、不同周期长度、不同样本长度和不同温度下的导热性能.模拟结果表明， Si/Ge超晶格的热导率随着界面倾斜角的增加而非单调变化.当周期长度为4—8原子层时，界面倾斜角为45°的热导率比其他界面倾斜角时热导率增大了一个数量级，且热导率随样本长度的增加而增加，随温度的增加而减小.然而当周期长度为20原子层时，由于声子局域化的存在，热导率对样本长度和温度的依赖性都较弱.

1 引言

热电材料是一种利用固体内部载流子运动实现热能和电能相互转换的功能材料，具有方便、清洁等优点，在热电制冷器、红外探测仪、超导电子仪、遥控导航系统等方面有着广阔的应用前景[1-4].目前热电材料的能量转换效率采用无量纲热电优值(ZT)来衡量，ZT=S2σT/κ，其中S，σ，κ和T分别为塞贝克系数、电导率、热导率和绝对温度.一般来说追求高的ZT值需要塞贝克系数S和电导率σ最大化的同时减小热导率κ.由于纳米结构材料在对电导率和塞贝克系数影响较小的前提下可显著降低材料热导率从而提高热电转换效率，因此近年来纳米结构材料的热输运问题引起了广泛的研究兴趣[5，6].超晶格是典型的由相同或者不同材料周期排列组成的纳米结构，是热电转换领域的优秀候选材料之一[7-10]，所以学者们针对超晶格周期长度[11-13]、样本长度[14-16]和界面粗糙度[17，18]等方面对其导热性能的影响做了广泛的研究，并取得了一系列显著进展.

界面角度是影响多晶和复合材料导热性能的重要因素，国内外研究者对其开展了许多研究，如Elapolu等[19]采用反非平衡分子动力学模拟方法研究了含对称倾斜晶界的二维六方氮化硼纳米带的热导率，研究发现在晶界角度θ＜ 27.79°时，热导率随着倾斜角的增加而减小.作者将其归因于随着倾斜角的增加，晶界缺陷密度增加，声子散射增加所致.我们[20]采用NEMD方法研究了双晶ZnO在更大倾斜角范围内的卡皮查热阻，结果发现卡皮查热阻具有明显的尺寸效应，在小角度 (θ＜36.86°) 区域受倾斜角的影响较大，其卡皮查热阻随着倾斜角呈线性增加，而在大角度 (θ＞ 36.86°)区域受倾斜角的影响较小，所以热导率变化也相对较小.

Fujii等[21]采用MD方法对MgO中对称倾斜晶界的晶格热传导进行了系统的研究，结果表明，小角度对称倾斜晶界主要降低低频声子热输运，而高角度对称倾斜晶界主要降低中、高频声子热输运，从而使得高角度对称倾斜晶界热导率总体小于低角度对称倾斜晶界热导率.Bagri等[22]研究了晶粒尺寸和倾斜角对石墨烯导热性能的影响，研究发现，界面角度对热导率的影响与晶粒尺寸相关，且存在一临界尺寸，当临界晶粒尺寸在0.1 μm以下时，倾斜边界对热导率的贡献可以与晶粒本身的贡献相当.而在晶粒尺寸大于0.1 μm时，石墨烯晶体热导率主要受声子在晶粒内部的散射决定，晶界散射对热导率的贡献随着晶粒尺寸的增加逐渐消失.研究人员还发现，某一界面倾斜角度下材料热导率可能与其他角度下有很大不同，例如Tan等[23]利用真空热蒸发技术制备了具有倾斜结构的(Sb， Bi)2Te3纳米线，发现倾斜角为45°的纳米线阵列在室温下具有很高的热导率，说明独特的倾斜结构在影响声子输运特性方面起关键作用.

受以上研究启发，可以知道界面角度也会对超晶格热导率产生重要影响，是调控热导率的一个重要参数.所以，为了揭示界面角度对热导率影响的规律，发现不同倾斜角下声子的输运机理，本文采用NEMD模拟方法研究含有不同倾斜角的Si/Ge超晶格的导热特性，同时综合考虑了周期长度、样本长度和温度对试样热导率的影响.

2 模拟方法及步骤

本文所采用的超晶格导热模型如图1所示.在X，Y方向上施加周期性边界条件，在结构完全弛豫后在Z方向施加固定边界条件.由于截面积过小会对计算结果产生较大的误差[24]，因此在本文中我们选取XOY截面积为4 UC × 4 UC (UC为晶胞长度).为了建立热流方向的温度梯度，在Z方向上布置恒温器.首先将两边固定10原子层厚度为固定边界，防止与外界产生热量交换，并将该区域内粒子的速度设置为0.然后分别将紧邻边界区域的3 UC厚度设置为热源和热汇，以建立温度梯度.我们将一对连续的Si层和Ge层的总厚度定义为周期长度，利用2Si× 2Ge、4Si× 4Ge和10Si×10Ge分别表示周期长度为4， 8和20原子层厚度的Si/Ge超晶格.本文中主要计算3种不同周期长度 (2Si× 2Ge， 4Si× 4Ge和10Si× 10Ge) 和5种不同的倾斜角 (30°， 45°， 60°， 75°， 90°) 的结构如图2所示.

图1 NEMD 模拟计算热性质的示意图Fig.1.Schematic diagram of the NEMD model for calculating the thermal properties.

图2 周期长度为20原子层厚度时不同界面倾斜角的示意图Fig.2.Schematic diagram of different tilted interface angles with the period length of 20 atomic layers.

所有的模拟过程都是基于 LAMMPS 软件进行的[25]，模拟计算中时间步长采用0.001 ps.在MD模拟中，采用Tersoff多体势函数描述Si-Si，Si-Ge， Ge-Ge原子间相互作用.首先，在零温度条件下对体系进行能量最小化，优化原子位置; 然后，在正则系综 (NVT) 下进行 1 ns的平衡态模拟来达到控制的温度和体积; 最后，使用朗之万热浴在微正则系综 (NVE) 下施加温度梯度并测量热流值.

在稳定状态下超晶格内部建立线性温度梯度，而在体系的两端存在非线性区域.在界面处温度不连续，具有非常明显的温度跳跃.卡皮查热阻R为界面处的温度跳跃ΔT与热流密度J之比，

Z方向热流密度JZ的计算公式为其中，E为能量，t为模拟时间，A为横截面积.然后通过傅里叶定律计算超晶格热导率:

声子态密度(phonon density of states， PDOS)是表征材料中声子活动的一种有效的方法.根据所有原子的速度自关联函数(velocity autocorrelation function， VACF)的傅里叶变换[26]计算得到PDOS:

其中速度自关联函数(VACF)由(5)式定义:

其中，N为原子总数，vi(t) 为t时刻i原子的速度矢量，vi(0) 为原子的初速度，〈〉为系综平均.

声子参与率(phonon participation ratio， PPR)是了解声子输运的另一种有效方法[27-29]，尤其是定量描述声子局域化效应[30]:

其中 P DOSi(ω) 的表示频率为ω的第i个原子的PDOS.

3 结果与讨论

3.1 界面角对热导率的影响

我们计算了样本长度为32.6 nm、温度为300 K下不同倾斜角时2Si× 2Ge， 4Si× 4Ge和10Si×10Ge周期超晶格热导率，如图3所示.从图中明显地看到10Si× 10Ge周期超晶格热导率对倾斜角具有较弱的敏感度.当倾斜角从30°增加到75°时，由于卡皮查热阻的增加，热导率呈现略微降低的趋势; 而倾斜角从75°继续增加至90°时，由于在此范围内Si/Ge界面粗糙度减小，界面变得相对光滑，部分长波声子会以相干模式输运，从而热导率会略微升高.

与10Si× 10Ge周期超晶格不同的是，在2Si×2Ge和4Si× 4Ge周期超晶格中，可以发现在倾斜角为45°时，热导率出现了一个明显的峰值，其热导率约为10Si× 10Ge周期Si/Ge 超晶格的3倍.这是因为在倾斜角为45°时， Si/Ge界面比较光滑，大部分声子在Si/Ge超晶格中进行相干性输运.为了验证此结论的合理性，我们计算了界面角度为45°时4Si× 4Ge和10Si× 10Ge两种周期Si/Ge超晶格的PDOS，如图4所示.从图中可以看出，在4Si× 4Ge周期Si/Ge 超晶格中，低频声子的峰值较多，这说明有更多的相干性声子参与热传导.而在10Si× 10Ge周期Si/Ge 超晶格的PDOS中，相对于4Si× 4Ge周期Si/Ge 超晶格的峰值，其明显的峰值发生了“左移”，而较窄的声子频带表明热导率较低[31]，这与图3中结果保持一致.为了进一步解释在倾斜角度为45°时短周期结构的热导率明显高于10Si× 10Ge结构的热导率，我们计算了通过4Si× 4Ge和10Si× 10Ge结构的界面声子频谱热流[32].Si-Ge界面的单向光谱热流q(ω) 定义为[33，34]:

图3 不同周期长度下热导率与倾斜角的关系Fig.3.The relationship between thermal conductivity and tilted angle as the different period length.

图4 不同周期长度Si/Ge超晶格的声子态密度Fig.4.The PDOS of Si/Ge superlattices with different period lengths.

其中，A为截面面积，vi为Si界面上的原子速度，Fij为原子通过Si-Ge界面时的力.光谱界面热导可以计算为G(ω)=|q(ω)|/ΔT，其中 ΔT为通过界面时的温度降.最后计算界面的总热阻为[33，34]

其中， Δf为离散的频率区间.图5为SiGe界面的NEMD模拟的频谱界面热导G(ω)的计算结果.根据图5可以发现，在中低频范围(ω＜ 16 THz)内，4Si× 4Ge结构的超晶格的频谱热导高于10Si×10Ge结构超晶格的频谱热导.这是由于中低频声子具有较长的波长，通过短周期长度的超晶格时不易发生散射和局域化现象，对热导贡献较大.对10Si× 10Ge结构的超晶格，声子通过界面时容易发生散射和局域化效应，频谱热导较低，因此4Si×4Ge结构的超晶格的热导率明显高于10Si× 10Ge结构的超晶格的热导率.随着倾斜角继续增加至60°的过程中， Si/Ge界面逐渐变得粗糙，破坏了声子的相干性，因此在界面角度大于45°时热导率急剧下降.为了解释这一现象，计算4Si× 4Ge结构的30°， 45°， 60°时的PDOS，如图6所示.从图6可以看出，倾斜角为45°时的PDOS在5， 7， 10 THz时比其他两个角度的PDOS多了三个明显的峰值，这说明粗糙度极大地恶化了声子相干性输运.而当倾斜角从60°增加到90°，由于界面变的相对平滑，导热系数又开始上升，这与10Si× 10Ge周期Si/Ge超晶格是一致的.

图5 倾斜角度为45°时4Si × 4Ge和10Si × 10Ge超晶格的频谱热导Fig.5.Spectral thermal conductance of 4Si × 4Ge and 10Si ×10Ge superlattices at the tilt angle of 45°.

图6 不同倾斜角下Si/Ge超晶格的声子态密度Fig.6.The PDOS of Si/Ge superlattices with different tilted angle.

3.2 声子局域化

图7 为2Si× 2Ge、4Si× 4Ge和10Si× 10Ge周期Si/Ge 超晶格在不同样本长度下的热导率的变化，其倾斜角和温度设置为45°和300 K.从图中可以看到，由于声子的相干性输运， 2Si× 2Ge和4Si× 4Ge周期Si/Ge 超晶格的热导率随样本长度的增加呈线性增加[35，36].由于声子的平均自由程大于或等于样本长度，声子通过界面时保持其相位信息不变，只与样本边界发生散射，由动力学理论和声子-边界散射所引起的弛豫时间τ-1=bsv/L可知[37]，热导率与样本总长度成正比.然而， 10Si× 10Ge周期Si/Ge 超晶格的热导率随着样本长度的增加呈现非线性趋势，其值先增大后减小.这是因为尽管10Si× 10Ge周期Si/Ge 超晶格的周期长度处于非相干性区域，但是仍有一少部分声子参与相干性传输.而随着样本长度继续增加，相干声子发生局域化，不再对热传导起作用，因此热导率随着样本长度的增加而呈现减小的趋势.

图7 Si/Ge超晶格热导率与样本长度的关系Fig.7.Thermal conductivity of Si/Ge superlattice vs.sample total length.

图8 显示了倾斜角为45°、样本长度为32.6 nm时不同周期长度结构超晶格的热导率与温度的关系.我们发现2Si× 2Ge和4Si× 4Ge周期Si/Ge 超晶格的热导率随温度的升高呈现下降的趋势，这是由于随着温度的升高，原子振动增强，导致Umklapp过程的散射增强，热导率逐渐降低[20，35，36].与2Si×2Ge和4Si× 4Ge周期Si/Ge 超晶格的热导率趋势相反的是， 10Si× 10Ge周期Si/Ge 超晶格的热导率对温度不太敏感.这主要由于声子的散射过程和声子局域化两种竞争效应引起的.在10Si× 10Ge周期Si/Ge 超晶格中，随着温度的升高，声子的相干性的恶化将削弱局域化机制，从而导致热导率的增加; 而同时由于超晶格中的非相干性声子的Umklapp散射，使热导率呈1/T变化趋势，这将导致10Si× 10Ge周期Si/Ge 超晶格热导率随温度的升高而下降.这两种机制共同作用使得10Si×10Ge周期Si/Ge 超晶格热导率基本保持不变.

图8 Si/Ge超晶格热导率与温度的关系Fig.8.Temperature dependent thermal conductivity of Si/Ge superlattice.

为了说明局域化现象，我们根据(6)式计算了不同周期长度结构中声子参与率，如图9所示.从图9(a)中可以看出无论是低频声子还是高频声子，其声子参与率都随着周期长度的增加而减小.即随着周期长度的增加，越来越多的声子表现出局域化的特征.正如我们计算的那样，由于声子局域化效应， 10Si× 10Ge周期Si/Ge 超晶格的热导率通常小于2Si× 2Ge和4Si× 4Ge周期Si/Ge 超晶格的热导率.从图9(b)中发现声子参与率随着样本长度的增加而减小.在12.5—16.0 THz频率下，样本长度为108.6 nm时的声子参与率明显小于54.3712 nm时的值，表现出更明显的局域化效应.这些结果也更好地支持了图7中声子局域化的定性讨论.

图9 (a)不同周期长度时的声子参与率; (b)不同样本长度的声子参与率Fig.9.(a)The participation ratio of superlattices with different period length; (b) the participation ratio of superlattices with different sample length.