Metal-insulator transition and the anomalous Hall effect in the layered magnetic materials VS2 and VSe2

We investigated the electronic structure of the layered transition-metal dichalcogenides VS2 and VSe2 by first-principles calculations. Both compounds exhibit metal-insulator transitions when crossing over from the bulk to the two-dimensional monolayer. In the monolayer limit, the Coulomb interaction is enhanced due to the dimension reduction, leading to the insulating state. Moreover, these monolayers are found to be ferromagnetic, supplying excellent candidates for ferromagnetic insulators. When increasing the thickness, the few-layer structure turns metallic and presents large anomalous Hall conductivity (∼100 S/cm), which oscillates with respect to the thickness due to the size effect. Our findings presents profound materials, such as ferromagnetic insulators and anomalous Hall ferromagnets, for the spintronic application.


Introduction
Semiconductor-based spintronics has attracted a lot of attention in recent years because of the spin current could transport without the presence of a net charge current [1]. For instance, the spin Seebeck effect [2,3] or the spin pumping [4] have successfully created pure spin currents through the thermal gradients across a ferromagnetic layer.
Moreover, in a magnetic insulator Y 5 Fe 3 O 12 (YIG) and Pt functioning as a spin current detector, [5,6] the spin current is transformed into an observable transverse voltage by the inverse spin Hall effect [7,8]. Nevertheless, due to the conflicting requirements in the crystal and electronic structure of semiconductors and ferromagnetic materials [9], the intrinsic semiconducting materials with strong ferromagnetism are very rare. Most widespread magentic semiconductor materials are diluted magnetic semiconductor (DMS) based on traditional semiconductors, but doped with transition metals. This kind of materials have the problem of low solubility of transition metals and of low mobility (10cm 2 V −1 s −1 ) of current carriers [10]. This is the reason why most of spintroics related works rely on YIG.
Very recently, a new 2D TMDs material of a few layer Vanadium disulfide (VS 2 ), has been successfully synthesized [11][12][13]. The synthesis procedures used are also applicable for other 2D TMD monolayers, such as the VSe 2 monolayer. The intrinsic magnetism and the potential applications as a high-performance functional nanomaterial attract particular interest [11,14,23]. Besides the bulk VX 2 [24], recent phonon dispersion calculations also reveal that monolayer VS 2 and VSe 2 with space group P6M2 are stable [25]. Moreover, the DFT calculation results also show that this type of materials is an intrinsic ferromagnetic semiconductor with a Curie temperature close to room temperature [20,21]. This property makes it an interesting material to study. Furthermore, the anomalous Hall effect (AHE) is an archetypal spin-related transport phenomenon in magentic materials, and hence has received renewed interests in recent years [26]. The ordinary Hall effect is driven by the Lorentz force in nonmagnetic conductors when an external magnetic field is applied, while the AHE depends on magnetization in ferromagnetic materials. The intrinsic mechanism can be interpreted in terms of the Berry curvature of the occupied Bloch states [27]. Recent first-principles studies based on the Berry phase formalism show that the intrinsic AHE is important in various materials, such as bulk Fe, Co, and Ni [28,29]. Also, first principles density functional calculations have been found to agree rather well with the experimental observations of AHC. Nonetheless, the physical origin of the AHE in 2D layered ferromagnets of VS 2 and VSe 2 is still not fully understood.
In this work, we therefore carry out first principles calculations not only to explain the metal-insulator transition mechanism, but also to interpret the intrinsic anomalous conductivities.

