$\textit{Ab Initio}$ Mismatched Interface Theory of Graphene on $\alpha$-RuCl$_3$: Doping and Magnetism

TL;DR

This paper introduces a new ab initio method called mismatched interface theory (MINT) for studying heterobilayers with no periodicity, applied to graphene on $ extit{α}$-RuCl$_3$, predicting doping and magnetic shifts.

Contribution

The paper develops MINT, a novel ab initio framework for mismatched heterobilayers, and applies it to predict doping and magnetic properties in graphene/$ extit{α}$-RuCl$_3$.

Findings

Predicts 4.77% doping from graphene to $ extit{α}$-RuCl$_3$

Shows magnetic interactions shift towards the Kitaev point

Demonstrates MINT's potential for guiding material design

Abstract

Recent developments in twisted and lattice-mismatched bilayers have revealed a rich phase space of van der Waals systems and generated excitement. Among these systems are heterobilayers which can offer new opportunities to control van der Waals systems with strong in plane correlations such as spin-orbit-assisted Mott insulator $\alpha$-RuCl$_3$. Nevertheless, a theoretical $\textit{ab initio}$ framework for mismatched heterobilayers without even approximate periodicity is sorely lacking. We propose a general strategy for calculating electronic properties of such systems, mismatched interface theory (MINT), and apply it to the graphene/$\alpha$-RuCl$_{3}$ (GR/$\alpha$-RuCl$_{3}$) heterostructure. Using MINT, we predict uniform doping of 4.77$\%$ from graphene to $\alpha$-RuCl$_3$ and magnetic interactions in $\alpha$-RuCl$_3$ to shift the system toward the Kitaev point. Hence we demonstrate that MINT can guide targeted materialization of desired model systems and discuss recent experiments on GR/$\alpha$-RuCl$_{3}$ heterostructures.

Equations

$$H_{JK\Gamma}=\sum_{\langle ij\rangle\in\alpha\beta(\gamma)}\big{[}KS_{i}^{\gamma}S_{j}^{\gamma}+J\mathbf{S}_{i}\cdot\mathbf{S}_{j}+\Gamma(S_{i}^{\alpha}S_{j}^{\beta}+S_{i}^{\beta}S_{j}^{\alpha})\big{]},$$

Full text

New capabilities for synthesizing atomic scale heterostructures with lattice-mismatched van der Waals materials have opened the floodgates to an infinite array of possibilities. Among them are twisted structures of identical monolayers, such as multi-layer grapheneCao et al. (2018); Yankowitz et al. (2018) and transition metal dichalcogenidesZhou et al. (2016) (TMDs), as well as structures involving two distinct monolayers, such as TMD heterobilayersChiu et al. (2014); Hong et al. (2014) and the GR/$\alpha$-RuCl3 heterostructureZhou et al. (2018); Mashhadi et al. (2019). These capabilities offer a new control parameter to design new systems. Unfortunately, traditional ab initio techniques for calculating the electronic structure of materials are powerless when the lattice mismatch between two crystals leads to the absence of periodicity Terrones and Terrones (2014); Tritsaris et al. (2016).

For twisted graphene bilayers, moiré materials offer a superlattice and the community poured on theoretical efforts to construct effective tight-binding models Bistritzer and MacDonald (2011); Fang and Kaxiras (2016); Zou et al. (2018); Yuan and Fu (2018); Massatt et al. (2017, 2018) and develop specialized techniques for solving those models in extended, aperiodic systemsFang and Kaxiras (2016); Carr et al. (2017). Pioneering efforts have also been made to develop perturbation theory for interlayer coupling that is homogenized in the in plane directionsTritsaris et al. (2016), though actual calculations within this latter approach also have required the use of non-self-consistent tight-binding models. Without self-consistency, such tight-binding-based approaches are limited to homobilayers in that they cannot account for effects such as screening and charge transfer.

