Quantum transport in graphene in presence of strain-induced pseudo-Landau levels

Inhomogeneous lattice deformations in graphene generate an effective gauge field modulating the electronic spectrum [1–3]. However, compared to a real magnetic field, the formation of a pseudomagnetic field preserves time reversal symmetry, having an opposite sign in the two inequivalent K and K' valleys [4–6]. This leads to different behavior than for real magnetic field especially when introducing disorder. Experimentally, scanning-tunneling measurements on graphene nanobubbles have revealed an electronic spectrum consisting of pseudo-Landau levels (pLL), including a zero-energy peak, showing that moderate spatial deformations can introduce pseudomagnetic field values reaching hundreds of Tesla [7–9]. Pseudomagnetic fields (PMF) have also been analyzed in deformed crystals by an atomically controlled arrangement of CO molecules on a gold surface [10], graphene on Ir with intercalated Pb monolayer islands [11], or have been harnessed for designing innovative electronic and photonic graphene devices [12–17].

However, the intrinsic quantum transport fingerprints of graphene in presence of PMF still needs investigation, especially in disordered systems. In particular, random strain fluctuations are believed to be the dominating disorder source in high-quality on-substrate graphene devices [18, 19]. A signature of the pseudomagnetic n = 0 pLL state has been predicted for stretched graphene ribbons in the form of a quadruplet low-energy conductance resonance split by edge-induced valley mixing [16], but its experimental confirmation remains challenging, and the variable range of transport features in deformed graphene still require further in-depth exploration.

Here we study the effect of pLLs, generated by an ordered network of graphene nanobubbles (or pseudomagnetic dots), using an efficient real space Kubo quantum transport methodology. We consider samples containing electron–hole puddles caused by substrate interactions where the presence of pLLs leads to several anomalous transport features resulting from an intertwined contribution of pseudomagnetic field and disorder effects. The formation of the zero-energy pLL reduces the mean free path by orders of magnitudes in the limit of low defect density (c, modeled by a random distribution of Gaussian impurities), but scales as ${{\ell }}_{e}\sim c$ up to $c\sim 1 \%$, in contradiction with the usual Fermi golden rule (FGR) argument [20] predicting ${{\ell }}_{e}\sim 1/c$ which we obtain by considering the unstrained structure. Additionally, our simulations show that ${{\ell }}_{e}\sim 1/\sqrt{{B}_{{s}}}$ which evidences a superimposed scaling with the pseudomagnetic length. Finally, we demonstrate further unconventional behavior of the conductivity close to zero-energy where it becomes largely independent of energy for large enough disorder (above $1 \%$). Although in the limit of sample length $L\to \infty$, one expects an insulating behavior at zero temperature. These features manifest the strong influence that strain induced PMF have on the quantum transport properties and also suggest possibilities for metal-insulator transition driven by strain fields.

To describe the electronic properties of graphene, we use the common tight-binding (TB) model of graphene

$$\begin{eqnarray}&&\hat{H}=\displaystyle \sum _{i}{\epsilon }_{i}{c}_{i}^{\dagger }{c}_{i}+\displaystyle \sum _{\langle i,j\rangle }{\gamma }_{0}{c}_{i}^{\dagger }{c}_{j},\end{eqnarray} \tag{ 1 }$$

where ${\epsilon }_{i}$ is the onsite energy and the sum over $\langle i,j\rangle$ runs over nearest neighbor sites with ${\gamma }_{0}=-2.7\,{\rm{eV}}$. To generate a PMF, we consider a dot with radius R subjected to a planar triaxial displacement, which in polar coordinates are given by [22, 23]

$$\begin{eqnarray}&&{\bf{u}}=({u}_{r},{u}_{\theta })=({u}_{0}{r}^{2}{\rm{\sin }}(3\theta ),{u}_{0}{r}^{2}{\rm{\cos }}(\theta )),\end{eqnarray} \tag{ 2 }$$

where r and θ are the polar coordinates and u₀ determines the strength of the strain field. This displacement field gives rise to a constant PMF [23] given by B_s. Here we ignore out-of-plane and curvature [21, 24, 25] components of the strain field in order to get a simple connection between the strain and the magnitude of the PMF. This is justified when in-plane strain dominates (to induce pLLs) and no sharp bends are present [26–28]. The PMF is inherently related to the first order expansion of the TB model, but we emphasize that the calculations reported below do not rely on such approximations, instead we use the modified atomic positions to modify the TB parameters. Changing the atomic positions according to the displacement field ${\bf{u}}$ alters the bond lengths d_ij, and thereby leads to renormalized TB hopping parameters

