Quantum Hall Effect in Graphene with Interface-Induced Spin-Orbit Coupling

We consider an effective model for graphene with interface-induced spin-orbit coupling and calculate the quantum Hall effect in the low-energy limit. We perform a systematic analysis of the contribution of the different terms of the effective Hamiltonian to the quantum Hall effect (QHE). By analysing the spin-splitting of the quantum Hall states as a function of magnetic field and gate-voltage, we obtain different scaling laws that can be used to characterise the spin-orbit coupling in experiments. Furthermore, we employ a real-space quantum transport approach to calculate the quantum Hall conductivity and investigate the robustness of the QHE to disorder introduced by hydrogen impurities. For that purpose, we combine first-principles calculations and a genetic algorithm strategy to obtain a graphene-only Hamiltonian that models the impurity.


I. INTRODUCTION
The interaction between a graphene sheet and different substrates has attracted great attention in recent years due to the appearance of interesting effects in the graphene layer [2][3][4][5][6] . This is particularly relevant in context of van der Walls heterostructures 7 , a vertical stack of bidimensional materials that promise to lead to novel electronics devices and applications. One possibility in this direction is the use of graphenebased van der Walls heterostructures for spintronics.
Graphene is a non-magnetic material with very weak spinorbit coupling (SOC) [8][9][10] , due to lightness of carbon atoms. Since its discovery, there were several proposals to introduce and control spin-dependent properties in graphene by inducing spin-orbit coupling [11][12][13][14][15][16] . Recently, there has been a significant progress in engineering those properties by proximity effect, which allows introducing spin-dependent features while preserving graphene's high electronic mobility. Graphene has been proximity-coupled with a magnetic thin layer of YIG for transport 17 and spin-pumping 18 measurements. SOC was also induced by proximity effect in graphene on top of different transition metal Dichalcogenites (TMDC) [19][20][21][22] and graphene decorated with Gold 23 . More recently, spins were optically injected in graphene/TMDC systems 24,25 , also indicating the presence of SOC in the graphene layer. However none of these measurements give clear indications of the type of underlying spin-orbit coupling mechanism that generates the observed phenomena.
In this article, we exploit the possibility of using quantum Hall measurements to extract the characteristics of the spinorbit coupling in graphene. For that purpose, we employ the effective model for graphene with interface-induced spin-orbit coupling from Ref. 1. We then use Landau operators to calculate the quantum Hall effect in the low-energy limit. By performing a systematic analysis of the spin-splitting of the Landau levels (LL) in the quantum Hall regime as a function of magnetic field and gate-voltage, we obtain characteristic scaling laws for the splittings produced each type of spin-orbit coupling. The same type of analysis can be used to charac-terise and estimate the SOC in transport experiments. Furthermore, we employ a real-space quantum transport approach based on a Chebyshev polynomial expansion of disordered Green functions to calculate the quantum Hall conductivity 40 . We use this numerical approach to investigate the robustness of the QHE to disorder introduced by hydrogen impurities. The impurities are modelled by an ab initio-derived tightbinding model. For the extraction of the tight-binding parameters, it is necessary to perform a multiparametric fit. Deterministic approaches for the fit can be quite difficult, because of the occurrence of large number of extrema. Therefore, here we proposed the use of an heuristic algorithm to perform this task efficiently.
The article is organised as follows: in section II, we introduce the tight-binding and low-energy models for graphene with interface-induced SOC. We also present the Hamiltonian of the system under an external magnetic field, written in terms of Landau operators. In section III, we discuss our analytical results for the energy spectra and Hall conductivity and introduce the scaling laws that can be used to discriminate the different SOC. In section IV we present the heuristic algorithm that was used to extract the tight-binding parameters from density functional theory spectra and our numerical approach for the conductivity calculations. We then discuss the results for the effect of hydrogenation on the QHE for graphene with interface-induced SOC. Finally, we present our conclusions in section V.