Computational methods
The electronic structure calculations of bulk and monolayer VX 2 are performed using the projector augmented wave (PAW) method with the Perdew-Burke-Ernzerhof (PBE) generalized gradient approximation (GGA) [30] as implemented in the VASP package [31,32]. The energy cutoff of 400eV is used for the plane-wave basis expansion with the total energy convergence criteria of 1×10 −5 eV per unit cell. Gamma-centered k-points grids 16×16×1 were sampled over the 2D Brillouin Zone. Optimized monolayer structures are obtained with the residual force and stress less than 0.01eV/Å and 1.0kBar, respectively. In calculations involving few layered VS 2 and VSe 2 , the van der Waals corrections (vdW-DF) [33] are adopted to optimize the lattice structural parameters and bondlengths. The van der Waals (vdW) correction barely changes the intra-layer bond lengths, while the inter-layer bondlengths are significantly reduced, confirming that the interactions between the VS 2 and VSe 2 layers is of a weak and non-local vdW type force.
The intrinsic anomalous Hall effect is caused by the spin orbital interaction. First-principles calculations must be based on a relativistic band theory. The Berry curvature and the intrinsic dc anomalous Hall conductivity (AHC) can be evaluated by using the Wannier functions by the maximally localized algorithm, as implemented in the Wannier90 [34][35][36] code. The construction of maximally localized Wannier functions is a non-self-consistent process on a uniform 16×16×1 grid ofk-points with convergent self-consistent charge potential. The chemical compounds VS 2 and VSe 2 are mainly formed by V d-and S (Se) p-orbitals. We chose ten d orbitals on atom V and six p orbitals on each atom S (Se) as the initial guess of the number of the Wannier functions. On other hand, there are 22, 44, 66, 88, 110 bands in the case of a monolayer, bilayer, threelayer, fourlayer, and five-layer, respectively. Also, in bulk systems, the unit cell contains two vanadium and four sulphide atoms, and there are 44 bands. After less than 100 iterative steps, the total Wannier spread was well converged down to 10 −11 Bohr 2 . To obtain accurate AHC, a fine Berry mesh was used. We set a convergence Berry k-mesh of 250×250×1 k-mesh in the first Brillouin zone, this mesh was used in monolayer to five-layer in Berry curvature calculation. Also, a fine 140×140×35 k-mesh in the first Brillouin zone was adapted to bulk.

The metal insulator transition mechanism
The bulk VS 2 and VSe 2 can be formed into 2H (D 3h ) and 1T (D 3d ) polymorphs [11,14,24]. A multilayer VS 2 in the 2H structure has been synthesized in recent years [14]. The bulk and muiltlayer VSe 2 can only be synthesized in the 1T structure to date [15][16][17][18]. The bulk 2H-TMDs in the unit cell with an AbA stacking sequence in each layer, such that the V ion is sandwiched by two S(Se) ions shown in figure 1(a) and also contains two hexagonal monolayers. The top view of the bulk VS 2 and VSe 2 is shown in figure 1(b). The bulk 1T-TMDs contain only one trigonal monolayer with an AbC stacking sequence. However, the 1T-VX 2 (X = S and Se) monolayer series is metallic, with higher formation energies than the 2H-VX 2 [14] therefore they are not considered in this work.
The structural parameters of bulk and monolayer VS 2 and VSe 2 are shown in table 1. Due to the lack of experimental structure parameters for the 2H-bulk VS 2 and VSe 2 , we first check the lattice parameters a and c of bulk 1T-VS 2 which are estimated to be 3.185 and 5.919 Å. This is in reasonably good agreement with the available experimental values 3.22 and 5.76 Å [19] and the theoretical values 3.17 and 5.91 Å [21]. In case of bulk 2H-VS 2 , the c value becomes almost doubled because of the AbA stacking as compared to its AbC stacking sequence in bulk 1T-VS 2 . The estimated values for a and c are 3.169 and 12.16 Å, in reasonably good agreement with the available theoretical values 3.17 and 12.13 Å [21]. Moreover, in our calculation, the lattice parameters a and c of bulk 1T-VSe 2 are estimated to be 3.330 and 6.271 Å also in agreement with the available experimental values of 3.356 Å and 6.104 Å [22]. The a and c lattice constants of bulk 2H-VSe 2 are 3.319 and 12.817 Å. The structure parameters of 2H monolayers of VS 2 and VSe 2 are smaller than the 1T VS 2 and VSe 2 monolayers. The lattice parameter a of monolayer is large than bulk except of 1T-VS 2 .
Electronic structure calculated for the 2H-VS 2 muiltlayer and the corresponding bulk materials is shown in figures 4(a) to (g). We only plot the electronic structure of the 2H-VS 2 muiltlayers because the 2H-VSe 2 are close to the 2H-VS 2 structures. The absence of interlayer interaction leads to a strong modification of the electronic properties in the VS 2 and VSe 2 monolayer systems. They are found to be indirect semiconductors with the VBM at the Γ point and the CBM nearly at the M point of the Brillouin zone. The high DOS of VS 2 and VSe 2 monolayer at the Fermi level cause the Stoner instability. Leading to the exchange splits of two pure spin polarization bands at E F and induces the intrinsic 2D ferromagnetic ordering [20].
In contrast to the semimetal of bilayer or bulk materials, the energy gap decrease is accompanied by a semiconductor metal transition from the monolayer to the bilayer. From the bilayer to the bulk, the negative energy gap increases and the semimetal evolves into metal. There is no experimental information available for the VS 2 and VSe 2 monolayers. Because of similar crystal structure as well as chemical composition of the 2H-MoS 2 monolayer are good in agreement between theory and experiment, we expect our results to be accurate also for these systems.
In figure 2, we give a physical explanation of the metal insulator phase transition which occurs when the number of layers is increased. The U value means the quantity of the Coulomb energy between atoms. The VS 2 monolayer has the largest Coulomb energy due to the lack of the interlayer atoms coupling van der Waals force. An increase in the number of layers induces a slow decrease of the U values. The most significant decrease of the U value takes place on transition from the monolayer to the bilayer VS 2 , which amounts to a 15.39% decrease. Otherwise, there is an only 0.93% decrease form the five-layer to the bulk. Moreover, the U value of the VSe 2 has the similar behavior as in the VS 2 . There is a 7.81% decrease on transition the monolayer to the bilayer, and only a 0.52% decrease on transition from the five-layer to bulk.
The d value signifies the intra band coupling between the atoms. The d value increases when increasing the layer thickness. The VS 2 monolayer has the smallest d value because it has larger lattice parameters than bulk. The largest lattice parameter induces the rare overlap between the atoms, therefore resulting in the smallest d value and in the smoothest bandstructure of the VS 2 monolayer near the Fermi level. The d value increases by about 37.12% on transition from the monolayer to the bilayer. On the other hand, there is only a 0.93% increase of the d value form the five-layer to the bulk. We also show the V-S distance, the V-S bond length, and the S-S distance on transition form the monolayer to the bulk. When increasing the number of layers all of the above lattice parameters decrease and display a tendency toward bulk properties. The bulk has the smallest lattice  Table 1. Structural parameters of bulk and monolayer VS 2 and VSe 2 : inplane (a) and out-of-plane (c) lattice parameters, vertical distance of the V and S(Se) layers (d). For the monolayer and bulk systems, both a and c have been optimized.