On the other hand, heterobilayers without superlattice structure are lacking a theoretical framework outside of making supercellsLazic (2015), which are computationally costly and introduce significant strain. Of our particular interest is the GR/$\alpha$-RuCl3 system. $\alpha$-RuCl3 (RuCl3) is a layered, spin-orbit-assisted Mott insulator that lies very close to forming the exotic quantum spin liquid ground statePlumb et al. (2014); Sandilands et al. (2015); Baek et al. (2017); Do et al. (2017); Zheng et al. (2017); Banerjee et al. (2017); Janša et al. (2018); Kasahara et al. (2018). Hence heterostructuring may offer a tantalizing possibility of exploring the phase diagram and doping the quantum spin liquid. However the large mismatch (see Fig.1(a)) rules out meaningful superlattice formation and the work-function difference ($\phi_{r}=6.1$ eV for $\alpha$-RuCl3Pollini (1996) and $\phi_{g}=4.6$ eV for grapheneYu et al. (2009)) suggests charge transfer.

We introduce a new framework for fully self-consistent electronic-structure studies of lattice-mismatched atomic heterostructures called mismatched interface theory (MINT). We then apply this approach to carry out full, direct density-functional theory studies of the GR/RuCl3 heterostructure and predict that electrons from the graphene layer dope RuCl3 while moving the system closer to the Kitaev pointKitaev (2006) [see Fig. 1(b)]. RuCl3 has been hypothesized to lie in the spin liquid regime at this phase pointKim et al. (2015), where spin correlation functions beyond nearest-neighbor (NN) vanish and, on a given NN bond, only the components of spins matching the bond type remain correlatedBaskaran et al. (2007).

Pure RuCl3 has been intensely studied since it was recognized to be a candidate system to materialize the honeycomb lattice Kitaev modelKitaev (2006) with extremely anisotropic spin-spin interaction among $j_{eff}=1/2$ pseudospin moments on the Ruthenium sites. In the bulk crystal the edge-sharing RuCl6 octahedra form two-dimensional RuCl3 layers with weak interlayer van der Waals coupling. Although there are signs of Kitaev physics in the bulk system, it orders into a zigzag antiferromagnet at $T_{Neel}$Sandilands et al. (2015); Baek et al. (2017); Do et al. (2017); Zheng et al. (2017); Banerjee et al. (2017); Janša et al. (2018) and evidence of Kitaev quantum spin liquid physics is only seen at temperatures above the ordering temperatureDo et al. (2017); Banerjee et al. (2017) or under a magnetic field, which suppresses orderingBaek et al. (2017); Janša et al. (2018); Zheng et al. (2017); Janša et al. (2018); Kasahara et al. (2018). Kim et al. (2015) obtained an effective model that captures competing interactions leading to the zigzag and other nearby orders from ab initio studies: 111Here we are keeping only the nearest neighbor coupling as further neighbor coupling constants were found to be suppressed by over 200%.

[TABLE]

where $i,$ $j$ designate the Ru3+ sites and the $S_{i}^{\alpha}$ are components of the $j_{eff}=1/2$ pseudospin operator $\mathbf{S}_{i}$, $\alpha\beta(\gamma)$ labels a bond on which the spin direction $\gamma$ is fixed. They further placed the model parameters relevant for RuCl3 in the zigzag-ordered phase close to the ferromagnetic- and $120^{\circ}$-ordered state in a classical phase diagram based on Luttinger-Tisza analysisLuttinger and Tisza (1946) [see the reproduced phase diagram in Fig. 1(b)]. However, little is known about how to move the system closer to the Kitaev point at the origin. Here we will use MINT to extract the tight-binding parameters for the graphene-RuCl3 heterostructure to obtain the Kitaev ($K$), Heisenberg ($J$) and symmetric off-diagonal exchange coupling ($\Gamma$) constants.

