Valley-polarised tunnelling currents in bilayer graphene tunnelling transistors

We study theoretically the electron current across a monolayer graphene/hexagonal boron nitride/bilayer graphene tunnelling junction in an external magnetic field perpendicular to the layers. We show that change in effective tunnelling barrier width for electrons on different graphene layers of bilayer graphene, coupled with the fact that its Landau level wave functions are not equally distributed amongst the layers with a distribution that is reversed between the two valleys, lead to valley polarisation of the tunnelling current. We estimate that valley polarisation $\sim 70\%$ can be achieved in high quality devices at $B=1$ T. Moreover, we demonstrate that strong valley polarisation can be obtained both in the limit of strong-momentum conserving tunnelling and in lower quality devices where this constraint is lifted.


I. INTRODUCTION
Electron tunnelling through a potential barrier is one of the most widely known physical consequences of quantum mechanics, responsible for effects as varied as nuclear fusion in stars, radioactive decay or spontaneous DNA mutation [1]. In particular, the probability of successful tunnelling decays exponentially with the width of the barrier, an effect best visualized in scanning transmission microscopy where moving the conducting tip away from the sample leads to rapidly decaying tunnelling currents, hence allowing for imaging of the corrugation of the sample surface [2]. Recently, the limit of singleatomic-layer barrier thickness has been achieved in van der Waals (vdW) heterostructures using, first, atomically thin graphene and hexagonal boron nitride (hBN) [3,4], as the electrode and barrier respectively, and later other two-dimensional atomic crystals [5]. The resulting tunnelling transistors offer a solution to the graphene 'band gap problem' -the lack of a band gap in the conical electronic dispersion of the material [4]. Moreover, in devices with ultra-high quality interfaces, momentum-conserving tunnelling was demonstrated, leading to negative differential resistance [6][7][8][9][10][11] and valley polarization due to an in-plane magnetic field [12]. It was also shown that electron tunnelling in vdW heterostructures can be accompanied by excitation of various quasiparticles, for example, phonons [13] or magnons [14] and influenced by defects in the tunnel barrier [15,16]. Finally, moiré superlattice effects can be used to engineer the electronic densities of states of the electrodes [17][18][19][20].
Here, we study theoretically the tunnelling current flowing between bilayer graphene (BLG) and monolayer graphene (MLG) electrodes through a hBN barrier, in the presence of a magnetic field perpendicular to the atomic layers. The impact of the applied magnetic field is twofold: firstly, the electronic density of states is modified * J.J.P.Thompson@bath.ac.uk due to Landau quantisation and second, layer polarisation of the low-energy Landau levels in BLG [21] leads to the efficient generation of valley polarisation [22,23]. We show here that valley polarisation of order unity is possible, in magnetic fields as low as ∼ 1 T, and that choice of the valley quantum number of the tunnelling current can be made electronically without reversing or changing the magnitude of the magnetic field. While the largest valley polarisation can be achieved in high quality devices in which tunnelling electrons conserve both energy and momentum, our results suggest that even in the absence of momentum conservation, polarisation ∼ 70% can be achieved at B = 1 T.

