Gyrotropic Zener tunneling and nonlinear IV curves in the zero-energy Landau level of graphene in a strong magnetic field

We have investigated tunneling current through a suspended graphene Corbino disk in high magnetic fields at the Dirac point, i.e. at filling factor ν = 0. At the onset of the dielectric breakdown the current through the disk grows exponentially before ohmic behaviour, but in a manner distinct from thermal activation. We find that Zener tunneling between Landau sublevels dominates, facilitated by tilting of the source-drain bias potential. According to our analytic modelling, the Zener tunneling is strongly affected by the gyrotropic force (Lorentz force) due to the high magnetic field.

The zero-energy Landau level is a unique feature of a graphene monolayer in a strong magnetic field 1 . Without interaction, the state is four-fold degenerate with respect to two values of true electron spin and two values of pseudospin, which describes distribution of electrons between two valleys in the graphene Brillouin zone. The whole state is half-filled, i.e. one half of possible Landau states are occupied and another half of them are empty (holes). However, Coulomb interaction lifts the degeneracy and the original four-fold degenerate zero-energy states split into two states with a gap between them. All states with the lower energy are now occupied, while all states with higher energy are empty. This leads to the integer quantum Hall effect, and electrons at the lower energy sublevel form an insulator state separated from the higher energy sublevel by an energy gap. It was expected that the insulator state is ferromagnetic 1 , although the experiment by Young et al. 2 and Giesbers et al. 3 did not confirm this. Other options (pseudospin "ferromagnetism", e.g.) have also been considered for the ground state 4 .
Independently from the character of the gap, the zero-energy Landau state is an insulator with the gap Δ of the order of the characteristic Coulomb energy equal to  e / B 2  in graphene. Here ε is the dielectric constant, is the magnetic length, B is the magnetic field, and Φ 0 = hc/e is the single-electron flux quantum. But any insulator in a high voltage, which we note by V cr , becomes a conductor (dielectric breakdown). At bias voltages V exceeding V cr , i.e. after the dielectric breakdown, the IV curve is close to linear (ohmic regime). On the other hand, at electric fields essentially less than V cr the conductance is very small and strongly nonlinear. At low temperatures, an exponentially small current I emerges due to Zener tunneling between the two bands 5 , creating an electron in the empty upper band and a hole in the full lower band 6 . Without a magnetic field Z b 3/2 must be close to V cr . Here a is the lattice constant, d is the distance between contacts, and E b is the band width.
In the quantum Hall state the two bands are flat in the bulk and hence the tunneling current vanishes. In a strong magnetic field the derivation, which leads to Eq. (1) becomes invalid. Motion of electrons in a strong magnetic field is not determined by the Newton's second law with the inertial force proportional to the electron mass, but by the equation of motion of guided centra (centra of Larmor circles around which electrons move). In this equation the external force on the electron is balanced by the Lorentz force and the inertia can be neglected. The equation of motion of guided centra is similar to the equation of motion of quantized vortices in superfluids and clean superconductors with the Magnus force (analog of the Lorentz force on an electron) balancing the external force.
Vortex-like equations of motion essentially modify the semiclassical theory of quantum tunneling. For vortices it was first demonstrated by Volovik 7 , who considered nucleation of a circular vortex half-loop near a plane boundary. Probability of quantum tunneling derived from the semiclassical theory based on the equations of vortex motion was confirmed by the results of the many-body approach based on the Gross-Pitaevskii equation 8 It is possible to demonstrate 9 how the theory of usual quantum tunneling for a massive electron governed by Newton's law transforms to the theory based on the equation of vortex motion when the inertial force proportional to the electron mass becomes much weaker than the gyrotropic force (Magnus force on the vortex, or the Lorentz force on electrons in a strong magnetic field).
Strong magnetic field dramatically affects probability of tunnelling through potential barriers in the configurational space, see e.g. Jain and Kivelson 10 or Dykman et al. 11 . Here we study Zener tunnelling in which a particle crosses a classically forbidden interval in the energy spectrum (energy gap). The important role of Zener tunneling on electron transport in the quantum Hall regime of a 2D electron gas has already been discussed in the context of quasi-elastic-inter-Landau-level scattering (QUILLS) 12,13 and in the connection with magnetoresistance oscillations in the ohmic regime 14,15 . Here we present direct measurement of a current produced by Zener tunneling and calculation of its probability in the subohmic regime, before the breakdown. We analyze Zener tunneling between two Landau sublevels emerging from the zero-energy Landau level (around filling factor ν = = 0 hn eB , where n is the charge carrier density) using the semiclassical theory of tunneling. Also, we measured nonlinear IV curves experimentally at low voltages (weak electric fields) and argue that they provide evidence of Zener tunneling governed by vortex-like equations of motion. This type of quantum tunneling was called Hall tunneling for vortices in superconductors 16 and Magnus tunneling for vortices in superfluids 9 . Here we call it gyrotropic tunneling because of the crucial role the gyroscopic force (Lorentz force on electrons or Magnus force on vortices) in the process of tunneling.
The IV measurements were done on undoped suspended graphene in the Corbino disk geometry. Corbino geometry has a distinct advantage over Hall bar geometry: no edge states persist and the breakdown takes place through the states in the bulk. Consequently, this geometry has been successfully employed for studying phenomena related to transitions between Landau levels in 2D electron gases: microwave-induced resistance oscillations & zero-resistance states 17 , and phonon-induced resistance oscillations 18 . It has also been employed in studies of edge channel tunneling 19 , regular Zener tunneling between different Landau-levels 20 , as well as the bootstrap electron heating (BSEH) model 21 leading to breakdown of the quantum Hall effect [22][23][24][25] .
By fitting theoretical formulas to the experimentally obtained IV curves, we could reveal the interval of voltages where the measured current is reasonably well described by the exponential law which follows from the gyrotropic Zener tunneling theory in a strong magnetic field. This provides evidence of quantum tunneling processes governed not by the particle mass but by the gyrotropic force on a particle.

