Interfacial Magnon-Mediated Superconductivity in Twisted Bilayer Graphene

Heterostructures of ferromagnets (FM) and conductors are currently attracting considerable attention in spintronics. The interfacial coupling between the localized spins and itinerant electrons gives rise to intriguing phenomena such as RKKY interactions, spin-transfer, and spin-pumping [1, 2]. Furthermore, the coupling between electrons and magnons can mediate an attractive electron-electron interaction [3, 4, 5]. This effect is analogous to the electron-phonon coupling in conventional Bardeen-Cooper-Schrieffer (BCS) superconductivity [6]. Superconductivity mediated by magnons has been studied experimentally in different materials [7, 8, 9]. Furthermore, superconductivity mediated by (antiferromagnetic) magnons might also appear in certain high-$T_{c}$ superconductors [10, 11].

Superconductivity induced by interfacial coupling to magnons could exist in various material combinations. Examples are normal metals coupled to ferro- and antiferromagnets [3, 5], as well as ferromagnets and antiferromagnets coupled to the surface of topological insulators [12, 4, 13, 14]. The ferromagnetic case has also been experimentally studied [15, 16], showing a superconducting state with a critical temperature $T_{c}$ significantly higher than the intrinsic superconductivity of two materials. These studies consider either surface effects or monolayer conductors.

Systems designed for interfacial magnon-mediated superconductivity require specific properties. In general, superconducting critical temperatures $T_{c}$ are exponentially sensitive to interaction strength. In the present case, this is the magnon-mediated electron-electron interaction. Therefore, the electron states should be localized at the interface. Thus, the conducting layer should be as thin as possible yet stable. Furthermore, the electron density of states (DOS) should be large at the Fermi level. In fulfilling these conditions, twisted bilayer graphene (TBG) stands out as an ideal candidate [17].

Twisted bilayer graphene is a two-dimensional material. This renders the electron-magnon-induced effects in TBG more robust than in 3D normal metals, where interactions are constrained to the surface. The relative twist angle of the two graphene layers creates a long-period moiré pattern that, in turn, gives rise to flat electronic bands at certain magic angles [18, 19, 20]. Flat bands at the magic angles greatly enhance the electron DOS. TBG is, therefore, a laboratory for studying the transition from weak- to strong coupling superconductivity by tuning the twist angle. Graphene can couple to conventional ferromagnets [21, 22, 23]. TBG is also an interesting component of van der Waals heterostructures such as spin valves [24, 25]. Moreover, it is an intrinsic superconductor with a critical temperature of about $T_{c}=1.7$ K at half filling [26]. The underlying mechanism is still under debate, and the explanations range from phonon-mediated superconductivity [27, 28, 29, 30] to non-BCS type mechanisms [31, 32, 33].

In this Letter, we consider another path to superconductivity in TBG via magnons in adjacent layers. Fig. 1 shows twisted bilayer graphene sandwiched between two identical ferromagnets with oppositely aligned magnetization. The interfacial coupling to magnons gives rise to an effective electron-electron interaction. In combination, TBG´s valley degree of freedom causes a multicomponent superconductor. We find that two of the interlayer coupling channels are suitable for Cooper pair formation. Using a BCS model, we find a $p$-wave superconducting state with a critical temperature of the same order of magnitude as that of the intrinsic mechanism.

To describe the heterostructure, we use a Hamiltonian $H=H_{\mathrm{TBG}}+H_{\mathrm{FM}}+H_{\mathrm{int}}$, where the first term describes the electrons in the twisted bilayer of graphene sheets, and the second term describes the magnons in the top and bottom layer ferromagnets. The last term describes the interfacial coupling between the ferromagnets and the graphene layers.

We consider a continuum model for the ferromagnets given by

where $\mathbf{m}$ is the magnetization, $J>0$ is the exchange coupling, and the easy-axis anisotropy is parametrized by $K_{z}>0$. A Holstein-Primakoff transformation to second order in magnon operators and a subsequent Fourier transform yields the magnon Hamiltonian

where $a_{\bm{q}}^{(\dagger)}$ and $b_{\bm{q}}^{(\dagger)}$ are the magnon annihilation (creation) operators of momentum $\bm{q}$ in the top and bottom layer, respectively. The magnon dispersion is $\omega_{\bm{q}}=2JMq^{2}+2MK_{z}$, where $M$ is the ground state magnetization.

