应变诱导单层NbSi2N4 材料磁转变的第一性原理研究*

姜楠 李奥林 蘧水仙 勾思 欧阳方平†

1) (中南大学物理与电子学院,长沙 410012)

2) (新疆大学物理科学与技术学院,乌鲁木齐 830046)

3) (中南大学粉末冶金研究院,粉末冶金国家重点实验室,长沙 410083)

二维材料磁性的有效调控属于国内外的前沿研究领域.本文运用基于密度泛函理论的第一性原理方法,研究了双轴拉伸应变对单层NbSi2N4 磁性的影响.声子谱和分子动力学的计算结果表明,单层NbSi2N4 结构具有良好的动力学与热力学稳定性.研究发现单层NbSi2N4为无磁金属,1.5%的双轴拉伸应变可使其转变为铁磁金属.对单层NbSi2N4 材料电子结构的分析表明,拉伸应变诱导的铁磁性具有巡游电子起源: 当不考虑自旋极化时,单层NbSi2N4 在费米能级处存在一条半满的能带,其主要由Nb 原子的dz2 轨道贡献,拉伸应变可使其更局域化,进而引起斯通纳不稳定性,导致铁磁性的产生.此外,对磁各向异性能的计算表明,应变可使单层NbSi2N4 的易磁化轴方向发生垂直-面内-垂直方向的翻转.基于海森伯模型的蒙特卡罗模拟结果表明,拉伸应变可显著提高单层NbSi2N4 的居里温度.单层NbSi2N4 的居里温度在2%应变时为18 K,在6%应变时提高到87.5 K,比2%应变时提高了386%.本研究为应变调控二维层状材料的磁性提供了理论参考,在力学传感器设计和低温磁制冷领域有着潜在的应用前景.

1 引言

以硅为主导的传统晶体材料已经无法满足电子器件的发展需求.二维磁性材料因独特的电子特性,在高速、高集成密度、低功耗逻辑运算与信息存储等方面具有强大的优势,是纳米级自旋电子的理论研究与实际应用的新平台[1-4].二维材料主要通过磁性掺杂[5]、近邻效应[6]或缺陷调控[7,8]等方法来引入磁性,由这些传统手段产生的磁性微弱且难以控制.近年来,具有本征磁性的二维材料如CrI3[9],Cr2Ge2Te6[10]等陆续在实验中成功合成,引起了研究者广泛研究兴趣.

探索新的二维材料可以创造新的性能和进一步的潜在应用.MA2Z4(M表示过渡金属元素,A表示Si 或Ge,Z表示VA 族元素)是一种新的人造二维材料,没有任何已知的三维层状近亲,通过过渡金属和层终端基团的适当组合可以实现任何期望的电子、磁性或催化性能,MA2Z4中元素的多样性使得它们的间隙和磁性能具有宽的可调性,这对于电子学、光电子学和自旋电子学的应用是必不可少的[11-16].随着Hong等[12]首次利用化学气相淀积(chemical vapor deposition,CVD)成功生长MoSi2N4和WSi2N4,掀起了对一系列具有独特七原子层Z-A-Z-M-Z-A-Z结构的二维MA2Z4家族的研究热潮.Chen等[16]通过第一性原理计算预测NbSi2N4为铁磁性(ferromagnetic,FM)金属,但该结论仅基于FM 和反铁磁(antiferromagnetic,AFM)两种磁序间的能量差,对单层NbSi2N4的磁基态有待进一步验证.NbSi2N4作为最近备受关注的MA2Z4体系成员之一,确定其正确的磁基态对后续的理论和实验研究均具有重要的参考意义.

应变可以有效调控二维材料的电子结构和磁性能[17-25].例如,单层CrI3在6%压缩应变下发生从FM 构型向AFM 构型的转变,在拉伸应变下磁取向由平面外翻转到平面内[19,20];当含单空位缺陷VSe和多空位缺陷VGaInSeTe的GaInSeTe 单层膜分别在拉伸应变超过8%和4%时,其从无磁性(nonmagnetic,NM)材料变为固有磁矩大于1.25µB和0.27µB的磁性材料[21];MnBi2Te4/CrI3异质结的磁各向异性能(magnetocrystalline anisotropy energy,MAE)在—5%—5%的应变范围内有约10 倍的幅度提升[22].此外,MnBi2Te4/CrI3异质结在应变为5%时,居里温度(Curie temperature,TC)达到91 K,比无应变时高13.8%[22];实验测得单层Cr2Ge2Te6的TC为22 K,在10%应变下,TC提高了191%[10,24];当施加应变达到10%时,二维过渡金属二硫化合物MnSe2的TC从300 K 提高到530 K[25],这些结果表明应变也是提高材料TC的一种可行方法.

