1.   引　言

二维材料以其独特的晶体结构、高载流子迁移率、高热导率和高机械强度等优良特性在光电、热电、柔性等微型器件领域受到了广泛关注^[1–3]. 作为二维材料家族的重要分支, 过渡金属二硫族化合物(transition metal dichalcogenides, TMDs)不仅和石墨烯一样能达到稳定的原子层厚度, 而且展现出了包括导体、半导体和金属在内的一系列丰富的电子性质, 因此被当作研究多种低维奇异物理现象的理想平台^[4–9]. 例如, 在ZrSe₂中, 用少量的Ti替代Zr可以有效调控ZrSe₂的带隙, 实现半导体到金属的转变^[10]; 将石墨烯和MoSe₂结合形成的异质结可以通过氧掺杂调控肖特基势垒, 有助于设计和制造高性能纳米场效应晶体管^[11]; 此外, 以能谷自由度为信息载体的能谷电子学受到了人们的极大关注, 最近也报道了通过异质结可以有效操控TMDs的谷极化^[12]. TMDs性质的多样性不仅得益于其家族的相对庞大, 同时也来自于其晶体结构的多样性^[13]. 单层TMDs最常见的结构是三棱柱相(H相)和八面体相(T相). 研究表明, 除了这两种高对称相, 单层TMDs还存在其他畸变相, 例如: 金属原子二聚化对应的$T^\prime$相, 该相在V, Mo, W, Ta等的二硫族化合物中被发现, 并且实验上已经可以制备^[14–18]; 通过偏振拉曼光谱确定的ReS₂与ReSe₂的四聚化$T^{\prime\prime\prime}$相^[19,20]; 理论计算预测了MoS₂含有非中心对称结构, 即Mo原子发生三聚化的$T^{\prime\prime\prime}$相, 在后来的实验中也得到证实^[21,22]. 畸变相往往具有与对称相不同的物理性质, 而且在某些性质上表现得更为优异^[23–25]. 例如, MoTe₂在畸变结构下, 载流子迁移率远高于对称结构, 可达$4000\, \rm cm^2{{\cdot}} {V}^{-1}{{\cdot}} s^{-1}$^[26]. 此外, 即使同为畸变相, 微小的结构差别常常会带来性质上的较大差异^[27–30]. 研究发现, Cr的二硫族化合物单层在不同畸变相中, 会因两Cr原子间的距离差异而展现出截然不同的磁学性质^[31]. 由于不同畸变相间的结构差异往往很小, 目前的实验表征手段还很难判断样品中是否同时存在着多种具有微小差异的畸变相, 因此从理论上探究不同畸变相的结构稳定性和物理性质是有意义的. 同时我们注意到, 最近的文献报道了压力可以诱导二维材料的相变^[32], 而不同的畸变相差距非常微小, 两者之间相变的可能性也值得探究.

本文预测了二维RuSe₂的一种$T^{\prime\prime\prime}$相结构, 这是一种不同于目前实验上获得的二维六方H相和T相RuSe₂^[33]的结构, 和已报道的理论预测的$T^\prime$相^[34]结构也不一样. 借助第一性原理计算, 通过对比4种不同结构的能量, 并结合声子谱和分子动力学模拟结果预测了$T^{\prime\prime\prime}$相结构在室温下的稳定性. $T^{\prime\prime\prime}$相与$T^\prime$相结构的差别导致了材料性质的差异, $T^\prime$相的RuSe₂单层具有间接带隙, 而$T^{\prime\prime\prime}$相则为直接带隙半导体. 此外, 计算了$T^{\prime\prime\prime}$相单层RuSe₂的载流子迁移率, 发现$T^{\prime\prime\prime}$相的迁移率远高于$T^\prime$相. 最后, 发现通过施加zigzag方向的单轴应力, 可以显著降低$T^\prime$与$T^{\prime\prime\prime}$相之间的转变势垒, 在室温下实现两个畸变相之间的转变.

2.   计算方法