Mismatched interface theory (MINT) — Electronic structure theory offers two broad sets of approaches for treatment of either isolated or periodic systems, respectively. Standard praxis for treating periodic structures within isolated-system methods is to construct large clusters of periodic material. Conversely, to treat aperiodic structures within periodic methods, one constructs large, periodic “supercells” containing the aperiodic structure. Both methods depend on the nearsightedness of electronic matter (NEM) Kohn (1996); Prodan and Kohn (2005) to ensure convergence toward exact behavior as the size of the calculation increases to infinity. This well-established principle is reflected in the mathematics underlying the recently developed tight-binding-based methods for twisted bilayersMassatt et al. (2017); Carr et al. (2017).

We here demonstrate for the first time that a simple combination of the supercell and cluster approaches allows treatment of incommensurate interfaces directly with standard density-functional theory software without the need for specialized techniques or reduction to non-self-consistent tight-binding models. We find, moreover, that nearsightedness ensures sufficiently rapid convergence that the calculations for our system of interest are quite practical. The basic approach, illustrated in Fig. 2, begins with a large, periodic supercell of the material system of primary interest $S$ (e.g., single-layer $\alpha$-RuCl3.) Next, we place into this supercell clusters $C$ of the material from the second subsystem, with terminating groups such that no bonds are left free and adjacent cluster images do not interact (e.g., hydrogen-terminated graphene flakes).

Finally, we study convergence as the cluster size is increased. Nearsightedness then ensures that sufficiently far from the boundary of the cluster, both materials behave just as they would for a truly infinite interface. Moreover, as the cluster grows, it eventually samples all possible registries with the other material. Appropriate finite-size scaling of thermodynamically intrinsic quantities then enables extraction of the behavior of the infinite interface. As a matter of practice, following past tight-binding workCarr et al. (2017), convergence with respect to sampling over registries (and even local rotational disorder) can be accelerated in systems with large moiré patterns by employing smaller clusters and averaging over different relative positions and orientations. Below, we present the first fully self-consistent density-functional theory (DFT) calculations carried out within this approach.

The final step of MINT is to use the “MINT representation” to make predictions for the electronic structure and magnetic interactions of the heterostructure. We include not only charge transfer but also the effects of strain by extrapolating the changes in bond length from the flake calculations. The benefit of the MINT representation is that it is a system that effectively models mismatched interfaces that is, nevertheless, well suited to all standard ab initio methods with periodic boundary conditions. For example, one can carry out the ab initio total-energy calculations for different magnetically ordered states to probe magnetism. We can also calculate the electronic structure of the MINT representation to obtain effective models for analyses suited to correlated electron methods.

Computational methods — * All ab initio calculations were carried out within the total-energy plane wave density-functional pseudopotential approach, using Perdew-Burke-Ernzerhof generalized gradient approximation functionalsPerdew et al. (1996) and optimized norm-conserving Vanderbilt pseudopotentials in the SG15 familySchlipf and Gygi (2015). Plane wave basis sets with energy cutoffs of 30 hartree were used to expand the electronic wave functions. We used fully periodic boundary conditions and a single unit cell of RuCl3* with a $6\times 4\times 1$ $k$-point mesh to sample the Brillouin zone. Electronic minimizations were carried out using the analytically continued functional approach starting with a LCAO initial guess within the DFT++ formalismFreysoldt et al. (2009), as implemented in the open-source code JDFTxSundararaman et al. (2017) using direct minimization via the conjugate gradients algorithmPayne et al. (1992). All unit cells were constructed to be inversion symmetric about $z=0$ with a distance of $\approx 60$ bohr between periodic images of the RuCl3 surface, using coulomb truncation to prevent image interaction.

Application of MINT to graphene/$\alpha$-RuCl3 — Here, we consider $\alpha$-RuCl3 as the system of interest $S$ and employ hydrogen-terminated graphene clusters $C$ [planar CH10H8, C14H10, C16H10, C24H12, and C30H16 as in Fig. 3(a)]. To calculate the expected charge transfer in the macroscopic system, we first determine charge transfer for each element of the convergence sequence and then scale to the transfer expected for a full graphene layer $L$ by multiplying by $N(L)/N(C)$, the ratio of the (incommensurate) number of carbon atoms expected for a full graphene layer $N(L)$ and the number in each cluster. Figure 3(b) shows that the intrinsic quantity $\delta$ (scaled charge transfer per Ru atom in $S$) converges reliably and rapidly to a value of about 4.77% $e$/Ru (electrons per ruthenium).

