Strain-induced effects on the electronic properties of 2D materials

Abstract

Thanks to the ultrahigh flexibility of 2D materials and to their extreme sensitivity to applied strain, there is currently a strong interest in studying and understanding how their electronic properties can be modulated by applying a uniform or nonuniform strain. In this work, using density functional theory (DFT) calculations, we discuss how uniform biaxial strain affects the electronic properties, such as ionization potential, electron affinity, electronic gap, and work function, of different classes of 2D materials from X-enes to nitrides and transition metal dichalcogenides. The analysis of the states in terms of atomic orbitals allows to explain the observed trends and to highlight similarities and differences among the various materials. Moreover, the role of many-body effects on the predicted electronic properties is discussed in one of the studied systems. We show that the trends with strain, calculated at the GW level of approximation, are qualitatively similar to the DFT ones solely when there is no change in the character of the valence and conduction states near the gap.

Introduction

Much of the emphases on 2D materials ^{1 –3} was born with the discovery of graphene, for which the Nobel Prize in physics was awarded to Novoselov and Geim in 2010. ⁴ Graphene is a 2D crystal made up of carbon atoms arranged in a hexagonal honeycomb form; it is one million times thinner than paper, almost transparent, and, at the same time, is the strongest material in the world. ^5,6 Its electronic structure can be easily derived from a simple tight-binding model, which explains the presence of bands with conical dispersion intersecting at the Fermi level, thus making graphene a semimetal. ⁷ Massless Dirac fermions move in graphene as fast as $v_F{=10}^6\text{\hspace{0.17em}m/s}$ at the Fermi level, and twisted bilayers of graphene have even been shown to be superconductors for very small twisting angles. ⁸ This peculiar linear dispersion of low-energy carriers, which can be mapped to an effective 2D Dirac Hamiltonian, is very different from the usual parabolic dispersion of bulk semiconductors. This has stimulated a lot of work to theoretically predict and experimentally observe novel physical effects in this 2D material. Nevertheless, despite the enormous interest both at the fundamental and applicative level, the lack of an electronic gap limits graphene use in applications like digital electronics, field-effect transistors, and optoelectronics at visible frequencies. For this reason, the last years have witnessed growing efforts to find a way to open its gap, leaving unmodified its peculiar electronic and mechanical properties, but also to grow and characterize other novel metallic and semiconducting two-dimensional materials beyond graphene.

Following the route of graphene, broad families of 2D materials are hence continuously developed and studied in view of their interesting physical properties and of a large number of their envisaged device-oriented applications. ^9,10 Because of their atomic-scale thickness, they are characterized by weak dielectric screening, strong light–matter interaction, and highly bound excitons. Moreover, atomic and molecular doping, ^{11 –13} external fields, ¹⁴ and also strain ^{15 –17} may have a very deep impact on their electronic and optical properties. In particular, strain engineering is very exciting since, differently from 3-D traditional materials, 2D materials can endure remarkably large mechanical strain (up to 10%), hence creating opportunities to modulate their physical properties for interesting device applications.

Silicene, ^{18 –20} the silicon-based counterpart of graphene, represents the first exciting material merging the exceptional physical properties of graphene with the simplicity of easily integrating it in already existing and largely developed silicon-based technology. ^21,20 A field-effect transistor has been reported at room temperature. ²² Further interest in silicene arises from its predicted nontrivial topological properties. ²³ Freestanding ideal silicene presents a buckled honeycomb structure. ¹⁹ Like graphene, silicene has a semimetallic behavior and possesses (in the absence of spin–orbit corrections) massless fermions at the k-point of the brillouin zone (BZ). ^24,25

Other members of the so-called X-enes family (borophene, germanene, stanene, phosphorene, arsenene, antimonene, bismuthene, and tellurene) are also of particular interest for their excellent physical, chemical, electronic, and optical properties. ²⁶

The evidence that bulk group-III nitrides are among the most important materials for solid-state lighting, as witnessed by the Nobel prize awarded in 2014 to Akasaki, Amano, and Nakamura, stimulated in the last years several theoretical ^27,28 and experimental ^{29 –31} studies on 2D honeycomb III-N sheets. From the experimental side, their growth is very challenging because, similar to the case of silicene, no simple route to mechanical or chemical exfoliation can be used, due to their 3-D wurtzite structures (only BN crystallizes in the hexagonal layered form in bulk). However, several promising experimental attempts to realize 2D III-nitrides have been reported so far.

