Correlated Insulating States in Twisted Double Bilayer Graphene

TL;DR

This study combines experiments and theory to explore correlated insulating states in twisted double bilayer graphene, revealing how electric and magnetic fields influence these states at various band fillings.

Contribution

It provides new insights into how moiré flat bands and external fields induce correlated insulators in twisted double bilayer graphene, highlighting the roles of spin and valley polarization.

Findings

Insulators at integer moiré band fillings 1 and 3, but not 2.

Electric field enhances insulating states at fractional fillings 1/4, 1/2, 3/4.

Magnetic field dependence indicates spin and valley polarization effects.

Abstract

We present a combined experimental and theoretical study of twisted double bilayer graphene with twist angles between 1{\deg} and 1.35{\deg}. Consistent with moir\'e band structure calculations, we observe insulators at integer moir\'e band fillings one and three, but not two. An applied transverse electric field separates the first moir\'e conduction band from neighbouring bands, and favors the appearance of correlated insulators at 1/4, 1/2, and 3/4 band filling. Insulating states at 1/4 and 3/4 band filling emerge only in a parallel magnetic field (B_{||}), whereas the resistance at half band filling is weakly dependent on B_{||}. These findings suggest that correlated insulators are favored when a moir\'e flat band is spectrally isolated, with spin polarization at 1/4 and 3/4 band filling and valley polarization at 1/2 band filling.

Figures

Full text

When two graphene sheets are stacked and twisted to an angle near $1.1\degree{}$, hybridization of the lowest lying energy bands gives rise to moiré superlattice bands with a very flat dispersion Bistritzer and MacDonald (2011); Suárez Morell et al. (2010), which greatly enhances the local density of states and can induce strong electron-electron interactions. Recent advances in fabrication techniques Kim et al. (2016) have enabled van der Waals heterostructures in which the relative layer orientation is controlled to sub-degree precision. By using the relative twist between layers as a new design parameter, many interesting phenomena have been observed, including moiré bands and Hofstadter’s butterfly spectra in twisted bilayer graphene Cao et al. (2016); Kim et al. (2017) and, most notably, correlated insulators Cao et al. (2018a) and superconductivity Cao et al. (2018b); Yankowitz et al. (2019); Lu et al. (2019) in ’magic angle’ twisted bilayer graphene (MATBG).

Bernal stacked bilayer graphene has parabolic bands touching at low energies, four-fold spin and valley degeneracy, and a band structure that is strongly modified by a transverse electric field McCann and Fal’ko (2006); Zhang et al. (2009). Twisted double bilayer graphene (TDBG) is an attractive platform to probe electron-electron interactions in flat bands because its band widths and band gaps can be controlled by electrostatic gating. We present an electrical transport study of TDBG at angles between 1$\degree$ and 1.35$\degree$. We observe resistivity maxima corresponding to the single particle band gaps of the moiré band structure, as well as correlated insulators at half filling of the first conduction band at finite transverse electric fields. Measurements as a function of in-plane magnetic fields show insulators developing at $1/4$ and $3/4$ moiré band (MB) filling factor, suggesting spin polarization at MB quarter filling. By comparing the dependence of experimental data on twist angle and transverse electric field with band structure calculations, we conclude that correlated insulators are most likely to appear when a moiré band is spectrally separated from neighboring bands, with the moiré band flatness playing a secondary role.

Figure 1(a) shows a schematic of the moiré pattern formed by two Bernal stacked bilayer graphene sheets stacked with a relative twist ($\theta$). The moiré pattern retains the hexagonal structure of the underlying layers, and is characterized by a wave-length $\lambda=(a/2)/\sin(\theta/2)$, where $a=$2.46\text{,}\mathrm{\SIUnitSymbolAngstrom}$$ is the graphene lattice constant. In reciprocal space, the Brillouin zone of the superlattice forms across the displaced $K$-points of the two graphene bilayers [Fig 1(b)].