Theoretical Results
Since we address the semiclassical theory of quantum tunneling, we need the classical equations of gyrotropic motion for guided centers (in cgs units): where  is the energy of the electron-hole pair, which depends on the position vector r connecting the hole in the lower subband with the electron in the upper subband. The right-hand side of the equation is the external force and the left-hand side is the Lorentz force. This vector equation is equivalent to two Hamiltonian equations, for the pair of canonically conjugate variables x − P x , where the conjugate momentum is connected with the second coordinate y of the guided center: We consider tunneling between two bands separated by the gap Δ. In an electric field the bands are tilted (Fig. 1a), and the energy of the electron in the upper band is equal to the electron energy in the lower band at the distance An electron tunneling from the lower band to the upper band leaves in the valence band a hole. The electric field E = V/d drives the electron and the hole in opposite directions but this is resisted by the Coulomb attraction between the electron and the hole. The energy of the electron is The position vector r connects the positions of the electron and the hole. For the sake of simplicity we assume that the hole is at rest, and only the electron participates in the tunneling event. We consider a weak electric field parallel to the x-axis. The shortest path across the interband barrier is at At a weak electric field, x f B , and two zeroes of the energy are very close to the points x = 0, y = 0 and x = x f , y = 0:   One should remember that we treat electrons as point-like objects, although any electron is distributed in space over the scale  B . Thus, an electron position can be determined only with uncertainty B  . Therefore differences between x 1 and x = 0 and between x 2 and x = x f are within the error bar of our analysis.
In the semiclassical theory of quantum tunneling one must find a trajectory of the electron at constant energy (equal to 0 in our case). Such a trajectory is possible only for imaginary y determined from Eq. (8) at 0  = : Analytical continuation of the energy as a function of coordinates from the real plane y = 0 to the complex plane = 0  is shown in Fig. 1c. Since x f B , it is accurate enough to approximate the trajectory in the complex plane by a simple linear function y ≈ −ix.
After the trajectory has been found we can calculate the exponent in the exponential law for a small current: Here ImS is the imaginary part of the classical action variation along the trajectory: Thus, the exponent is The exponent of the law is of the order of the number of flux quanta in the area x f 2 ∼ , while for vortex tunneling the exponent is the number of bosons in the same area 9 . Remarkably, the exponent depends only on the thickness of the barrier but not on its height. Qualitatively it is of the order x / f B 2 2  , i.e. of the order of the exponent of the overlapping integral for wave functions of two Landau electron states at the distance x f between them.
Equation (14) is exact if we know the gap Δ. However, we know it only by the order of magnitude: This yields the exponential law (3) for the current with the characteristic electric field: the exponential growth of the current saturates, and a crossover to the ohmic regime on the IV curve takes place. Our theory addresses only the exponential law of the tunneling probability ignoring the prefactor. The calculation of the prefactor is much more difficult and would require a lot of information on scattering processes and on spin (pseudospin) structure of electron-hole pairs.