Among layered materials, for which exfoliation is possible, transition metal dichalcogenide (TMD) occupies a prominent place in recent worldwide research. ^{32 –35} They are of type MX ₂, where one layer of transition metal atoms (M) is sandwiched between two layers of chalcogen (X) atoms, crystallizing mainly in the hexagonal or rhombohedral forms with metal atoms having octahedral or trigonal prismatic coordination. Among this broad family, group-VI TMD ( $M{X}_2$ with M = Mo, W, and X = S, selenium, tellurium (Te)) has received a lot of attention due to their tunable bandgap, strong light–matter interaction, ^36,37 large spin–orbit effects, strong influence of doping, functionalization, external field, and strain. In the hexagonal form, they are direct gap semiconductors when monolayer (ML) while exhibiting indirect gaps when BL and thicker multilayers, being promising materials for flexible electronics, light emission, energy storage, solar energy conversion, as well as electrochemical catalysis and biosensors. ^{38 –41}

In this article, we use first-principle DFT calculations to study the role of biaxial uniform strain on the structural and electronic properties of 2D materials beyond graphene: silicene, group III-nitrides, and TMD. In particular, we focus on the behavior of their band structure, electronic gap, workfunction (WF), ionization potential (IP), and electron affinity (EA) and perform a comparative analysis for nonpolar (graphene and silicene), polar (aluminum nitride (AlN) and gallium nitride (GaN)) honeycomb lattices, and TMD (molybdenum disulfide (MoS₂) and molybdenum ditelluride (MoTe₂)) in the monolayer and bilayer form.

Methods and computational details

All the DFT calculations to obtain the structural and electronic properties have been performed using the Quantum-Espresso code ⁴² within the local density approximation (LDA), ⁴³ perdew becke ernzerhrof (PBE), ⁴⁴ XC functional for nitrides (TMD). For the case of the two X-enes (graphene and silicene), we performed calculations using both LDA and PBE with the specific goal to look at the effect of local or semilocal XC on the calculated values of the work functions. The van der waals (vdW) correction ⁴⁵ has been applied to take into account the interaction between the layers in TMD BL.

For X-enes and nitrides (TMD), we used scalar (fully)-relativistic optimized norm-conserving ⁴⁶ pseudopotentials from the quantum-espresso (QE) repository.

The self-consistent density for graphene and silicene was computed using a 15 × 15 × 1 k-point sampling for the ground state with plane wave cutoff of 100 Ry. For nitrides, we used a nonshifted 18 × 18 × 1 k-point sampling for the ground state and a plane wave cutoff of 100 (200) Ry for the structural optimization and band structure of AlN (GaN). For TMD, we sampled the BZ using a uniform 18 × 18 × 1 mesh and used a plane wave cutoff of 120 Ry.

To simulate isolated layers, vacuum sizes of 15 and 20 Å have been chosen for X-enes and nitrides, respectively, while 17.8 (16.5) Å has been used for TMD-ML (BL). The geometry of each system has been relaxed at each value of strain until the forces on the atoms were less than ${10}^{-2}$ eV ${\AA }^{-1}$ . The values of applied strain are reported in the form $\epsilon =\frac{a-a_0}{a_0}$ , where a ₀ is the equilibrium lattice constant. Positive values correspond to tensile strain, while negative values correspond to compressive strain. The equilibrium geometries at zero strain are in very good agreement with the existing literature ^47,48 and the main equilibrium structural parameters are reported in Table 1.

Table 1.

Equilibrium lattice constant a ₀, buckling for X-enes and nitrides, and X-M-X angle for TMD.

	a 0 (Å)	Buckling (Å)/angle (°)	d 1 (Å)	d 2 (Å)
Graphene	2.47	0	—	—
Silicene	3.87	0.45 Å	—	—
GaN	3.17	0	—	—
AlN	3.03	0	—	—
MoS2-ML	3.18	81.0°	3.14	—
MoTe2-ML	3.56	82.6°	3.62	—
MoS2-BL	3.18	80.7°	3.12	3.05
MoTe2-BL	3.54	83.1°	3.62	3.46

d ₁: the vertical distance chalcogen–chalcogen; d ₂: vertical distance between layers; TMD: transition metal dichalcogenide; GaN: gallium nitride; AlN: aluminum nitride; MoTe₂: molybdenum ditelluride; MoS₂: molybdenum disulfide.