Figure 1(c) shows the moiré band structure calculated for TDBG with $\theta=1.01\degree$, assuming that the middle two graphene layers are coupled by the same sub-lattice and position dependent hybridization as in twisted bilayer graphene Bistritzer and MacDonald (2011, 2010). The $\mathbf{k}$-dependent Hamiltonian is then constructed in a plane-wave representation. To be consistent with atomistic calculations that account for strain and out-of-plane relaxation effectsvan Wijk et al. (2015); Jain et al. (2016), we set $|T_{AA}|/|T_{AB}|=0.8$, where $T_{AA}$ and $T_{AB}$ are respectively the interlayer hopping amplitudes between A and A and between A and B sublattices. The outer two graphene layers are described by the minimal model of Bernal stacked bilayer graphene, considering only the hopping ($\gamma_{1}$) between dimer sites. We use $\gamma_{1}=0.33$ eV, consistent with infrared spectroscopy studies Zhang et al. (2008); Li et al. (2009). For twists in the 1-1.3° range, the lowest lying conduction and valence bands become relatively flat, which promotes electron-electron interactions and can lead to correlated transport phenomena Cao et al. (2018a, b); Yankowitz et al. (2019); Sharpe et al. (2019); Lu et al. (2019). Figure 1(c) exhibits gaps at neutrality between the lowest conduction and valence bands ($\Delta_{\mathrm{0}}$), and also gaps within the conduction ($c$) and valence ($v$) band structures between the first and second ($\Delta_{\mathrm{12}}^{c,v}$), and third and fourth MBs ($\Delta_{\mathrm{34}}^{c,v}$). The mechanism responsible for band flatness in TDBG is quite distinct from that in twisted single-layer systems, as explained more fully in the supplemental material.

Our TDBG samples are encapsulated in 20-50 nm thick boron nitride with graphite top and bottom gates [Fig 1(d)]. All layers are mechanically exfoliated, and the TDBG is realized using techniques similar to those described in Kim et al. (2016), with both graphene bilayers originating from the same single crystal and subsequently twisted to a precise angle during the transfer process. The TDBG samples are shaped into Hall bars using reactive ion etching, and independent edge contacts Wang et al. (2013) are made to the active area and gates. We discuss three TDBG samples, with twist angles $\theta=$ $1.01\degree{}$, $1.10\degree$, and $1.33\degree{}$. The $1.33\degree{}$ sample has a uniform twist angle to within $\mkern 1.5mu{<}\mkern 1.5mu{}0.01\degree{}$ over a channel length of $5.5\text{\,}\mu{}\mathrm{m}$, while the $1.10\degree$ sample shows a variation in the twist angle of $\mkern 1.5mu{\sim}\mkern 1.5mu{}0.015\degree{}/\mu{}\mathrm{m}$ along the channel. The $1.01\degree{}$ sample did not have a sufficient number of contacts to quantify the angle uniformity.

The longitudinal ($\rho_{\mathrm{xx}}$) and Hall ($\rho_{\mathrm{xy}}$) resistivities of the TDBG samples are probed using low frequency lock-in techniques. Our dual gated structure allows top ($V_{\mathrm{TG}}$) and bottom ($V_{\mathrm{BG}}$) gate voltages to independently tune the carrier density, $n=(V_{\mathrm{TG}}C_{\mathrm{TG}}+V_{\mathrm{BG}}C_{\mathrm{BG}})/e$, and the transverse electric field, $E=(V_{\mathrm{TG}}C_{\mathrm{TG}}-V_{\mathrm{BG}}C_{\mathrm{BG}})/2\varepsilon_{\mathrm{0}}$, where $C_{\mathrm{TG}}$ and $C_{\mathrm{BG}}$ are the top and bottom gate capacitances, $e$ is the electron charge, and $\varepsilon_{\mathrm{0}}$ is the vacuum permittivity. Figure 2(a) shows $\rho_{\mathrm{xx}}$ vs $V_{\mathrm{TG}}$ and $V_{\mathrm{BG}}$ in the 1.01$\degree$ sample at a temperature $T=$1.5\text{,}\mathrm{K}$$. The data show resistivity maxima at the charge neutrality point (CNP), and at a fixed density $n_{\mathrm{s}}$ and 3$n_{\mathrm{s}}$. Based on Fig. 1(c), which show gaps at neutrality as well as between the the first and second, and third and fourth MBs, we associate $n_{\mathrm{s}}$ with filling of the first MB. Using the four-fold spin and valley degeneracy for each MB state, the wave-length, and therefore the twist angle, can be determined using $\lambda^{2}=4(2/\sqrt{3}n_{\mathrm{s}})$.

