Angle dependence of the local electronic properties of the graphene/MoS2 interface determined by ab initio calculations

We present a full theoretical study of the graphene/MoS2 interface, using density functional theory (DFT) calculations and scanning tunneling microscopy (STM) simulations. In particular, we show that contrary to previous theoretical predictions, the rotation angle between the layers has no influence on the global electronic properties of the interface, providing a careful choice of lattice vectors and supercells is made, in order to avoid artificial modifications in the electronic structure. However, small modifications of the local electronic properties do appear, as revealed by the calculated STM images. This result might be exploited in nanoelectronic devices by specific local contacting.

S Supplementary material for this article is available online (Some figures may appear in colour only in the online journal)

Letter
Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
does not seem to be clear yet [14]. As such, several theoretical and experimental papers have exhibited the specific modification of graphene and MoS 2 electronic properties, according to the specific rotation angles between the layers [15][16][17][18]. From a theoretical point of view-except in the case of some specific rotation angles, for which the mismatch between the single layer lattices is compensated exactly-the optimization of the supercell using periodic boundary conditions results in an equilibrium structure characterized by an artificial strain on both the graphene and MoS 2 . The strain on MoS 2 or graphene has an important effect on its electronic properties, in particular affecting the value and type of energy band gap [19]. For example, in their optimized structure, Yandon Ma et al [20] found a small band gap opening in the Dirac cone of graphene of some meV that increases as the interlayer distance is reduced. More recently Ebnonnasir et al [21] found a dependence of the MoS 2 thickness with respect to the orientation between the constituent layer, affecting the value and the type of the semiconductor energy band gap. In this work, we consider the different rotation angles between graphene and MoS 2 , in order to elucidate the mutual influence between the two monolayers. We first discuss the structural aspects of the graphene/MoS 2 interfaces, namely the unit cells for the density functional theory (DFT) calculations considered for the different rotation angles and the corresponding interaction energies. Then we analyze the electronic band structures and density of states (DOS), and finally we present the simulated scanning tunneling microscopy (STM) image calculations, which exhibit different Moiré patterns for the different structures.
Our calculations have been performed using the localized orbital DFT-Fireball code [22,23], which includes a specific treatment of vdW interactions [24]. In addition, a scissor operator has been used for electronic level alignment corrections [25,26]. This operator is introduced in the Hamiltonian in order to correct for the misalignment between graphene and MoS 2 . This misalignment is a result of the small size of the basis set used in our calculations. By virtue of this operator, we can rigidly shift the band of one system with respect to the other and obtain the correct energy level alignment of both systems using our very efficient basis set. In particular, if we want to shift the band ( ) ε α k by the value ( ) ∆ α k , we can write the scissor operator as: where | ( )⟩ α k is the eigenorbital with energy ( ) ε α k . A detailed mathematical derivation of the scissor operator matrix elements on a numerical atomic orbital basis set can be found in the supplementary information (see stacks.iop.org/JPhysD/50/17LT02/mmedia). The STM images are calculated within a Keldysh-Green function formalism using the Hamiltonian obtained in the Fireball simulation [27][28][29]. Full details of the method are provided in the supplementary information. As is well known in the theoretical design of material interfaces, due to the periodic conditions imposed in DFT to reproduce an infinite interface, a common basis of lattice vectors for the new superstructure composed by the rotation of one material with respect to the other has to be found [30]. Therefore, considering an × n n unit cell of graphene and an × m m unit cell of MoS 2 , the indices n and m are related by the following equation: where → a 1 and → b 1 are the lattice vectors of the isolated graphene and the MoS 2 unit cells respectively, and θ is the rotation angle of one layer with respect to the other. However, the perfect matching of both structures is almost impossible to obtain, making DFT calculations impractical. As a consequence, the relative matching of the two structures in the new superstructure can be obtained, inducing a necessarily small error in the optimization process. Indeed, depending on the lattice vector used, namely na 1 or ( ) θ mb cos 1 , the structures relax differently in most cases. This induces artificial strain and corrugation that lead to changes in the electronic properties of the deformed material at the interface. As an attempt to overcome this difficulty, we have designed supercells with reasonable numbers of atoms for the DFT calculations, favoring optimization along one material lattice vector or the other. Our criterium for designing the cell is that the lateral distance between the matching atoms of the different planes should not exceed 0.5 Å. In this way, we can observe the differences in the electronic properties induced by the calculation conditions, and remove them from the possible physical effects at the interface. According to their corresponding lattice vectors, these different structures have then been optimized and the equilibrium distance between the two layers is determined, making use of the LCAO-S 2 + vdW formalism implemented in the Fireball code [24]. The corresponding interaction energy curves as a function of the average carbon-sulfur distance between the graphene and MoS 2 are represented in figures 1 (a)-(d). The geometries of the unit cell designed for each angle and for the graphene or MoS 2 lattice vectors are represented in the insets. Except for the 0 degree structure, there is no noticeable difference between the interaction energies corresponding to the optimization, either with the graphene or the MoS 2 lattice vectors. In all the structures, the interaction energy is around 22 meV − Å 2 . This value is a bit lower than the one found for the AB stacking of graphene, which is around 40 meV − Å 2 , calculated using the same formalism. This smaller interaction energy can be explained by the honeycomb structure of MoS 2 , composed of alternating sulfur and molybdenium atoms, with the last kind in a lower plane. Consequently, the molybdenium atoms are located farther from the graphene plane, which reduces the overall interaction energy with the graphene sheet. Regarding the 0 degree graphene/MoS 2 interface, which is composed of a × 3 3 MoS 2 and a × 4 4 graphene unit cell, this superstructure presents the most important strain after optim ization, either for graphene with a 3.1 % extension in the MoS 2 lattice vectors, or for MoS 2 , with a 2.8 % compression in the graphene lattice vectors. As a consequence of this important strain, the graphene plane presents a larger corrugation which is responsible for the energy difference between the two lattice vector optimizations. Since there is no significative difference in the other rotation angles, we can deduce that the interaction energy mainly depends on the strain and the graphene corrugation. Regarding the equilibrium distance, defined as an average distance for the corrugated system, in Table 1. Evolution of the C-C distance in graphene, the MoS 2 lattice parameter, the strain and the corrugation for 0, 10, 20 and 30 degrees according to the corresponding lattice vector optimization, either for graphene (Gr) or MoS 2 . A positive strain corresponds to a compression whereas a negative strain corresponds to an extension of the layer.  The calculated electronic structure of the previous configurations shows that the optimization with respect to different lattice vectors leads to rather different electronic properties in the equivalent unit cells. At 0 degrees, in the lattice vector optimization of graphene, the MoS 2 band gap is substantially reduced from 1.83 eV to 1.1 eV, and the bottom of the conduction band is located right at the Fermi level (see figure S2 in the supplementary information). We can also observe a small amount of p-doping in graphene, as the Dirac point is shifted above the Fermi level. This is related to the very important artificial strain on the MoS 2 layer, favoring an important charge transfer from graphene to MoS 2 . Since the same strain on graphene does not show any important effect on the electronic structure, we will consider the MoS 2 lattice optimization for the rest of the study, in order to avoid any artificial strain effects.
The band structure and the DOS for 0, 10, 20, 30 degrees are represented in figure 2. For each angle, the Dirac point is now located at the Fermi level, as expected for the isolated graphene layer. Notice that the Dirac cone is mapped at the Γ point for the 10 and 30 degree structures, due to the particular symmetry of the corresponding supercells. Indeed, as is wellknown [31,32], when an × n n supercell is considered, the corresponding Brillouin zone (Bz) is reduced by a factor of × n n. As a consequence, the k-points are re-mapped into the shrinked Bz by projection. This is the so-called Bz folding effect. While Γ is always located at the center of the Bz, the other symmetry points (K and M in the hexagonal structure) can be projected in different points depending on the periodicity and the crystal structure of the supercell. In the specific case of a hexagonal crystal structure, when n is a multiple of 3, the K symmetry point is projected on Γ, as happens in our case for the 10 and 30 degree configurations.
Also, the bottom of the MoS 2 conduction band is now located at 0.3 eV above the Fermi level and the MoS 2 band gap presents the same value for all the considered angles. Since the 20 and 30 degree configurations are made with the same number of MoS 2 unit cells, a MoS 2 band to band comparison does not reveal any band structure dependence on the orientation. The similarity of the DOS for each configuration also demonstrates that the orientation does not affect the global electronic properties, which result in the simple superposition of the electronic properties of each single system, as suggested by Geim et al [8]. However, in these initial calculations, there is an electronic level misalignment that should be corrected. Using the 0 degree cell and comparing the level positions with respect to the isolated layers, we have defined a unique scissor potential that will be applied to all the angles considered here. For example, the scissor-corrected DOS and band structures for 0 and 20 degrees are represented in figure 3. It is important to notice that the MoS 2 conduction band is now shifted to 0.7 eV, above the Fermi level.
We can thus conclude that the global electronic structure of the graphene/MoS 2 interfaces is unaffected by the rotation angle. However, small modifications of the local electronic structures appear as the STM calculations reveal. The simulated images of graphene/MoS 2 interfaces for 0 and 20 degrees within the scissor-correction approach are represented in figures 4(a) and (b) respectively. The applied voltage is +0.7 V in order to include the conduction band of MoS 2 .
As a first remark, even though the global electronic properties remain the same, different Moiré patterns can be obtained for the STM images for each rotation angle. It is important to notice that the bright spots are placed on the hollow sites of the graphene sheet due to the large contribution of the non-directional d-orbitals in the tip (see the DOS of the W   tip in [27]) and the C-C distance reduction due to the graphene compression. These Moiré patterns have already been observed in ARPES experiments for example, as a mini-gap opening in the global band-structure of the interface [33]. As the electronic structure of the C atoms is almost unaltered, the difference in the brightness of the spots is directly linked to the corrugation in the graphene sheet. For example, the Moiré pattern for the 0 degree unit cell exhibits a large black area around the C-atoms in the lower position-i.e. the hollow site with a S atom below-while the brightest sites correspond to the coincident points where a C atom of the graphene sheet falls over a S atom of the MoS 2 layer. In order to check the effect of the MoS 2 layer in the current, in figures 4(c) and (d) we represent the STM image for an isolated graphene plane in the configuration of the graphene/MoS 2 interface for 0 and 20 degrees, respectively. The main features of the graphene/MoS 2 STM image are caught by the image of the isolated graphene layer, but less defined spots are obtained. This result means that there is a modulation effect in the image due to the inclusion of the MoS 2 layer. The STM images of the other cells show a less pronounced contrast between the brightest and the darkest areas due to the lower corrugation.
As illustrated in figure S3, the corrugation decreases as the rotation angle is increased (up to 30 degrees). Also notice that the corrugation is mainly due to the interaction between graphene and MoS 2 , since we have checked that the isolated graphene monolayer remains flat even in the most corrugated case (see supplementary information for more details on corrugation). On the other hand, the graphene-MoS 2 interaction is reflected in the STM image change when the MoS 2 underlayer is removed (compare images in figures 4(a), (c) and (b), (d)). Finally, the images calculated in the gap of MoS 2 , i.e. at −0.1 V, (see figure S4), present the bright spots over the C atoms of graphene, exhibiting a lower contrast. In that respect, graphene acts as a grid for the MoS 2 electronic structure. We have also calculated the charge transfer between the two layers, which is very small, in agreement with the unaffected electronic structure with the rotation angle. As a consequence of these results, the electronic device seems to be unaffected by the inclusion of graphene/MoS 2 interfaces with different rotation angles. However, their local electronic property modifications can be exploited through specific connections to the area with high or low electronic density. To illustrate this feature, we have calculated the conductance variation when approaching an STM tip to a black or white area in the 0 degree graphene/MoS 2 interface. As represented in figure S5 (supplementary information), we can observe significative differences in the corresponding conductance. These results, however, present small differences from some experimental results, where the stretch of one of the 2D materials [34] or band bending [35] have been observed at such interfaces, modifying the global electronic properties of the system. Here, we propose a tentative explanation for this slight discrepancy; however, without any calculation data for the moment. This will be the goal of a full future work on lateral interfaces between 2D materials. Theoretically, the interface is composed of two infinite planes, which implies full vdW interaction between the two structures. However, experimentally, a graphene/MoS 2 interface is usually composed of MoS 2 triangles deposited on graphene, which therefore present connections at the border, whose nature is slightly different from vdW interactions. Indeed, the triangle edges present dangling bonds which are much more reactive than the π orbitals involved in the weak vertical interaction. Consequently, we think that the main difference between experiment and theory in the structural and electronic behavior of vdW heterostructures might arise from those dangling bonds forming a lateral heterostructure between MoS 2 and graphene, and not weak horizontal heterostructures. A future work will then be devoted to the study of such heterojunctions, in order to estimate their weight in the experimentally measured vdW heterostructures between 2D materials.
To summarize, we have presented a full study of the influence of the rotation angle on the electronic properties of the graphene/MoS 2 interface. As a result, the global electronic structure remains unaffected by the rotation, even though calcul ations have to be conducted with a careful choice of lattice vectors and supercells to avoid artificial effects. For example, Yandong Ma et al [20] have theoretically demonstrated the small gap opening in the graphene band-structure for a specific rotation angle at the graphene/MoS 2 interface, which might be due to the strain induced by the choice of the supercell and the corresponding lattice vectors. The rotation angle, however, has an influence on the local electronic properties through the different Moiré patterns observed in the calculated STM images. These findings can be important in the design of future nanoelectronic devices employing the local contacting of vdW heterostructures.