Giant Spin Lifetime Anisotropy in Graphene Induced by Proximity Effects

We report on fundamental aspects of spin dynamics in graphene interfaced with transition metal dichalcogenides (TMDCs). By using realistic models derived from first principles we compute the spin lifetime anisotropy, defined as the ratio of lifetimes for spins pointing out of the graphene plane to those pointing in the plane. In the presence of strong intervalley scattering the anisotropy can reach unprecedented values of tens to hundreds, while it reduces to 1/2 for weak disorder. This behavior is mediated by spin-valley locking, which is strong in TMDCs and is imprinted onto graphene. Such giant spin transport anisotropy, driven by proximity effects, provides an exciting paradigm for designing novel spin device functionalities.

Introduction. Following the discovery of graphene in 2004 [1], a host of other two-dimensional (2D) materials have been synthesized and studied, each demonstrating unique properties and showing promise for technological applications [2]. Currently, there is a great deal of interest in layered heterostructures of these materials [3,4], where the combined system might be engineered for specific applications [5] or might enable the exploration of new phenomena [6,7]. In the field of spintronics, graphene has exceptional charge transport properties but weak spin-orbit coupling (SOC) on the order of 10 µeV [8], which makes it ideal for long-distance spin transport [9][10][11] but ineffective for generating or manipulating spin currents. To advance towards spin manipulation, recent work has focused on heterostructures of graphene and magnetic insulators [12][13][14][15][16] or strong SOC materials such as transition metal dichalcogenides (TMDCs) and topological insulators [17][18][19]. The SOC induced in graphene by a TMDC could enable phenomena such as topological edge states [20] or the spin Hall effect [21][22][23].
To this end, a variety of recent experiments have explored spin transport in graphene/TMDC heterostructures [21,[24][25][26][27][28][29]. Magnetotransport measurements revealed that graphene in contact with WS 2 exhibits a large weak antilocalization (WAL) peak, revealing a strong SOC induced by proximity effects [24][25][26]30]. Fits to the magnetoconductance yielded spin lifetimes τ s ≈ 5 ps, which is two to three orders of magnitude lower than graphene on traditional substrates [10,31]. It was later asserted that after the removal of a temperatureindependent background, τ s becomes at most only a few hundred femtoseconds [26]. Nonlocal Hanle measurements, meanwhile, have revealed spin lifetimes up to a few tens of picoseconds [27][28][29] that can be tuned by a back gate [28,29]. Finally, charge transport measure-ments on a Hall bar demonstrated a large nonlocal signal that was related to the spin Hall effect [21]. Fits to experimental measurements have estimated the induced SOC in graphene to be 10-20 meV [21,26], while most density functional theory (DFT) and tight-binding (TB) calculations find values closer to 1 meV [18-20, 24, 25, 32]. While these studies have demonstrated that TMDCs induce strong SOC in graphene, the estimated values of τ s vary by three orders of magnitude and nothing is yet known about the mechanisms governing spin dynamics and relaxation in these systems.
In this Letter, we employ dissipative quantum spin dynamics arguments, and quantum mechanical numerical simulations, to elucidate the nature of spin relaxation in graphene/TMDC heterostructures. We find that spin relaxation follows the D'yakonov-Perel' (DP) mechanism, with τ s = 1 − 100 ps for realistic momentum relaxation rates and Fermi energies. Remarkably, the spin lifetime anisotropy, defined as the ratio of lifetimes for spins pointing out of the graphene plane to those pointing in the plane, can reach unprecedented values of tens to hundreds in the presence of intervalley scattering. This behavior is mediated by spin-valley locking induced in graphene by the TMDC, which ties the in-plane spin lifetime to the intervalley scattering time. In the absence of valley mixing this ratio reduces to 1/2, typical of systems dominated by Rashba SOC [33]. A giant spin lifetime anisotropy thus represents a qualitatively new regime of spin relaxation not typically seen in 2D systems, and its measurement [34,35] should be an experimental probe of systems with strong spin-valley coupling, which includes both graphene/TMDC heterostructures and TMDCs themselves. Furthermore, systems with giant spin lifetime anisotropy could serve as an exciting new platform for the manipulation of spin and the im- Dissipative spin dynamics model. To clarify the nature of spin relaxation in graphene/TMDC systems, we follow the approach in [33], which describes spin dynamics in a randomly fluctuating magnetic field. The low-energy (E F < 300 meV) Hamiltonian of graphene on a TMDC substrate is given by In Eq. (1), v F is the Fermi velocity, κ = 1(−1) for the K (K ) valley, σ i (s i ) are the sublattice (spin) Pauli matrices, k i are the wave vector components relative to K or K , and a = 0.246 nm is the graphene lattice constant. H 0 represents the graphene Dirac cone, and H ∆ is a staggered sublattice potential induced by the TMDC. H