II. THEORY
With a combination of density functional theory (DFT) and group theory analysis, Kochan et al. proposed a tight-binding Hamiltonian and its corresponding low-energy approximation for graphene stacked on transition metal dichalcogenides (TMDC) 1 are the creation and annihilation operators for an electron on a lattice site i and spin s belonging to sublattice A or B, respectively. The first term represents the hopping with amplitude t between π orbitals in a honeycomb lattice and the second term is an energy offset ∆ between sub-lattices A (η ai = 1) and B (η bi = −1) due to the super-lattice effect originated by the incommensurability of the two lattices. The third contribution is the typical Rashba SOC with strength λ R 26,27 . It arises because the inversion symmetry is broken when graphene is placed on top of a TMDC.The fourth contribution is a sublattice-dependent intrisic SOC with coupling intensities λ a and λ b . H ξ , the valleydependent low-energy limit of the Hamiltonian of Eq. 1, has four terms 1 , The first two terms are the spinless contributions where v F is Fermi velocity given by v F = 3ta 2 with lattice constant a, σ is a pseudospin Pauli matrix related to sublattices A and B and k x and k y are the components of the electronic moment relative to the Dirac points. ξ is related with the valley degree of freedom, ξ = + for valley K and ξ = − for valley K . The third and fourth terms are the Rashba and sublattice resolved intrinsic spin-orbit couplings respectively. The fourth term contains the well known Kane-mele term 28,29 with strength λ I = (λ a + λ b )/2, and a valley-Zeeman SOC with strength λ VZ = (λ a − λ b )/2 that couples spin and valley degrees of freedom.
In Figure 1 we show the energy spectrum of the Hamiltonian of Eq. 2 for different combinations of λ a , λ b and λ R to better understand the characteristics of the novel intrinsic spin-orbit coupling. Differently from the usual case of λ a = λ b , if only the valley-Zeeman contribution λ a = −λ b (i. e. λ VZ = λ a ) is present, the spectrum is gapless and the spin degeneracy is broken (see Figure 1(a)). If only Rashba SOC is present, the spectrum is also gapless 27 . However, any combination λ R = 0 and λ VZ = 0 opens a gap in the energy spectrum as shown in Figure 1, panels b and c. For comparison, Figure 1 (d) presents the energy spectrum for λ a = λ b (i. e. λ I = λ a ) and λ R = 0. It is important to note that valley-Zeeman combined with Rashba coupling preserves particlehole symmetry ( Figures1(b) and 1(c)), while, when Kane-Mele and Rashba are both present, the particle-hole symmetry is broken (Figure1(d)).
If a uniform perpendicular magnetic field is applied, p → π = p + e A, where A is the vector potential in the Landau gauge A = B(−y, 0, 0), and B is the intensity of the magnetic field. The low-energy Hamiltonian can be written as where Landau operators a ξ are creation and annihilation operators on valley K and K for ξ = ± : The Hamiltonians H ξ=± are block diagonal with each block indexed by an occupation number n. The two lowest blocks are 1 × 1 and 3 × 3 matrices and higher blocks are 4 × 4 matrices for both valleys. The energies E ξ n,i and eigenvectors |ψ ξ n,i are indexed by the valley ξ, the occupation number n and i, that labels the eigenvalues and eigenvectors of a given block n. For more details, see Appendix A. The transverse Hall conductivity σ xy can be calculated in the framework of the Kubo formula 30 where f (E) is the Fermi-Dirac distribution function. We will not consider thermal effects here and the Fermi-Dirac distribution has a Heaviside function profile, Expressions for the eigenvectors are computed in appendix A.