在当前对二维材料的研究中,关于应变诱导FM-AFM 转变的报道很多,但呈现出NM-FM 转变性质的材料则较少,单层NbSi2N4为研究磁性转变现象提供了一个新的平台,具有一定的科学意义.虽然单层NbSi2N4金属性使其在自旋电子学中的应用受到了一定的限制,但其固有磁矩和TC可受应变调控,理解应变对单层NbSi2N4电子结构与磁性以及TC的影响为其在力学传感器设计和低温磁制冷等领域的应用具有重要意义.本文采用基于密度泛函理论的第一性原理方法,研究面内双轴拉伸应变对单层NbSi2N4磁性的影响,通过计算确定NM-FM 转变的临界应变,并对磁性转变的可行性和内在物理机制进行探究,进一步通过蒙特卡罗(Monte Carlo,MC)模拟预测了拉伸应变下单层NbSi2N4的TC.

2 计算模型与方法

本文基于密度泛函理论(density functional theory,DFT)的第一性原理计算均通过VASP(Viennaab-initiosimulation package)[26,27]软件实现.交换关联势采用广义梯度近似(generalized gradient approximation,GGA)下的PBE (Perdew-Burke-Ernzerhof)泛函[28],平面波基组采用投影缀加平面波(projected augmented wave,PAW)方法[29],平面波截断能设置为600 eV.结构优化采用共轭梯度算法,对晶格常数和原子位置均进行优化,优化总能收敛标准设置为1×10—7eV/cell,力收敛标准设置为—0.01 eV/Å (1 Å=0.1 nm).金属材料的总能计算需要稠密的k点网络,经测试发现只有当原胞的k点网络不低于12×12×1 时,对单层NbSi2N4总能的计算结果才能给出正确的磁基态,因此,在计算中第一布里渊区采样使用以Γ为中心的Monkhorst-Pack 方法[30]生成20×20×1 的网格.测试发现自旋轨道耦合(spin-orbit coupling,SOC)作用对单层NbSi2N4电子结构的影响很小,因此本文的计算均不考虑SOC 的作用.声子谱通过密度泛函微扰理论(density functional perturbation theory,DFPT),利用开源软件包PHONOPY实现[31],使用3×3×1 的超胞;在第一性原理分子动力学(ab initiomolecular dynamics,AIMD)方法模拟中,使用正则系统(canonical ensemble,NVT)条件,模拟温度为300 K,电子步的能量收敛标准为1×10—4eV/cell,原子运动的时间步长为10 ps.对单层NbSi2N4居里温度的计算,采用基于海森伯模型和Metroplis 算法的蒙特卡罗方法,通过开源软件Mc_solver 实现[32，33].

3 计算结果与讨论

3.1 单层NbSi2N4 的电子结构及稳定性

结构稳定是二维材料可实验制备的前提条件之一,本文首先通过声子谱和AIMD 计算验证了单层NbSi2N4的动力学和热力学稳定性.由图1(c)可知,单层NbSi2N4的声子谱仅ZA 支在Γ点存在可忽略的虚频,这主要由计算精度的限制引起,证明了其动力学稳定性.图1(d)给出了单层NbSi2N4总能量随时间的波动,能量涨落均在20 meV/atom左右,表明其具有良好的热力学稳定性,能够在室温下稳定存在.

图1 (a)单层NbSi2N4 的俯视图(左)和侧视图(右),自旋密度分布用黄色表示;(b)不考虑自旋极化时单层NbSi2N4 的能带结构和态密度图;(c)单层NbSi2N4 的声子谱图;(d)温度300 K,时长10 ps 的分子动力学模拟下,系统总能量的变化;(e)考虑自旋极化时单层NbSi2N4 的能带结构和态密度图Fig.1.(a) Top view (left) and side view (right) of monolayer NbSi2N4,and the spin density distribution is represented in yellow;(b) energy band structure and density of states of monolayer NbSi2N4 without considering spin polarization;(c) phonon spectra of monolayer NbSi2N4;(d) the change of total energy of the system under the molecular dynamics simulation of temperature 300 K and duration 10 ps;(e) energy band structure and density of states of monolayer NbSi2N4 considering spin polarization.