II. DEVICE DESCRIPTION AND TUNNELLING MATRIX ELEMENT
The schematic of the device we study is shown in Fig.  1 with the assumed direction of the magnetic field, B. We consider BLG, two layers of carbon atoms arranged arXiv:1811.04687v1 [cond-mat.mes-hall] 12 Nov 2018 in regular hexagons and stacked in a Bernal (AB ) formation, adjacent to a few-layer hBN sheet with a monolayer graphene (MLG) sheet placed onto the opposite face, such that the hBN acts as a tunnelling barrier between the two graphene materials. The BLG and MLG electrodes can be rotationally misaligned so that their crystallographic directions are rotated by an angle θ. Following experimental device architectures [8,12,25], we assume that both MLG and BLG are encapsulated with hBN. Finally, a silicon back gate which enables doping of the tunnelling electrodes is located on the side of MLG, separated by SiO 2 , while a gold top gate is attached adjacent to the BLG electrode.
In order to find the tunnelling current through the device we use Bardeen's formalism [4,9,[26][27][28][29], which utilises the wave functions of the source and drain electrodes to model the tunnelling probability. The matrix element, M (ε), associated with the probability of an electron with energy ε tunnelling through the barrier (which we take to lie in the xy-plane), is calculated (up to some constant prefactor with dimension of energy × distance) as where Ψ S (Ψ D ) describes the wave function on the source (drain) electrode and the integration is over the volume of the tunnelling junction. In a multilayer electrode, such as bilayer graphene, the overall wave function can be described as a linear combination of the wave functions on the constituent layers. Therefore, in the case of an Nlayer source electrode we can write where ψ S,i is the electronic wave function of the i-th layer of the source electrode and |n i | 2 describes the relative occupation of the ith layer by an electron in state Ψ S . Assuming a clean sample, all ψ S,i (r, ε) are separable into the in-plane and perpendicular components, ψ S,i (r, ε) = ϕ S,i (x, y, ε)φ S,i (z, ε). This enables us to decompose the matrix element, Eq. (1), into transverse, z, and in-plane, x, y, components. We model the transverse component of the wave functions, φ S,i (z, ε), as exponentially decaying, taking care to change the decay constant in the regions corresponding to different materials. In our case, there are two materials to consider: hBN, comprising the tunnelling barrier, and graphene. We assign to them decay constants c(ε) and c (ε), respectively (the decay constants vary as a function of the tunnelling state energy ε [27]). With these assumptions, we integrate the expression in Eq. (1) in the direction transverse to the barrier. For our BLG source electrode which consists of two graphene layers, the electrons from the layer closest to the barrier tunnel directly to the drain, passing only through one material, hBN. However, the electrons in the layer further from the barrier have to travel an increased distance. Since there are no available states on the other graphene layer the electrons have to pass through, it can be effectively treated as an insulator and so the only mechanism for transport is tunnelling. Integrating over the total width of the barrier for these electrons, d + d 0 , we obtain an expression depending only on the in-plane components of the wave function and energy, where we have absorbed all normalisation factors into ϕ S,i (x, y, ε) and used the fact that in our device the drain electrode is built of only one layer. Note that the expression in Eq. (3) conserves the in-plane electron momentum in the tunnelling process, as is the case in experiments performed on the highest quality devices [6,8,9,12]. While the argument above can be extended to any number of layers in both the source and drain electrodes, the exponential dependence of tunnelling probability on the barrier width means that only tunnelling from/into the first few layers next to the barrier is measurable. For hBN, experimental works [3,4,9,12,30,31] suggest a value of the decay constant c(0) ≈ 5nm −1 . In the case of graphene, studies of its role as a barrier in magnetic tunnel junctions [32][33][34] and between metal contacts [35,36] showed that it behaves as a strong out-of-plane insulator. In fact, in experiments conducted in the absence of a magnetic field and in the presence of a field parallel to the graphene layers, the measured tunnelling current has been well described by assuming that all tunnelling from the further BLG layer is suppressed [10,12]. For this reason, here, we take the limit c = c, corresponding to the decay through graphene being significant and similar to that through hexagonal boron nitride. However our conclusions hold for notably smaller c (we discuss what happens for differing estimates of c in Appendix A).