本文使用了第一性原理计算软件包VASP (Vienna ab initio simulation package)^[35,36]完成了基于密度泛函理论的计算模拟工作, 赝势采用的是投影缀加波PAW (projector augmented wave)近似^[37]. 计算采用了广义梯度近似GGA (generalized gradient approximation)下的PBE (Perdew-Burke-Ernzerhof)交换关联泛函^[35–38], 波函数的截断能取为550 eV. 本文涉及的几种晶相结构都采用矩形晶格进行建模计算, 如图1所示. 通过Monkhorst-Pack方法^[39]对矩形晶胞的第一布里渊区进行采样, 以Γ为中心取$19\times 19\times 1$的K点网格进行晶体结构和电子结构优化. 每个原子上的力收敛阈值为$10^{-3}\; \text{eV/Å}$, 电子步收敛精度设置为$10^{-6}\; \text{eV}$. 为了考虑电子之间的强关联效应, 采用DFT+U的方法进行修正. 采用线性响应法对Ru原子的有效U值${\left(U_\text{eff}\right)}$进行测试^[40], 获得${U_\text{eff}} = 4.36\; \text{eV}$. 为了验证结构的稳定性, 通过密度泛函摄动理论DFPT (density functional perturbation theory)方法, 利用PHONOPY程序包计算声子色散谱^[41]. 本文还使用了基于第一性原理的分子动力学模拟AIMD (ab initio molecular dynamics)来验证材料的热力学稳定性, 计算基于正则系综NVT (canonical ensemble, particles N, volume V, temperature T) 并采用Nosé–Hoover方法进行^[42]. 载流子迁移率是通过基于有效质量近似的形变势理论计算的, 对于二维系统, 其计算公式如下^[43]:

图 1  (a) 单层RuSe₂ 4种结构的相对能量, 以图中矩形超胞的长度作为横坐标进行排列, 绿球代表Ru原子, 红球代表Se原子; (b) RuSe₂的$T^\prime$相晶体结构; (c) RuSe₂的$T^{\prime\prime\prime}$相晶体结构

Fig. 1.  (a) The relative energies of four structures of monolayer RuSe₂, arranged along the lengths of rectangular supercells, green spheres represent Ru atoms, while red spheres represent Se atoms; (b) crystal structure of RuSe₂ in the $T^\prime$ phase; (c) crystal structure of RuSe2 in the $T^{\prime\prime\prime}$ phase.

下载: 全尺寸图片 幻灯片

其中, $k_{\rm B}$为玻尔兹曼常数, $m_{\rm a}$和$m_{\rm b}$是图1中zigzag和armchair方向上的载流子有效质量, $E_{\rm 1 a}$和$E_{\rm 1 b}$是两个正交方向的形变势, $C_{\rm 11}$和$C_{\rm 22}$为弹性模量.

3.   结果与讨论

图1(a)展示了4种二维RuSe₂的晶体结构, 为了便于比较, 这里统一采用了大小为T相单胞$\sqrt{3 }$ × 3倍的矩形超胞. 超胞在zigzag方向上的长度(L)与相对能量(relative energy, RE)的关系表明$T^\prime$相结构的能量最低, H相的能量最高, 而T相与$T^{\prime\prime\prime}$相则介于两者之间. 对T相结构进行分子动力学模拟发现, 50 K下的结构会自动弛豫到非对称$T^\prime$相, 且其声子谱也存在较大的虚频 (见补充材料图S1(a)(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557)和S1(b)(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557)), 结果与先前理论研究结果一致^[34]. 这和实验上发现T 相并不稳定, 高温退火后结构就会被破坏的现象也是一致的^[33,34,44]. 但对于$T^{\prime\prime\prime}$相结构, 其自由能比T相低165 meV/f.u., 且获得的声子谱和分子动力学模拟结果都预测了它的稳定性, 因此很有可能在实际中存在.

为比较$T^\prime$相和$T^{\prime\prime\prime}$相结构的细微区别, 图1(b)和图1(c)分别给出了RuSe₂的这两种结构的俯视图和侧视图. 它们均可以被视为扭曲的T相, 关键的结构差异是相邻的Ru原子之间的距离. $T^\prime$相属于$P21/m$空间群, 其中的2个Ru原子和4个Se原子的Wyckoff位置都为2e, 晶格在b方向形成之字链结构, 并在a轴方向出现了明显的二聚化, 其最小重复单元是由两个Ru原子相靠近构成的2×1晶胞(相对于T相而言), 如图1(b)中的绿色矩形框所示. $T^{\prime\prime\prime}$相是$T^\prime$相在b 轴方向进一步发生周期性团聚形成的, 三个相邻的Ru原子靠近形成正三角形结构, 元胞大小约为$T^\prime$相的3倍. $T^{\prime\prime\prime}$同样属于P21/m空间群, 和$T^\prime$相不同的是, 其中有4个Ru原子和8个Se 原子的Wyckoff位置为4f, 2个Ru原子和4个Se原子的Wyckoff位置为2e, 如图1(c)中的红色矩形框所示, 详细的原子位置信息见补充材料表S1(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557).