如图1(a)所示,单层NbSi2N4由N-Si-N-Nb-N-Si-N 7 个原子层组成,该结构关于Nb 原子呈镜面对称.对于NM 和FM 两种磁结构,优化后单层NbSi2N4的晶格常数a0均为2.97 Å,与文献结果相近[16].图1(b),(e)分别给出单层NbSi2N4在NM态和FM 态下的电子结构.当不考虑自旋极化时,一条半满的能带穿过费米能级,使单层NbSi2N4表现出NM 金属性;考虑自旋极化后,这条半满能带的自旋简并解除,每个原胞具有0.40µB的净磁矩,从自旋密度分布图1(a)可知,其主要由Nb 原子贡献.为确定单层NbSi2N4正确的磁基态,分别考虑NM,FM 与AFM 三种自旋结构.单层NbSi2N4单位原胞内FM 态与NM 态间的能量差记为ΔEFM-NM,Zigzag 型AFM 结构平均到每个原胞的总能与NM态间的能量差记为ΔEAFM-NM.根据总能计算结果显示,ΔEFM-NM=2.012 meV/cell,ΔEAFM-NM=1.8745 meV/cell,这表明单层NbSi2N4的基态为NM 态,而非文献报道的FM态[16],这限制了其在自旋电子学中的应用.

由于ΔEFM-NM非常小,磁基态或受到温度的影响.在有限温度下,电子能级占据满足费米-狄拉克函数:

其中,σ=kBT表示能级展宽,kB为玻尔兹曼常数,T为电子温度.在DFT 计算中,通过改变展宽σ可定性评估温度对材料体系电子性质的影响.已有文献表明,利用该方法估算的电荷密度波相的转变温度与实验观测值定性一致[34],证明了该方案的有效性.图2(a)模拟了无应变时单层NbSi2N4中ΔEFM-NM与电子温度T间的函数关系,结果显示温度不超过500 K 时,ΔEFM-NM恒为正值,表明在室温或低温下单层NbSi2N4的基态始终为NM 态.

图2 (a)单层NbSi2N4 铁磁态与无磁态间的能量差随电子温度的变化关系;(b)单层NbSi2N4 在2%应变下的声子谱;(c) 单层NbSi2N4 6%应变下的声子谱Fig.2.(a) Variation of energy difference between ferromagnetic and nonmagnetic states of monolayer NbSi2N4 with electron temperature;(b) the phonon spectrum of monolayer NbSi2N4 at 2% strain;(c) the phonon spectrum of monolayer NbSi2N4 at 6% strain.

为使单层NbSi2N4转变为FM 态,最简单的方法是施加外磁场.外磁场可在FM 态和NM 态间引入额外的能量差:

其中,B为磁感应强度,M为每个Nb 原子的固有磁矩.为使EFM＜ENM,应有ΔEB＞ ΔEFM-NM.当温度为600 K 时,ΔEFM-NM＜ 0,此时单层NbSi2N4由NM 转变为FM 基态,取其所对应的能量差值—0.054 meV/cell 来预测将单层NbSi2N4磁性转变所需的最小临界磁场强度B,可得B=2.3 T,这接近于大型粒子加速器中所需磁场.因此,通过外磁场控制单层NbSi2N4由NM 转变为FM,虽理论上具有可行性,但并不适用于低功耗、高集成度的自旋电子器件应用.在本工作的后续部分,将演示通过应变可灵活调控单层NbSi2N4的NM-FM 相变.

3.2 应变对单层NbSi2N4 磁性的影响

面内双轴拉伸应变的定义为ε=(a—a0)/a0,其中a为施加应变后单层NbSi2N4的晶格常数.对结构施加0%—6%的拉伸应变,在拉伸范围内,单层NbSi2N4的结构未发生明显形变.图2(b),(c)给出单层NbSi2N4应变为2%和6%时的声子谱图,结果表明,在拉伸应变的作用下,单层NbSi2N4的声子谱也仅ZA 支在Γ点存在可忽略的虚频,因此施加应变不会影响体系的动力学稳定性.

如图3(a)所示,无外加应变时,ΔEFM-NM和ΔEAFM-NM都大于0,NM 态时体系的能量最低,说明单层NbSi2N4基态是NM 态;当应变增加到1.5%时,ΔEFM-NM＜ 0,此时单层NbSi2N4基态由NM转变为FM;在1.5%—6%的拉伸应变范围内,ΔEFM-NM始终小于ΔEAFM-NM,表明单层NbSi2N4在一定程度的应变下具有FM 基态.如图3(b)所示,1.5%应变下,Nb 原子的磁矩从0 变为0.4473µB,并随着应变的增大而增大,在6%应变时增至0.7548µB.由应变诱导Nb 原子产生固有磁矩,其大小不受热扰动影响,因此在室温下单层NbSi2N4顺磁态的磁化率也受应变调控,使其在传感器方面具有一定的应用潜力.