III. WAVE FUNCTIONS OF GRAPHENE ELECTRODES
In order to obtain the wave functions of electrons in our BLG and MLG electrodes, we use the low-energy description for electrons in these materials, applicable in the vicinity of the inequivalent Brillouin zone corners (valleys) K ξ = ξ(4π/3a, 0), where ξ = ±1 and the graphene lattice constant a = 2.46Å. A single graphene layer consists of two sublattices, A and B, and in the case of BLG, an effective model can be constructed using Bloch states, φ(A1) and φ(B2), formed from p z -orbitals on the nondimer sites (sites which do not have a neighbour directly above/below them) [21] which we refer to as A1 and B2 with the labels 1 and 2 corresponding to the layer closer and further from the barrier, respectively. For an electron with momentum p = (p x , p y ) measured from the centre of valley K ξ , the resulting Hamiltonian, written in the basis Above, the velocity, v ≈ 10 6 ms −1 , is related to the inplane nearest-neighbour hopping while γ 1 ≈ 0.38 eV is the vertical interlayer coupling. The termĤ u captures the effect of energy difference between sites on different layers, u (interlayer asymmetry), due to the electric field perpendicular to the BLG electrode induced by the applied voltages. Misalignment between the MLG and BLG electrodes, generated by a small clockwise rotation of the MLG sheet about the z-axis by angle θ , leads to an identical rotation between the corresponding Brillouin zones. As a result of this rotation, the position of MLG valley centres is offset from that of BLG by a vector ∆K ξ = ∆K x ξ , ∆K y ξ = (1 −R θ )K ξ . Taking into account this shift as well as the rotation between the two materials, electrons in the MLG electrode are described by a Hamiltonian which acts on the basis {φ(A), φ(B)} T for the K + valley and {φ(B), −φ(A)} T in the K − valley and π = (∆K x ξ + i∆K y ξ ). We include the magnetic field B applied perpendicular to the graphene planes in Eq. (4) and (5) by using the Peierls substitution, p → p + eA, and the Landau gauge A = (0, −Bx, 0). As a result, the operators π and π † become lowering and raising operators, respectively, for functions built of quantum harmonic oscillator states along the x-direction, φ m (x), and plane waves along the y-direction, where λ B = /eB is the magnetic length, H m (x) is a Hermite polynomial of order m, L is the dimension of the flake along the y-direction and X = λ 2 B k y is the position of the centre of cyclotron orbit of an electron with wave vector k y . Using Eq. (6), the energies of the BLG Landau levels, can be expressed as [21] ε 0,ξ = ξu 2 , where C = C 2 1 + C 2 2 is the normalisation constant. In rotated monolayer graphene, the displacement of the momentum-origin modifies the harmonic oscillator state, and shifts the cyclotron orbits toX = λ 2 B (k y + ∆K y ξ ). The resulting energy levels in MLG are where s = ±1 is the conduction/valence band index in MLG, defined when n ≥ 1. Therefore, the wave function corresponding to the n-th Landau level in rotated MLG can be written as In MLG, all the Landau levels have an additional fourfold degeneracy due to spin and valley. Moreover, the n = 0 Landau level is positioned at the Dirac point and its energy does not depend on the magnetic field. In BLG, for u = 0, in addition to the valley and spin degeneracies of each level, both the m = 0 and m = 1 Landau levels sit at the neutrality point, leading to an unusual eightfold degenerate zero-energy state. Non-zero u lifts both the m = 0, 1 and valley degeneracies.
In both MLG and BLG, the wave functions are distributed asymmetrically between the two sublattices concerned (A, B in MLG and A1, B2 in BLG). In particular, the electrons in the n = 0 MLG level and m = 0, 1 BLG states occupy only one of the sublattices [38]. For the case of BLG, this results in two states that in the K + valley are located only on layer 1 and in the K − valley only on layer 2. In the vertical tunnelling transistor geometry like in Fig. 1, we expect electrons from layer 2 to have a smaller chance of tunnelling through the barrier than electrons from layer 1, due to the additional effective barrier thickness. As a result, more electrons from BLG K + valley will tunnel through than from the K − valley, leading to valley-polarized current arriving in the MLG drain electrode.
In order to quantify valley polarisation of the tunnelling current, we use the wave functions from Eq. (8) and (11) to compute the tunnelling matrix element, Eq. (3) (see Appendix B for more details). The Landau level wave functions of BLG are already written so that their components correspond to different layers and we can identify the states in Eq. (8) with (ϕ S,1 (x, y, ε), ϕ S,2 (x, y, ε)) T in Eq. (3). However, although the MLG wave functions are also written as spinors in Eq. (11), both of their components correspond to wave function amplitudes on sublattices in the same layer. Hence, for a given Landau level state we take In clean samples, with a small misalignment angle between the source and drain electrodes the shift ∆K ξ is small and for all the dominant processes the valley quantum number is conserved in the tunnelling process [9,12]. Therefore, we use Fermi's golden rule to relate the tunnelling matrix element to current of electrons originating in the K ξ valley of BLG, where we have already taken the spin degeneracy into account. We measure the energy ε from the charge neutrality point of the BLG electrode while µ BLG and µ MLG represent the distance in energy between the charge neutrality point and the chemical potential in the BLG and MLG electrode, respectively. Finally, we define ∆ as the shift between the source and drain neutrality points such that, in the low temperature limit, the local chemical potentials in the source and drain electrodes, µ BLG and µ MLG + ∆ respectively, determine the energy window within which tunnelling processes can occur, while the number of the initial and final states at a given energy is provided by the densities of states D MLG and D BLG in the monolayer and bilayer graphene, respectively. In a device with high quality layers, free from defects and in a quantizing external magnetic field, these densities of states consist of a series of sharp peaks at the energies of the Landau levels. We model the latter using a Lorentzian shape with the same full width at half maximum for all Landau levels, 2 meV and 4 meV for B = 1 T and B = 4 T, respectively [9,39]. Finally, we define the valley polarisation, P , of the tunnelling current, Because our tunnelling matrix element, Eq. (1), is defined up to a proportionality constant, values of tunnelling current in this paper are given in arbitrary units. Polarisation, however, as a ratio of currents, does not depend on that constant itself.
We set the thickness of the hBN barrier (separation between the MLG and BLG) to d = 13Å and the interlayer distance in BLG as d 0 = 3.3Å. We also relate the energies ε, µ BLG and µ MLG to the applied voltages V t , V b and V g (see Fig. 1) through the electrostatic equations discussed in more detail in Appendix C. We define n MLG , n BLG and n Au as the carrier densities on the MLG, BLG and gold electrodes respectively. The distance between the gold top gate and BLG, d Top , is set as 30 nm in our numerical calculations. Furthermore, d hBN and d SiO2 represent the thicknesses of the hBN and SiO 2 substrates, which, following previous experimental works, we set as 30 nm and 300 nm respectively. Finally, 0 is the permittivity of free space while hBN ≈ 3 and SiO2 ≈ 3.9 are the relative permittivities of hBN and SiO 2 . We also take into account that the electric field between the graphene layers of BLG induces the interlayer asymmetry u which we compute self-consistently, The quantities n BLG,1 (u) and n BLG,2 (u) describe the carrier density on the two layers of BLG, which are calculated by taking into account the wave function distribution, for each Landau level, at a given u.

