1.   引　言

实现“双碳”的宏伟目标, 关键在于大力推广清洁能源和提高能源的利用效率. 在工业领域温度位于中温区(60—225 ℃)的低品质余热约占全社会总能耗的15%, 大量的工业余热资源未能被有效利用, 随介质直接排放到环境中, 造成严重的环境污染和能源浪费. 热电材料是可将热能与电能直接相互转换且对环境友好的新能源材料, 其温差发电效应在低品质工业余热的回收利用上具有极大的优势和应用前景, 其热电转化效率可以通过无量纲热电优值ZT来描述^[1,2]:

$$ZT = \frac{{{S^{\text{2}}}\sigma T}}{{{\kappa _{\rm{l}}} + {\kappa _{\rm{e}}}}} \text{, }$$

(1)

式中, 分子S²σ是功率因子(包括塞贝克系数S和电导率σ); T是绝对温度; 分母为晶格热导率κ_l和电子热导率κ_e之和. 一般来说, 热电材料的商业应用标准要求其ZT值要高于1, 而具有高ZT值的热电材料通常要有高的功率因子或低的晶格热导率, 或者两者同时具备^[3-5]. 传统的热电材料在中温区的转换效率不够高, 严重制约其大规模应用, 因此寻找在中温区具有优良热电性能的材料具有非常重要的科学意义. 近年来由于纳米技术的发展, 不断有新的热电材料被发现, 其中层状材料因其特殊的晶体结构和输运性质引起人们的广泛关注.

近年来层状材料在热电领域因其独特的晶体结构和优异的电子性能备受关注^[6-11]. 经研究发现将其降维后不仅会增强界面散射, 减小晶格热导率^[12,13], 同时还会导致价带(导带)变窄, 增大载流子有效质量和塞贝克系数, 进而提高了热电性能^[14-16]. 例如, 块体相的Bi基三元化合物Bi₂O₂S^[17]和Bi₂O₂Se^[18,19], 它们本身就具有很小的晶格热导率, 但二维单层结构却拥有更小的晶格热导率和更好的热电性能. Huang等^[20]预测的单层Mg₃Sb₂表现出了小的晶格热导率和大的ZT值(600 K时>2), 其热电性能也远高于对应的块体结构的. 因此, 低维层状材料的引入为调控热电性能提供了新的思路, 并极大地刺激了其他二维层状热电材料的探索.

近期, 由Zhu等^[21]理论预测的磷烯衍生物—单层Ge₂X₄S₂ (X = P, As), 具有类似于GeP₃的层状六方晶体结构^[22,23]. 目前已知的类似结构中, 如单层SnP₃^[24], InP₃^[25], SbP₃和GaP₃^[26]都表现出了较低的晶格热导率和优异的热电性能, 这大大地提高了研究单层Ge₂X₄S₂热电性质的兴趣. 基于第一性原理和玻尔兹曼输运理论, 本文系统地研究了这两种材料的声子和电子输运性质. 结果证明单, 层Ge₂P₄S₂和Ge₂As₄S₂在300 K时沿armchair方向具有较低的晶格热导率, 分别为3.93 W·m^–1·K^–1和3.19 W·m^–1·K^–1. 另外, 价带的双重简并使得这两种材料都具有较大的塞贝克系数. 单层Ge₂X₄S₂基于其弱的声子输运能力和良好的电子输运性质, 展现了优异的热电性能. 其中, 单层n型掺杂的Ge₂P₄S₂在500 K时最大ZT值达到了3.51(zigzag方向), 在中温区优异的热电性能表明它在工业余热发电方面有很好的应用前景. 本文的理论计算研究也将激发实验科学家对这类材料制备、表征以及热电性能调控的研究兴趣.

2.   计算方法