A/B I and H
A/B P IA are the intrinsic and the pseudospin inversion asymmetry (PIA) SOC, respectively, the latter of which is permitted by broken z/-z symmetry in graphene [36]. Due to the broken sublattice symmetry, these terms can have different strengths and signs on the A and B sublattices (λ H P IA = aλ P IA σ z (k x s y − k y s x ), where λ I = (λ A I + λ B I )/2, λ V Z = (λ A I − λ B I )/2, λ P IA = (λ A P IA + λ B P IA )/2, and ∆ P IA = (λ A P IA − λ B P IA )/2. In this form, H I is the usual intrinsic SOC in graphene, which opens a topological gap 2λ I at the Dirac point [37]. H V Z is a valley Zeeman term, which locks valley to spin and polarizes the bands out of the graphene plane with opposite orientation in the K and K valleys. H P IA renormalizes the Fermi velocity, while H ∆ P IA leads to a k-linear splitting of the bands, as in traditional 2D electron gases with Rashba SOC [38]. Except for the PIA terms, this Hamiltonian is the same as that considered in previous works [24][25][26]32].
The next step is to derive the effective spin-orbit field arising from the SOC terms. This is done by rewriting Eq. (2) in the basis of the eigenstates of H 0 and projecting onto the conduction and valence bands. At Fermi energies away from the Dirac point (E F 1 meV), this gives where k is the wave vector magnitude, θ is the direction of k with respect to k x , and ω is the spin precession frequency of the effective spin-orbit field. The in-plane components of ω give a Rashba-like spin texture, where +(−) is for the conduction (valence) band. Strong PIA SOC thus leads to electron-hole asymmetry, as will be seen for graphene on WS 2 . The out-of-plane component of ω is determined by λ V Z and changes sign between valleys. The overall texture of the effective spin-orbit field is depicted in Fig. 1.
Owing to momentum scattering, each component of ω will fluctuate in time. A simple model for the correlation of the fluctuating spin-orbit field is [33] where the correlation time of fluctuation τ c,α depends on the component of ω. The in-plane components ω x/y depend only on θ, and thus τ c,x = τ c,y = τ p , the momentum relaxation time. Meanwhile, the out-of-plane component ω z depends only on the valley index, giving τ c,z = τ iv , the intervalley scattering time. Assuming that τ c,α ω α 1, applying Eqs. (3) and (4) to the equation of motion for the density matrix [33] yields the final expressions for the spin relaxation rates In Eq. (5), the out-of-plane spin relaxation follows the usual DP relation, τ −1 s,⊥ ≡ τ −1 s,z = [2(ak∆ P IA ± λ R )/ ] 2 τ p , with the Rashba SOC augmented by the PIA term. However, the in-plane relaxation includes contributions from both the intervalley and the overall momentum scattering, and is given by The nature of the spin relaxation, with τ s, determined by τ iv and τ s,⊥ by τ p , is shown schematically in Fig. 1. Ignoring the PIA term, the spin lifetime anisotropy is Equation (6) is the main result of this work, and indicates that a giant spin lifetime anisotropy, with the inplane spins relaxing much faster than the out-of-plane spins, should be a defining characteristic of systems with strong spin-valley locking. Using DFT values of λ V Z = 1.2 meV and λ R = 0.56 meV for graphene on WSe 2 [20], and assuming relatively strong intervalley scattering (τ iv ∼ 5τ p ), we obtain a spin lifetime anisotropy of ∼20. This represents a qualitatively different regime of spin relaxation than the usual case of 2D Rashba systems, where without valley Zeeman SOC the anisotropy is 1/2, with the in-plane spins relaxing more slowly than the out-of-plane spins. Equation (5) assumes strong intervalley scattering, τ iv ω z 1, such that fast fluctuation of ω z results in motional narrowing of the in-plane spin precession and an inverse dependence of τ s, on τ iv . In contrast, when τ iv → ∞, electrons experience a constant out-of-plane spin-orbit field and only the in-plane components fluctuate with time. In this limit, the procedure above yields where τ * p = τ p /(ω 2 z τ 2 p +1). Without intervalley scattering the spin lifetime anistropy thus collapses to 1/2, as found in Rashba systems [33]. Interestingly, in this regime an external perpendicular magnetic field B z can induce an imbalance in the spin population of each valley by enhancing (canceling) the spin-orbit field at K (K ). The ratio of spin lifetimes in each valley thus becomes where g is the electron g-factor and µ B is the Bohr magneton. For graphene on WSe 2 with τ p = 100 fs, the difference in τ s can reach 10% for B z ≈ 4 T. Although this difference is too modest to achieve a complete valleyspin imbalance, it should be considered when observing spin relaxation in these structures in a magnetic field. Numerical simulations. To verify the above results, we perform numerical simulations of spin relaxation. The graphene/TMDC system is described by the TB form of Eq. (1), to which we add a disorder term is the creation (annihilation) operator at site i with spin s, and V dis ( r i ) is the potential at site i. We assume the disorder consists of Gaussian-shaped electron-hole puddles [39], such that , with the strength j of each scatterer randomly chosen within [− , ], and with a uniform width ξ = √ 3a. In the dilute limit, τ p and τ iv are inversely proportional to the number of scatterers N , while controls their relative magnitude, with larger giving stronger intervalley scattering [40,41].
To calculate charge and spin transport, we employ a real-space wavepacket propagation method that allows for efficient simulation of large-scale disordered graphene systems [42][43][44]. For charge transport we use the meansquare spreading of the wavepacket X 2 (E, t) to calculate the diffusion coefficient D(E, t) = ∂ X 2 (E, t) /∂t, which in turn gives the momentum relaxation time τ p (E) = max D(E, t)/2v 2 F . We simultaneously calculate the expectation value of the spin of the wavepacket s(E, t), from which the spin lifetime is evaluated by fitting to exp(−t/τ s,α ) or exp(−t/τ s,α ) cos(ω z t), as appropriate. The density of charge scatterers is characterized as a percentage of the number of carbon atoms, n = N/N C × 100%. We consider a 500 nm × 500 nm system with 9.2 million carbon atoms, and TB parameters are taken from Table I of Ref. [20]. Figures 2(a) and (b) show s and τ p for disorder profiles corresponding to strong and weak intervalley scattering, respectively. In the former we set n = 0.1% and = 2.8 eV, and in the latter n = 1% and = 0.5 eV. The τ p for these two cases are shown in the insets, with values typical of those found experimentally [24][25][26]. The dif- ferent energy dependence of τ p , with a minimum or maximum at the Dirac point, is indicative of the contribution of intervalley scattering [45]. In Fig. 2(a), where intervalley scattering is strong, the in-plane component of s decays much more quickly than the out-of-plane component, and spin precession is suppressed. Meanwhile, in Fig. 2(b) the in-plane spin precesses about the effective spin-orbit field with frequency ω z = 2λ V Z / , and relaxes more slowly than the out-of-plane spin. This behavior is consistent with Eqs. (5)- (7). Figure 3 shows the numerical spin lifetimes in the case of strong intervalley scattering for graphene on (a) WSe 2 and (b) WS 2 . The solid lines, for n = 0.1%, indicate a giant anisotropy with τ s,⊥ = 20 − 200 ps and τ s, ≈ 1 ps. There is also a significant electron-hole asymmetry in τ s,⊥ for graphene on WS 2 , arising from the larger PIA SOC in this system; λ R = 0.56 meV and ∆ P IA = 75 µeV for WSe 2 , while λ R = 0.36 meV and ∆ P IA = 1.4 meV for WS 2 [20]. The open circles are the values of τ s estimated from Eq. (5), showing good agreement between the numerical simulations and the spin dynamics model. To fit τ s, we assumed τ iv = 5τ p ; although our calculations do not permit an exact determination of τ iv , this ratio is consistent with prior numerical results [40]. As shown by the dashed lines, increasing the disorder density to n = 1% scales τ s by a factor of 10, confirming the inverse relationship between τ s and τ p,iv .
The numerical spin lifetimes in the absence of intervalley scattering are shown in Fig. 4, where τ s, is now larger than τ s,⊥ . The agreement with the predictions of Eq. (7), shown as the open circles, is very convincing. However, we note that the agreement worsens at low energies, as the effective spin-orbit field in Eq.  show the numerical values of the spin lifetime anistropy. As predicted by the semiclassical theory, the anisotropy is giant in the case of strong intervalley scattering, and collapses toward 1/2 otherwise.
Summary and conclusions. Using realistic quantum spin dynamics modeling and numerical simulations, we have presented a unified picture of the spin relaxation in graphene on TMDCs. We predict a giant spin relaxation anisotropy, which emerges in graphene due to proximity effects but should exist in any system with strong spinvalley locking, including TMDCs themselves. In the absence of spin-valley locking or intervalley scattering the anisotropy falls to 1/2, as expected for Rashba systems. This large variation indicates a qualitatively new regime of spin relaxation in graphene and other 2D materials.
It should be noted that the theory presented here is applicable when spin relaxation is dominated by SOC, but other spin relaxation mechanisms can take over when the SOC is small. This appears to be the case for graphene on SiO 2 , where measurements yielded no anisotropy, i.e., τ s,⊥ = τ s, [34,35]. In these systems the SOC is small and spin relaxation is likely dominated by paramagnetic impurities [46,47]. Meanwhile, very recent measurements have confirmed our prediction of giant spin lifetime anisotropy in graphene/TMDC heterostructures, with an anisotropy of ∼11 (40) for graphene on MoSe 2 (WSe 2 ) at a temperature of 75 K [48]. Another recent measurement found an anisotropy of ∼10 in graphene on WS 2 at room temperature [49], suggesting that temperaturedependent effects driven by electron-electron or electronphonon scattering should have a weak impact.
These results also have important implications for the WAL analysis of magnetotransport in graphene/TMDC heterostructures. Previous analyses have concluded that the spin relaxation is dominated by Rashba SOC [25,30], which is seemingly at odds with the presence of giant spin lifetime anisotropy. By reanalyzing the magnetoconductance measurements of Ref. [30], and introducing valley Zeeman SOC into the analysis, the experimental results can be shown to be consistent with our theory [50].
On the more applied side, the giant spin lifetime anisotropy in graphene/TMDC heterostructures might be utilized for practical purposes in spin logic devices [51,52] or in relation with opto-valleytronic spin injection in graphene/TMDC spin valves [53,54]. One possible application would be the design of a linear spin polarizer, where the in-plane components of an incoming spin current would be filtered out, leaving only the net out-of-plane polarization. In this Supplemental Material, we present a brief complementary analysis of the weak antilocalization (WAL) measurements of Ref.
[? ] (Ref. [30] in the main text). This analysis is intended to show that the theory of spin lifetime anisotropy presented in this Letter can be consistent with prior measurements of WAL in graphene/TMDC heterostructures. It also underscores the need for further theoretical and experimental study of the relationship between spin transport and quantum conductivity in these systems.
We analyze the magnetoconductance (MC) measured in a graphene/WSe 2 heterostructure, shown in Fig. 3 where ∆σ is the quantum correction to the conductivity in units of e 2 /h, F (z) = ln(z) + Ψ(1/2 + 1/z), Ψ is the digamma function, τ −1 B = 4DeB/ , D is the diffusion coefficient, τ φ is the dephasing time, and B is the external magnetic field. The spin-orbit time τ asy arises from Rashba SOC, while τ sym is typically assigned to intrinsic, or Kane-Mele, SOC.
The dashed line in Fig. S1 shows a fit to the experimental data using Eq. (1). For this fit we set D = µ v F √ πn/2e, where µ is the electron mobility, n is the charge density, and v F = 10 8 cm/s is the graphene Fermi velocity. Using µ = 10 4 cm 2 /V-s and n = 2 × 10 12 cm −2 gives D = 0.08 m 2 /s. Fitting the rest of the parameters then yields τ φ = 20 ps, τ asy = 13 ps, and τ sym = 0.6 ps.
The first thing to consider from this fit is the magnitude of τ asy . Assuming this arises solely from Rashba SOC, and assuming the D'yakonov-Perel' mechanism of spin relaxation, τ −1 asy = (2λ R / ) 2 τ p , we can estimate the Rasbha SOC strength as λ R = / 4τ asy τ p . For a momentum relaxation time τ p ≈ 100 fs (cf. Fig. 4(b) of Ref. [? ]), we get λ R ≈ 0.3 meV. This rough estimate is right in the range of the values predicted by DFT [? ].
[? ] predict λ I on the order of tens of µeV, and those of Ref.
[? ] found it to be vanishingly small, but even a value of λ I = 1 meV would give τ sym ≈ 10 ns, four orders of magnitude larger than what is found in the above fit. This analysis thus shows that the small value of τ sym does not arise from intrinsic SOC.
[? ] it was argued that this term does not relax the spin, but we have shown that it does relax the in-plane spin. Indeed, it should be noted that the EY mechanism induced by intrinsic SOC is also an in-plane spin relaxation process, as shown in the Supplementary Information of Ref.
[? ]. Thus, the valley Zeeman SOC can contribute to the MC through τ sym . If we choose τ iv ≈ 10τ p ≈ 1 ps then λ V Z = / 4τ sym τ iv ≈ 0.4 meV, which is a reasonable value in line with DFT simulations [? ]. Finally, assuming that τ asy is driven by Rashba SOC and τ sym by valley Zeeman SOC, the spin relaxation anisotropy would be τ asy /τ sym = 13/0.6 = 22. Therefore, the MC data appear to support the presence of giant spin lifetime anisotropy.
[? ]. Thus, in this instance, ignoring the valley Zeeman term could lead to an overstimate of the Rasbha SOC by a factor of ∼5.
In summary, the original analysis of Ref.
[? ] concluded that spin relaxation in a graphene/WSe 2 het-erostructure was dominated by Rashba SOC. This conclusion was supported by a large WAL peak in the magnetoconductance, and by the linear scaling of τ −1 SOC with τ p . A Rashba-dominated spin relaxation implies a small spin lifetime anisotropy, which is seemingly in contradiction with our theory as well as recent measurements of large spin lifetime anisotropy in graphene/TMDC heterostructures [? ? ]. However, our revised analysis of the WAL results, shown above, indicates that large spin lifetime anisotropy actually could be present in the devices measured in Ref. [? ]. The discrepancy could lie in the fits of Eq. (1) to the magnetoconductance. These fits have a degree of parametric freedom, and thus can lead to different conclusions depending on the extracted spin-orbit times. We also propose that the valley Zeeman SOC should be considered in any WAL analysis of graphene/TMDC heterostructures, as it is responsible for the in-plane spin relaxation. More theoretical and experimental work is clearly needed to fully reconcile the spin lifetimes extracted from WAL analysis to those from Hanle measurements.