Monolayer
Bulk parameters, but the largest d value which indicates that here exists the largest overlap between atoms. The variation percentage of the V-S distance, the V-S bond length, the S-S distance are 0.37%, 0.21%, and 0.37% on transition from monolayer to bilayer. On the contrary, the lattice parameter in the transition from the five-layer to the bulk is one tenth smaller than on transition from the monolayer to the bilayer. The variation percentages of the V-S distance, the V-S bond length, and the S-S distance are 0.04%, 0.02%, 0.03% on the transition from the five-layer to the bulk.
In the VSe 2 case the d value has the same behavior as in the VS 2 . There is a 44.24% increase on transition from the monolayer to the bilayer, and an only a 0.52% decrease on transition from the five-layer to the bulk. The percentage variation of the V-Se distance, the V-Se bond length, and the Se-Se distance is 0.20%, 0.25%, and 1.10% on transition from the monolayer to the bilayer. The percentages variation of the V-Se distance, the V-Se bond length, and the Se-Se distance are 0.031%, 0.02%, and 0.1% on transition from the five-layer to the bulk. This is one tenth smaller than in the monolayer to the bilayer transition.
Moreover, the correlation effects are included beyond GGA shown in figure 3. The on-site Coulomb energy U=1 and 2 eV, J=0.87 eV for V 3d electrons [38] are taken into account for the correlation effects of the localized d orbital in the GGA+U [39]. The VS 2 and VSe 2 with U=1 eV cases show the same phase diagram with GGA calculation but the U value is upward of 0.3 and 0.4 eV shown in figures 3(a) and (c). The semiconductor to metal transition from the monolayer to the bilayer also occurs in VS 2 and VSe 2 . For the  U=2 eV case, the semiconductor to metal transition from the monolayer to the bilayer also occurs in VS 2 but the VSe 2 exhibits a gap of about 0.1 eV from monolayr to bulk shown in figures 3(b) and (d).