基于第一性原理和密度泛函理论, 结合Vienna Ab-initio Simulation Package (VASP)软件^[27] 计算了单层Ge₂X₄S₂ (X=P, As)的热电性质. 为了避免层间的耦合作用, 沿Z方向设置了一个厚度为30 Å的真空层, 并通过vdW-DF2方法^[28]修正了范德瓦尔斯相互作用. 计算中, 本文采用Generalized Gradient Approximation (GGA)中的Perdew-Burke-Ernzerhof(PBE)交换关联泛函^[29,30]. 在结构弛豫时, 使用了500 eV的平面波截断能和11×11×1的K点网格, 总能和每个原子间相互作用力的收敛标准分别设置为10^–7 eV·Å^–1和0.01 eV·Å^–1.

通过求解玻耳兹曼输运方程计算了单层Ge₂X₄S₂ (X = P, As)的声子输运性质. 其中, 在3×3×1的超胞中使用3×3×1的K点网格, 使用VASP结合Phonopy软件^[31]获得二阶原子间力常数矩阵(2^nd IFCs), 并计算得到了声子色散曲线. 基于有限差分的方法, 使用thirdorder.py程序^[32]求得三阶力常数(3^rd IFCs). 本文考虑到第七近邻原子的相互作用. 将2^nd和3^rd IFCs作为输入, 通过ShengBTE软件包求解并获得了材料的晶格热导率^[33]. 为了获得更精确的晶格热导率, 本文引入了50×50×1的致密K点网格.

本文还通过玻耳兹曼输运理论和刚性带近似的方法对其电子输运性质进行了研究. 以塞贝克系数和电导率为代表的电子输运性质是在BoltzTrap软件包^[34]中实现的, 这里的计算选取了45×45×1的K点网格. 电子能带结构是基于Heyd-Scuseria-Ernzerhof(HSE06)杂化泛函^[35]计算获得. 另外, 基于形变势理论^[36,37], 本文还采用弛豫时间近似的方法计算了这两种材料的载流子迁移率和弛豫时间. 这种方法在能量尺度的影响下变化较小, 并且已经准确地验证了许多材料的热电性质^[38,39].

3.   结果与讨论

3.1   晶体结构与声子色散

图1(a)和图1(b)呈现了单层Ge₂X₄S₂的几何结构(空间群为$P\bar3m1$, No.164)和对应的布里渊区. 从图中可以看到单层Ge₂P₄S₂和Ge₂As₄S₂都表现出了类似于单层InP₃的褶皱结构, 这能够增加声子散射并降低晶格热导率. 这两种材料的每个原胞都由两个Ge原子、四个X原子和两个S原子组成, 其中每个Ge原子同时与两个X原子和一个S原子相连. 完全弛豫后的结构参数如表1所列, 这与Zhu等^[21]的理论预测值是一致的. 为了了解这两种材料的成键特性, 计算了它们的电子局域函数(ELF). 如图1(c)所示, 当ELF = 0时代表没有电子, ELF = 0.5对应于自由电子气, ELF = 1表示完全局域化. 可以发现所有原子键区ELF值都在0.5以上, 尤其是P—P(As—As)键的ELF值在0.75以上, 说明它们是由强的共价键组成, 大大提高了晶体稳定性^[40-42]. 为了比较准确地计算Ge₂X₄S₂的热导率和电导率, 并与其他单层材料进行比较, 应用文献[24]的方法计算了单层Ge₂X₄S₂的有效厚度. 如图(1)所示, Ge₂X₄S₂最上层和最下层的原子都是Ge原子, 这两层Ge原子的垂直距离加上2倍的Ge原子平均范德瓦耳斯半径, 定义为单层Ge₂X₄S₂的有效厚度. Ge原子的平均范德瓦耳斯半径为2.32 Å^[43], 可计算出单层Ge₂P₄S₂和Ge₂As₄S₂的有效厚度分别为6.57 Å和6.49 Å.