Interestingly, $\rho_{\mathrm{xx}}$ at the CNP increases as a function of $E$-field, while the $\rho_{\mathrm{xx}}$ values at $\pm 1n_{\mathrm{s}}$ are large at $E=0$, but decrease with increasing $E$. The resistivity at $\pm 3n_{\mathrm{s}}$ is relatively constant within the accessible gate voltage range. In addition, we observe resistivity peaks at fractional MB fillings, most notably at half filling of the first conduction band within two small $E$-field windows. The appearance of a gap in a half filled band is a hallmark of strong correlation physics, as demonstrated explicitly in ultra-cold atoms trapped in optical lattices Greiner et al. (2002); Stöferle et al. (2004); Jördens et al. (2008).

To help elucidate the TDBG transport characteristics, we examine the band structure evolution as a function of $E$-field by introducing an on-site energy difference between adjacent layers ($\Delta_{\mathrm{V}}$), which should be viewed as the external field potential difference $\Delta_{\mathrm{V}}=eEd$ corrected for screening; $d$ is the graphene layer separation. The screening effect due to induced carriers is non-negligible (see supplemental material), and results in a significant reduction in $\Delta_{\mathrm{V}}$ compared to $eEd$. Figure 2(b) shows how the gaps introduced in Fig. 1(c) vary with $\Delta_{\mathrm{V}}$ when $\theta=1.01\degree$. $\Delta_{\mathrm{0}}$ behaves like the gap of a single graphene bilayer and increases with $E$ Zhang et al. (2009). On the other hand, $\Delta_{\mathrm{12}}^{c,v}$ decreases and vanishes at a finite $E$, while $\Delta_{\mathrm{34}}^{c,v}$ remains relatively constant. For both $\Delta_{\mathrm{12}}^{c,v}$ and $\Delta_{\mathrm{34}}^{c,v}$, the gaps between valence sub-bands are smaller and tend to zero at lower $E$-fields than for the corresponding conduction sub-bands. A comparison between Fig. 2(a) and 2(b) data shows good qualitative agreement between experiment and calculations, and validates the assignment of $n_{\mathrm{s}}$ to full MB filling, and therefore also the angle extraction from the transport characteristics. Interestingly, the $\rho_{\mathrm{xx}}$ maxima at $n_{\mathrm{s}}/2$ is observed at $E$-field values that yield insulators at both charge neutrality and $+1n_{\mathrm{s}}$, i.e. when $\Delta_{\mathrm{0}}\simeq\Delta_{\mathrm{12}}^{c}$. Under this condition, the first conduction band is maximally separated from both of its neighboring bands, suggesting that isolated bands favor the emergence of correlated insulators, and localized mainly in the outside layer that has the highest on-site energy (Fig. S1). This observation is also consistent with the absence of features at half filling in the valence band of Fig. 2(a) data, since $\Delta_{\mathrm{0}}$ and $\Delta_{\mathrm{12}}^{v}$ are never large enough at the same $E$-field to sufficiently separate the first valence band from surrounding bands.