TBG consists of two graphene monolayers with a relative twist angle $\theta$, as shown in Fig. 2. In the decoupled limit, the top layer has two Dirac cones $K_{1}$ and $K^{\prime}_{1}$. Similarly, the bottom layer has two Dirac cones $K_{2}$ and $K^{\prime}_{2}$. Because of the relative twist, the Dirac points are related by the three vectors $\bm{q}_{1}=k_{\theta}(0,-1)^{T}$, $\bm{q}_{2}=k_{\theta}\left(\sqrt{3}/2,1/2\right)^{T}$ and $\bm{q}_{3}=k_{\theta}\left(-\sqrt{3}/2,1/2\right)^{T}$. Here, $k_{\theta}=2k_{D}\sin{\theta/2}$ is the magnitude of $\bm{q}_{j}$. The coefficient $k_{D}$ is the magnitude of the Dirac cone momenta in the monolayers. For decoupled layers, the electrons at crystal momentum $\bm{k}$ near $K_{1}$ or $K_{2}$ are governed by the single layer Dirac Hamiltonian $h^{K}(\bm{k})=v\bm{\sigma}^{*}\cdot\bm{k}$, whereas the electrons near the cones $K^{\prime}_{1}$ and $K^{\prime}_{2}$ are governed by $h^{K^{\prime}}(\bm{k})=-v\bm{\sigma}\cdot\bm{k}$. Here, $\bm{\sigma}$ is the Pauli matrix vector, and $v$ is the graphene Fermi velocity.

To model the low-energy electron bands of TBG, we use the Bistritzer-MacDonald Hamiltonian [20]. In the following, we present the main features required to apply the model for magnon-induced superconductivity. We start with an effective spin-degenerate $4\times 4$ Hamiltonian describing the electrons at cone $K_{1}$ and interlayer hopping in momentum space to cone $K_{2}$. It is given by

The three other Hamiltonians describing electrons similarly at the cones $(K_{1}^{\prime},K_{2},K_{2}^{\prime})$ must also be considered (see Supplemental Material). In Eq. (3), each element is a $2\times 2$ matrix in sublattice space. The Hamiltonian acts on four two-component spinors $\left(\psi_{K_{1}}(\bm{k}),\psi_{K_{2}}(\bm{k}_{1}),\psi_{K_{2}}(\bm{k}_{2}),\psi_{K_{2}}(\bm{k}_{3})\right)^{T}$. Here, the first spinor component describes electrons near $K_{1}$, whereas the other three components describe the electrons near $K_{2}$. The crystal momentum $\bm{k}_{j}\equiv\bm{k}-\bm{q}_{j}$, where $\bm{k}$ is measured relative to the top layer‘s Dirac cone $\mathrm{K}_{1}$.

The Hamiltonian at $K_{1}$ couples to the Dirac cone $K_{2}$ through three hopping processes. The interlayer hopping is captured by the hopping matrices $T_{1}=w(\sigma_{0}+\sigma_{x})$, $T_{2}=w\left(\sigma_{0}-1/2\sigma_{x}-\sqrt{3}/2\sigma_{y}\right)$ and $T_{3}=w\left(\sigma_{0}-1/2\sigma_{x}-\sqrt{3}/2\sigma_{y}\right)$. The hopping strength $w\approx 113$ meV [34].

The Hamiltonian in Eq. (3) exhibits low-energy eigenstates

where we used the shorthand notation $h_{j}=h^{K}(-q_{j})$. Projecting the Hamiltonian on these low-energy states, we get an effective $2\times 2$ sublattice space Hamiltonian on the form

This effective Hamiltonian has the form of a single-layer Dirac Hamiltonian with a renormalized velocity

where $\alpha=w/vk_{\theta}$. Note how $\alpha^{2}=1/3$ yields a vanishing Fermi velocity and flat bands. This value corresponds to the largest ”magic angle.”

The Hamiltonian in Eq. (3) models electrons at cone $K_{1}$ in the top layer and three allowed interlayer hoppings to $K_{2}$. As noted above, the related Hamiltonians for the $K_{1}^{\prime}$ point in the top layer and $K_{2}$ and $K_{2}^{\prime}$ in the bottom layer must also be considered (see Supplemental Material). Diagonalizing all four effective $2\times 2$ Dirac Hamiltonians, we find the resulting Hamiltonian

with a linear electron dispersion $\varepsilon_{n\bm{k}}=nv^{\star}\absolutevalue{\bm{k}}$. The $c$-operators are creation- and annihilation operators for the upper $n=1$ and lower $n=-1$ bands. The eigenstates are superpositions of states at Dirac cones in both layers. Each state is ”based” at one of the four cones $I\in\{K_{1},K_{1}^{\prime},K_{2},K^{\prime}_{2}$} with a threefold contribution from the opposing layer. The relative weight of the contributions depends on the twist angle $\theta$ and interlayer hopping strength $w$. We index the eigenstates according to the cone at which it is based. Furthermore, $s$ is the spin index, and the crystal momentum $\bm{k}$ is taken with respect to cone $I$.

