Interlayer quantum transport in Dirac semimetal BaGa2

The quantum limit is quite easy to achieve once the band crossing exists exactly at the Fermi level (EF) in topological semimetals. In multilayered Dirac fermion systems, the density of Dirac fermions on the zeroth Landau levels (LLs) increases in proportion to the magnetic field, resulting in intriguing angle- and field-dependent interlayer tunneling conductivity near the quantum limit. BaGa2 is an example of a multilayered Dirac semimetal with its quasi-2D Dirac cone located at EF, providing a good platform to study its interlayer transport properties. In this paper, we report the negative interlayer magnetoresistance induced by the tunneling of Dirac fermions between the zeroth LLs of neighboring Ga layers in BaGa2. When the field deviates from the c-axis, the interlayer resistivity ρzz(θ) increases and finally results in a peak with the applied field perpendicular to the c-axis. These unusual interlayer transport properties are observed together in the Dirac semimetal under ambient pressure and are well explained by the model of tunneling between Dirac fermions in the quantum limit.

BaGa 2 is recently predicted to be a candidate of Dirac semimetals. The crystal structure of BaGa 2 consists of Ga honeycombnet layers and Ba layers. Its quasi-2D Dirac cone originates from Ga p z orbitals with the Dirac point located at E F 37 , while Ba layer between two Ga honeycomb-net layers act as the interlayer transport barrier. Thus, BaGa 2 is believed as an ideal candidate of quasi-2D Dirac material and it provides a good platform to investigate the unusual interlayer transport properties induced by the tunneling of Dirac fermions on the zeroth LLs. Figure 1a presents the schematic illustration of the Dirac fermion interlayer tunneling on the zeroth LLs in BaGa 2 .
In this work, the quasi-2D Dirac cone at K with the Dirac point located at E F is confirmed by the first-principles calculations, angle-dependent dHvA oscillations and angle-resolved photoemission spectroscopy (ARPES) measurements. The abovementioned NIMR and peak in angle-dependent interlayer resistivity are both observed clearly in BaGa 2 for the first time and well explained by the model of Dirac fermion tunneling between the zeroth LLs 34 .