IV. MOMENTUM-CONSERVING TUNNELLING
A. Total Tunnelling Current at B = 1 T Our simulation of the total tunnelling current, I = I + + I − , between the BLG and MLG electrodes, produced using Eq. (12), is shown in Fig. 2. For momentumconserving tunnelling, the strength of the coupling is dictated by the magnitude of the applied magnetic field, relative orientation of the electrodes and Landau level indices of the involved electronic states.
In panel (a), we show the current for V t = 0 V and ideally aligned electrodes, θ = 0 • . For V g = V b = V t = 0 V, the chemical potentials in the BLG and MLG electrodes, µ BLG and µ MLG + ∆, are located at their respective neutrality points, which are at the same energy, resulting in zero tunnelling current. As the bias voltage is increased, a shift between the local chemical potentials in the two electrodes is induced. This opens an energy window, within which electrons occupying states in one electrode can tunnel into empty states at the same energy in the other electrode, thus leading to a non-zero tunnelling current. The coupling strength between the initial and final states in the tunnelling process is set by |M s,s n,m,ξ (ε)| 2 which, because of the spinorial nature of the MLG and BLG wave functions, is expressed as a sum of four terms, each of which contains an integral φ * n φ m dA of two oscillator statesφ n and φ m . For ideal alignment of the electrodes, θ = 0 • , the set {φ n } is equivalent to {φ m } and the integrals express orthonormality of functions with different indices, φ * n φ m dA = δ n,m . As a result, tunnelling only occurs if one of the four conditions is fulfilled: i) n = m, ii) n − 1 = m, iii) n = m − 2 and iv) n − 1 = m − 2. Hence, the central region of large current in Fig. 2(a), for V g = 0 and non-zero V b , corresponds to the coupling between m = 0 and n = 0 Landau levels in BLG and MLG respectively. Although at low voltages the m = 0 and m = 1 Landau levels in BLG are degenerate, transitions between the m = 1 and n = 0 level in MLG are forbidden. Although increasing V b increases the size of the tunnelling energy window to include higher Landau levels, due to the selection rules for θ = 0, these do not contribute to the tunnelling current.
Setting non-zero V g , at constant V b dopes the graphene electrodes, shifting the two chemical potentials together such that the difference between them remains unchanged. At small V b and zero V g in Fig. 2(a) clear current is observed. However, as V g is increased (decreased), the electrodes become hole-doped (electron-doped) and the tunnelling energy window, set by V b , moves away from the positions of the m = 0 and n = 0 Landau levels and the current decreases. This behaviour results in the horizontal, finger-like regions of large current visible Fig. 2 (a). Additionally, V b and V g induce an electric field between the graphene layers in BLG which leads to non-zero interlayer asymmetry u. This opens a band Finally, the current diagram as shown in panel (a) seems to have inversion antisymmetry with respect to the point (V b , V g ) = (0, 0), I(V b , V g ) = −I(−V b , −V g ). In fact, within the voltage window presented in Fig. 2, this antisymmetry is only weakly broken by the energy dependence of the decay coefficient c(ε) [see Eq. (3)]a feature also observed experimentally [3,4]. We investigate this symmetry in more detail in panel (e) where we present current plots for changing V b and constant V g from -9 V to 9 V in the steps of 3 V corresponding to solid/dashed lines marked in (a). We show in solid lines current for negative V g and V b and in dashed lines current for positive V b and V g . The same colour is used for curves of the same magnitude of V g and, for all V g , we have I(V g , V b ) almost equal to −I(−V g , −V b ).
In Fig. 2(b), we show the tunnelling current as a function of V b and V g for V t = 0.5 V. Non-zero V t induces interlayer asymmetry u even for V b = V g = 0 V while also introducing a shift between the MLG and BLG neutrality points, ∆. The latter leads to energy misalignment of the m = n = 0. However, because in BLG the position of the m = 0 Landau level depends linearly on u [see Eq. (7)], the impact of V t can be counterbalanced by choosing V g such that the overall u shifts the m = 0 BLG Landau level in the K + valley back into alignment with the n = 0 Landau level in MLG. This restores the finger-like feature, in panel (b) visible for some positive V g .
In Fig. 2(c) and (d), we show the impact of misalignment θ = 0.25 • and θ = 0.5 • , respectively, between the two electrodes on the tunnelling current. For non-zero θ, the oscillator functionsφ n and φ m are no longer orthonormal and transitions between any pair of states are allowed. Moreover, as shown in Appendix B, the coupling strength between states also depends on the misalignment angle. For these reasons, in panels (c) and (d) the finger-like feature present in (a) and (b) becomes increasingly smeared out with increasing θ and the tunnelling current also decreases as compared to (a) and (b). Additionally, misalignment between the electrodes breaks the approximate inversion antisymmetry of the current diagram in (a). In the presence of θ = 0, each of the four terms of the kind φ n φ m dA appearing in the calculation of the matrix element M s,s n,m,ξ (ε) (see also Appendix B) comes with a prefactor that depends on the MLG and BLG band indices s and s. Upon inversion from I(V g , V b ) to I(−V g , −V b ), the interference between these terms leads to different results depending on whether the initial and final states originate in the conduction or valence band. As a consequence, the approximate inversion antisymmetry about (V g , V b ) = (0, 0), present in Fig. 2(a) is strongly broken in both (c) and (d). This is demonstrated in more detail in Fig. 2(f), where we show similar current curves as in panel (e) (changing V b for constant V g from -9 V to 9 V in the steps of 3 V with the same colour scheme) produced for θ = 0.25 • (the cuts are also indicated in panel (c)). In particular, the magnitude of the current for V g = 0 V is much larger for As the Landau level wave functions of BLG are not distributed equally between its two constituent graphene layers and this distribution is reversed between the valleys, tunnelling in the device shown in Fig. 1 can be used to produce unequal electron occupations in the MLG drain electrode. Such an effect can be characterised by the valley polarisation, P , of the tunnelling current, introduced in Eq. (13), which we plot in Fig. 3 as a function of the gate voltages V g and V b for B = 1 T. Panels (a)-(d) in this figure correspond to the same parameters for which we presented total tunnelling current in Fig.  2(a)-(d).
In panel (a), we show the case of θ = 0 • and V t = 0 V. The polarisation diagram has inversion symmetry with respect to (V b , V g ) = (0, 0). A bright red cross-like feature corresponds to P ∼ 50% with a region of P ∼ 80% in the centre of the diagram. This high valley polarisation is due to the tunnelling between the m = 0 BLG and n = 0 MLG Landau levels. In the BLG K + valley, electrons in the m = 0 state occupy exclusively the layer closer to the barrier, whereas in the K − valley all of them sit on the layer further from the barrier. Consequently, the current in the K + valley is significantly larger than in the K − valley.
As V g is increased, the m = 0 and n = 0 Landau levels move out of alignment and the dominant source of tunnelling current becomes the m = 2 to n = 0 transition. From Eq. (8), the BLG m = 2 wave function is ψ ξ 2,s ∝ (C 1 φ 2 , C 2 φ 0 ) T , where C 1 and C 2 are complex numbers, and the first and second components of ψ ξ 2,s are located respectively on layer 1 (layer 2) and layer 2 (layer 1) in valley K + (K − ). For the MLG n = 0, ψ = (φ 0 , 0), so that for θ = 0 tunnelling in K + is only possible for BLG electrons from layer 2, further from the barrier, while in K − it is the electrons from layer 1 that can tunnel into MLG. As a result, the overall current has negative (K − ) polarization as shown by dark blue regions above and below the central red cross in Fig. 3(a). Similar arguments can be used to explain other regions of the polarisation diagram.
For non-zero top gate voltage, as shown in Fig. 3(b) for V t = 0.5 V, the polarisation map is modified as a result of the shift between the neutrality points, ∆, as well as non-zero interlayer asymmetry, u, at (V b , V g )=(0,0). The latter lifts the valley degeneracy of the BLG m = 0 Landau level. Alignment of the K + BLG m = 0 and MLG n = 0 states, responsible for the red cross-like feature in (a), now requires compensating with positive gate voltage. However, for negative V b it is not possible to achieve both alignment of these two states and position the BLG and MLG chemical potentials such that the aligned states contribute to the current. As a result, the left arm of the red cross disappears and the m = 2 to n = 0 transitions lead to negative polarisation in this region.
In Fig. 3(c) and (d), we show valley polarisation as a function of V b and V g for increasing misalignment between the electrodes, θ = 0.25 • and θ = 0.5 • corresponding to total current plots in Fig. 2(c) and (d). Similarly to the current features, when the graphene electrodes are misaligned, individual polarisation features become smeared out and the variation of polarisation throughout the (V b , V g )-space becomes more gradual. The oscillator states φ m andφ n are not orthonormal for θ = 0 • so that many different transitions contribute to the overall polarisation for given (V b , V g ). Importantly, interference of electronic states tunnelling between any of the BLG layers and any of the MLG sublattices which leads to different outcomes for conduction band-conduction band and valence band-valence band transitions, strongly breaks the inversion symmetry of polarisation present in (a) for θ = 0 • . This symmetry breaking grows with increasing θ.