$$\begin{eqnarray}&&{\gamma }_{{ij}}={\gamma }_{0}{{\rm{e}}}^{[-\beta ({a}_{{ij}}/{a}_{0}-1)]},\end{eqnarray} \tag{ 3 }$$

where $\beta =\partial {\rm{log}}(\gamma )/\partial {\rm{log}}(a){| }_{a={a}_{0}}\approx 3.37$ [1]. Here ${a}_{{ij}}=({a}_{0}^{2}+{\epsilon }_{{xx}}{x}_{{ij}}^{2}+{\epsilon }_{{yy}}{y}_{{ij}}^{2}+2{\epsilon }_{{xy}}{x}_{{ij}}{y}_{{ij}})/{a}_{0}$ denotes the modified bond length caused by the deformation, where ${a}_{0}=0.142\,{\rm{nm}}$ is the equilibrium bond length and the strain tensor is ${\epsilon }_{\nu \mu }=({\partial }_{\mu }{u}_{\nu }+{\partial }_{\nu }{u}_{\mu })/2$ with $\nu ,\mu =x,y$ [29]. Using these definitions the deformation field in equation (2) gives rise to a magnetic field of ${B}_{{s}}=8{u}_{0}{\rm{\hslash }}\beta /2{{ea}}_{0}$. At last, we note that equation (3) can easily incorporate out-of-plane components of the strain as it only depends on the change in bond lengths. The connection between strain tensor and bond length deformation is only approximate [30] but sufficient for the present analysis of generic quantum transport fingerprints of strain-induced pLLs. Especially, we note that the effect of the PMF on the LDOS is qualitatively unchanged by relaxation using molecular dynamics methods [22, 27, 31–35] and that second nearest neighbor terms only contribute with a scalar potential without generating pseudomagnetic effects [22, 26].

Outside the central dot region, we apply a smoothing to the strain tensor to assure a soft transition to the strain-free pristine regions. This is accomplished by applying a Gaussian transformation $\epsilon \to \epsilon ^{\prime} =\epsilon \ {\rm{\exp }}[-{(r-R)}^{2}/2{\sigma }^{2}]$ for $r\gt R$, where R is the radius of the PMF region as shown in the left inset of figure 1. The triaxial strain for $r\lt R$ gives rise to a constant PMF, whereas for $r\gt R$ an r-dependent PMF of opposite sign develops [21]. The field of opposite sign within the smoothing region arises because we apply the smoothing to the physical strain, and not to the derived PMF. We repeat this pseudomagnetic dot deformation in a periodic array with lattice constant L (see left inset of figure 1). This system is an idealized representation of the array of self-formed bubble deformations giving rise to the pLL structure envisioned experimentally by Lu et al [8]. We note that similar calculations have been performed for other array symmetries yielding qualitatively the same conclusions.

Zoom In Zoom Out Reset image size

We study the quantum transport using an order-N, real space implementation of the Kubo approach for the conductivity ${\sigma }_{{xx}}(E,t)$ [20, 36, 37]. The scaling properties of ${\sigma }_{{xx}}$ are followed through the dynamics of electronic wavepackets using ${\sigma }_{{xx}}(E,t)={e}^{2}\rho (E){D}_{x}(E,t)/2$, where $\rho (E)\ =\mathrm{Tr}[\delta (E-\hat{H})]$ is the density of states (DOS), ${D}_{x}(E,t)={\rm{\Delta }}{X}^{2}(E,t)/t$ is the diffusion coefficient and the time- and energy-dependent mean square displacement of the wavepacket is

$$\begin{eqnarray}&&{\rm{\Delta }}{X}^{2}(E,t)=\mathrm{Tr}[\delta (E-\hat{H}){| \hat{X}(t)-\hat{X}(0)| }^{2}]/\rho (E),\end{eqnarray} \tag{ 4 }$$

where $\hat{X}(t)$ is the position operator in Heisenberg representation. The semi-classical conductivity can be calculated using ${\sigma }_{{xx}}^{{\rm{SC}}}(E)={e}^{2}\rho (E){{\rm{\max }}}_{t}({D}_{x}(E,t))/2$. Calculations are performed on systems containing approximately $5\times {10}^{6}$ atoms. The state propagation is followed through a total simulation time of 14 ps and split into three intervals with different time steps ${\rm{\Delta }}{t}_{1}=1$ fs, ${\rm{\Delta }}{t}_{2}=15$ fs and ${\rm{\Delta }}{t}_{3}=60$ fs. The expansion of the time evolution operator uses Chebyshev polynomials corresponding to coefficients larger than ${10}^{-12}$. Finally, we approximate the traces using random phase states and the Lanczos method with 1500 iterations and a broadening of $\eta =5$ meV [38].