$T^{\prime\prime\prime}$相单层RuSe₂结构的弹性模量计算为 $C_{11} = 97.80\; {\rm N}{\cdot}{\rm m}^{-1}$, $C_{22} = 99.31\; {\rm N}{\cdot}{\rm m}^{-1}$, $C_{12} = C_{21} = 23.66\; {\rm N}{\cdot}{\rm m}^{-1}$, $C_{66} = 33.21\; {\rm N}{\cdot}{\rm m}^{-1}$. 这些数据满足稳定性标准($C_{11}C_{22}-C_{12}C_{21}>0$, $C_{66}> 0$)^[45], 确保了这些材料的力学稳定性. 晶体结构的AIMD模拟结果如图2(a)所示, 在300 K下对体系进行了10 ps的动力学模拟后, 晶体结构没有明显的变形, 预示了这一三聚相畸变结构具有很好的热力学稳定性. 此外, 还计算了该晶相的声子谱, 如图2(b)所示, 仅在Γ点附近具有很小的虚频, 这类虚频一般认为是因计算采用有限晶胞造成的, 可以通过进一步扩胞来消除^[46], 因此也说明了这一$T^{\prime\prime\prime}$相RuSe₂单层结构的动力学稳定性. 同时也对$T^\prime$相RuSe₂单层结构进行了AIMD模拟和声子谱的计算(见补充材料图S2(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557)), 结果和之前已报道的结果一致^[34], 预示了该结构在热力学和动力学上具有稳定性.

图 2  (a) AIMD模拟300 K时$T^{\prime\prime\prime}$相单层RuSe₂总能量随时间的变化, 插图为始末状态的晶体结构图; (b) $T^{\prime\prime\prime}$相单层RuSe₂的声子谱

Fig. 2.  (a) Variation of the total energy of monolayer RuSe₂ in the $T^{\prime\prime\prime}$ phase during AIMD simulation at 300 K. The crystal structures of the initial and final states shown in the inset; (b) phonon spectrum of monolayer RuSe₂ in the $T^{\prime\prime\prime}$ phase.

下载: 全尺寸图片 幻灯片