Results
Materials characteristic and in-plane transport properties. BaGa 2 single crystals are grown by self-flux method (see the "Methods" section). The crystal structure consists of Ba planes and the Ga honeycomb-net planes which shows 2D characteristic. The (00l) lattice plane and lattice parameters are determined by the Xray diffraction (XRD) measurements on BaGa 2 single crystal and powder as shown in Fig. 1b and c. Both the in-plane resistivity (ρ xx ) and the interlayer resistivity (ρ zz ) show metallic behavior (Fig. 1d). The ratio of ρ zz /ρ xx is about 26 at 2.5 K, indicating a moderate electronic structure anisotropy. The measurements of Hall resistivity ρ xy reveal that the in-plane dominant carrier is hole-like with the carrier concentration and mobility estimated to be n h = 4.28 × 10 21 cm −3 and μ h = 3277 cm 2 V −1 s −1 at 2 K. The MR of BaGa 2 reaches 400% and shows no sign of SdH oscillation (Fig. 1e).
dHvA quantum oscillations. dHvA oscillations at various temperatures with the magnetic field along [00l] direction is shown in The dHvA oscillations can be described by the Lifshitz-Kosevich (LK) formula 38 where R T = (λTμ/B)/sinh(λTμ/B), R D = exp(−λT D μ/B) and λ = (2π 2 k B m 0 )/(ℏe). μ is the ratio of effective mass m * to free electron mass m 0 , and T D is the Dingle temperature. γ = 0, δ = 0 for a 2D system, and γ = 1/2, δ = ±/8 for a 3D system. β = ϕ B /2π and ϕ B is the Berry phase. The inset of Fig. 2c displays the temperaturedependent FFT amplitudes and fittings using the thermal factor R T in LK formula. The obtained effective cyclotron masses are quite small, comparable with those in topological semimetals, due to the almost linear dispersion though they all originate from trivial bands lately determined by the magnetotransport analysis and first-principles calculations. In order to extract the Berry phase conveniently, we applied four-band LK formula fitting directly. The fitting result is displayed in Fig. 2d with the violet dots being experimental results and the red line being the fitting curve. All of the Berry phase of these four pockets greatly deviate from the non-trivial value π implying that these bands may be trivial. Detailed data is exhibited in Table 1. The angle-dependent dHvA oscillation measurements are applied to further reveal the Fermi surface of BaGa 2 as exhibited in Fig. 3a. Figure 3b and c shows the FFT spectra of the corresponding dHvA oscillations with the B rotating from B//c to B//ab. The fundamental frequencies increase with the angle θ and vanish at θ = 90°, indicating the moderate electronic structure anisotropy of these trivial bands. If these trivial bands are highly 2D-like, these frequencies should vanish once the magnetic field rotates by a small angle, which is inconsistent with the experimental measurements that F α , F β , and F γ are still observable at θ = 56°. Besides, the twofold symmetry of polar plot MR ( Fig. 3d) with B always vertical to I at 2.5 K and the ratio of ρ xx (14 T, 0°)/ ρ xx (14 T, 90°) ≈ 5 also indicate a moderate electronic structure anisotropy in BaGa 2 .
First-principles calculations. The band structure and Fermi surfaces of BaGa 2 calculated with spin-orbit coupling (SOC) effect included are shown in Fig. 4. As exhibited in Fig. 4a, the Dirac cone locates exactly at K point (E F ) in the Brillouin zone (BZ), which is consistent with previous work 37    of k z dispersion, which is consistent with the angle-dependent dHvA oscillations (Fig. 3). According to the Onsager relation ARPES measurement. The electronic structure of BaGa 2 is measured by ARPES and the quasi-2D Dirac cone at the K point is revealed. Figure 5a shows a schematic drawing of the 3D BZ and its projected surface BZ onto the (00l) surface measured. The FS is shown in Fig. 5b. The most obvious feature is a nearly circular pocket around the K point. Besides this pocket, there are several much weaker features around the Γ and K points where several bands cross the Fermi energy. This can be seen more clearly in the ARPES spectra along the high-symmetry directions as shown in Fig. 5c. The experimental spectra show an overall agreement with the calculated dispersions at different k z values.
There also exist a few extra band dispersions near the Fermi energy (indicated by red circles in Fig. 5c) originating from surface states which is consistent with the first-principles calculations, as shown in Fig. 5d. Besides, the experimental dispersions which are away from Fermi energy show a mixture of calculated dispersions ranging from k z = 0 and k z = 0.5c * (the dispersions at the k z boundaries are colored black in Fig. 5d), suggesting a very strong k z -broadening 39 . We zoom in the dispersion near E F to further reveal the Dirac cone at the K point. The dispersion at k y = 0 is shown in Fig. 5f. Due to the matrix element effect 40 , the left branch is much weaker than the right branch. To enhance the visualization of the dispersion on both sides, we show in Fig. 5g the normalized spectra by integrating each energy distribution curve (EDC). In addition to the Dirac cone at the K point, there is an extra dispersion originating from the surface state, which is consistent with the calculations. The evolution of Dirac cone along the k y direction is shown in Fig. 5f. Moving away from the K point, the dispersion shifts down in energy and becomes more parabolic, in agreement with a Dirac cone at the K point as shown in Fig. 5a. The k z -dependence of the Dirac cone is shown in Fig. 5g. Similar behavior is also observed, with the dispersion shifting down in energy when moving away from K point (k z = 2c * ). However, the k z dispersion of the Dirac cone is much weaker than k y dispersion and the Fermi velocity along k z is almost one order of magnitude less than k y , suggesting that the Dirac cone is quasi-2D like centered at K point. Interlayer transport properties. BaGa 2 provides an ideal platform to study the unusual interlayer transport properties caused by the tunneling of Dirac fermions between the zeroth LLs since the Ga honeycomb-net layer contributes the quasi-2D Dirac cone and the Ba layer between adjacent two Ga honeycomb-net layers provides the barriers. Figure 6b presents the angle-dependent interlayer resistivity ρ zz (B, θ) at 2.5 K under different fields. It shows the interlayer resistivity peak at θ = 90°. θ is the angle between B and I as defined in Fig. 6a. Besides, the interlayer resistivity decreases with the increasing field, resulting in the NIMR at θ = 0°, as shown in Fig. 6b, d and e. The NIMR decreases quickly when temperature changes from 2.5 to 70 K. Both of these anomalous interlayer transport properties can be understood based on the tunneling of Dirac fermions between the zeroth LLs. According to previous works 34, 36 , the interlayer Dirac fermion tunneling conductivity σ LL0 t ðB; θÞ can be described as where A is a field-independent parameter, θ is the angle between B and I, and c is the distance of the adjacent two Ga layers. When the field is applied parallel to I along c-axis (θ = 0°), the tunneling conductivity is simplified as σ LL0 t ðB; 0 Þ = A B j j. In this case, the tunneling conductivity increases proportional to the magnetic field, leading to the NIMR as exhibited in Fig. 6d. When θ increases from 0°to 90°, the interlayer tunneling conductivity decreases gradually because of the suppression of 2D LL quantization, resulting in the interlayer resistivity peak as shown in Fig. 6b.
In an ideal 2D Dirac system, the interlayer transport is only contributed by the tunneling of Dirac fermions between the zeroth LLs (tunneling channel) and the resistivity exhibits the equivalent peak values at θ = 90°under different fields 34,35 . However, as shown in Fig. 6b, the interlayer resistivity peak increases with the increasing field. The inconsistency is reasonable since there also exist trivial bands contributing the positive interlayer MR except for the Dirac band. As discussed in the angle-dependent dHvA oscillations, the results from firstprinciples calculations, the ratio of ρ xx (14 T, 0°)/ρ xx (14 T, 90°) ≈ 5, and ρ zz /ρ xx ≈ 26, the electronic anisotropy of these trivial bands are analyzed as moderate. Based on these discussion, we assume the trivial bands contribute to the interlayer transport through a momentum relaxation mechanism (i.e. coherent band transport) and describe it with Drude model. Considering both the tunneling mechanism and momentum relaxation mechanism, the total interlayer resistivity is described as where σ c (B, θ) is the conductivity from the coherent trivial bands. As shown in the Supplementary Table 1, all of the trivial bands do not reach the quantum limit, so the conductivity σ c (B, θ) is taken as σ 0 =ð1 þ k Á B 2 xy Þ for low field approximation (σ 0 is the Drude conductivity, k is a constant, B xy = |B|sin θ). B 0 is a fitting parameter 34 . At θ = 90°, the tunneling channel is suppressed [σ LL0 t ðB; 90 Þ = 0] and the interlayer resistivity mainly originates from the trivial bands. In such case, the interlayer resistivity is described as ρ zz (B, 90°) = 1/σ zz (B, 90°) = 1/(σ c (B, 90°) + B 0 ). As shown in Fig. 6e the fitting curve matches well with the experimental data. Therefore, to a certain extent, the contribution from trivial bands to interlayer transport can be described by Drude model. The total interlayer transport properties can be considered as the result of the competition between the tunneling mechanism and momentum relaxation mechanism. When θ is close to 0°, the trivial band-induced positive MR is relatively small because the field-induced scattering is very weak (absence of Lorentz force). However, the density of states (DOS) at Dirac point increases dramatically with the increasing field since the Dirac cone is located exactly at E F . In this case, the tunneling mechanism dominants the interlayer transport, resulting in the NIMR as shown in Fig. 6b, d and e. The angle-dependent and fielddependent unusual interlayer resistivity can be well fitted by Eq. (3), as shown in Fig. 6c and e. When θ increases from 0°to 90°, the tunneling of Dirac fermions between the zeroth LLs is suppressed while the trivial band induced positive MR increases due to the enhancement of the Lorentz force-induced scattering. It is consistent with the observed interlayer MR, where the NIMR vanishes gradually and positive MR increases with θ varying from 0 ∘ to 90 ∘ . These resistivity curves are well fitted by Eq. (3) as shown in Fig. 6c and e. These anomalous interlayer transport properties induced by the tunneling of Dirac fermions between the zeroth LLs have also been observed in (BEDT-TTF) 2 I 3 , where the Dirac cone only exists under pressure 34,35 . YbMnBi 2 is another quasi-2D topological material that exhibits the same anomalous interlayer transport properties 36 , in which the NIMR was not observed possibly due to the large positive MR induced by other bands overwhelming the signature of NIMR. To our knowledge, BaGa 2 is the only material that exhibits both NIMR and interlayer resistivity peak at ambient pressure.