We model the interfacial coupling with a conventional s-d Hamiltonian

where $\bm{s}$ is the spin operator of the itinerant electrons and $J_{\mathrm{s-d}}$ is the coupling strength.

As explained above, we take into account spin fluctuations by means of a Holstein-Primakoff transformation. We only consider the interfacial coupling between the ferromagnets and their nearest graphene layer. For instance, the first component of the electron state in Eq. (4) couples to the top ferromagnetic layer, whereas the three next components couple to the bottom layer. The interfacial coupling gives rise to a layer-dependent spin splitting

where the spin splitting is given by

The layer index $l$ can be extracted directly from $I$. It takes the value $l=1$ for $I\in\{K_{1},K^{\prime}_{1}\}$ and $l=-1$ for $I\in\{K_{2},K^{\prime}_{2}\}$. The spin index $s=1$ for spin up and $s=-1$ for spin down. The spin-split electron dispersion in the moirè Brillouin zone is shown in Fig. 3. In the case of a monolayer conductor sandwiched between two oppositely aligned ferromagnets, the spin splitting vanishes exactly. This is not necessarily the case for a bilayer, as seen from Eq. (9). Each of the electron eigenstates is asymmetrically delocalized in the two layers. Hence, the net magnetic field shifts the energy of the electron states. The sign of the shift depends on spin and the layer in which the state is ”based.” At the special twist angle $\alpha^{2}=1/6$, the states are symmetrically distributed in the two layers such that the spin shift cancels.

We next consider the electron-magnon coupling

The coupling parameters $g_{I,\bm{k},\bm{q}}^{a}$ and $g_{I,\bm{k},-\bm{q}}^{b}$ depend on the index $I$, in addition to the ferromagnetic layers $a$ and $b$ to which it couples. Their full form is given in the Supplemental material.

We derive an effective electron-electron interaction between electrons in state $I$ and $I^{\prime}$ via a Schrieffer-Wolff transformation [35], obtaining

where the interaction strength $V_{\bm{k}\bm{k}^{\prime}II^{\prime}}$ is given by

The rich $I,I^{\prime}$ structure of Eq. (13) yields, in principle, a large number of possible pairing channels. However, not all of them are suitable for the formation of Cooper pairs. For instance, states based in the same layer with opposite spins are spin split and thus not suitable for spin-unpolarized Cooper-pairing either in the spin-singlet or spin-triplet channels. Spin-unpolarized Cooper-pairing is the only possibility when the electron-magnon coupling originates from collinear spin ground states in the magnetic insulator [36]. Hence, we focus on the coupling between electrons in different layers. This excludes half of the pairing channel candidates. Furthermore, requiring the Cooper pairs to have a zero net momentum with respect to the moiré Brillouin zone leaves two possible pairing channels. These are inter-layer intra-valley pairs denoted by

Within these pairing channels, the interaction can be further decomposed in terms of pairing symmetry. To that end, we approximate the magnon frequency to be constant $\omega_{\bm{k}^{\prime}+\bm{k}}=\omega_{M}$. In this approximation, the full angular dependence originates from the coupling constants $g$ in Eq. (13). The effective interaction decomposes into an $s$-, a $p$-, and a $d$-wave component. In the low-frequency limit $(\varepsilon_{\bm{k}^{\prime}}-\varepsilon_{\bm{k}})^{2}\rightarrow 0$, the $s$- and $d$-wave components are repulsive. In contrast, the $p$-wave symmetry component is attractive and enables Cooper pair formation. Hence, we expect TBG to exhibit an interlayer magnon-mediated superconducting state with $p$-wave symmetry.

We now give an estimate of the critical temperature from the conventional BCS expression

where $\omega_{M}$ is the characteristic magnon frequency. The coupling $\lambda=VN^{\prime}_{D}$ depends on the effective interaction $V$ and the DOS $N^{\prime}_{D}$ per valley per layer per spin. The electron DOS is enhanced near the magic angle due to the flat energy bands. Ref. [37] reports a total DOS $N_{D}\gtrapprox 100\;\mathrm{eV}^{-1}\mathrm{nm}^{-2}$ close to the magic angle. This suggests a DOS $N_{D}^{\prime}=12.5\;\mathrm{eV}^{-1}\mathrm{nm}^{-2}$. Although the BCS theory does not predict the critical temperature, it may be used to obtain estimates of $T_{c}$. An important feature of Eq. (15) is the non-perturbative renormalization of the magnon energy scale, in that $T_{c}$ depends exponentially on the inverse of the DOS and the interfacial s-d coupling.