To mimick the effects of electron–hole puddles induced by the substrate [39], the Hamiltonian incorporates long range impurities defined by onsite potential on the n'th site ${V}_{n}={\sum }_{i}{\epsilon }_{i}{{\rm{e}}}^{-| {{\bf{r}}}_{n}-{{\bf{r}}}_{i}{| }^{2}/(2{\xi }^{2})}$, where ${{\bf{r}}}_{i}$ is the center of the ith impurity, ${\epsilon }_{i}\in [-W/2,W/2]$ ($W=2| {\gamma }_{0}|$) is the maximum onsite energy and the range is $\xi =3{a}_{0}$. A low concentration of impurities ($c=0.05 \%$, $c=0.1 \%$ and $c=0.2 \%$) is distributed randomly throughout the sample allowing us to reach the diffusive regime [20, 37, 39, 40].

We first consider the DOS of a pseudomagnetic dot array with $L=200a$, $R=40a$, $\sigma =10a$ ($a=\sqrt{3}{a}_{0}$), and a maximum strain of approximately $15 \%$. The parameters lead to a PMF of 450 T (see figure 1), corresponding to a pseudomagnetic length ${{\ell }}_{B}=\sqrt{{\rm{\hslash }}/{{eB}}_{{s}}}\sim 1.22$ nm. The pseudomagnetic length is smaller than the dot size which is necessary for the appearance of Landau quantization [6]. Similar conclusions are obtained for other geometries by varying L and R (not shown). The DOS of the PMF array in figure 1 shows the formation of pLLs following the expected $\sim \sqrt{| n| }$ behavior (figure 1, right inset), where n is the pLL index. In particular, a strong zero-energy peak is formed, as for real magnetic fields. The peak features are superimposed with the linear dispersion characteristic of the unstrained graphene calculated (black, dashed).

To characterize the impact of pLL states on the quantum transport, we first consider the mean free path (see figure 2(a)), extracted in the diffusive regime as ${{\ell }}_{e}(E)={{\rm{\max }}}_{t}(D(E,t))/2{{v}}_{{\rm{F}}}$ where $D(E,t)={D}_{x}(E,t)+{D}_{y}(E,t)$. We use ${v}_{{\rm{F}}}=8.6\times {10}^{5}$ m s⁻¹ [20] even though the Fermi velocity v_F is weakly strain dependent [15, 41–46]. In the absence of a strain field (dashed lines), the mean free path follows the FGR and scales with the impurity density as ${{\ell }}_{e}(E)\sim 1/c$, with values and energy dependence dictated by intervalley scattering [37]. When superimposing the strain field, we observe a large difference between ${{\ell }}_{e}(E)$ with (full lines) and without (dashed lines) the deformation field (see figure 2(a)). For a given impurity concentration c, we find a systematic decrease of ${{\ell }}_{e}(E)$ compared to the unstrained case. The enhanced scattering induced by strain is particularly significant at the pLL energies where the DOS is large, as revealed by the dips in ${{\ell }}_{e}(E)$ in figure 2(a). At the centers of pLLs, the mean free path scaling is driven by the pseudomagnetic length ${{\ell }}_{B}$ i.e. ${{\ell }}_{e}(E)\sim 1/\sqrt{{B}_{{s}}}$ (see figure 2(b)).

Zoom In Zoom Out Reset image size

Next, we focus on the low energy regime in the inset of figure 2(a). Here, we observe a surprising crossover in the mean free path near $E=\pm 50$ meV (indicated by vertical dashed lines in the inset of figure 2(a)) where the scaling of ${{\ell }}_{e}(E)$ with impurity density is reversed going from higher to lower energies. While the high-energy ($| E| \gt 50$ meV) behavior follows the FGR, ${{\ell }}_{e}(E)\sim 1/c$ (filled triangles in figure 2(c)), the opposite behavior is observed for $| E| \lt 50$ meV (filled circles in figure 2(c)). Consequently, we observe an anomalous regime (filled circles) at low energies for which ${{\ell }}_{e}$ increases with disorder. This suggests that the pseudomagnetic field counteracts disorder effects. Increasing the disorder beyond a certain concentration disrupts this mechanism and we recover the more conventional decay of the mean free path with increasing disorder.