Discussions
To further understand the observed NIMR in BaGa 2 , we discuss on the possibly different origins of the negative magnetoresistance (NMR) below. BaGa 2 is a non-magnetic material, spin flip may not be a suitable reason for the NIMR in BaGa 2 41 . Secondly, the NMR has also been observed in topologically trivial materials at high field such as graphite, where the NMR is induced by the ellipsoidal Fermi surfaces approaching the quantum limit 42 . As displayed in the Supplementary Table 1, all of the trivial bands in BaGa 2 are not under the quantum limit with the field of 14 T. In addition, the NIMR in BaGa 2 appears once a quite small field is applied. Thus, the scenario of trivial bands under quantum limit does not apply in BaGa 2 . Thirdly, the NMR induced by chiral anomaly has been widely observed in 3D Dirac and Weyl semimetals 5,9,12,21,22 , which usually has the form of Δσ xx ∝ B 2 or B depending on the Fermi energy 43 . However, the Dirac cone in BaGa 2 opens a small energy gap (~10 meV) when SOC is considered and it will not evolve into Weyl cones when time reversal symmetry is broken (applying the magnetic field). So the chiral anomaly is absent in BaGa 2 when I//B. Besides, the in-plane NMR is not observed when the field is parallel to the current in BaGa 2 (Supplementary Fig. 1). Thus the scenario of chiral anomaly induced NMR does not apply either. Fourthly, current jetting can often cause the NMR 44 , which has been explained in detail in ref. 45 . Specifically, it is caused by the inhomogeneous spatial distribution of the current in the sample 22,45 . In order to eliminate the current jetting effect, the four-probe AC transport method is applied in the interlayer transport measurement and the current probe with silver paste almost contact fully on the samples' upper and lower surfaces ( Supplementary Fig. 2) to make the current flow cross the sample homogeneously. Besides, the NIMR in BaGa 2 decreases quickly with the increasing temperature, inconsistent with the current jetting effect 44 . Finally, the NIMR can arise in high-purity layered metal PdCoO 2 with the residual resistivity ρ xx ranging from about 10 to 40 nΩ cm in most single crystals 46  ℏω c (ε n = (n + 1/2)ℏω c + 2t c cos(k z c), t c is interlayer transfer integral) describing the Landau levels for quasi-2D Fermi surface under the magnetic field in the quasi-2D metal material PdCoO 2 may not be fulfilled here. Besides, the interlayer MR in PdCoO 2 increases quickly at low field and decreases at high field when field deviates from c-axis which is quite different from that of BaGa 2 . From the above, we conclude that the interlayer NMR and the interlayer resistivity peak result from the tunneling of Dirac fermions between the zeroth LLs. In summary, BaGa 2 is a multilayered Dirac semimetal with quasi-2D Dirac cone located at E F , which is confirmed by the first-principles calculations, angle-dependent dHvA oscillations and ARPES measurements. The unusual interlayer transport properties (NIMR and peak in angle-dependent interlayer resistivity) originating from the tunneling of Dirac fermions on the zeroth LLs are both observed in BaGa 2 under ambient pressure for the first time. Our findings highlight the unusual role of the zeroth LLs in interlayer magnetotransport and enrich the novel transport properties of topological materials at quantum limit.