To assess the role of many-body effects, for one of the studied systems, that is, the MoTe₂ monolayer, we went beyond the single-particle approach. On top of the DFT simulations, we performed one-shot perturbative $GW$ calculations at $2\%,0,-2\%$ strain using the Yambo code. ^49,50 A cutoff of 240 Ry has been used for the exchange part of the quasiparticle corrections $<{\Sigma }_x-V_{xc}>$ , while 8 Ry and 100 bands were employed for the correlation part of the self-energy $<{\Sigma }_c>$ . To speed up the convergence with respect to the empty states, we adopted the technique described by Bruneval F and Gonze. ⁵¹ A plasmon–pole approximation for the inverse dielectric matrix has been applied. ⁵² Moreover, to guarantee the simulation of isolated layers, a cutoff in the bare Coulomb potential has also been used. ^50,53 The k-point sampling was selected to be 42 × 42 × 1 in the BZ.

Results

The main goal of this work is to describe how the electronic properties of different classes of metallic and semiconducting 2D materials can be tuned via the application of a uniform biaxial strain with the aim to highlight similarities and differences among the different systems.

For all the 2D materials, we have calculated the WF obtained as the energy difference between the Fermi energy E_F and the vacuum potential $E_{\text{vac}}$ . For 2D semiconductors, we report the electronic band structures and bandgap values at several biaxial uniform compressive and tensile strains, the IP ( $\text{IP}=E_{\text{vac}}-E_{\text{VBM}}$ ) and the EA ( $\text{EA}=E_{\text{vac}}-E_{\text{CBM}}$ ) to show their specific dependence on strain and to understand which are the different effects leading to the peculiar behaviors observed. For nitrides (TMD ML and their homobilayers), the DFT band structures have been calculated for applied strains ranging from −10% to 10% (−7% to 7% for MoS₂, −10% to 10% for MoTe₂) but, here, we report only those at zero strain and −4%, 4% (−5%, 5% for MoS₂ and MoTe₂). For each atomic structure at a given strain, all the quantities are determined self-consistently in a ground-state electronic structure calculation. The vacuum level is obtained as the asymptotic value of the electrostatic potential in the direction perpendicular to the layer at a far distance, 10 Å or more, from the system.

Graphene and silicene

For these two X-enes, we observe that for increasing compressive strain, the WF decreases, whereas it increases for increasing tensile strain (see Figure 1). For graphene, this behavior is in substantial agreement with the existing literature ^54,55 and it has been explained in terms of a strain-induced enhancement of the density of states (DOS) close to the Fermi level. ^56,57 Results for silicene are qualitatively similar to the ones obtained for graphene. In the tensile region, the WF variation due to a 10% of strain is found to be 0.15 eV, in good agreement with the results of Qin et al., ^58,59 where a change of about 0.2 eV has been obtained.

Figure 1.

Dependence of the WF by the strain value using PBE and LDA XC potentials. Top panel: results for graphene. Bottom panel: results for silicene.

We observe that, for both the systems, a simple argument of plausibility of this behavior is based on the fact that as long as the materials are stretched, the interaction among the ions of the lattice decreases, thus approaching the behavior of isolated atoms in the limit of infinite tensile deformation. The ionization potentials for C and Si atoms are 11.2 and 8.1 eV, respectively, and these values are much higher than the WF of their correspondent 2D forms at equilibrium, 4.24 eV for graphene, and 4.35 eV for silicene. For this reason, it is expected that the WF of these materials, characterized by fully covalent bonds, should grow with increasing uniform tensile deformation.

To test the possible dependence of our results on the choice of the XC potential used to describe the various systems, we have analyzed the strain dependence of the WF of graphene and silicene both within the LDA and PBE scheme for the XC term. The results of this analysis are shown in Figure 1 in different colors. It is noted that the choice of the functional affects the value of the WF by adding a constant shift of 0.24 eV for graphene and of 0.41 eV for silicene. In particular, the zero strain WF for graphene ranges from 4.24 eV in PBE to 4.48 eV in LDA, thus producing a spread of results, which is fully consistent with values in the range from 4.28 eV to 4.5 eV found in the literature. ^{54,55,60 –62} But, interestingly, we observe that the functional induced variation turns out to be a constant shift, which is independent of the strain value. Thus, we can conclude that the strain-induced variations of the WF do not depend on the choice of functional and can be assessed unambiguously within the DFT framework.