A. Weak spin-orbit coupling
Here, we show and discuss our results on the quantum Hall effect in under the effect of Rashba and valley-Zeeman SOCs, that are present in graphene-TMDC heterostructures. Experiments and first-principle calculations point out to values of the coupling constants in the range of 0.1 − 10 meV 21,22 . We express our results in terms of a gate voltage which controls the Fermi energy. This voltage is related to the electrons' Fermi momentum in graphene by where α depends on the substrate. Here, we use α = 7.2 × 10 10 V −1 cm −2 , which is an appropriate value for either silicon oxide or TMDC substrates.
To be compatible with the analysis of experimental results, instead of showing our energy spectra, we present fan diagrams of the density of states as a function of gate voltage V G and magnetic field B. We begin by presenting the well known spin-splitting generated by λ R . In this case, the n = 0 Landau level is spin-degenerate, while all other levels are spinsplit, as illustrated in Fig. 2. It is also clear, from the Landau level splittings of Fig2(a) and the quantum Hall conductivity of Fig.2(b) that, as expected, the splitting increases with the level number n.
We proceed to analyse the spin-splitting generated by λ VZ : similarly to the previous case, the n = 0 Landau is spindegenerate while all other levels are spin-split, as illustrated in Fig.3. Again, the Landau level splittings of Fig. 3(a) and the Quantum Hall conductivity of Fig. 3(b) show the increase of the splitting with the level number n and it is considerably larger than the spin-splitting produced by Rashba, even for weak couplings and small n.
To compare the spin-splitting produced by Rashba with the one by the valley-Zeeman SOC, we select the difference in energy of the n = 3 levels ∆V G as a function of the spinorbit strength for various values of the external magnetic field . From Fig.4, we see a very different dependence of ∆V G as a function of λ R (panel (a)) and λ VZ (panel (b)). ∆V G varies linearly with λ VZ and quadratically with λ R . For the case of valley-Zeeman SOC, ∆V G is 1-5 V even for very weak couplings of the order of λ VZ = 1 − 5 meV and can be resolved  ure 5, we consider the case where λ VZ > λ R while in Figure 6 we consider the case where λ VZ < λ R . Because of the different functional laws presented in Fig.4, if both couplings have similar strengths the splitting is dominated by λ VZ . However, as can be seen in these figures, it is difficult to extract information on the nature and strength of the spin-orbit coupling in a graphene heterostructure from fan diagrams and quantum Hall measurements. The spin-splitting of the Landau Levels and the quantum Hall plateaus for graphene with valley-Zeeman SOC are very similar to the ones for graphene with Rashba SOC.
To address this issue, we found two different scaling laws for the spin-split states that can be used to distinguish between systems with valley-Zeeman SOC and ones with Rashba SOC. For the case of valley-Zeeman SOC, the scaled quantity is ∆V G , the voltage difference between spin-split states with same quantum number n. If we consider ∆V G /n as a function of the magnetic field, all curves for different values of n collapse perfectly into a single one, as can be seen in Figure  7(a). If λ R = 0, the scaling starts to fail whenever the main contribution for the spin-splitting is λ R (see Figure 7(b)-(c)), indicating that this scaling law is a characteristic of valley-Zeeman SOC.
On the other hand, if we consider λ R = 0 and λ VZ = 0, we need to use a different scaling law to collapse all curves. In this case, we use the quantity 4(V 2 G,s − V 2 G,s )/n 2 where V G,s and V G,s are the voltage of the two spin-split states (represented by s and s ) for a given n. This simple analysis, that can be used in transport measurements, is able to determinate what is the main SOC in graphene and the deviations from   Let us now discuss the case where the spin-orbit coupling is of the order of tens of meV, as estimated via weak anti-weak localisation measurements on graphene-TMDC heterostructures 21 . Rashba combined with valley-Zeeman SOCs produce particle-hole symmetric Landau levels spectra in both strong and weak coupling regimes. This is a conse-   quence of particle-hole symmetric spectra in absence of mag-netic field (see Figure 1 (b) and (c) ). However, as we can seen in Figure 9, the structure of spin-splittings and plateaus is more complex in the strong coupling regime, due to the presence of several level crossings, that make it difficult to perform quantitative analysis and estimations based solely on spectra and conductivity profiles. On the other hand, Rashba imprints a clear signature in the fan diagram: for strong λ R , the fan diagram acquires two extra lateral fans and the separation between the main fan and the satellite ones is proportional to λ R . Furthermore, if both λ R and λ VZ are strong, there is a splitting of the n = 0 level that does not occur in the case of pure Rashba SOC. Instead of considering the interplay between Rashba and valley-Zeeman SOCs, we can look at the structure of the Landau fan diagram of graphene with interface-induced Rashba and intrinsic SOCs (λ R and λ I ,) as in the case of graphene intercalated with gold. An example of this intercalated structure is a recent experiment on van der Waals heterostructures of graphene-gold-hBN 23 . We show in figure 7, the spectra for (a) λ I > λ R and (b) λ I < λ R . In both cases, the particle-hole symmetry is broken, a characteristic of the interplay between Kane-Mele SOC and Rashba SOC. For λ I > λ R , it presents a topological gap at V = 0 while the gap is closed if λ I < λ R .