C. Tunnelling at B = 4 T
The Landau level structures in BLG and MLG depend on the strength of the magnetic field differently, hence the tunnelling current and polarisation features in the (V b , V g ) diagrams depend on B. For this reason, to contrast our results for B = 1 T presented in Fig. 2 and 3 with the case of stronger magnetic field, in Fig. 4 we show the tunnelling current and its valley polarisation for B = 4 T, θ = 0 • and V t = 0 V. Due to the increased electron density per Landau level at B = 4 T, it is necessary to increase the voltage range in order to compare features arising from similar electronic tunnelling transitions. Similarly to the B = 1 T case, a finger-like structure is present in Fig. 4(a). In fact, it is more pronounced because the separation between the m = n = 0 Landau levels and the rest of the electronic spectra in the corresponding materials increased. Consequently, the central cross-like region of K + -polarised current is also sharper, including polarisation of P ∼ 90% in the vicinity of (V b , V g ) = (0, 0). The maximum K − polarisation in the blue region dominated by m = 2 to n = 0 tunnelling is also increased.

V. TUNNELLING WITH STRONG MOMENTUM SCATTERING
In the presence of a poor interface between the electrodes and the hBN barrier, the scattering length-scale becomes very small such that the momentum resolution of the tunnelling electron becomes lost. In this limit, the momentum-nature of the initial and final state has no effect on the magnitude and valley polarisation of the current across the device. Instead, the tunnelling current depends only on the density of states of the source and drain electrode. As a consequence, we expect the valley polarisation of the tunnelling current to arise purely due to differences in valley occupations of the two BLG layers. We model this regime by setting each harmonic oscillator integral, φ n φ m dA, equal to 1 for all the transitions independently of their initial and final states [43] and present our results for θ = 0 • and V t = 0 V in Fig.  5.
In panels (a) and (b), we show the tunnelling current for B = 1 T and B = 4 T, respectively. Because all of the transitions are now allowed, the graphs look similar to that in Fig. 2(d), corresponding to θ = 0.5 • . In (b), the increased magnitude of the magnetic field leads to larger spacing between the Landau levels so that, as compared to (a), a larger voltage window is necessary to capture features due to transitions between the same pair of Landau levels.
In Fig. 5(c) and (d), we present valley polarisation of the currents shown in Fig. 5(a) and (b), respectively. The relaxation of the selection rules discussed in the previous section results in polarisation maps which are heavily weighted in favour K + valley. In particular, the maximum valley polarisation occurs at low voltages, where the participating Landau levels are those with low index (in particular m = 0 and m = 1). This is because, for the m = 0 and m = 1 BLG Landau levels, the valley and layer degrees of freedom are coupled. Furthermore, the interlayer asymmetry, u, generates a layer population difference in the m ≥ 2 Landau levels in BLG, which is opposite in the two valleys. This induced interlayer asymmetry is responsible for small regions of minor K − polarisation which occur at higher voltages. These two principles are responsible for all polarisation features observed in Fig. 5(c) and (d).
The relative misalignment of the graphene electrodes has no effect on the tunnelling probability in this limit. Similarly to the case of momentum-conserving tunnelling between misaligned electrodes, θ = 0 • , lack of restrictions on allowed transitions leads to chiral interference. This interference results in the asymmetry in inversion about the origin that is observed in Fig. 5. We expect momentum non-conserving tunnelling to be the dominant mechanism when the misalignment angle between the graphene electrodes is large. This is because, while increasing misalignment angle decreases the magnitude of the momentum-conserving current, it should have no effect on the transitions involving scattering.