To test the sensitivity of the calculation with respect to the relative displacements, we have integrated over different relative displacements for a heterogeneous incommensurate bilayer for the first time in a full density-functional theory context. As a matter of expediency, we considered the CH6 cluster and sampled the primitive $d=2$-dimensional surface cell of the $\alpha$-RuCl3 with 12 CH6 configurations: six positions in the plane covering the unit cell and two rotations of the CH6 at each position (see Supplemental Material). We found the averaged doping to be 4.41%$e$/Ru compared to our original estimate for CH6 of 4.57%$e$/Ru. This demonstrates that convergence with respect to relative displacement is extremely rapid in this system. Consequently, for all results below, we sample a single relative displacement only.

We next explore the effect of uniaxial pressure on charge transfer. Holding the $C$ cluster at fixed distance from the $\alpha-$RuCl3 system and computing pressure from the resulting force per unit area, we find that $\delta$ increases monotonically with compression and is much more responsive to positive than negative pressure (Figure 4). (For these data we used CH10H8, as it already exhibits good convergence.) We find perpendicular pressure to be an effective tool for controlling doping of GR/RuCl3.

Effects on magnetism — To understand the effect of heterostructuring on magnetism, we consider the charge transfer predicted by MINT added to a single layer of $\alpha$-RuCl3 placed in an effective medium of dielectric constant $1$ that models screening with a debye length 3 Å, as implemented within JDFTxGunceler et al. (2013); Sundararaman et al. (2017). This allows us to simulate the doping due to graphene while suppressing unwanted interaction between $\alpha$-RuCl3 planes. We first perform the ab initio total-energy calculations for the different magnetic ground states.

For pure $\alpha$-RuCl3, various ab initio studies taking spin-orbit coupling into account have found that the two lowest energy states are ferromagnetic- and zigzag-ordered states that are extremely close in energy. Indeed we find the energy difference between the ferromagnetic state and the zigzag state to be far less than our energy resolution for pure $\alpha$-RuCl3 (see Fig. 5(a)). On the other hand, the MINT representation for GR/RuCl3 displays a dramatic change in this energy hierarchy. First, the ferromagnetic state experiences a large increase in relative energy. Second, the antiferromagnetic state comes closer to the zigzag state, which remains the lowest energy state (see Fig. 5(b)). These results indicate that GR/RuCl3 should be closer to the AFM state in the phase diagram of the effective model in Eq. (1).

The full description of GR/$\alpha$-RuCl3 requires understanding of how charge transfer affects the interatomic overlaps that enter the strong coupling expansion of the Kanamori HamiltonianKim et al. (2015) that results in the magnetic Hamiltonian of Eq.1. In addition, one should consider how to describe the doped magnetic system with now more nontrivial Hamiltonian. This second step is beyond the scope of this Letter. Nevertheless we will investigate how the parameters of the magnetic interactions $J$, $K$, and $\Gamma$ are expected to change.

To accomplish this, we extract the the tight-binding parameters (intra-$t_{2g}$ and nearest-neighbor $t_{2g}-e_{g}$ orbital overlaps) from our ab initio calculations on the MINT representation using the maximally localized Wannier orbital methodMarzari and Vanderbilt (1997). We estimate the on-site Coulomb interaction $U$ following Ref.Sundararaman and Arias (2013). We then use the expressions for the coupling constants in terms of these parameters given in Ref. Rau et al. (2014). The resulting estimates of the NN exchanges give $J/K\simeq-0.3$ and $\Gamma/K\simeq 0.3$ for the doped system, corresponding to the green star on the Luttinger-Tisza phase diagram in Figure 1 (b)Kim et al. (2015). When compared to the previously obtained values of $J/K\simeq-0.7$ and $\Gamma/K\simeq 0.7$ for plain RuCl3, this clearly indicates that the charge transfer from graphene to RuCl3 has moved the system closer to the Kitaev point.222This calculation was repeated at an intermediate doping of $\delta=$2.30% e/Ru to illustrate the trend in magnetism as doping increases. At this value of $\delta$, the parameters extracted via the maximally-localized Wannier-orbital method give $\Gamma/K\simeq 0.4$ and $J/K\simeq-0.5$.