IV. NUMERICAL RESULTS
A. Graphene-only Hamiltonian for hydrogen on graphene.
Hydrogen adsorption is one of the main sources of contamination when manufacturing graphene. Hydrogen hybridizes with p z orbitals in graphene, modifying the local symmetry from sp 2 to sp 3 , and create midgap-states 32,33 . To analyse the effect of hydrogenation in the quantum Hall effect, we employ a real-space quantum transport approach 40 . For an optimized use of computational resources, that allow us to study large systems, we need to find a single-orbit tightbinding Hamiltonian for hydrogen on graphene that reproduces a band-structure obtained by density functional theory (DFT).
We propose the following graphene-only Hamiltonian for hydrogen adatoms incorporated in a graphene layer where c † I and c I are the creation and annihilation operators for an electron on the adsorption site I, and I, j and I, j represent the sum over nearest and next-nearest neighbours of the adsorption site respectively. t I and t (2) I are parameters that have to be adjusted to properly fit the DFT band-structure.
In Fig. 12.a we show the band-structure of a 5 × 5 graphene supercell with one hydrogen adatom. The first-principles calculations are based on DFT 34,35 . as implemented in the SIESTA code 36 , using GGA functional approximation following the PBE approach 37 . The pseudopotential were obtained through Troullier-Martins scheme 38 and a double-ζ polarised basis set was used to described the electronic orbitals. The self-consistent cycle was performed using 16 × 16 × 1 ksampling of the Brillouin zone. The structural relaxation was performed using conjugate gradient minimization until the forces were smaller than 0.01 eV/Å. The gap in the band structure artificially originates from the broken sub-lattice symmetry due to the arbitrary choice of an adsorbtion site, which in this case belongs to the sub-lattice A, and disappears when the adatoms are placed randomly, because in average the symmetry will be restored.
Multiparametric fits such as the present one are usually difficult to be performed by deterministic approaches due to the occurrence of large number of extremas. Therefore, we propose the use of an heuristic algorithm to efficiently perform this task. A simple algorithm that is inspired in the natural evolutions to find an optimal solution for a given problem is called genetic Algorithm, the scheme is outlined in Fig. 11 and follows the logic 1. Primordial Generation: A population N t individual with random genome vector x p is generated, within a hypercube of allowed genomes.

2.
Breeding: 2N p < N t pairs of individuals are choicen randomly and is combined to give raise to N o < N p ofsprings.
3. Ranking: All individuals parents are offspings are avaliated and ranked through the fittness function.

Reinsetion:
The least fitned N o individuals are terminated while keeping the best N t individuals.

5.
Mutation: An stochastic alteration of the i-th genome vector occurs with probability p M . Then the breeding phase occurs again For the primordial generation phase we consider each component of the genome vector to lie within a range of x ∈ [−x 0 , x 0 ], with x 0 =1eV. The selection of pairs is performed by following stochastic universal sampling with N p = N t /2, and it is combined using an intermediated recombination (8) with p = −0.25, producing a set of N o = N t /4 offspings, the fitness function is where Y is the band structure shown in Fig. 12(a) Figure 11. Schematics of a genetic algorithm.
In Fig. 12, we show a comparison between (a) the band structure and (b) density of states of our fitted tight-binding Hamiltonian and the first-principle results. In the inset of Fig.  12 we also compare our results with a fit that considers an independent orbital for hydrogen 39 . The fit obtained with a genetic algorithm for a graphene-only Hamiltonian is comparable with the one with independent orbitals and has the advantage (if compared with other multiparametric approaches) of being a semi-automatic procedure and using a reduced Hilbert space.