Nitrides monolayers: AlN and GaN

The 2D nitrides have the same flat honeycomb structure as graphene but with polar bonds due to the different electronegativity of the lattice atoms. To understand the trends of the IP, WF, and EA under strain, it is important to analyze their band structures. These are shown in Figure 2 for GaN and AlN at 0, −4%, and 4% uniform strain. In agreement with existing literature, ²⁷ we have found that at zero strain, both AlN and GaN have an indirect bandgap ${\Gamma }_c-K_v$ and the gap (increases) decreases for (compressive) tensile uniform strain, as shown in Figure 2. Interestingly, an important difference occurs between AlN and GaN.

Figure 2.

Band structrure of AlN (top) and GaN (bottom) monolayers for biaxial uniform strain of −4% (blue), 0% (red), and 4% (green). The dashed gray line indicates the VBM for the unstrained material, set as zero energy.

While AlN remains an indirect gap material for all the values of considered strain, GaN undergoes an indirect to direct gap transition at the compressive strain of about −4%. It is worth to mention that unstrained GaN has an indirect gap only at the DFT level, which switches to direct when the GW corrections are introduced. ^27,63 In the case of strain, the VBM switches from K to Γ, because the state at Γ goes up in energy under compressive strain, while the corresponding eigenvalue at K does not show a big variation. This behavior can be explained by analyzing the character of the KS wavefunction at K_v and at ${\Gamma }_v$ (see bottom panel of Figure 3). While the main contribution of K_v is due to p_z N orbitals, ${\Gamma }_v$ has a large contribution from $p_{x,y}$ (VBM and VBM-1 are degenerated) orbitals of N atoms. The analysis of the corresponding states in AlN (see Figure 3, top panel) reveals similar contributions for K_v and ${\Gamma }_v$ , but the more localized nature of the orbitals implies a smaller sensitivity of the bonds to the applied strain and, in particular, a minor energetic variation of ${\Gamma }_v$ in AlN under compressive strain, with respect to the case of GaN. We report in Figure 4, for 2D AlN and GaN, the dependence of valence band maximum (VBM), conduction band minimum (CBM), vacuum potential (left panel); IP, EA, WF (center panel), and electronic gaps (right panel). In both materials, the vacuum potential, as well as CBM and VBM, decreases (increases) when increasing tensile (compressive) strain (left panel). Looking at the central panel, we note that the IP of GaN increases (decreases) with increasing tensile (compressive) strain, showing similar behavior to what observed for the WF of X-enes. The two different slopes between $-10\%$ to $-4\%$ and between $-4\%$ to 10% are related to the change of VBM from ${\Gamma }_v$ to K_v , as discussed before.

Figure 3.

Plot of $|\psi {|}^2$ for ${\Gamma }_v$ (left) and K_v (right) of AlN (top panel) and GaN (bottom panel) monolayers at zero strain. In green, the isosurfaces at 1%. The represented atoms are Al (orange), Ga (magenta), and N (blue).

Figure 4.

AlN (solid lines) and GaN (dashed lines). Left panel: vacuum level (black), CBM (blue), and VBM (red) dependence on the strain. Center panel: IP (black), EA (blue), and WF (red) dependence on the strain. Right panel: gap dependence on the strain.

It is worth to notice that the corresponding IP of AlN shows an opposite trend, increasing (decreasing) when the compressive (tensile) strain increases. This can be explained, again, in terms of the stronger localization of the electronic charge near the nitrogen atoms that occurs in AlN, with respect to the case of GaN, due to the larger difference in the electronegativity between anion and cation. The EA is given by the position of the conduction bands bottom, and both materials show a similar behavior: EA increases (decreases) for increasing tensile (compressive) strain. The CBM is always at Γ and from the analysis of the KS wave functions, it is due to a hybrid of s states of the cation (Al/Ga) and anion (N), where the former gives a larger contribution. Also, for GaN, a small contribution from the d states appears. This analysis explains why the slope of CBM and then of EA is steeper in GaN with respect to AlN. The dependence of WF and gaps from strain is then a direct consequence of the IP and EA behavior.

TMD: MoTe₂ and MoS₂