In Figures 2(c-d) we examine the temperature dependence of the $\theta=1.01\degree{}$ TDBG. Figure 2(c) shows $\rho_{\mathrm{xx}}$ vs $n$ along a horizontal line-cut of Fig. 2(a) data, at $V_{\mathrm{TG}}=$-0.8\text{,}\mathrm{V}, from $1.5\text{\,}\mathrm{K}$ to $200\text{\,}\mathrm{K}$. The upper inset shows a zoomed view of the line-cut at $T=$1.5\text{\,}\mathrm{K} in which several developing $\rho_{\mathrm{xx}}$ maxima are observed at half and quarter MB fillings. For densities between the $\pm 1n_{\mathrm{s}}$ insulators, the temperature dependence is metallic at lower $T$’s, with $\rho_{\mathrm{xx}}$ increasing with $T$ up to $\mkern 1.5mu{\sim}\mkern 1.5mu{}$50\text{,}\mathrm{K}$$, followed by a decrease with increasing $T$. Figure 2(d) shows the $T$-dependence at select densities, indicated by the arrows in Fig 2(c). In the metallic regions the resistivity has a nearly linear dependence on temperature, similar to experimental observations in MATBG Cao et al. (2019a); Polshyn et al. (2019), which has been theoretically attributed either to acoustic phonon scattering in flat bands Wu et al. (2019); Li et al. (2019) or to strong correlation effects Cao et al. (2019a). The decrease in $\rho_{\mathrm{xx}}$ at higher temperatures is attributed to thermal activation to higher bands which are more dispersive and therefore have higher electron velocities Polshyn et al. (2019). We note that exceptionally low $\rho_{\mathrm{xx}}$ values can be observed at certain MB fillings [blue arrow and trace in Fig. 2(c) and 2(d), respectively], which might signal an emerging superconducting state. Figure 2(d) inset shows that the slope of $\rho_{\mathrm{xx}}$ vs. $T$ in the linear regions is relatively agnostic to $n$, apart from slight dips indicating emergent insulators at fractional band fillings.

Figure 3(a) shows $\rho_{\mathrm{xx}}$ vs $n$ and $E$ measured in a TDBG sample with $\theta=1.33\degree$, at $T=1.5$ K. Similar to the Fig. 2(a) data, single particle gaps appear at charge neutrality and $\pm 1n_{\mathrm{s}}$, and $n_{\mathrm{s}}/2$ correlated insulators arise in a range of positive and negative $E$-fields. The $\rho_{xx}$ vs. $T$ data of the $\theta=1.33\degree$ sample is similar to Fig. 2(d) data, showing a linear dependence below $50\text{\,}\mathrm{K}$, with density dependent $d\rho_{xx}/dT$ values of $\mkern 1.5mu{\approx}\mkern 1.5mu{}$100\text{,}\Omega\mathrm{/}\mathrm{K} (Fig. S2). Figure 3(b) shows $\rho_{\mathrm{xx}}$ and $\rho_{\mathrm{xy}}$ vs. $n$ measured in a perpendicular magnetic field $B=$1\text{\,}\mathrm{T} at $E\approx$0.4\text{,}\mathrm{V}\mathrm{/}\mathrm{n}\mathrm{m}$$ [dashed line in Fig 3(a)]. The $\rho_{\mathrm{xy}}$ vs. $n$ data changes sign at the single particle and correlated insulators, consistent with a transition between hole-like and electron-like bands when the Fermi level crosses an energy gap, and signaling that the four-fold band degeneracy is lifted in the first conduction MB. As in the Fig. 2 data, the insulators at $n_{\mathrm{s}}/2$ are present at $E$-field values concomitant with insulating states at both CNP and $+n_{\mathrm{s}}$. However the correlated insulators at $n_{\mathrm{s}}/2$ are now more prominent, with $\rho_{\mathrm{xx}}$ values similar to those of the single particle insulators.