The repulsive Coulomb interaction can be detrimental to the superconducting state. TBG´s Coulomb interaction is largely screened at the magic angle for long wavelengths due to a large twist-angle dependent dielectric constant $\epsilon>250$ [38]. The Coulomb coupling strength is $\mu=V_{C}(k_{\theta})N_{D}\approx 1$. It slightly exceeds the attractive magnon-mediated interaction. However, for two reasons, the superconducting state is robust in the presence of this Coulomb repulsion. First, the superconducting gap function is of interlayer $p$-wave symmetry. Hence, the Cooper pairs circumvent the significant on-site $s$-wave contribution of the Coulomb interaction. However, we will not consider the decomposition of the Coulomb interaction as the approximation is already at a crude level. Second, the Coulomb interaction is frequency independent at the scale of the magnon cut-off frequencies. To account for this, we adopt the Morel-Anderson model [39] to find an effective coupling strength

Here, we used the observed interband plasmon frequency $\omega_{p}\approx 200\;\mathrm{meV}$ as the Coulomb interaction cut off [40]. The effective interaction strength takes the value $\lambda^{*}=\lambda-\mu^{*}$.

Interfacial coupling between graphene and ferromagnets has been studied theoretically and experimentally for numerous materials [41, 42, 22, 43, 44, 45, 46]. Hence, there are several candidates. Here, we consider two specific ferromagnets.

EuO is a ferromagnetic semiconductor with Curie temperature $T_{c}=69$ K. It has an fcc unit cell with lattice constant $5.1$ Å. Hence, two magnetic Eu²⁺ ions per unit cell, each with spin $S=7/2$, are located at the interface and thus accessible for interfacial s-d coupling [47]. EuO thin films can be deposited on graphene epitaxially [41, 42]. The induced exchange splitting is found to be $\Delta=36$ meV [48]. At the wave number $k_{\theta}$, $\omega_{M}=0.15$ meV is an appropriate frequency cut-off. These parameters suggest an effective coupling strength $\lambda^{*}\approx 0.8$ and a critical temperature $T_{C}\approx 0.5$ K.

CrI₃ is a van der Waals ferromagnet down to the monolayer limit [49]. The crystal has two magnetic ions per unit cell. Each magnetic ion carries a magnetic moment $S=3\mu_{B}$ [50]. CrI₃ hosts two magnonic modes accessible for electron-magnon coupling. Their respective energies at momentum $\bm{q}=0$ are $0.3$ meV and $17$ meV [51]. In a graphene-CrI₃ heterostructure, CrI₃ is theoretically found to induce an exchange splitting of $20$ meV [46]. Considering CrI₃ as the ferromagnet, we find a coupling constant $\lambda^{*}\approx 0.4$ and critical temperature $T_{c}\approx 0.3$ K.

Van der Waals ferromagnets, such as CrI₃, are particularly interesting candidates. This is because the interfacial exchange splitting increases significantly under compression. For CrI₃, a slight decrease in the interlayer gap can enhance the exchange splitting. A moderate reduction of the interlayer distance $\Delta d=-0.5$ Å leads to an exchange splitting of 80 meV. The splitting can reach values up to 150 meV [45]. The enhanced interfacial interaction renders the higher energy magnon branch of CrI₃ accessible for electron-magnon coupling. We expect this to increase the critical temperature significantly. In this way, the van der Waals spin valve exhibits a tunable compression-controlled superconducting state that connects the weak- and strong-coupling regimes. Due to the limited validity of the weak-coupling BCS model, we do not estimate the critical temperature for the compressed heterostructure.

In conclusion, we have demonstrated that interfacial magnons can induce superconductivity in twisted bilayer graphene. The magnons yield a multicomponent superconducting state due to the valley structure of TBG. We find an attractive interlayer channel suitable for Cooper pairs with $p$-wave symmetry. Moreover, we have considered two promising candidate ferromagnets, EuO and CrI₃. Both exhibit a critical temperature of the same order of magnitude as the intrinsic superconducting mechanism. Using van der Waals magnets is particularly interesting because the interlayer interaction strength is tunable through compression. For this reason, the superconducting state is tunable both via the twist angle and external compression. Twisted bilayer graphene sandwiched between ferromagnets is, therefore, a promising platform in which to explore magnon-mediated superconductivity.

The Research Council of Norway (RCN) supported this work through its Centres of Excellence funding scheme, project number 262633, ”QuSpin”, as well as RCN project number 323766.

Supplemental material to ”Interfacial Magnon-Mediated Superconductivity in Twisted Bilayer Graphene”