Upon the application of strain, the electronic properties of TMD show a more complicate behavior with respect to those of X-enes and nitrides. Moreover, despite the similarity in their atomic structures, the electronic properties of MoS₂ and MoTe₂ monolayers and bilayers have similar but not completely equivalent trends. ^64,65

MoTe₂ and MoS₂ monolayer

In agreement with the existing literature, the two $Mo{X}_2$ -ML studied here exhibits, at zero strain, a direct bandgap with VBM and CBM located at the six K points corners of the BZ (see the central panels of Figure 5). The bandgap values of 1.61 eV for MoS₂ and 0.93 eV for MoTe₂ are in agreement with previous DFT calculations performed at the same level of theoretical approximation. ⁶⁶ Beyond the K points also Γ and Q (where a minimum of the CB is visible in the band structure, along the $\Gamma -K$ direction) critical points have to be mentioned, because they are energetically close to K_v and K_c and, under strain or in multilayers, Q_c may become the global CBM and ${\Gamma }_v$ the global VBM. For all the group-VI TMD, it is well-known that the edges of the conduction and valence bands are composed predominantly by d metal (M) and p chalcogen (X) orbitals. ⁶⁷ However, the types (symmetries) and contributions of these orbitals vary with the chosen k-point and with the constituent atoms. As shown in Figure 6 at K, the VB (CB) is mainly described by $d_{xy},d_{x^2-y^2}$ ( $d_{z^2}$ ) M-orbitals with a small contribution of $p_{x,y}$ X-orbitals. At Γ, the VB has mainly a M- $d_{z^2}$ character with a small X-p_z contribution.

Figure 5.

Electronic band structure of monolayer MoTe₂ (top) and MoS₂ (bottom): under compressive strain (left/blue), under tensile strain (right/green), and with zero strain (center/red).

Figure 6.

Orbital projected band structure of MoTe₂ ML: p-orbitals of Te (top panel) and d-orbitals of Mo (bottom panel).

The analysis of Q_c shows a similar composition of K_v with minor contribution from $d_{z^2}$ and $p_{x,y,z}$ orbitals. ⁶⁸ Uniform tensile strain causes a general reduction of the direct gap (at K), which can be understood, using a simple tight-binding picture, as due to a minor overlap among orbitals. Nevertheless, because of the different orbital composition, the application of strain acts in a different way at the various high-symmetry k-points. In particular, due to the dominant $d_{xy},d_{x^2-y^2}$ contribution, a downward shift of K_v is observed, while ${\Gamma }_v$ , for its $d_{z^2}$ main composition, remains essentially unshifted. ⁶⁹ This effect induces a direct-to-indirect bandgap crossover ( ${\Gamma }_v$ − K_c ) when tensile strain is applied, which occurs at $\epsilon <1(\sim 5)\%$ for MoS₂ (MoTe₂) monolayer. Indeed, the extreme delocalized nature (see, for instance, Figure 7) of the p-orbitals of Te with respect to those of sulfur makes MoTe₂ less sensitive to strain. This means that smaller amounts of tensile strain are needed to shift down K_v in MoS₂, as evident in the right panels of Figure 5.

Figure 7.

Plot of $|\psi {|}^2$ for Q _c of MoS₂ (left) and MoTe₂ (right) bilayers at zero strain green isosurface at 5%. Mo (red), S (pink), and Te (yellow).

Interestingly, the behavior of MoS₂ and MoTe₂ under uniform compressive strain is very different. In MoS₂-ML, the CBM moves from K to Q for strain above $-2\%$ , hence inducing the formation of an indirect gap $K_v-Q_c$ . On the other hand, in MoTe₂-ML, when compressive strain increases, the two CB minima, along the Λ and Σ directions decrease their energies, while the valence band at M increases in energy. These two facts lead to the first transition from a direct gap $K_c-K_v$ to an indirect gap $Q_c-K_v$ at about $-2\%$ of strain and to a second transition to another indirect gap, $Q_c-M_v$ , for strain of about $-5\%$ . From the analysis of the wavefunctions at M_v , it results that contributions come from $d_{x^2-y^2}$ , $d_{z^2}$ orbitals of the metal atoms and more dominant contributions come from $p_{x,y,z}$ orbitals of chalcogen atoms. Again, due to the more delocalized nature of the p orbitals of Te with respect to S, M_v is much more sensitive to compressive strain in MoTe₂ than in MoS₂ (left panel of Figure 8). The calculated WF, IP, EA, and electronic gap trends with strain for TMD monolayers are shown in Figure 8. VBM, CBM, and vacuum potential curves (left panel) are qualitatively similar to those of 2D nitrides, but we note that the slopes are very different. In particular, the vacuum potential decreases more rapidly in MoTe₂ with respect to the other cases. For MoS₂-ML, the IP (central panel in Figure 8) decreases when compressive strain is applied with an almost linear behavior due to the fact that the VBM remains K_v . Also, in MoTe₂, the IP decreases by applying a compressive strain, but at approximately −5%, there is a change in the slope, caused by the change of the VBM from K_v to M_v .