To further characterize the transport fingerprints caused by the strain-induced pLLs, we analyze the diffusion coefficients $D(E,t)$ in figure 3 at energies marked with symbols in figures 1 and 2, allowing us to distinguish various transport regimes. A diffusive regime is observed at the energy of the first pLL (figure 3(c), square), with a constant asymptotic behavior of $D(E,t)$ both with (full lines) and without a strain field (dashed lines). The additional scattering due to the strain field results in a lower diffusion coefficient leading to the shorter mean free path at the pLL energies as discussed above.

Zoom In Zoom Out Reset image size

The low energy regime, however, is qualitatively different (figures 3(a)–(b)). The DOS is dominated by the zeroth pLL, which induces a strong sublattice polarization [6, 22, 47]. Here the decrease of $D(E,t)$ with time reveals the onset of localization effects. In figures 3(a)–(b), we also show the time evolution of the diffusion coefficient in the absence of the strain field (dashed lines). Comparing the curves with and without strain, we conclude that localization effects are generated by the presence of the pseudomagnetic field. This is especially clear in figure 3(d) showing the variation of $D(E=0.1\ {\rm{eV}},t)$ while decreasing the pseudomagnetic length. The enhancement of the localization effects are evident through the reduced value and sharper decay of $D(t)$ for higher pseudomagnetic field strengths.

The anomalous transport regime around zero-energy is also manifested in the diffusion coefficients $D(E,t)$. Indeed, $D(E,t)$ exhibits a qualitatively different behavior for $| E| \lt 50$ meV (figure 3(a)) compared to $| E| \gt 50$ meV (figure 3(b)) even though the strain-induced pLL causes localization in both regimes. In agreement with the discussion of the mean free path, the absolute value of the diffusion coefficient at low energy increases with impurity concentration.

Finally, we consider the semi-classical (${\sigma }^{{\rm{SC}}}(E)$) and quantum ($\sigma (E)$) conductivities at low energies for different impurity concentrations (figure 4(a)). As seen in the figure, ${\sigma }^{{\rm{SC}}}(E)$ always remains larger than or equal to ${\sigma }_{0}=4{e}^{2}/\pi h$ (horizontal dashed line) over the whole energy spectrum, in agreement with the lower bound of the semi-classical value [48, 49]. The larger value of ${\sigma }^{{\rm{SC}}}(E\sim 0)$ may be related to a low energy behavior, similar to the effect caused by zero-energy modes around lattice monovacancies [50–52].

Zoom In Zoom Out Reset image size

The value of ${\sigma }_{0}$ actually separates two different transport regimes, which are also identified by scrutinizing the long time behavior of $D(E,t)$ [51–56]. As long as $\sigma (E)\gt {\sigma }_{0}$, the system remains diffusive in a metallic (i.e. non insulating) regime. Correspondingly, the diffusion coefficient saturates at long times. On the other hand, the condition $\sigma (E)\lt {\sigma }_{0}$ implies the existence of localization effects, whose strength (weak or strong localization) is dictated by the impurity density and the length-scale (or time-scale) at which the quantum conductivity is evaluated.

The anomalous scaling of the mean free path with disorder is further translated to similar intriguing scaling for the quantum conductivity. In presence of strain, the zero-energy quantum conductivity (taken at the maximum calculation time, t_max) remains very close to the fundamental semi-classical limit ${\sigma }_{0}$ (see figure 4(b)-inset). The dynamical conductivity, $\sigma (E,{t}_{{\rm{\max }}})$ also shows a complex energy-dependent profile evolution with defect concentration. According to figure 4(b), the conductivity for $| E| \geqslant 0.1\,{\rm{eV}}$ decays with defect density, i.e. $\sigma (| E| \geqslant 0.1\ {\rm{eV}},{t}_{{\rm{\max }}})\sim 1/c$. However, at lower energies such behavior is interrupted and $\sigma (| E| \leqslant 0.1\ {\rm{eV}},{t}_{{\rm{\max }}})$ instead increases with defect density. A special situation is observed in figure 4(b)-inset where $\sigma (| E| \sim 0\ {\rm{eV}},{t}_{{\rm{\max }}})$ is almost constant for a wide range of disorder.

The asymmetric disorder effect between low and high energy is in sharp contrast to many types of disordered graphene systems with similar mean free paths (namely on the order of a few nanometers) at zero magnetic field [51–56], and resembles the percolation mechanisms [56–67], suggesting that the zeroth order pLL plays a critical role. This requires additional insight, provided in the next section.