VI. SUMMARY
We have explored the tunnelling characteristics of a vertical field effect transistor comprising monolayer and bilayer graphene electrodes, in the presence of a perpendicular magnetic field. The coupled layer and valley polarisation in the Landau levels of bilayer graphene gives rise to a valley-polarised tunnelling current through the device resulting in unequal valley populations in monolayer graphene. Importantly, this valley polarisation can be tuned solely by electrostatic means without the need to reverse the direction of the magnetic field. Our modelling suggests that P ∼ ±70% is possible in high quality devices in homogeneous fields of B = 1 T. Fields of such magnitude could be in principle generated by placing ferromagnets on top of the device [44]. While the homogeneity of the field distribution across the device would then depend on the size of the ferromagnet, thickness of the tunnelling junction and distance between the two, valley polarisation might still be possible in such a setup. In order to detect the valley polarisation produced using the proposed device, two of these could be connected in series: the first one to produce unequal valley populations and the second to act as a detector. Alternatively, the produced valley polarisation can be measured using optical means [45]. In the main text, we take the decay constant c characterising tunnelling through monolayer graphene to be equivalent to the decay constant c corresponding to tunnelling through hBN. Here we discuss the effect of significantly reducing the decay constant, c , on our results. Following previous work [4,41], we relate the decay constant to the height of the tunnelling barrier. For the case of hBN, we treat it as an isotropic potential step with barrier height Φ 0 = −1.5 eV, corresponding to experimental measurements of the valence band maximum (VBM) of hBN [4,5]. The hBN energy dispersion around the VBM is roughly parabolic in k z and this allows us to write [41] c(ε) = Im where m * ≈ 0.5m 0 is the effective mass [4,5]. Notably, this relation predicts weak electron-hole asymmetry in tunnelling current as observed in experiment [4,5]. In our work, we use the expression in Eq. (A1) to obtain the decay constant c for the tunnelling of BLG electrons from the layer further from the barrier across the layer closer to the barrier. While both the barrier height and the effective mass would be different for graphene as compared to hBN, our main conclusions are quite insensitive to the details of the functional form of c (ε). For example, in Fig. (6), we demonstrate the valley polarisation at B = 1 T and θ = 0 • for barrier height decreased to (a) 0.5Φ 0 and (b) 0.25Φ 0 (while using the same effective mass . While the asymmetry between P (V b , V g ) and P (−V b , −V g ) increases for smaller barrier height, qualitative features of the valley polarisation graphs remain the same. Also, the maximum valley polarisations are still significant, 58% and 48% respectively, compared to 70% in Fig. 3(a).

Appendix B: Landau Level Couplings and Matrix Element in the Momentum-Conserving Limit
The matrix element determining the tunnelling between the BLG and MLG electrodes depends on the Landau level indices, n, m, as well as the magnetic field, B, and misalignment angle between the two sheets, θ.
In both cases we define 6(a) , increasing misalignment redistributes the coupling strength amongst other transitions (in particular, for any non-zero θ transitions between any two oscillator states are in principle allowed). Interestingly, changing the misalignment angle also changes the preferred transition (the one with the largest coupling strength). However, because the matrix element depends on products of the type ∆K 2 ξ λ 2 B , a change of angle (which determines ∆K ξ ) can be to some extent counterbalanced by changing the magnetic field (and hence λ B ).

Appendix C: Electrostatics
The bias voltage V b and gate voltages V g and V t control the local Fermi levels in BLG and MLG as well as the shift between the neutrality points and the interlayer asymmetry in BLG. The Fermi levels and vertical shift of the energy bands would be controlled via a bias V b and gate voltages V g,t in an experimental setup. We use a four-plate capacitor model to express the electric fields between the gates and consecutive graphene layers (we treat hBN and SiO 2 as homogeneous insulators with dielectric constants hBN and SiO2 , respectively). The carrier densities per graphene layer on the BLG source (j =BLG) or MLG drain (j =MLG) electrodes can be expressed as Through charge conservation, for each combination of µ BLG , µ MLG and u, we obtain corresponding bias, bottom and top gate potentials. Furthermore, the charge build up on the bilayer sheet acts as a capacitance leading to a difference in neutrality points of the two spectra For simplicity we set E BLG 0 = 0 and therefore the position of the charge neutrality point in monolayer is related to the number of excess charge carriers in the bilayer graphene by relation The voltages of the system can therefore be expressed as