Band structure calculations offer an explanation for the differences between the $\theta=1.01\degree$ and $\theta=1.33\degree$ samples (Fig. S3). As a function of $E$-field, the $\Delta_{\mathrm{0}}=\Delta_{\mathrm{12}}^{c}$ condition provides the maximum separation of the first MB from the neighbouring bands. This separation varies with $\theta$, and reaches a maximum at $\theta\approx 1.3-1.4\degree$. We note that the width of the first conduction band at $\Delta_{\mathrm{0}}=\Delta_{\mathrm{12}}^{c}$ also increases with $\theta$. Together the experimental data and calculations suggest that correlated insulators at fractional band fillings are most likely to emerge if the band is maximally separated from neighbouring bands, and that this criterion is more important than extreme band flatness. If only band flatness is considered, correlated insulators should be more prominent at $\theta\simeq 1\degree$ where the lowest bands are narrower. We emphasize again that the mechanism for band flatness in TDBG is different than in the twisted single layer case, and that the band width in TDBG is a less sensitive function of twist angle.

Further insight into the nature of the correlated insulators is provided by measurements in an in-plane magnetic field ($B_{\mathrm{\parallel}}$), which couples only to the electron spin. Figure 3(c) shows $\rho_{\mathrm{xx}}$ vs $E$ and $n$ at $B_{\mathrm{\parallel}}=$14\text{,}\mathrm{T}, and $T=1.5$ K. Correlated insulators emerge at $n_{\mathrm{s}}/4$ and $3n_{\mathrm{s}}/4$ in narrower $E$-field windows of $\mkern 1.5mu{\approx}\mkern 1.5mu{}$0.05\text{\,}\mathrm{V}\mathrm{/}\mathrm{n}\mathrm{m}. In addition, the insulator at $n_{\mathrm{s}}/2$ extends over a slightly larger range of $E$-fields. The presence of gaps at each MB quarter filling suggests a full lifting of the four-fold spin and valley degeneracy, in which each electron added per moiré unit cell is polarized in both spin and valley. Figure 3(d) shows $\rho_{\mathrm{xx}}$ normalized to its zero field value ($\rho_{\mathrm{xx,0T}}$) vs. $B_{\mathrm{\parallel}}$ at different MB fillings. At $n_{\mathrm{s}}/4$ and $3n_{\mathrm{s}}/4$, $\rho_{\mathrm{xx}}$ grows rapidly with $B_{\parallel}$, whereas the growth is more gradual at $n_{\mathrm{s}}/2$. Figure 3(c-d) data suggest the insulators at $n_{\mathrm{s}}/4$ and $3n_{\mathrm{s}}/4$ are spin polarized, with gaps enhanced by the applied $B_{\mathrm{\parallel}}$, and correspondingly that the insulator at $n_{\mathrm{s}}/2$ is valley polarized.