Figure 8.

Dependence on the strain of vacuum level (black), CBM (blue), VBM (red) (left panel); IP (black), EA (blue), WF (red) (central panel); and (c) electronic gap (right panel) for monolayer MoS₂ (solid lines) and MoTe₂ (dashed lines).

For what concerns tensile strain, we see that the IP of MoS₂ decreases by increasing tensile strain. This can be explained in the following way: from very small strain values (less than 1%), the VBM moves from K_v to ${\Gamma }_v$ and it remains almost unaffected by the strain (see the blue dots in the left panel of Figure 8); in this way, we can attribute the IP variation to the decreasing behavior of the vacuum potential. For MoTe₂ (black triangles, center panel in Figure 8), IP first increases until it reaches a maximum around +5% and then begins to decrease due to the fact that the VBM shift from K to Γ at higher tensile strains.

For both TMD-ML, EA decreases when reducing the lattice parameter up to the crossover of the CBM from K_c to Q_c . Then, EA begins to increase with a very small rate when Q_c becomes the global conduction band minimum (central panel of Figure 8).

From the analysis of the data shown in the central panel of Figure 8, and similar to the case of nitrides, it is clear that the behavior of the WF depends on both IP and EA curves and shows an intermediate behavior. We can say that the WF and EA values for the two TMD-ML have a qualitatively similar behavior: an increase with a change in the slope that occurs when the direct to indirect bandgap crossover happens for low compressive strain. The linearity of EA is due to the smoothness of the vacuum potential and CBM (left panel), while the WF values are affected by the more complex behavior of the VBM (left panel) with strain. Finally, the bandgaps (Figure 8, right panel) have a totally different trend for MoS₂-ML and MoTe₂-ML, showing a much more symmetrical behavior of MoTe₂ for tensile and compressive strain.

Beyond the monolayer

TMD bilayers and multilayers are known to be indirect bandgap semiconductors. MoS₂ (MoTe₂) shows an indirect gap ${\Gamma }_v-K_c$ of 1.1 eV ( $K_v-Q_c$ of 0.93 eV) in partial agreement with previous results ⁴⁷ (see the central panels of Figure 9).

Figure 9.

Band structure of bilayer MoTe₂ (top) and MoS₂ (bottom): under compressive strain (left/blue), under tensile strain (right/green), and with zero strain (center/red).

The p X-orbital composition at Q_c and ${\Gamma }_v$ plays an important role in the crossover from the direct to indirect bandgap from monolayer to bulk. Indeed, the close distance between X-p orbitals from neighboring layers leads to large hopping, which changes the energy of Q_c and ${\Gamma }_v$ substantially. ⁶⁷ For this reason, the bandgap is very sensitive to small changes in the interlayer distance, thus dependent on the approximated vdW functional used.

In this regard, it is worth to mention that for MoTe₂ bulk, we specifically tested how our simulated structural parameters (a = 3.54 Å, Mo–Mo vertical distance between two different layers 7.33 Å) are obtained using the selected vdW functional, reproduce the experimental data (a = 3.54 Å, Mo–Mo distance between two different layers 7.33 Å) reported by Knop and MacDonald, ⁷⁰ finding a very good agreement.

Moreover, always for MoTe₂-BL, we also tested other two vdW functionals taken from the study by Grimme ⁷¹ and Barone et al., ⁷² respectively. In both cases, the final optimized structures are very similar to the one reported in Table 1, with lattice parameter a and Te–Te intralayer vertical distance, which changes less than 0.3%. Only using the second functional, we observe a reduction of 2.6% of the interlayer vertical distance. The electronic band structure does not show relevant changes near the gap region with the largest downshift of the VBM at Γ of 0.05 eV.