为确定磁矩的空间取向,定义MAE=E100—E001,其中E100和E001分别为自旋空间取向沿面内及面外方向时FM 的总能.MAE为正表示材料具有垂直磁各向异性(perpendicular magnetic anisotropy,PMA),反之则具有面内磁各向异性(inplane magnetic anisotropy,IMA).单层NbSi2N4在无应变时具有NM 基态,因此这里只考虑应变大于1.5%的情况.如图3(c)所示,单层NbSi2N4的MAE 在应变下呈“M”型非单调变化: 当应变1.5% ≤ε≤ 2.7%时,单层NbSi2N4具有PMA,其中,当ε=1.5%时MAE为0.081 meV,而当ε=2.5%时MAE 增至0.439 meV,较1.5%应变时增大了5 倍左右.当应变ε继续增大到3%时,单层NbSi2N4由PMA 转变为IMA;而当应变ε≥ 4%时,单层NbSi2N4又从IMA 变为PMA,这表明应变可作为调控单层NbSi2N4磁各向异性的有效手段.

图3 (a) 单层NbSi2N4 反铁磁态和铁磁态与无磁态间能量差随应变变化;(b)单层NbSi2N4 磁矩随应变变化;(c) MAE 随应变变化;(d)不同应变下,Nb 原子轨道磁矩与自旋磁矩对MAE 的贡献Fig.3.(a) Energy difference between antiferromagnetic state,ferromagnetic state and non-magnetic state of monolayer NbSi2N4 with strain;(b) magnetic moment of monolayer NbSi2N4 with strain;(c) MAE with strain;(d) contribution of orbital magnetic moment and spin magnetic moment of Nb atom to MAE under different strains.

为理解应变下单层NbSi2N4磁各向异性变化的原因,考虑其SOC 的哈密顿量[35]:

其中,ξ为Nb 原子的自旋轨道耦合常数,L和S分别表示原子的轨道磁矩和自旋磁矩.由于晶体场的作用使Nb 原子4d 轨道电子的轨道角动量冻结,为此必须考虑(3)式的微扰项,其对单层NbSi2N4磁各向异性能的贡献为

为便于后续分析,定义ΔL和ΔS分别为Nb 原子的磁矩沿[100]和[001]方向时轨道磁矩和自旋磁矩的变化量,即:

图3(d)给出了ΔL和ΔS随应变的变化规律,其中ΔS的相对变化量小于0.01%,可将(4)式进一步改写为

随应变增大,ΔL呈现“W”形变化,与图3(c)中MAE 的变化在总体上呈相反趋势,这表明单层NbSi2N4磁各向异性能的改变主要由ΔL的变化引起.

3.3 单层NbSi2N4 铁磁性的起源

为了理解单层NbSi2N4在拉伸应变下FM 产生的物理机制,首先分析单层NbSi2N4的态密度,如图4(a)—(e)所示.不考虑自旋极化时,单层NbSi2N4的费米能级处态密度的峰主要由Nb 原子的轨道贡献.在0—6%应变范围内,费米能级处的态密度随着拉伸应变增大呈增大的趋势,说明拉伸应变使费米能级处的能带更加局域化,局域化导致其态密度峰值升高.根据斯通纳不稳定性,当费米能级处的态密度达到一定程度时,材料便转换成FM态,即当过渡金属在费米能级处的态密度达到一定程度时,巡游电子在同一格点上的库伦作用会导致FM.当考虑自旋极化时,Nb-轨道的自旋简并解除,使费米面附近的态密度减小,并随着应变的增大而不断减小,自旋极化导致费米能级处态密度减小,这有利于降低系统能量,从而增强体系的稳定性.这表明单层NbSi2N4在拉伸应变下FM的产生可通过斯托纳不稳定性来解释.

为进一步验证单层NbSi2N4的FM 是否具有巡游电子起源,考虑不同应变下磁矩M与总能ΔE间的关系,如图4(f)所示,并通过如下公式进行拟合[36]:

图4 (a) 单层NbSi2N4 分波态密度;(b)—(e) 单层NbSi2N4 在不同应变下的态密度;(f)铁磁态与无磁态间的能量差随磁矩变化的DFT 计算曲线,将无磁态的能量设为0 meVFig.4.(a) Fractional density of states of monolayer NbSi2N4.(b)—(e) The density of state of monolayer NbSi2N4 under different strains.(f) DFT calculation curve of the energy difference between the FM state and the NM state change with the magnetic moment,the energy of the NM state is set to 0 meV.