In this section, we provide a tentative interpretation for the opposite effect of the Gaussian disorder on the low and high-energy behavior, in the classical transport regime (dictated by the mean free path). To this end, we scrutinize the role of real-space distribution of electronic states subjected to a pseudomagnetic field. It has been shown that low-energy states are preferentially located inside the deformed region, while high-energy states correspond to states outside the bubbles [68]. Furthermore, the zeroth pLL induces sublattice polarization [22], where the preferred sublattice depends on the deformation direction [6]. This sublattice polarization of the first Landau level is a unique feature of PMFs. The strain direction considered here causes the n = 0 pLL states in the central region (red color in left inset of figure 1) to localize on the B-sublattice. On the other hand, the states in the outer ring (blue color) with an opposite pseudomagnetic field sign, is A-sublattice polarized. Consequently, the low-energy states inside the bubbles tend to be confined by this opposite sublattice polarization from the normal conducting states outside the bubbles in the undeformed regions. Now, because we introduce long-range disorder randomly on both A and B sublattices everywhere in the sample, the increased disorder concentration disturbs the perfect low-energy sublattice polarization. Thus, the spatial sublattice-confinement of the states is broken, leading to an increased mean free path and conductivity. We note that this scenario is compatible with the fact that the quantum localization effects increase with disorder, as expected (see for instance figure 3(a)).

To further confirm the sublattice-polarization scenario, we introduce a sublattice-selective (vacancy) disorder, which has been shown to provoke asymmetric transport phenomena [69, 70]. We model such unreconstructed vacancies by removing carbon atoms randomly everywhere in the sample. In figure 5, the diffusion coefficients for sublattice-selective vacancy disorder, located either on the A- or B-sublattice, are shown, both with and without strain. The inequivalent effect of A- and B-vacancies exists only when the pseudomagnetic field is present. Such vacancies increase the number of states in the opposite sublattice from where they reside. As such, the A-vacancies increase the number of states in the B-sublattice. Because the low-energy states are mostly polarized on the B-sublattice, they are strongly affected by short-range valley-mixing effects induced by A-vacancies, leading to stronger localization effects (blue curve on figure 5). Indeed, unlike the zeroth order Landau level caused by a real magnetic field where valleys are located on opposite sublattices, the valleys are concentrated on the same sublattice in the presence of a pseudomagnetic field [6]. Finally, even if the A-vacancies yield stronger localization effects, they concurrently lead to longer mean-free paths, when compared to B-vacancies. This observation further supports our interpretation since A-vacancies (increased B-sublattice polarization) increase the number of states to percolate through the blocking A-polarized ring around the bubbles.

Zoom In Zoom Out Reset image size

The presence of strain-induced PMF forming pLL has been shown to trigger unconventional quantum transport features compared to other types of disordered graphene systems. The transport physics close to half-filling is particularly remarkable, with a huge drop of the transport mean free path upon switching on a pseudomagnetic field, indicating the possibility of a mechanically induced metal-insulator transition.

We note, however, that despite of the weak localization caused by the deformation field, the strength of quantum interferences, even in the large disorder limit, remains anomalously weak for the sublattice polarized zeroth pseudo Landau level (compared to other types of disorder). This anomalous scaling with disorder gives a peculiar transport fingerprint of the combination of pseudomagnetic field and disorder effects, especially for the sublattice-polarization unique to the pseudomagnetic field compared to real magnetic fields.

Superlattice networks of deformation fields have been experimentally demonstrated [71], and could be tunable under pressure [27, 72, 73] or by temperature. In this way, varying the mean free path upon deformation, could reversibly switch the system between metallic and insulating states (for $\xi /L\ll 1$, with ξ the localization length and L the sample size). Finally, the role of pseudomagnetic field in spin transport deserves further consideration, in the context of spin manipulation and graphene spin-based devices [74].

The work by MS is funded by the Danish Council for Independent Research (DFF5051-00011). Funding from the European Union Seventh Framework Programme under grant agreement 604391 Graphene Flagship is acknowledged. SR acknowledges the Spanish Ministry of Economy and Competitiveness for funding (MAT2012-33911), the Secretaria de Universidades e Investigacion del Departamento de Economia y Conocimiento de la Generalidad de Cataluña and the Severo Ochoa Program (MINECO SEV-2013-0295). The Center for Nanostructured Graphene (CNG) is sponsored by the Danish National Research Foundation (DNRF103).