图2(a)和图2(b)为单层Ge₂X₄S₂的声子色散曲线和声子态密度(PhDOS). 从中可以清楚地看到这两种材料的声子谱都没有虚频, 这说明它们的结构是动态稳定的. 通过使用ab initio molecular dynamic(AIMD)模拟表明单层Ge₂X₄S₂可以在800 K下保持热力学稳定^[21]. 图2(a)和图2(b)中还可以看到声子色散曲线由3条声学支(其中包括纵向模LA、横向模TA和面外模ZA)和21条光学支组成. 由于As原子质量大于P原子, 导致单层Ge₂As₄S₂的声子色散曲线的频率都低于单层Ge₂P₄S₂的, 这意味着单层Ge₂As₄S₂具有更小的声子群速度. 进一步分析发现, 单层Ge₂As₄S₂声学支与最低光学支在M点简并在一起, 这意味着强的光声耦合作用和小的晶格热导率^[44-47]. 图2(c)为单层Ge₂X₄S₂的ZA, TA, LA模和最低光学分支的振动模式. 从图2(c)中可以看出, 低频下ZA模的所有原子沿面外方向运动, TA模与LA模分支的所有原子都是在面内做同向运动. 而在最低的光学模Opt₁中, Ge与S原子在中间层中沿相反方向移动, P原子以及As原子在平面内无序运动, 这种多原子的混合运动将大大地抑制声子的传输并降低晶格热导率.

3.2   声子输运性质

材料的晶格热导率是所有声子模贡献的总和^[48,49], 是影响其热电性能的重要因素之一, 可通过求解玻尔兹曼输运方程获得^[50]

$${\kappa _{\rm{l}}} = \dfrac{1}{V}\displaystyle\sum\limits_\lambda {{C_\lambda }{v_\lambda }{\tau _\lambda }} \text{, }$$

(2)

其中C为热容; ν为声子群速度; τ为弛豫时间; λ为声子模式. 基于玻耳兹曼输运方程并通过ShengBTE软件包本文计算得到了单层Ge₂X₄S₂的晶格热导率与温度的函数图, 如图3(a)所示. 其中, 在300到800 K的温度范围内, 晶格热导率与温度之间存在明显的反比关系, 其值随温度的增加而减小. 这种现象主要是由于声子散射随温度的升高而增强造成的, 这在许多材料中都得到了证实^[51,52]. 在300 K时单层Ge₂P₄S₂沿armchair方向和zigzag方向的晶格热导率分别为3.93 W·m^–1·K^–1和4.38 W·m^–1·K^–1, 而单层Ge₂As₄S₂的相应值分别为3.19 W·m^–1·K^–1和3.79 W·m^–1·K^–1, 这远低于一些经典的层状材料, 如MoS₂(26.2 W·m^–1·K^–1)^[53]和phosphorene(83.5 W·m^–1·K^–1)^[54]. 图3(b)是室温下累积导热系数与声子频率的关系图, 它反映了不同频率的声子对总晶格热导率的贡献. 结果表明, 单层Ge₂P₄S₂和Ge₂As₄S₂的声子分支(<2 THz的低频区域)在两个方向上对总晶格热导率的都贡献超过了90%. 为了进一步研究各声子分支对晶格热导率的贡献, 本文计算了如图3(c)所示的各声子模对总晶格热导率的相对贡献. 对单层Ge₂P₄S₂而言, 当热流沿armchair方向时, 各分支对总晶格热导率贡献的顺序为TA>LA>ZA>Opt, 而沿zigzag方向的贡献顺序为LA>TA>Opt>ZA; 单层Ge₂As₄S₂在armchair方向的顺序依次为LA>TA>Opt>ZA, 在zigzag方向的贡献依次为TA>LA>Opt>ZA. 从图3(c)中可以看出, 这两种材料的热传输主要由TA模和LA模贡献. 此外, 本文还研究了尺寸效应对单层Ge₂X₄S₂热电性能的影响. 图3(d)是300 K下累积晶格热导率与平均自由程(MFP)的关系图, 可以看到随着MFP的增加, 累积晶格热导率显著增加, 直到出现最大值. 值得注意的是, 自由程100 nm以上MFP的声子对热导率的贡献超过50%, 这表明可以通过尺寸效应和纳米化进一步降低晶格热导率.

为了进一步分析晶格热导率, 使用以下公式计算了单层Ge₂X₄S₂在室温下的声子群速度:

$${v_{\lambda ,q}} = \frac{{\partial {\omega _{\lambda ,q}}}}{{\partial q}} \text{, }$$

(3)

其中ω_{λ, q}是声子频率; q是声子波失的大小. 观察图4(a)和图4(b)发现这两种材料在低频区域都具有较小的声子群速度, 其中各声子模的群速度大小顺序依次为LA>TA>ZA, 这与图2中各声子模不同的声子色散曲线斜率有关. 我们还可以看到其中单层Ge₂As₄S₂的声子群速度更小, 这是它具有较小晶格热导率的重要原因.

另一个影响晶格热导率的关键参数是格林艾森参数(γ), 它常用来描述非简谐作用, 而非简谐作用又决定了声子—声子相互作用的强度, 格林艾森参数绝对值越大, 声子输运能力越弱. 格林艾森参数可以通过下列公式计算得出:

$${\gamma _{\lambda ,q}} = - \frac{V}{{{\omega _{\lambda ,q}}}}\frac{{{\partial _{\lambda ,q}}}}{{\partial V}} .$$

(4)

如图4(c)和图4(d)所示, 单层Ge₂X₄S₂在低频区域格林艾森参数的绝对值非常大, 尤其是ZA模, 这大大地抑制了声子输运能力, 对获得低的晶格热导率是非常有利的.

还通过求解各散射过程计算出声子弛豫时间, 如图4(e)和图4(f)所示, 发现单层Ge₂X₄S₂的低频声子分支具有较小的弛豫时间, 可与单层SnP₃相比较^[24], 这与其三声子散射相空间(P₃)有很大关系. P₃是指三声子散射过程的可用空间^[55], 是决定声子弛豫时间的重要参数. 较高的P₃值意味着有更多的空间用于三声子散射过程, 这通常是小的声子弛豫时间和低的热导率的标志^[56]. 从图4(e)和图4(f)的插图来看, 发现三个声学分支都具有较高的P₃值. 通过对声子群速度、格林艾森参数和声子弛豫时间的深入分析, 揭示了单层Ge₂X₄S₂具有较小晶格热导率的原因.

3.3   电子输运性质

计算了这两种材料的电子能带结构和对应的电子态密度(DOS). 如图5所示, 可以清晰地看到单层Ge₂P₄S₂的导带最小值与价带最大值分别位于Г点和K点, 带隙值为1.13 eV, 这说明它是一种间接带隙半导体. 而单层Ge₂As₄S₂是具有带隙值为1.21 eV的直接带隙半导体. 能带简并是材料具有优异的热电性能的重要特征, 目前常通过能带工程等方式实现^[57]. 这两种材料的价带在费米能级附近出现了明显的双重简并, 这意味着它们具有大的塞贝克系数, 有利于获得优异的热电性能. 观察它们的电子态密度可以发现价带对应区域主要由Ge原子控制, 而导带对应的部分Ge原子贡献较小, 主要由其他两种原子做贡献, 图5(c)也证明了这点. 从最低导带所对应的电荷密度中可以观察到电荷集中在P(As)原子和S原子附近, 同时最高价带的电荷在各原子附近都可观察到, 但Ge原子附近相对更多一些. 此外, 态密度在费米能级附近出现了陡峭的峰, 这是Seebeck系数在此处急剧增大的重要标志^[58], 有利于实现高的功率因子.

准确的载流子迁移率和弛豫时间是评估电子输运性质的重要基础. 基于形变势理论^[36], 通过下列公式计算了这两种材料的载流子迁移率^[59,60]:

$${\mu _{{\text{2D}}}} = \frac{{2e\hbar {C^{2{\rm{D}}}}}}{{3{k_{\rm{B}}}T{{\left| {{m^*}} \right|}^2}E_{\text{l}}^2}} \text{, }$$

(5)