Methods
Crystal growth and magnetotransport measurements. The single crystal BaGa 2 was grown by self-flux method. The ingredient ratio with Ba:Ga = 41.7:58.3 were put into the alumina crucible and sealed into a quartz ampoule. Then the ampoule was heated to 950°C and kept for 50 h to make sure the materials have melted and mixed thoroughly. The process of crystal growth was carried up in cooling at 3°C/h to 700°C. The BaGa 2 single crystal was obtained after centrifuging the excess flux at that temperature. The atomic composition of BaGa 2 was checked by energydispersive x-ray spectroscopy (EDS, Oxford X-Max 50). The crystal structure was determined by XRD in Bruker D8 Advance x-ray diffractometer. TOPAS-4.2 was employed for refinement. The electric transport and magnetic properties were collected in Quantum Design PPMS-14 T and MPMS-3 SQUID VSM system.
Angle-resolved photoemission spectroscopy. ARPES measurements were taken at BL13U of Hefei National Synchrotron Radiation Laboratory and our home laboratory. The crystals were cleaved in situ and measured at a temperature of T ≈ 15 K in vacuum with a base pressure better than 1 × 10 −10 torr.
Band structure calculations. The electronic structures of BaGa 2 were studied by using first-principles calculations. The projector augmented wave (PAW) method 47,48 as implemented in the VASP package [49][50][51] was used to describe the core electrons. The SCAN type exchange-correlation potential was adopted 52 . The kinetic energy cutoff of the plane-wave basis was set to be 400 eV. A 20 × 20 × 18 kpoint mesh was utilized for the BZ sampling and the Fermi surface was broadened by the Gaussian smearing method with a width of 0.05 eV. Both cell parameters and internal atomic positions were fully relaxed until the forces on all atoms were smaller than 0.01 eV/A. Once the equilibrium structures were obtained, the electronic structures were calculated with the SOC effect. The Fermi surfaces were studied by using the maximally localized Wannier functions (MLWF) method 53,54 . To study the surface states, a 1 × 1 two-dimensional supercell with an BaGa 2 slab containing 59 atoms and a 20 Å vacuum was employed to simulate the BaGa 2 (001) surface.