Experimental results
The IV characteristics measured at B = 2 T using positive bias on the inner Corbino contact are presented in Fig. 2. The data are compared with two different theoretical models based on Zener tunneling: the conventional model ∝ − V V exp( ( / )) Z of Eq. 1, and the gyrotropic model of Eq. 3. The third curve is based on thermal activation, where the number of independent activated regions is a fitting parameter (see below). The data presented in this section is from our sample S2 (see Methods section); very similar results were found on sample S1.
In the case of Zener tunneling, both models can be fitted to the data but the gyrotropic Zener tunneling yields a better agreement. More importantly, the fitting parameter V Z finds more reasonable values when using the gyrotropic model. By reasonable we refer to the fact that the characteristic voltage V Z should be close to the critical voltage V cr , the point where significant current starts to flow through the bulk of the Corbino disk. We define the critical voltage V cr = V(I = 10 nA), which correspond to the voltage in the middle of the region where the current grows significantly. One can note that V cr ≈ 10 mV in Fig. 2. By comparing this to the values of V Z in Fig. 2, 21.5 meV for the gyrotropic tunneling and 99 meV for the conventional case, along with the better agreement of the fit, we conclude that the gyrotropic model is the best to account for our observations at B = 2 T. Moreover, the same conclusion is found at negative bias voltages.
In addition to the Zener tunneling models, we have also considered conduction by thermal activation across localized states. In Fig. 2, there is a fit based on such an activation model with tunneling current between impurity islands given by: where N is the number of islands in the percolation chain, Δ is the activation gap energy, and I 0 is a prefactor. By assuming that the tunneling happens through a series of equidistant jumps, the effective voltage is reduced by a factor of N, which in turn makes the dI/dV slope shallower. Without setting N = 5-10, the activation models become too steep to fit the data. Additionally, this model requires us to set the effective temperature to T = 1.7 K and the gap energy down to 1 3 Δ = . meV (compared to the measured eV cr ≈ 10 meV from the critical voltage). Even with these assumptions in the activation model the gyrotropic model nevertheless results in better fits in the region where the dielectric break down takes place V V ( ) cr ∼ .
IV curves were also investigated at different magnetic fields, B = 1-9 T, and fits such as the ones in Fig. 2 were made in order to extract the magnetic field dependence of V Z and V cr . V Z was obtained from the gyrotropic Zener tunneling fits, whereas the critical voltage V cr was taken as the voltage needed to drive I = 10 nA through the device, as defined before.
The extracted V Z and V cr are displayed in Fig. 3 for the positive bias data. It can be seen that an approximate correspondence ∼ V V cr Z is valid over the whole range as well as approximately linear magnetic field dependence in accordance with Eq. (15). There is a factor of two difference between V Z and V cr throughout the magnetic field range, which was already noted in the data at B = 2 T in Fig. 2. The conventional Zener tunneling yields V Z that is about ten times larger than V cr . A theoretical estimate for V Z can be calculated using Eq. (14), converted to SI units   tion), and ε = ε 0 ε r the permittivity of graphene. The effective value of relative permittivity ε r = 5 was used 4 . The model yields the observed linear behaviour, V Z ∝ B, but with 10  times higher theoretical threshold values calculated from the equation above for V Z . Since the measured V Z are reasonable in the sense that they correspond to V cr , there must be a substantial inaccuracy in the parameters employed in the theoretical estimation of V Z . The most likely culprit is the estimate for the gap  ∆ ∼ e B 2  that is poorly known and sample dependent 2,26 .