We also consider the effect of a magnetic field perpendicular to the sample plane. Electrons in a periodic potential and perpendicular magnetic field develop a self-similar energy spectrum known as Hofstadter’s butterfly Hofstadter (1976), which has been studied extensively in graphene/boron nitride moiré patterns Ponomarenko et al. (2013); Dean et al. (2013); Hunt et al. (2013); Yu et al. (2013), and in twisted bilayer graphene Cao et al. (2016); Kim et al. (2017); Cao et al. (2018a, b); Yankowitz et al. (2019); Lu et al. (2019). In such systems, quantum Hall states (QHSs), indexed by a Landau level filling factor ($\nu$) and a sub-band filling factor ($s$), emerge when $n=\nu\cdot(eB/h)+s\cdot n_{\mathrm{s}}$; $e$ is the electron charge and $h$ is Planck’s constant. Figure 4(a) shows a contour plot of $\rho_{\mathrm{xx}}$ vs. $n$ and $B$ measured in a TDBG sample with $\theta=1.10\degree{}$, at $T=$0.3\text{,}\mathrm{K}. The data shows QHS fans corresponding to $s=0,-1,-2,\pm 4$, with emerging fans at $s=-3,-3/2$, and $+2$. While the $\nu$ values are predominately multiples of four, the $s=0$ fan shows $\rho_{\mathrm{xx}}$ minima at single integer filling factors, which indicates a lifting of the Landau level spin and valley degeneracy, as well as unexpected states at $\nu=\pm 11$ and $\pm 13$ [Fig. 4(b)]. In Fig. 4(c) we show a line-cut of Fig. 4(a) along $B=$8.5\text{\,}\mathrm{T} where the various minima are labeled according to MB and $\nu$. At high fields, the $s=0$ fan exhibits $\rho_{\mathrm{xx}}$ minima corresponding to $\nu=\pm 1/2$ (Fig. S4). It is also interesting to examine the fan intersections. For example, where the $s=0$ and $s=-2$ fans meet, only the $\nu=-16$ minima from the $s=0$ fan persists while all others are suppressed. A similar behavior is observed for $\nu=+16$. This may point towards topological distinctions between the different moiré bands or Landau levels Lian et al. (2018).

We have presented a study of electrical transport in twisted double bilayer graphene, a system that exhibits single particle gaps in the moire band spectrum along with correlated insulators at half and quarter fillings of the first conduction band. A combination of theoretical and experimental data suggests that correlated insulators are most likely to occur in a band when it is maximally separated from neighbouring bands. Measurements in an in-plane magnetic field indicate that the correlated insulators are predominantly spin polarized at 1/4 and 3/4 fillings, and valley polarized at 1/2 filling.

Note added.—During the preparation of this manuscript we became aware of two related studies Liu et al. (2019); Cao et al. (2019b).

I Supplemental Material:

Correlated Insulating States in Twisted Double Bilayer Graphene

In TDBG, the mechanism of flat band formation is different from that of TBG, as noted in the main text. The nature of the low energy bands can be understood qualitatively by first considering the limit in which the inner two graphene layers, the layers which have a relative twist, are decoupled. As illustrated in Fig. S1(a), the low-energy states at the moiré Brillouin zone (MBZ) corners are then localized on individual layers from 1 to layer 4, and on the honeycomb sub-lattice in that layer which does not have a near-neighbor in the same bilayer. The energies at the MBZ corners are $3\Delta_{V}/2$, $\Delta_{V}/2$, $-\Delta_{V}/2$ and $-3\Delta_{V}/2$ respectively, where $\Delta_{V}$ is the on-site energy difference between adjacent layers. While in TBG, low-energy states localize on both sub-lattices of each layer. The color-coded bars ($y$-label on the right) in Fig S1 give the layer-dependent probability distributions for the corresponding high-symmetry $k$-points in the MBZ marked by colored points in the band structure. When the inner layers are then coupled, level repulsion tends to push the energies of states with high weight in these layers to high energies. The interplay between sublattice-coupling within layers and sublattice-dependent interlayer tunneling which yields the magic angle flat bands in bilayers does not take place because only one sublattice is available in each layer at low energies. This leaves low energy bands that are dominantly localized in the outer two layers, as illustrated in Fig. S1(b).

At $n=n_{s}/2$, the screened potential energy difference between two adjacent graphene layers is $\Delta^{\epsilon}_{V}=\frac{4\pi e^{2}nd}{\epsilon}=\frac{2\pi e^{2}n_{s}d}{\epsilon}\approx\frac{145\text{meV}}{\epsilon}$, where $d=0.34$nm, and the superlattice density $n_{s}\approx 4\times 10^{12}$ cm*-2* for $\theta=1.33^{\circ}$. Using a relative permittivity of twisted bilayer graphene of $\epsilon\sim 5$ /citenamefont Chung et!al. (2018) leads to a screened potential energy $\Delta^{\epsilon}_{V}\approx 30$meV.