其中,ΔE为单层NbSi2N4单位原胞内FM 态与NM 态间的能量差,系数a2与斯托纳系数IS存在对应关系:

其中,N(EF)为非自旋极化系统费米能级处的态密度.根据斯托纳理论,自发磁化的必要条件为[37]

由(9)式和(10)式可知,若材料具有巡游FM,则有

如图4(f)所示,在0%应变下,ΔE在M=0时取得最小值,表明单层NbSi2N4的基态为NM;当M=0.4µB时,ΔE达到局部极小 值,说 明FM 是作为体系的亚稳态而存在,这与自旋极化计算的结果一致.拉伸应变达到2%时,ΔE在M=0时转变为局部极小值;而当M=0.5µB时,ΔE转变为全局最小值,表示单层NbSi2N4在2%应变附近发生NM 到FM 转变.当应变为6%时,M=0不再是局部最小值,仅FM 在能量上是稳定的.此外,随着应变的增加总能的最小值持续降低,说明增加一定程度的应变可以有效增强FM 耦合的稳定性.

通过最小二乘法拟合式(8)式可得到在0%,2%和6%应变时,系数a2分别为4.029,—23.332 和—70.923 meV,对应的N(EF)×IS分别为0.99,1.06和1.27.这表明在0%应变时不足以诱导出巡游电子磁性;在2%应变时,N(EF)×IS略大于1,对应FM-NM 转变的临界点;在6%应变时,N(EF)×IS达到1.27,FM 变得更稳定.上述拟合结果与图3(a)中第一性原理计算相一致,证实了单层NbSi2N4的FM 具有巡游电子起源.

3.4 单层NbSi2N4 的居里温度

居里温度TC是二维磁性材料的关键性质之一.由于平均场理论容易高估TC[38],本文使用基于海森伯模型的MC 方法来计算不同应变下单层NbSi2N4的TC.MC 算法通过开源项目Mc_solver实现,该软件已成功应用于其他二维材料TC的预测,如单层FeTe2、单层FeSi2,Fe 与Co 掺杂的WS2单层等[33,39,40].海森伯哈密顿量可表示为

其中,Si为格点i处Nb 离子的自旋,J代表最近邻Nb 离子间的磁交换作用常数,D是单离子磁各向异性能参数.(10)式中的计算参数均可通过第一性原理方法得到.对于单层NbSi2N4,取S=1/2,J=(EFM-EAFM)/(4S2),D=MAE/S2.不同应变下对应的J和D如表1 所示,不考虑没有磁性的情况.

表1 在不同应变下的交换常数J 和各向异性参数DTable 1.Exchange constant J and anisotropy parameter D under different strains.

从图5 可以看出,在2%临界应变下,NbSi2N4的磁矩与磁化率在18 K 处同时达到峰值,所对应的温度即为单层NbSi2N4的TC;TC随着应变的增大而增大,当应变达到6%时,TC增加到87.5 K,相对于2%应变时增大了386%.这些结果表明,单层NbSi2N4的TC受应变显著调控,这一特点可作为磁热工质应用于低温磁制冷领域中,应变或可有效降低所需的磁场强度.

图5 (a)不同应变下居里温度的变化;(b),(c)在2%和6%应变下磁矩与磁化率随居里温度的变化Fig.5.(a) Variation of Curie temperature under different strains;(b),(c) variation of magnetic moment and susceptibility with Curie temperature under 2% and 6% strain.

4 结论

本文通过基于密度泛函理论的第一性原理计算和海森伯模型下的MC 模拟研究面内双轴拉伸应变对单层NbSi2N4电子结构、磁特性和TC的影响.单层NbSi2N4动力学稳定性以及热力学稳定性通过声子谱分析与分子动力学模拟得到证实.无外加应变时,理论预测单层NbSi2N4的基态为NM,FM 作为其亚稳态存在.在施加1.5%应变时,单层NbSi2N4由NM 向FM 转变,MAE 表现为PMA.随着应变的增大,MAE 发生PMA-IMA-PMA 翻转.当无外加应变时,单层NbSi2N4费米能级处存在一条由Nb 原子的dz2轨道贡献的自旋简并的半满能带,拉伸应变可使该能带更加局域化,引起斯通纳不稳定性,导致单层NbSi2N4巡游电子呈现FM 有序.增大一定程度的应变,可以有效增强FM耦合的稳定性.同时,在2%—6%应变范围内,TC由18 K 变至87.5 K.应变工程可有效调控单层NbSi2N4的磁基态和TC,研究结果有望促进MA2Z4材料在力学传感器件设计和低温磁制冷领域的发展.

感谢中南大学高性能计算中心.