式中, C ^2D, k_B, m^*和E_l分别代表二维弹性模量、玻尔兹曼常数、载流子有效质量和形变势常数. 这些参数和电子弛豫时间(${τ}=μ{{m}}^{{*}}{/}{e}$)的具体值如表2所列. 从表1中可以看到, 不管是沿armchair方向还是zigzag方向, 这两种材料的电子的迁移率和电子弛豫时间都是远高于空穴的, 其中单层Ge₂P₄S₂沿zigzag方向的值更是达到了3277.22 cm^–2·V^–1·s^–1, 这远高于单层δ-InP₃的2684 cm^–2·V^–1·s^–1[61], 这对获得好的热电性能是非常有利的.

在玻尔兹曼输运理论和弛豫时间近似的基础上, 通过BoltzTrap软件计算出了单层Ge₂X₄S₂的电子输运性质. 还通过刚性带近似理论研究了掺杂对电子输运性质的影响. 图6呈现了300至500 K下各电子输运参数与化学式μ的关系图. 其中, 当μ是负值时代表p型掺杂, 反之为n型掺杂. 塞贝克系数和电导率一般通过下列式子描述^[34]:

$$\begin{split} & {S_{\alpha \beta }}(T,\mu ) =\dfrac{1}{{eTV{\sigma _{\alpha \beta }}(T,\mu )}} \\& ~~\times \displaystyle\int \displaystyle\sum\limits_{\alpha \beta } {\left( \varepsilon \right)\left( {\varepsilon - \mu } \right)\left[ { - \dfrac{{\partial {f_\mu }\left( {T,\mu } \right)}}{{\partial \varepsilon }}} \right]} {\rm{d}}\varepsilon \text{, } \end{split}$$

(6)

$${\sigma _{\alpha \beta }}\left( {T,\mu } \right) = \dfrac{1}{V}\displaystyle\int {\displaystyle\sum\limits_{\alpha \beta } {\left( \varepsilon \right)\left[ { - \dfrac{{\partial {f_\mu }\left( {T,\mu } \right)}}{{\partial \varepsilon }}} \right]} } {\rm{d}}\varepsilon \text{, }$$

(7)

$$\sum\limits_{\alpha \beta } {\left( \varepsilon \right) = \frac{{{e^2}}}{{{N_0}}}\sum\limits_{i,k} {\tau {v_\alpha }\left( {i,k} \right){v_\beta }\left( {i,k} \right)} } \frac{{\delta \left( {\varepsilon - {\varepsilon _{i,k}}} \right)}}{{{\rm{d}}\varepsilon }} \text{, }$$

(8)

其中V是晶胞的体积; $\sum\nolimits_{\alpha \beta } {\left( \varepsilon \right)}$是输运分布函数; ${f_\mu }\left( {T, \mu } \right)$是费米狄拉克分布函数.

从图6(a)和图6(b)中可以看出, 塞贝克系数具有明显的温度依赖行为, 会随温度的升高而降低. 在室温下, 单层Ge₂P₄S₂沿armchair方向和zigzag方向在–0.06 eV位置具有最大的p型塞贝克系数, 分别为1800 μV·K^–1和1680 μV·K^–1, 小于单层Ge₂As₄S₂的2070 μV·K^–1和2040 μV·K^–1; 同时单层Ge₂P₄S₂和Ge₂As₄S₂在0.06 eV也获得了较大的n型塞贝克系数, 其中最大值分别为–1720 μV·K^–1(zigzag方向)和–1940 μV·K^–1 (armchair方向). 这两种材料如此大的塞贝克系数可与单层InP₃的2000 μV·K^–1相比拟, 有利于获得大的功率因子. 显然, p型塞贝克系数明显高于n型值, 这主要与态密度在费米能级附近出现的陡峭的峰有关. 另外塞贝克系数的最大值出现在低化学势区域与态密度的能量依赖行为有关, 它们之间遵循莫特关系^[62], 表现为

$$S = \frac{{{\pi ^2}k_{\rm{B}}^2T}}{{3e}}{\left[ {\frac{1}{n}\dfrac{{{\rm{d}}n\left( \varepsilon \right)}}{{{\rm{d}}\varepsilon }} + \dfrac{1}{\mu }\frac{{{\rm{d}}n\left( \varepsilon \right)}}{{{\rm{d}}\varepsilon }}} \right]_{\varepsilon = \mu }} .$$