Discussion
There have been several investigations on electric breakdown in QH systems. Various models have been developed, which in addition to Zener tunneling also include thermal instabilities. Thermal runaway, i.e. positive feedback where increase in input power increases conductance that in turn increases power dissipation, would be one option 21 . The critical field E c has been found to scale approximately linearly with the magnetic field and the largest electric field reported for E c amounts to 25 kV/m at 9 T 27 (see the analysis and compilation in ref. 25 for filling factors ν = 2, 4, 6). Our 9 T result at ν = 0, E c = 100 kV/m, is by a factor of four larger than the quoted result in GaAs heterostructures. On the basis of inelastic electron -acoustic phonon scattering, the BSEH model leads to E c ∝ B 3/2 , which has been verified experimentally 28 , but even in these experiments E c extrapolates only to 14 kV/m at 9 T. Taking also into account the weak coupling to acoustic phonons in similar graphene samples 29 , we conclude that in graphene we are not dealing with thermal runaway at the onset of the nonlinear IV characteristics. Further evidence of nonthermal origin of the IV characteristics is provided by the absence of hysteresis, which one would expect to exist for this type of instabilities. The inset in Fig. 2 displays IVs when ramping the bias voltage up (blue) and down (red), and subsequently the curves laying on top of each other indicates negligible hysteresis. One factor which may contribute to the absence of the hysteresis is the equilibration of the edge states at the metal/graphene contacts in Corbino samples 30 and the shortness of the bulk graphene: the avalanche type of excitation generation, typically assumed to be present in BSEH analysis, requires a substantial distance which can be even in excess of 100 μm 31,32 . Note that BSEH type of behavior with E c ∝ B 3/2 has been observed in epitaxial graphene with sample lengths 5-35 μm 33 .
In Fig. 2, we compared our data against two quantum tunneling models based on Zener tunneling: the conventional form and the gyrotropic model presented here. Evidently, the fit to the law of gyrotropic tunneling looks better over the full current variation by nearly three orders of magnitude. For the gyrotropic tunneling, moreover, the fitting parameter V Z is closer to the critical voltage V cr than for the usual Zener tunneling. According to the theory of Zener tunneling, the critical voltage V cr should be on the same order of magnitude as the voltage scale of tunneling V Z . This relation is clearly better valid for the gyrotropic tunneling than for the regular Zener tunneling. Additionally, we observe linear field dependence of V Z and V cr in the range B = 1-7 T, which is characteristic for the gyrotropic Zener tunneling but not for the conventional model. While our observations seem to agree well with the presented model for the gyrotropic Zener tunneling, the measured value of V Z is smaller by a factor of ten compared to that obtained from the order of magnitude estimation in Eq. 15 using an energy gap of  e / B 2  Δ ∼ . We assign this discrepancy between the gyrotropic theory and our experiment to the uncertainty in the effective gap value in the transport experiment at zero filling factor. In similar context, observation by the Singh and Deshmukh et al. 34 shows discrepancy in the observed activation gap energy and the breakdown Hall voltage. Moreover, if we compare the transport result 10  Δ meV of ref. 2 with the result Δ = 58 meV of the scanning SET measurement of ref. 35 (both at at B = 9 T), we note a difference by a factor of six. Note also that our Arrhenius data in Fig. 3 inset scaled up to B = 9 T (with linear field dependence) corresponds to 7  Δ meV. Taking this effective gap reduction factor into account, we find that our theory and experiment agree within a factor of two when using our measured Δ in Eq. 14 instead of e B 2 ∆ ∼  . Reduction of the gap might be induced by the nonuniformity of the charge density caused by the bias voltage. In our experiment at low bias, the charge density is not altered by the source-drain biasing, and the measured state will settle to the ν = 0 quantum Hall state with a robust gap 2 . The relatively high asymmetric source-drain bias offsets the chemical potential in a manner dependent on the transport path and contact resistances (possible with Schottky barriers, see below). The asymmetric bias is due to the use of a low-noise transconductance amplifier, which acts as a virtual ground for the terminal connected to it. As a consequence, there is an increase in the charge density near the biased lead, which will reduce the energy separation with respect to the next Landau sublevel. Since the other end is connected to virtual ground, there is going to be a region where the gap is unmodified. The voltage distribution of the sample across the disk follows resistance R(r), where r is the radial distance from the virtual ground. Since dR(r)/dr grows towards the unbiased terminal, the electric field is the largest in that regime. Hence, we expect the largest gap regime to dominate at the onset of charge transport, and our results correspond to the regime close to the Dirac point with the filling factor ν = 0. This conclusion was checked at B = 2 T by measuring the same IV characteristics in a symmetrized biasing configuration, where the sample was biased on both sides using current biases that resulted in chemical potential offset Δμ = ±eV/2 on the opposite ends. This way, the chemical potential remained at μ = 0 in the middle of the graphene ring, and the tunneling must have happened through the strongly gapped region.
At much higher fields, the cyclotron radius of charged particles is greatly reduced and the carriers can be treated as point-like objects. The breakdown mechanism will be dominated by the electric-field-induced local breakdown within the localized states in the bulk, leading to an avalanche type of breakdown 25 . Large Fano factors, observed in recent noise measurements on Corbino disks 28,36 , have yielded support to the avalanche picture of electron transport at the breakdown point of the quantum Hall effect in GaAs heterostructures. It seems that our IV results also include similar avalanche-induced transport at large magnetic fields 37 where the current increase in the IV curves is actually faster than predicted by our gyrotropic tunneling model. Such enhanced current by avalanches could explain the deviation of the observed V Z values in excess of the linear behavior in Fig. 3. Fermi level pinning of metallic contacts may also play a role in the gapped transport phenomena. In semiconducting CNTs e.g., this leads to FET operation that is based on Schottky barrier modulation 38 . In our devices this effect may come into play because the charge density in graphene is always nonuniform due to the contacts inducing n-type doping. The charge doping by contacts leads to pn interfaces at negative bias voltages 39 . In the presence of the magnetic-field-induced energy gap, the nonuniform doping may lead to Schottky-barrier-like structures that would have influence on the Zener tunneling phenomena. Consequently, we have concentrated on the Zener tunneling results on the positive bias where the sample is unipolar.
In conclusion, we have measured tunneling current through suspended graphene Corbino rings in the middle of the zero-energy Landau level corresponding to the filling factor ν = 0. We found that the tunneling current near the dielectric breakdown is consistent with a nonstandard quantum tunneling process, namely the gyrotropic Zener tunneling where the quantum tunneling takes place along a curved path due to the presence of a large magnetic field. The gyrotropic tunneling model also accounts for the observed linear increase of V Z and V cr as a function of magnetic field. The presented model fits the qualitative features of our data well by reproducing the exponential growth of current with a characteristic voltage V Z of expected magnitude in relation to the onset of current. This agreement in the magnitudes of the measured activation gap energy and V Z , calls attention to the IV characteristics as an alternative approach for Arrhenius type measurements. We also note that the ν = 1 state supports skyrmionic excitations 40,41 ; this gyrotropic Zener tunneling model could be a promising tool to investigate such excitations.