进一步研究了上述几种不同晶相的电子结构, 从计算得到的能带结构来看, H相的单层RuSe₂展现出金属性质, 如补充材料图S3(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557), 而$T^\prime$相和$T^{\prime\prime\prime}$相结构则都是半导体, 并且两种结构都没有磁性. 首先对比了单层RuSe₂两种不同相使用PBE泛函在${U_\text{eff}} = 0.00 \; \text{eV}$, ${U_\text{eff}} = 4.36 \; \text{eV}$以及使用HSE杂化泛函计算的能带结构图, 如补充材料图S4(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557)和图S5(https://wulixb.iphy.ac.cn/article/doi/10.7498/aps.73.20240557)所示, 可以看到, 对于不考虑${\left(U_\text{eff}\right)}$的情况, 两种相的带隙都远小于HSE的结果, 对于${U_\text{eff}} = 4.36 \; \text{eV}$的情况, 带隙相比于${U_\text{eff}} = 0 \; \text{eV}$明显增大, 并且能带的色散分布和HSE的结果更为匹配, 因此基于计算效率和计算精度的综合考虑, 后面的计算都使用${U_\text{eff}} = 4.36 \; \text{eV}$的结果. 图3(a)和图3(b)分别为单层RuSe₂在两种不同相的轨道分辨能带投影图. 从图中可以看出, $T^\prime$相为间接带隙半导体, 带隙为1.11 eV; 而$T^{\prime\prime\prime}$相则为直接带隙半导体, 带隙为0.71 eV. 对于这两种结构, 费米面附近都是由Ru原子的d_xy, d_yz, d_xz贡献的, 而费米面以上主要是由Ru原子的${\text{d}}_{z^2}$和$\text{d}_{x^2 - y^2}$轨道贡献的. 进一步分析能带的轨道成分可以发现, Ru-d轨道主要还是按八面体晶体场的能级分裂规律进行排布的, 分成了能量较高的${{\mathrm{e}}}_{\rm g}$和较低的${{\mathrm{t}}}_{\rm 2 g}$两组. 八面体的示意图和计算所选取的笛卡尔坐标系如图3(c)和图3(d)所示. 由于两种晶相下的Ru-Se八面体都发生了或多或少的形变, 导致${{\mathrm{e}}}_{\rm g}$和${{\mathrm{t}}}_{\rm 2 g}$轨道进一步分裂. 对于$T^\prime$相结构, 由于八面体中靠近xy 平面的四根键并非严格的垂直于z轴, 所以d_yz和${{\mathrm{d}}}_{xz}$的能量被提升, 和d_xy轨道的能量不再简并, 如图3(e)所示. 而靠近z轴方向上的两Ru-Se键本身也略微偏离了z轴, 因此$\text{d}_{z^2}$和$\text{d}_{x^2-y^2}$轨道的简并也被打开了. 对于$T^{\prime\prime\prime}$相结构, 三聚化的结构使其元胞扩大为$T^\prime$相的3倍, 而且胞内的Ru-Se八面体发生了更大程度的扭曲变形, 再加上元胞扩大造成的能带折叠, 导致各轨道能级出现了明显的展宽(如图3(e)所示). 这种展宽使得$T^{\prime\prime\prime}$相的电子结构发生了较大变化, 带隙减小了约30%, 并由间接带隙转变成了直接带隙.

图 3  (a) 单层RuSe₂在$T^\prime$相的能带结构; (b) 单层RuSe₂ 在$T^{\prime\prime\prime}$相的能带结构; (c) $T^\prime$相中的Ru-Se八面体结构; (d) $T^{\prime\prime\prime}$相中的Ru-Se八面体结构; (e) Ru-Se 八面体的Ru-d轨道在不同晶相结构中的分裂情况示意图

Fig. 3.  (a) Band structure of monolayer RuSe₂ in the $T^\prime$ phase; (b) band structure of monolayer RuSe₂ in the $T^{\prime\prime\prime}$ phase; (c) octahedral structure of Ru-Se in the $T^\prime$ phase; (d) octahedral structure of Ru-Se in the $T^{\prime\prime\prime}$ phase; (e) schematic illustration of the splitting of Ru-d orbitals in the Ru-Se octahedron in different crystal phases.

下载: 全尺寸图片 幻灯片

载流子迁移率(µ)是半导体材料的一个重要特性, 它表征了电子或空穴的移动速度^[36]. 从图3可以看到, $T^\prime$相和$T^{\prime\prime\prime}$相的能带结构有着显著的区别, 这很可能会使单层RuSe₂在两种相下的迁移率有很大的差别. 为此, 我们基于形变势理论计算了两者的电子和空穴在armchair和zigzag方向 上的迁移率, 结果如表1所列. 无论是$T^\prime$相还是$T^{\prime\prime\prime}$相, 空穴的迁移率都远高于电子, 达到10倍以上. 而且, 单层RuSe₂还体现了较为显著的载流 子迁移率的各向异性, 在zigzag方向上的迁移率要高于armchair方向. 对比两种相, 可以发现 $T^{\prime\prime\prime}$相的迁移率相较于$T^\prime$相有大幅提高, 在zigzag和armchair方向分别为$3.22 \times 10^3\, {\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}$和$1.18 \times 10^3\, {\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}$, 和磷烯材料的空穴迁移率很接近^[43,47,48]. 分析其原因可以发现, 载流子迁移率的提高主要是因为该体系相变后带边处电子和空穴有效质量$\left({m}^*\right)$的减小造成的, 而弹性模量$\left(C_{\rm 11} , C_{\rm 22}\right)$和形变势$\left(E_{\rm 1} \right)$的变化则并不显著. 价带顶附近的空穴有效质量降低了一半以上, 再加上迁移率是和有效质量的二次方呈反比的, 因此相变后的空穴迁移率有了大幅的提升.

表 1  计算得到的单层RuSe₂在armchair和zigzag方向上的载流子有效质量$m^*\left({m}_{0}\right)$, ${{m}}_0$为单个电子的质量; 弹性模量$C\;\left( {\rm N}{\cdot}{\rm m}^{-1}\right)$; 形变势$E_1 \;(\text{eV})$和迁移率$\mathit{µ}\;\left({\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}\right)$.

Table 1.  Calculated effective mass of carriers $m^*\left({m}_{0}\right)$, here the ${m}_{0}$ represents the mass of a single electron; elastic modulus $C\;\left({\rm N}{\cdot}{\rm m}^{-1}\right)$; deformation potential $E_1 \; \left(\text{eV}\right)$; and mobility $\mathit{µ}\; \left({\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}\right)$ of monolayer RuSe₂ in the armchair and zigzag directions.

Phase/Carrier	$m^*_{{\mathrm{arm}}}$	$m^*_{{\mathrm{zig}}}$		$C_{11}$	$C_{22}$		$E_{1 {\mathrm{arm}}}$	$E_{1 {\mathrm{zig}}}$		$\mu_{\rm arm}$	$\mu_{\rm zig}$
	${m}_{0}$			$/({\rm N}{\cdot}{\rm m}^{-1})$			/eV			$/\left({\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}\right)$
$T^\prime \text{/e}$	37.59	9.64		95.73	92.01		0.74	–2.87		2.44	3.53
$T^\prime \text{/h}$	3.21	1.22					1.46	–0.71		224.42	823.25
$T^{\prime\prime\prime}\text{/e}$	3.48	1.16		97.80	99.31		–1.46	–2.98		70.53	150.57
$T^{\prime\prime\prime} \text{/h}$	1.10	0.60					1.78	0.50		1182.65	3219.20

下载: 导出CSV 

| 显示表格

从上述结果可以看到, $T^{\prime\prime\prime}$相的电子结构相较于$T^\prime$相发生了较大的变化, 载流子迁移率更是有了大幅的提高, 为此本文进一步研究单层RuSe₂的$T^{\prime\prime\prime}$相的形成性问题. 图4(a)为计算得到的两种相的能量随zigzag方向晶格长度的变化. 可以看到, 当晶格由$T^\prime$相的最低能量位置$\left( {L}=10.72 \;\text{Å} \right)$压缩到10.40 Å时, $T^{\prime\prime\prime}$相的能量开始低于$T^\prime$相. 但此时计算得到的两种结构发生转变的NEB势垒相对较高, 达到141 meV/f.u.. 两相的能量差随着晶格的压缩会不断增加, 而转变势垒则不断减小. 因此, 我们猜测, 在一定的温度下, 压缩应变可能使得单层RuSe₂越过转变势垒, 由$T^\prime$相转变到$T^{\prime\prime\prime}$. 我们在室温(300 K)下对多个位置的结构进行了AIMD模拟, 发现在图4(a) 所示的$L=10.00\;\text{Å}$位置处(压缩应变约为6.5%), 单层RuSe₂有可能从$T^\prime$相自发转变为$T^{\prime\prime\prime}$相, 对应的转变势垒为96 meV. AIMD模拟的结果如图4(b)所示, 经过6 ps之后, 总能量有一个明显的跃变降低, 结构在此时发生了相变, 由初始的$T^\prime$相转变到了$T^{\prime\prime\prime}$相. 由此可见, 单层RuSe₂在实际情况下很有可能会因为局部的应力作用而发生结构相变, 出现两相共存的现象. 如果能在实验上实现对这类材料施加均匀应力, 则能引起材料发生整体的结构相变, 从而实现对其性能的调控.

图 4  (a) 单层RuSe₂在$T^\prime$与$T^{\prime\prime\prime}$相下能量随zigzag方向晶格长度变化的关系曲线, 插图为部分位置处发生相变的NEB势垒图; (b) 300 K下, $T^\prime$ 相在zigzag方向施加6.5%压缩应变时的分子动力学模拟, 其中的插图分别是0和15 ps的结构图

Fig. 4.  (a) The energy of monolayer RuSe₂ in the $T^\prime$ and $T^{\prime\prime\prime}$ phases variation with the lattice along the zigzag direction. The insets are the NEB barriers showing phase transitions at certain positions; (b) AIMD of the $T^\prime$ phase under 6.5% compressive strain along the zigzag direction at 300 K, with insets showing crystal structures at 0 and 15 ps.

下载: 全尺寸图片 幻灯片

4.   总　结

通过第一性原理计算, 对单层RuSe₂的不同晶相结构进行了深入研究, 发现除了报道的H相、T相和$T^\prime$相之外, 还存在着$T^{\prime\prime\prime}$相的稳定态. 通过能量比较以及机械稳定性、热力学稳定性和动力学稳定性的分析, 预测了$T^{\prime\prime\prime}$相结构实际存在的可能性. 对单层RuSe₂的电子结构进行了分析, 发现$T^{\prime\prime\prime}$相的能带相较于$T^\prime$相发生了较大的变化, 由间接带隙半导体转变为了直接带隙半导体, 且载流子迁移率有了大幅的提升. 特别是$T^{\prime\prime\prime}$相的空穴迁移率达到了$3.22 \times 10^3\; {\rm cm}^{2}{\cdot}{\rm V}^{-1}{\cdot}{\rm s}^{-1}$, 远高于$T^\prime$相, 和磷烯材料的空穴迁移率相接近. 此外, 研究发现当在zigzag方向发生应变时, $T^\prime$相向$T^{\prime\prime\prime}$相的转变势垒将不断减小, 在300 K时进行的AIMD模拟预测了单层RuSe₂在应力作用下会出现$T^\prime$相向$T^{\prime\prime\prime}$相转变的现象. 本工作不仅拓展了二维RuSe₂的研究和应用范围, 也为分析其他二维TMDs材料的结构相变和多相共存等现象提供了理论参考.