(9)

因此, 根据公式可知通过调节载流子浓度能够进一步提高塞贝克系数.

图6(c)和图6(d)为单层Ge₂X₄S₂电导率与化学势的关系图. 从图6(c)和图6(d)中可以看出, 由于服从费米狄拉克分布, 使电导率最大值出现在较高的化学势区域. 与塞贝克系数明显不同的是, 电导率的n型值远大于其p型值, 这主要归因于其高的电子迁移率和弛豫时间. 本文还研究了电子热导率在不同温度下对热电性能的影响. 电子热导率遵循Wiedemann-Franz 定律^[63], 可通过下式计算获得

$${\kappa _{\rm{e}}} = L\sigma T \text{, }$$

(10)

其中L为洛伦兹数($L = {\pi ^2}k_{\rm{B}}^2/3{e^2}$). 因为电子热导率与电导率成正比关系, 所以可以观察到如图6(e)和图6(f)所示的电子热导率具有与电导率非常相似的函数曲线, 这也使得它们有相同的掺杂类型. 两种输运参数的不同之处在于电子热导率对温度更为敏感, 会随温度升高而增大. 从图6(e)和图6(f)中还可以看到, 单层Ge₂P₄S₂在zigzag方向表现出一个较大的n型电子热导率, 这主要与其较大的n型电导率有关.

通过已获得塞贝克系数和电导率计算出了单层Ge₂X₄S₂的功率因子(S ²σ). 如图6(g)和图6(h)所示, 功率因子也表现出类似于塞贝克系数的温度依赖行为, 但其值会随温度升高而增大. 在大的n型电导率的影响下, 功率因子在500 K时获得了一个较大的n型值, 这意味着它们可能是潜在的n型热电材料.

3.4   热电优值

利用前面得到的声子输运参数和电子输运参数, 计算了Ge₂X₄S₂的热电优值ZT, 评价其热电性能. 在图7(a)—(d)中, 可以看到单层Ge₂X₄S₂在armchair和zigzag两个方向的ZT值具有相似的曲线, 都会随着温度的升高而增大, 并且它们的最大值都表现出了n型大于p型的现象. 在最优的n型掺杂下, 单层Ge₂P₄S₂和Ge₂As₄S₂在500 K时沿armchair方向的ZT值分别为3.06和3.21, 室温下其值分别为1.75和1.86, 这超过了相同温度下磷烯的ZT值; 而这两种材料在zigzag方向的最大值分别为3.51和2.54. 在此对比发现, 由于单层Ge₂P₄S₂其在zigzag方向拥有更大的功率因子, 因此它表现出了一个更大的ZT值. 载流子浓度的变化对塞贝克系数和电导率的作用趋势是完全相反的, 表现为浓度的增加会导致塞贝克系数减小, 电导率增大, 反之亦是如此. 因此, 合适的载流子浓度可以使功率因子最大化. 图7(e)—(h)为一定温度下ZT值与载流子浓度的关系图. 结果显示, 最佳的空穴(电子)浓度为1×10¹⁰至1×10¹² cm^–2之间, 在此区间可以获得最大的n型(p型)ZT值.

4.   总　结

本文通过第一性原理和玻尔兹曼输运理论系统地研究了单层Ge₂X₄S₂的热电性能. 结果表明, 这两种材料具有较小的声子群速度、大的格林艾森参数和低的声子弛豫时间, 导致了非常低的晶格热导率. 电子能带结构表明单层Ge₂P₄S₂和Ge₂As₄S₂是分别具有带隙值为1.13 eV和1.21 eV的半导体, 并且价带的双重简并为这两种材料提供了很大的p型塞贝克系数. 基于弱的声子输运能力和很好的电子输运性质, 在最优的n型掺杂中获得了最佳的热电性能. 在500 K时单层Ge₂P₄S₂和Ge₂As₄S₂的最大ZT值分别为3.51(zigzag方向)和3.21(armchair方向), 这意味着它们在中温区具有优异的热电性能, 在工业余热发电领域具有应用前景.