Methods
The Corbino disk geometry has the unique trait that there are no edge states connecting the measurement ports of the device at high magnetic fields, unlike in Hall bar devices. This allows direct probing of the bulk and dielectric break down of the system in the Corbino disk. We study this break down in current annealed suspended graphene devices.
The samples used in this work were fabricated in a manner described in ref. 42 . The fabrication was based on a sacrificial layer of lift-off-resist (LOR) which acted as support for our leads. The LOR was spun on a SiO 2 covered p + + doped silicon substrate that also served as a back gate. The graphene part, and its immediate surroundings, were suspended by exposing the LOR in that area to e-beam and dissolving the exposed parts in ethyl lactate, and then rinsing in hexane where low surface tension allowed us to simply lift the chip out of the liquid without destroying the delicate suspended graphene membrane.
Two suspended graphene Corbino samples of similar characteristics were measured: S1 and S2, both depicted in Fig. 4. S1 had inner and outer diameters of D i = 0.8 μm and D o = 3.2 μm. In turn, S2 had inner and outer diameters of D i = 1.5 μm and D o = 3.8 μm. Before the actual measurements, the samples were current annealed using voltages 1.6 V (S1) and 2 V (S2), which resulted in high mobility in both samples. The Dirac point resided at V D ≈ −2 V and the electron side mobility exceeded 10 5 cm 2 /(Vs) at 8 · 10 10 cm −2 for both samples. The residual charge density was estimated to be around n 0 = 2 × 10 9 cm −2 , as illustrated in loglog-plot of conductance in Fig. 5.
The samples were mounted in a BlueFors LD-400 dilution refrigerator with base temperature of 10 mK. The cryostat was equipped with a superconducting magnet producing up to 9 T magnetic field perpendicular to the graphene surface. The DC-conductance measurements were conducted using RC-filtered (R = 450 Ω, C = 30 nF) twisted pair lines. At high fields, the sample resistance in series with RC-filters remained high even at high bias voltages, which resulted in long RC time constants in the system (e.g. τ = R sample C = 3 s, when R sample = 100 MΩ). Consequently, a pause time of 30 s was used for each bias point in the IV-sweeps.
The IV-measurements were carried out by applying a voltage bias at the top of the cryostat and measuring the current through the device using a transimpedance amplifier (Stanford Research Systems SR570). Thus, the IV data consisted of sets DC-voltages and corresponding currents through the Corbino disk. The current noise level in the measurements was 10 −13 A. The charge carrier density was tuned by applying a voltage V g to the back gate. The inverse conductance for sample S1 as a function of the charge carrier density n and the perpendicular magnetic field B is shown in Fig. 6 displaying the expected fan of quantum hall states where ν = 0 is situated at n = 0. For the measurements reported here, the gate voltage was set to V D = −2 V in order to stay in the middle of the 0th Landau level at the filling factor ν = 0. Data availability. The datasets generated and analysed during the current study are available from the corresponding author on reasonable request.  . Fan diagram of sample S1 where the logarithm of the inverse zero-bias conductance measured using a lock-in amplifier operated at f = 3.333 Hz is displayed as a function of the gate-induced charge density n and the perpendicular magnetic field B. The quantum Hall states at specific filling factors are marked for the positive field.