Anomalous Hall effect
The main equation for the intrinsic AHC in terms of the Berry curvature, The indices α and β denote Cartesian directions, W c is the cell volume, N k is the number of k-points used for sampling the Brillouin zone, and f nk is the Fermi-Dirac distribution function. W ab ( ) k n, is the Berry curvature in the first Brillouin zone [37], The normalized AHC values of the multilayer structure is defined as The s ab ( ) 0 is the orinigal AHC value from the first principles calculation, the l z is the height along z direction which includes the thickness of the vaccum layer, and l z 0 is the lattice constant c of bulk 2H-VS 2 and 2H-VSe 2 . In our calculations predicted that the AHC values of VS 2 and VSe 2 are −99.34 and −116.42 S/cm as shown in figure 5. In figures 6(a) and (b) we have plotted the total Berry phase from equation (1) and the corresponding bandstructure of the bulk 2H-VS 2 . The sharp two peaks derive from the high symmetry path from A to L, but we find that the dominant peak is opposite in sign, to the total AHC. On the other hand, there are very small positive peaks along the high symmetry path from H to L. They are less than 2% in magnitude, and give only a very small positive contribution to the total AHC. Other paths do not contribute the total AHC.
Haldane has pointed out that the AHC quantity may alternatively be expressed as a Fermi-surface property [40] and Ivo Souza's group [41] present an ab initio approach for computing the AHC of bulk Fe, Co, and Ni who convert the integral over the Fermi sea. Thus, we plot the four Fermi surface sheets associated with the bands 1-4 of bulk 2H-VS 2 which help us to understand the original contribution of the AHC. In figures 7(a)-(d), the Fermi surface of band 1 can be seen to be spherically symmetric and the band 2 can be regard as cylindrically symmetric around the first Brillouin zone. Both of then have a smooth curvature, and do not touch along the H to L edge. This implies that they do not affect the Berry phase results. The Fermi surface of bands 3 and 4 can be seen to reflect the six-fold crystal symmetry. Band 3 opens two gaps around the A to L and L to M, and changes the Berry phase contribution. The Fermi surface of the band 4 gives rise to a very small hole pocket in the A to L points direction. The hole pocket gives rise to a sharp negative Berry phase distribution. Surprisingly, the AHC oscillates within a certain range when increasing the number of layers as is shown in figure 5. First, there are no bands across the Fermi surface, and induce zero AHC in monolayer VS 2 and VSe 2 . Furthermore, the Fermi surface of bilayer to five-layer can be regard as the superposition of bulk Fermi surface, and are very similar to each other. In figures 7(e) and (f), we plot the two essential Fermi surfaces of VS 2 threelayer, which are the cause of the AHC oscillation from bilayer to five-layer. In our calculation, the three-layer and five-layer have large AHC values of −58.2 and −75.1 S/cm and also contain the figures 7(e) and (f) shape Fermi surfaces.
If we take the bilayer and four-layer as an example, then there is only one Fermi surface of figure 7(e), but there is no Fermi surface of the type shown in figure 7(f). The bilayer and four-layer AHC are −77.7 and −101.4 S/cm. The AHC of the bilayer is smaller than that of the three-layer, same as the AHC of a four-layer is smaller than that of a five-layer and displays oscillation. The oscillation amplitude becomes smaller as the numbers of layer increase and approach the bulk AHC.
The 1T-VS 2 and 1T-VSe 2 monolayer are metallic with higher formation energies than the 2H-VS 2 and 2H-VSe 2 [14]. We only present the AHC of 1layer and bulk 1T-VS 2 and 1T-VSe 2 . The AHC of 1T-VS 2 and VSe 2 monolayer of are 32.95 and 59.43 S/cm, respectively. The AHC of bulk 1T-VS 2 and 1T-VSe 2 are −12.60 and 264.62 S/cm. The AHC of monolayer and bulk 1T-VS 2 and 1T-VSe 2 are not significance different.

Conclusions
In summary, we have calculated the electronic structure and intrinsic anomalous Hall conductivity of layered and bulk VS 2 and VSe 2 via a GGA scheme. We find that a metal insulator transition occurs due to the competition between the intra band coupling (d value) between the atoms and the Coulomb energy (U value) between the layers. The monolayer of VS 2 and VSe 2 have the largest lattice constant smallest, and as a result of having the smallest intra band coupling. Also, they have the largest Coulomb repulsion force due to the lack of the interlayer van der Waals force. We conclude that there is an indirect energy gap in the monolayer VS 2 and VSe 2 when the U values are larger than the d values. On the contrary, there is no gap in the bilayer to five-layer VS 2 and VSe 2 where the U values are smaller than the d values.
The present calculations indicate that theoretical anomalous Hall conductivity (s xy ) obtained from the GGA calculations is about zero in the monolayer VS 2 and VSe 2 . When increasing the number of layers the anomalous Hall conductivity begins to oscillate. The most significant effect of the oscillation is that at all even number of layers of the small M hole pocket disappears. The M hole pocket is now pushed completely above the E F , and causes smaller oscillations of the anomalous Hall conductivity. Finally, the oscillation amplitude becomes smaller as the numbers of layers increases and approaches the bulk AHC values.