In Figure 7, we show the charge density plot of the Q_c state for MoS₂ (left) and MoTe₂ (right) bilayers, where the important contribution coming from p-orbitals of chalcogen atoms is evident. Similar to what discussed above, the very delocalized nature of these orbitals for the case of Te is evident from the figure. MoS₂ and MoTe₂ bilayers show, already at zero strain, different values of VBM and CBM due to the different hopping between p of the two adjacent layers, as shown in Figure 7. When tensile strain is applied in both BL, K_c and K_v go down in energy in agreement with the monolayer behavior (right panels of Figure 9). For compressive strain, similar to the monolayers, Q_c goes down in energy and, for MoTe₂, M_v becomes the VBM when $-4\%$ strain is applied (left panel).

In other words, trends previously discussed for monolayers remain qualitative similar for bilayers (see Figure 9), but different values of strain are needed to obtain changes in VBM and CBM with respect to the MLs. Moreover, MoTe₂-BL results are more sensitive to strain with respect to MoS₂-BL. This is also highlighted in Figure10, where the IP curves of the monolayer and bilayer of MoS₂ and MoTe₂ are reported, respectively, in the top and bottom panels.

Figure 10.

IP versus strain for monolayer (red) and bilayer (blue) of MoS₂ (top panel) and MoTe₂ (bottom panel).

Role of many-body effects

Despite the DFT electronic properties of 2D materials that are affected by the well-known bandgap problem, previous works have suggested that their behavior under strain reproduces, at least qualitatively, the one obtained with more refined excited state methods, such as $GW$ . ⁶⁵ In particular, for MoTe₂, it has been shown ⁷³ that the response to strain of CBM and VBM studied within DFT plus hybrid functional HSE06 qualitatively agree. To confirm these results, we have calculated, for MoTe₂-ML, how CBM and VBM change with strain, using the perturbative $G_0{W}_0$ method. The obtained G ₀ W ₀ bandgap is 1.73 eV for the unstrained structure. This value is in good agreement with that one reported by Robert et al. ⁷⁴ and Rasmussen and Thygesen ⁷⁵ of 1.72 eV.

Because of the heaviness of GW calculations, we have limited our analysis to zero strain and to one value of compressive ( $-2\%$ ) and tensile (2%) strain. In Figure 11, we report the energy values of the CBM (red) and VBM (blue) at both DFT and $G_0{W}_0$ levels of approximation. The VBM and the CBM are at K for 0 and +2% tensile strain, but the CBM moves at Q for $-2\%$ strain. While the KS eigenvalues (dashed lines) for K_c and K_v follow a trend with strain qualitatively similar to the correct quasiparticle one (solid lines), some deviation occurs for Q_c . This is clearly due to the fact that the quasiparticle correction of the state at Q_c is quite different from that of the state at K_c , due to the different contributions of atomic orbitals composing them. Our analysis of many-body effects is clearly limited, having considered only three values of strain, but it is clear that some caution has to be taken when trends with strain of the electronic gaps and work functions are extracted at the DFT level if a change in the character of the valence and conduction states near the gap occurs.

Figure 11.

Comparison of GW (solid) and DFT (dashed) energy of CBM (red) and VBM (blue). Values for the local minimum of the conduction band at K (green) are also reported.

Conclusions

We have presented a systematic study of the effect of uniform biaxial compressive and tensile strain on the electronic properties of several 2D materials ranging from semimetallic ones, like the X-enes, to semiconductors like nitrides and TMD.

Each material shows its own behavior strictly linked to the changes in the band structures and to the orbital character of the CBM and VBM. In all monolayers, the WF always increases with the lattice parameter, while the IP and EA show material-dependent trends.

The slope of the IP versus strain is usually positive, but exceptions arise for AlN (very small negative slope) and MoS₂ (where the right and left derivatives are different). The bilayers appear to be much more sensitive to strain than the isolated layers, hence, the curves are similar only in some regions. The difference in behaviors for nitrides and TMD is explained in terms of the localization of the VBM and CBM wave functions. We have also shown that the choice of the XC potential while affecting the absolute value of the WF does not change the general behavior of the electronic properties with respect to the strain. Finally, we have demonstrated that special care has to be used when the character of the VBM or CBM changes with the strain, as the DFT and $G_0{W}_0$ trends may not be equivalent.