Summary and outlook – In summary, we have introduced MINT: a new framework for studying lattice-mismatched atomic heterostructures ab initio. It is a two step process of (1) constructing the MINT representation by combining the cluster and supercell methods to exploit the principle of nearsightedness of electronic matter and then (2) using this representation to study the electronic structure of mismatched interfaces with full self-consistency in the description of charge transfer across a heterogeneous incommensurate bilayer. We then applied MINT to the GR/$\alpha$-RuCl3 system that has recently been realized, finding results that quickly converge to the overall doping of $4.77\%$ electron per Ru atom. A rough estimate based on the work-function difference yields a doping level that is about half our prediction333This estimate of $2.5\%$e$/\textnormal{Ru}$ comes from integrating the linear density of states of graphene up to the workfunction difference of $\Delta\phi=1.5$ eV.. This rapid convergence provides an internal check of how well the members of the MINT sequence converge to the full GR/$\alpha$-RuCl3 system. We also predict this doping to increase readily under positive perpendicular pressure. Finally we predict the doping to bring the GR/$\alpha$-RuCl3 system much closer to the Kitaev point in the phase diagram in terms of effective exchange parameters. Interestingly the enhancement in conduction observed in two recent experiments Zhou et al. (2018); Mashhadi et al. (2019) is consistent with our prediction.

The implications of our results are twofold. First, we presented the first framework for studying mismatched interfaces in a systematic yet efficient manner. Although the MINT is a new framework, it is based on a simple principle and it uses established and widely available standard ab initio methods in each of its steps. Hence MINT is versatile and accessible and we anticipate the application of this approach to produce many more interesting results in mismatched interface systems previously out of reach of ab initio studies. Second, using MINT we found the GR/$\alpha$-RuCl3 system to accomplish two sought-after controls: bringing the RuCl3 closer to the Kitaev point and doping. To the best of our knowledge, this is the only known case of control that can make the elusive quantum spin liquid physics accessible to RuCl3 without a magnetic field. Our results lay the field wide open to future experiments on GR/$\alpha$-RuCl3 to be compared to MINT predictions. Moreover, it will be interesting to study other heterostructures involving $\alpha$-RuCl3 partnered with different van der Waals systems and continue to explore this uncharted territory.

Supplemental Material

.1 Graphene/hBN

We applied MINT to graphene/hexagonal Boron Nitride (GR/hBN, shown schematically in Figure 6 (a)) and compared the result with a large, commensurate graphene/hBN supercell in which the graphene lattice constant is initially set to match that of hBN and allowed to relax. The red dashed line in Figure 6 (b) indicates the extrapolated value of charge transferred from graphene to hBN to be 0.073% electron/B atom (which is minimal, as expected, since hBN is relatively inert). The value obtained for the commensurate supercell is 0.071% electron/B atom, thus the two values agree to within 2.77%, an impressive agreement for such a subtle effect.

.2 Displacement Averaging

In sampling the primitive $\alpha$-RuCl3 cell, we averaged charge transfer over 12 spatial configurations, shown in Figure 7. We considered six CH6 positions in the plane covering the unit cell (three along the long edge of the rectangular unit cell as in Figure 7 (a) and (c), and three down the middle of the unit cell as in Figure 7 (b) and (d)) and two rotations of the CH6 at each position. We found the averaged result to be 4.41%e/Ru compared to our original estimate from CH6 of 4.57%e/Ru. This close agreement indicates rapid convergence with respect to relative displacement.