B. Effects of hydrogenation in the quantum Hall effect.
The transport coefficients were computed using a real-space O(N ) method based on the Chebyshev expansion of a variant of the Kubo formula, the Kubo-Bastin formula 40 We proceed to analise the effect of disorder in the QHE in graphene with SOC. Figure 13(a) shows de density of states for λ R = 20 meV, λ VZ = 30 meV -which are the same values of the SOC of Figure 9(a)-for B = 13T and different concentrations of hydrogen. The spin-split LL states are still visible with 0.1% of hydrogen, that is the concentration that can be expected from contamination of the samples in experiments. Still, for 1% of hydrogen, the gap at E = 0 is closed and the LL spectrum is destroyed. 13(b) presents longitudinal and transverse conductivities for 0.1% of hydrogen. Although not all quantum Hall plateaus are visible in the presence of disorder, the Laudau levels, including the spin-split states, are still visible in the longitudinal conductivity. This indicates that our analysis, based on the scaling of spin-split states from longitudinal conductivity data, should be effective even in the presence of hydrogenation that produces both intra-valley and inter-valley scattering. Inferring the type and strength of spin-orbit coupling in a graphene heterostructure is a difficult experimental task. Quantum Hall measurements provide a useful tool to understand the physics that govern the charge carriers in 2D materials. The Landau spectrum offers hints about the low energy physics of electrons (or holes) and the possible SOCs. Each Landau level in pristine graphene is eightfold degenerate due to spin, pseudospin and isospin quantum numbers. In the presence of SOC, this degenerescence is lifted in a specific way depending on the type of coupling induced in graphene. Motivated by experimental results on magnetotransport in graphene heterostructures, we studied QHE on graphene with different SOCs and provided a simple way to characterize SOCs induced by proximity effect by analyzing the spin-splitting of the Landau Levels. By using two different scaling laws, we were able to determinate what was the dominant contribution to the SOC in a graphene layer. Additionally, we used an efficient genetic algorithm strategy together with ab-initio calculations to obtain a realistic allgraphene tight-binding Hamiltonian that models hydrogenation in graphene. With this novel Hamiltonian and a quantum transport approach, we analysed the effect of hydrogenation on the QHE in graphene with interface-induced SOC. The numerical results indicate that the scaling laws can in principle be applied even with 0.1% of hydrogenation.
All results presented here are based on a free particle approximation. This approximation is appropriate in most situations but with evolution of the manufacturing of graphene, cleaner samples have been obtained and the effects of Coulomb interactions may become important. Due to the Dirac-Weyl nature of electrons on graphene, the role of interactions are different from typical materials, and still under intensive discussion 44 . The signature of Coulomb interaction for the quantum Hall physics of graphene is the appearance of extra four-fold splittings in the energy spectra, implying intermediate plateaus in the Hall conductivity. For low Landau levels, this splitting has already been observed experimentally with high magnetic fields 45,46 . The importance of Coulomb interaction increases exponentially with magnetic field 47 and the symmetry breaking is more pronounced for lower Landau levels, high magnetic fields and high mobility 48 . The quality of sample is a important factor to observe Coulomb splitting as disorder decreases the mobility of electrons. Approaches based on a free particle picture cannot capture the physics of Coulomb interactions, what limits our calculations for samples with moderate electron mobilities and moderate magnetic fields (typically smaller than 15 Tesla). The study of effect of Coulomb interaction on the Landau spectra is beyond scope of present work. However, it is important to mention that for experimental transport measurements, the analysis presented here shows that a two-fold splitting of LLs can be associated to the presence of spin-orbit coupling while a four-fold splitting is a signature of Coulomb interaction.