Valley splitter and transverse valley focusing in twisted bilayer graphene

We study transport through electrostatic barriers in twisted bilayer graphene and show that for certain configurations, electrons from the $K$ ($K'$) valley are transmitted only to the top (bottom) layer, leading to valley-layer locked bulk currents. We show that such a valley splitter is obtained when the potential varies slowly on the moir\'e scale and the Fermi energy in the barrier exceeds the kinematic barrier between Dirac electrons from the top and bottom layer. Furthermore, we show that for a given valley the current is transversely deflected, as time-reversal symmetry is broken in each valley separately, resulting in valley-selective transverse focusing at zero magnetic field.

Introduction -Twisted bilayer graphene (TBG) hosts a rich phenomenology ranging from its single-particle properties  to the recent discoveries of correlated insulating phases [24,25], superconductivity [26,27], nematicity [28,29], and ferromagnetism [30]. In TBG, two graphene layers are stacked at an angle that differs from Bernal (AB) or hexagonal (AA) stacking, giving rise to a triangular moiré pattern of alternating stacking regions. For small twist angles, the moiré pattern varies slowly relative to the interatomic scale and the low-energy electronic properties are significantly modified [1,4]. Its hallmarks are the narrowing of the bands near charge neutrality with decreasing twist angles and the concomitant reduction of the Fermi velocity [10,31], vanishing at the so-called magic angle, and emergence of low-energy Van Hove singularities (VHSs) [3,11,13,14].
In this paper, we argue that TBG for small twist angles above the magic angle is also an ideal platform for valleytronics [32] where the two twisted layers provide a natural setting to separate valleys. Contrary to valleytronics proposals in monolayer graphene which rely on specific edge termination or defects [33][34][35], we show that properly oriented electrostatic barriers in bulk TBG lead to valley-layer locked currents when the Fermi energy in the barrier exceeds the kinematic barrier between Dirac electrons localized on different layers. Consequently, one can design a valley valve [33,34] by putting two such valley splitters in series. We also show that the current is transversely deflected and in opposite directions for the two valleys, leading to transverse valley currents as time-reversal symmetry is effectively broken within a single valley of TBG [36]. The magnitude of deflection can be controlled by the slope of the potential and hence this effect, which we call transverse valley focusing, can be thought of as an electronic counterpart of transverse magnetic focusing in graphene [37].
We first present the basic idea of the valley splitter using semiclassical arguments. We then solve the scattering problem with a two-band lattice model for TBG [38][39][40] and discuss signatures of the valley splitter in the two-terminal conductance. Finally, we consider a multi- terminal setup to demonstrate the transverse deflection of the valley-layer locked current.
Twisted bilayer graphene -Twisted bilayer graphene can be obtained by rotating the top (bottom) layer of two layers of graphene in perfect registry, over an angle θ/2 (−θ/2), about a common hexagon center, where θ is the twist angle [1]. This results in a moiré pattern of alternating stacking regions with moiré lattice constant L m = a/[2 sin(θ/2)] where a ≈ 2.46Å is the graphene lattice constant, as shown in Fig. 1(a). For small θ, i.e. L m a, the low-energy physics becomes independent of twist center or commensuration due to approximate symmetries of the moiré pattern that are effectively exact above µeV energy scales [4,36,40,41]. Hence, it is sufficient to consider the emergent moiré lattice in the limit of small twist angles, which we do in this work.
In Fig. 1(a), the primitive moiré vectors are given by L 1,2 = L m √ 3/2, ±1/2 and the corresponding reciprocal basis G 1,2 is defined by G i · L j = 2πδ ij . Since the moiré pattern varies slowly (∼ a/θ) relative to a for small twists, interlayer scattering between valleys is suppressed. Above the first magic angle, the low-energy physics is then well described by the continuum model for TBG that treats the valleys, which are related to each other by time-reversal symmetry, independently [1,4,15,36]. In momentum space, shown in Fig. 1(b), the Dirac points of the rotated layers are given by K ( ) l = R ±θ/2 K ( ) for layer l = t, b, where R θ is the standard rotation matrix, K ( ) = τ |K|e x are the unrotated Dirac points with |K| = 4π/(3a), and τ = ± for the K and K valley, respectively. We choose the moiré Brillouin zone (MBZ) such that K ( ) l sit at the corners, and refer to these low-energy regions as minivalleys K/K = K t /K b or K b /K t for the other valley, where the bar indicates points in the MBZ, and which are separated in momentum space by k θ = 2|K| sin(θ/2). In Fig. 2, we show the low-energy bands calculated with the continuum model for θ = 2 • (details in Supplemental Material (SM) [42]).
Valley splitter and transverse valley focusing -Consider scattering in TBG at an n-n -n barrier, as illustrated in the top panel of Fig. 3. Here, the Fermi energy is close to charge neutrality in the outer regions, while inside the barrier it is above EM , the energy of the lowest conduction band at theM point, i.e. the midway point between Dirac points of the top and bottom layer. Here, we consider a potential U (y) in the zigzag direction of the honeycomb lattice dual to the moiré lattice, given by the y-direction in Fig. 1(a), such that inter-minivalley scattering is allowed. Now suppose there is a small bias voltage, such that electrons incident on the barrier come from the left. When the potential varies slowly on the scale of λ F ∝ E −1 F , scattering can be treated semiclassically as adiabatic transport of modes in momentum space. In this case, we expect that incident modes from K b are always reflected back to K t as they pass through a local extremum, i.e. a classical turning point. Moreover, because the minivalleys correspond to the layers at low energies, the transmitted current runs mostly in a single layer. Up to this point, we only considered the K valley, which is allowed since the valleys are decoupled for small twist angles and we assume the potential varies slowly on the interatomic scale. Since the valleys are related by time reversal (Fig. 1(b)), the total transmitted current consists of a K-polarized current in the top layer and a K -polarized current in the bottom layer. The valley polarization of the layers can be reversed by reversing the bias voltage. Furthermore, we find that incident modes at K b that are reflected to K t are transversely deflected, by inspection of the transverse velocity, shown in Fig.  2(b), and the classical trajectories in Fig. 4(a). It can be understood as a consequence of broken time-reversal symmetry in a single valley of TBG, leading to transverse valley currents and valley-dependent deflection.
The range of twist angles for which the valley splitter works is determined by the required carrier density nM ≈ where m is a band index, k is the crystal momentum constrained to the first MBZ, V is the sample area, G = n 1 G 1 + n 2 G 2 (n 1,2 ∈ Z) is a reciprocal lattice vector of the moiré lattice, and ψ l are two-component spinors containing the amplitudes on the A/B sublattices of layer l = t, b. The layer polarization is then given by and is shown in Fig. 2(a) for θ = 2 • . At the Dirac points, P L (K/K ) = ±0.32 for the K valley, so that (1±P L )/2 = 66% of the density remains localized in its layer, which increases up to 82% for θ = 3 • . Hence, the valley splitter regime lies roughly in the range 2 • < θ < 10 • , where the lower (upper) limit is determined by P L (nM ). Devices with twist angles in this range have already been realized in several experiments [3,11,[17][18][19].
Model -To study electronic transport, we use an effective two-band lattice model for TBG [38,39], defined on a honeycomb lattice with sublatticesĀ/B corresponding to AB/BA stacking regions ( Fig. 1(a)). This model describes the two (spin degenerate) bands near charge neutrality, for each valley independently, that are separated from the other bands by a minigap, as shown in Fig. 2(a). We then consider transport in the zigzag direction of the effective honeycomb lattice and take periodic boundary conditions in the armchair direction. We obtain a one-dimensional (1D) chain for each transverse mo- withĥ − (k ) =ĥ + (−k ) * by time-reversal symmetry.
Here,ĉ † nσ (ĉ nσ ) creates (destroys) a particle on sublattice σ of the nth cell of the chain, U n is the onsite energy, σ = ± forĀ/B, andσ =B,Ā for σ =Ā,B, respectively. The nearest-neighbor (nn) hopping amplitude t 1 is chosen real and positive: for θ = 2 • , we find t 1 ≈ 33 meV using the method of Ref. 39. Here, we only show the nn hopping term explicitly (see SM [42]). The potential barrier is modeled , where f (y) = tanh(y), d determines the variation of the potential in units L m /2, and the length of the barrier is given by L = (L m /2) (n r − n l ). Henceforth, we focus on θ = 2 • , for which L m ≈ 7.05 nm, and the K valley.
Results and discussion -The 1D scattering problem was solved with Kwant [43] and the resulting zero-bias differential conductance is shown in the middle panel of Fig. 3 as a function of the Fermi energy in the barrier E F − U 0 . When the barrier varies slowly on the moiré scale (d > 1) and E F − U 0 > EM , the conductance reaches a plateau of G 0 /2, where G 0 is the ballistic conductance per valley. This case corresponds to the schematic in the top panel of Fig. 3: the barrier is transparent for one minivalley and opaque for the other, which is corroborated by the minivalley polarization where T l←l is the transmission function for scattering from minivalley l to l, and which is shown in the bottom panel of Fig. 3. The beating pattern in the conductance is due to a mismatch of oscillations for modes at ±k due to time-reversal symmetry breaking within a single valley, which displaces the VHS away fromM , and hence observation of such beatings is a signature of a type-II VHS in TBG [44][45][46]. In general, interference oscillations are more prominent in the bipolar regime (E F − U 0 < 0) due to reflections at n-p junctions on both sides of the potential barrier [47,48]. For the same reason, the minivalley filter does not work perfectly for E F − U 0 < −EM and its efficiency decreases for the slow-varying potential barrier (d = 7, blue line). This can be understood as follows. When the barrier is also smooth on the scale of λ F ≈ 3.6L m , transmission through the n-p junctions is suppressed away from normal incidence [49]. Hence, the filter in the bipolar regime works better for barriers with intermediate smoothness (d ≈ 1, red line). Since a smooth barrier in the filter regime only transmits electrons from one minivalley for a given valley, one can try to design a valley valve [33,34] by placing two such barriers in series. When the barriers are gated similarly, one minivalley can go through. For opposite polarities, however, electrons that pass the first barrier are reflected and change minivalley in the second barrier, at which point they can pass the first barrier from the other side. More detailed discussions of the minivalley filter and valley valve are given in the SM [42].
Finally, to demonstrate the transverse deflection of the reflected current, we consider a multi-terminal setup of length L = 1 µm and width W = 3 µm, as shown in the center inset of Fig. 4(c), for a smooth potential step modeled by a linear carrier density profile. The corresponding on-site energy profile follows from the energy-density relation shown in the upper-right inset of Fig. 4 clarity, all transmission functions are normalized by the number of incoming modes in lead 1. In Fig. 4(b), the transmission function T 4←1 is shown as a function of the density n R and n L on the right and left side of the potential step, respectively. The upper part of the side panel of Fig. 4(b) zooms in near charge neutrality, and the lower part shows T 2←1 for comparison. When n L is near charge neutrality and n R corresponds to E > EM , we have T 4←1 = T 2←1 = 1/2 as before, which is better seen in Fig. 4(c) where we show the transmission (T ) for n L = 10 11 cm −2 , corresponding to the dashed lines in the insets of Fig. 4(b). The plateau of 1/2 for n R > 5 × 10 12 cm −2 is imperfect, as few modes escape through other leads. Nevertheless, almost all reflected modes end up in lead 2, whereas T 6←1 is negligible. The transverse deflection is best visualized by the current density; see, for example, Refs. 50 and 51 for technical details. The magnitude of the current density is shown in Fig. 4(d) for n L = −n R = 10 11 cm −2 and shows a Klein-collimated electron beam [51], similar as in graphene n-p junctions, due to the smooth junction [49]. The current density in the valley splitter regime is shown in Figs. 4(e) and (f), where in (f) the width of the injector lead is slightly reduced and the classical trajectories of the reflected current are superimposed. More exemplary current density images are given in the SM [42]. The magnitude of the deflection can be tuned by changing the slope of the potential which controls the classical turning point. Importantly, the deflection is in the opposite direction for the other valley (lead 2 ↔ lead 6) as the valleys are related by time reversal.

(c). For
Disorder on the moiré scale can induce scattering between minivalleys K t and K b (or K b and K t ) as they are separated by k θ ∝ L −1 m . However, since the minivalleys are localized on different layers, such disorder has to couple the layers. The valley splitter is also robust against an interlayer bias as long as the bias is smaller than the VHS energy, since it only depends on the change in topology of the Fermi surface in the barrier, such that two disjoint Fermi circles connect. Lastly, we want to address the orientation of the barrier. So far, we assumed that the barrier is oriented along the zigzag direction of the dual honeycomb lattice. In principle this could be achieved by visualizing the moiré pattern with scanning tunneling microscopy [3] and rotating the sample accordingly. If the orientation differs by an angle χ relative to the zigzag direction, there is always a finite part of electrons at the Fermi surface demonstrating the same physics provided that tan χ < 2k F /k θ = √ 2πn/k θ . Conclusion -We have shown how electrostatic potentials in TBG give rise to valley-layer locked currents and how these can be used to design a valley valve. Such currents are characterized by a conductance plateau in each layer given by one fourth of the total ballistic conductance, which could be measured by contacting a single layer. Furthermore, we demonstrated the valley-selective transverse focusing of the current. This could be observed by injecting electrons from a narrow contact and measuring the response at nearby contacts. We believe this work opens up new means of manipulating the valley degree of freedom in graphene systems.
Acknowledgments In the limit of small twist angles, the low-energy physics of twisted bilayer graphene (TBG) becomes independent of the geometric details of the lattice, and is well described by the continuum model of TBG. In this model, the two valleys are treated separately as intervalley processes are suppressed due to the slow variation of the moiré pattern on the atomic scale. Furthermore, only the dominant terms in the interlayer tunneling term are retained. Following Refs. 4, 15, and 36, the continuum model for TBG for a single valley can be written in real space aŝ whereĤ l (l = t, b) is the Dirac Hamiltonian corresponding to the top and bottom layer, respectively, and T is the interlayer moiré coupling: where K l = R ±θ/2 K and σ θ = e iσzθ/2 σe −iσzθ/2 with σ = (σ x , σ y ) the Pauli matrices. Here, we also have G 1,2 = √ 3k θ ±1/2, √ 3/2 (note that G 2 is defined with an opposite sign in the main text) for K = Ke x , which is a reciprocal basis of the moiré pattern and k θ = |K t − K b | = 2K sin(θ/2) = 4π/(3L m ) with K = 4π/(3a). In general, the interlayer coupling matrices T n (n = 0, 1, 2) are given by with u AA and u AB real parameters and w = e i2π/3 . The phase factor w arises since the AB/BA regions are located away from the AA region by ± (L 1 + L 2 ) /3, respectively, and u AA < u AB accounts for the corrugation of the twisted bilayer [9,20]. The form of the T n also follows from symmetries. Assuming T 0 is real [9], one finds that C 2x symmetry, which flips the layers but leaves the valleys invariant, together with C 2z T symmetry (with the T the spinless timereversal operator) which also leaves the valleys invariant, require T 0 = u AA σ 0 + u AB σ 1 . The other two are determined by C 3z symmetry, which gives T 1,2 = exp (±iσ z 2π/3) T 0 exp (∓iσ z 2π/3). The Hamiltonian is then diagonalized by plane-wave expansion. To this end, we consider a general Bloch wave, where G = n 1 G 1 + n 2 G 2 with n 1,2 integers and ψ l are 2-spinors containing amplitudes of A l and B l sublattices. In momentum space representation, the Hamiltonian can then be written aŝ whereĉ † l,k−G (ĉ l,k−G ) are two-component fermion operators that create (destroy) a Dirac particle in a single valley in layer l with momentum k − G − K l , H l (q) = e −iq·rĤ l e iq·r , and G 0 = 0. Here, k is constrained to the moiré Brillouin zone (MBZ) and defined with respect to theΓ point. Diagonalization of the Bloch Hamiltonian requires solving an infinite system of equations. In practice, to achieve convergence of the low-energy bands in the first MBZ, we only include G for which |K l + G| < 6k θ . In this way, we obtain the energy bands E mk , shown in Fig. 2 of the main text and Fig. S3  In Fig. S1, we show the reciprocal space of the moiré lattice of TBG in the extended zone scheme where the equivalent Dirac points of the top (bottom) layer are indicated by blue (red) dots. The form of the interlayer coupling in (2) corresponds to scattering between Dirac cones of different layers in momentum space with the smallest momentum transfer. Higher-order Fourier components (higher moiré harmonics) lead to longer-range "hopping" terms in momentum space but are exponentially suppressed [36].

S2. TWO-BAND LATTICE MODEL
In this section, we follow Ref. 39 and construct a Wannier basis from the two low-energy bands for a single valley and spin. We use a slightly different ansatz based on their results: where R = n 1 L 1 + n 2 L 2 (n 1,2 ∈ Z) is a moiré lattice vector, N is the number of moiré cells, k runs over all momenta in the first MBZ, and |Ψ mk are the Bloch states (4) obtained from the continuum model. Unlike Ref. 39, we take a gauge where the B t component of |Ψ mk is real at the center of the BA region of the moiré lattice ( Fig. 1 (a) in the main text), and the phase φ mk is chosen so that the B b component of e iφ mk |Ψ mk is real at the center of the AB region. This is still consistent with the symmetry analysis of Refs. 38 and 39. In this way, Bloch functions in |R, 1 interfere constructively at the BA region, and we expect it to be localized around R + r BA . Similarly, |R, 2 is expected to be localized around R + r AB . We do not further optimize these Wannier states, so our Wannier basis is not maximally localized. However, we find that it yields well-localized Wannier functions away from the magic-angle, as shown in Fig. S2 (a), which is sufficient for the regime of twist angles that we are interested in. The hopping amplitudes between the Wannier orbitals are then obtained as follows: where σ, σ =Ā,B, and we used Ψ mk |Ĥ|Ψ m k = E mk δ mm δ kk , and In the Wannier basis, the effective Hamiltonian for the two low-energy bands, can then be written aŝ where n, m run over the sites of the dual honeycomb lattice with sublatticesĀ = AB andB = BA, as shown in Fig.  1 (a) in the main text, on which the Wannier orbitals are centered. We also defined fermion operatorsĉ nσ (ĉ † nσ ) which destroy (create) particles in the Wannier state |R n , σ and hopping matrices t σσ nn = t σσ (R n , R n ). The hopping amplitudes are illustrated in Fig. S2 (b) for θ = 2 • . Note that (10) lacks time-reversal symmetry within a single valley because the hopping amplitudes are complex in general. The hopping amplitudes for the other valley are obtained through complex conjugation, restoring time-reversal symmetry when both valleys are accounted for.

Bulk spectrum
The bulk spectrum of the lattice model (10), taking into account only the largest hoppings: t 1 , t 2 , t 4 , t 5 , and t 6 , as defined in Fig. S2 (b), is given by where with L 1,2 = L m ( √ 3/2, ±1/2). The energy bands of the two-band lattice model for θ = 2 • is shown in Figs. S3 (a) and (c) where we compare it with the energy bands calculated with the continuum model.
The origin of the transverse deflection of the reflected current can be traced back to the imaginary intrasublattice hopping t 2 ≈ it 2 [38,39]. In k space, this term becomes ∝ τ g(k)σ 0 where τ = ± corresponds to the valley and g(k) is an odd function with g(k)|M = −4 √ 3k x L m + O(k 3 L 3 m ). This pushes the VHS fromM along theΓM line [23]. Hence, the transverse velocity nearM is asymmetric, giving rise to transverse valley currents.  (8) Table S1 and defined in Fig. S2  top panel of Fig. S4 shows almost perfect conductance G 0 /2 with very small interference oscillations. In this case all transverse modes (different k ) from one minivallley are reflected. However, in the energy range approximately given by where v θ is the renormalized Fermi velocity, the number of transverse modes that are reflected via the semiclassical mechanism, gradually diminishes and the conductance increases almost linearly from G 0 /2 to G 0 . When U 0 increases further for E F − U 0 > 0, the conductance is close to G 0 with small oscillations, until the barrier height becomes comparable to E F . In the bipolar regime (E F − U 0 < 0) one might naively expect also a minivalley filtering effect for E F − U 0 < −EM (with reversed sign) since electrons in the barrier which were perfectly reflected via the semiclassical mechanism are now fully transmitted and vice versa. However, in order to get inside the barrier electrons need to pass through the n-p-junction at the left side of the barrier (and then leave the barrier through the p-n junction at the right side). Due to Klein tunneling the n-p junction is mostly transparent when E F is sufficiently small [49] and we still observe a conductance that is reasonably close to G 0 /2 for smooth barriers (blue and red lines) as well as a high minivalley polarization of the current. Interestingly, both the conductance plateau and minivalley polarization in the bipolar region holds better not for the smoothest potential (d = 7, blue) but for the intermediate one (d = 1, red). To understand this, we need to consider another length scale present in our system. In order for electrons to have a well-defined minivalley that corresponds well to the layers (left panel of Fig. 2 of the main text), we need to consider incoming/outgoing electrons with E F close to the Dirac point. The validity of the semiclassical mechanism requires the potential to be smooth at the scale given by the VHS energy λ VHS ∼ v θ /E VHS ∼ L m . However, the electrons at the n-p-junction have much larger wavelength λ F ∼ v θ /E F . If the length scale of the barrier potential is large compared to λ F , tunneling through the junction is suppressed [49]. Hence the reduced conductance and minivalley polarization for d = 7 (blue lines). On the other hand, if the potential scale is small compared to λ F (still large compared to λ VHS ), the barrier at the n-p junction becomes effectively sharp. The total transmission probability for a sharp junction is independent of E F and turns out to be close to unity. This is because the transmission probability reaches a maximum (unity) at normal incidence due to Klein tunneling and for sharp junctions the width of this maximum effectively spans the entire range of transverse modes. Scattering at an n-p junction in TBG (above the magic angle) for E F close to the Dirac point is described by the Hamiltonian which describes a generic crossing of two bands in the presence of the potential U (x). Here, we take |U (x)| |ε| for x < 0 and U (x) > 0 with U (x) |ε| for x > 0. For concreteness, we take ε > 0 from now on. In this case, the right-moving scattering state for x > 0 can be approximately written as while at x < 0, the wave function is given by where we take p x > 0 and the spinors u + , u − may be chosen as The coefficients a and b are found from matching the two solutions at x = 0: Current conservation is assured as We thus find the reflection probability where both ε > 0 and p x > 0. For ε < 0, the sign of p x is also reversed as electrons in the valence band propagate in a direction opposite to their momentum and hence the results holds in general. The normalized conductance is then given by which is valid for a single n-p junctions. For two n-p junctions on either side of the barrier, we find where we neglected interference effects. We believe that this is the reason for the red curve (d = 1) for E F − U 0 < −EM showing better filtering than the blue curve (d = 7) in the bipolar regime. On the other hand, if in the unipolar filtering regime (E F − U 0 > EM ) the potential barrier is not sufficiently smooth with respect to λ F , electrons can be reflected at the "sharp" edge of the barrier. That is why the minivalley filter for d = 1 is similarly imperfect for both the bipolar E F − U 0 < −EM and unipolar E F − U 0 > EM regime.
The oscillations of the conductance in the bipolar regime, as mentioned in the main text, appear due to the interference with electrons bouncing inside the barrier [48]. The interference mechanism in the bipolar filtering regime is illustrated in Fig. S5 (a). Here, we sketch the trajectories of an electron, which in the absence of the n-p junction (e.g. for E F < 0) would be reflected via the semiclassical mechanism illustrated in the top panel of Fig. 3 of the main text. As we see, the wave function is first split into transmitted and reflected parts at the left n-p junction. The transmitted part is then reflected inside the barrier due to the semiclassical mechanism which changes the minivalley (blue to red). This trajectory goes back to the left junction and is split again. However, the transmitted (red) wave does not interfere with the previously reflected (blue) wave as they correspond to different minivalleys. Electrons bouncing inside the barrier start emitting interfering waves only after a full period consisting of two reflections that flip the minivalley and four reflections at two n-p junctions. This may be the reason why interference oscillations are stronger in the bipolar Dirac regime (−EM < E F − U 0 < 0) than in the bipolar filtering regime (E F − U 0 < −EM ) because of the smaller length of interfering trajectories in the former.
Finally, if |E F − U 0 | exceeds the bandwidth of the two lowest bands, the conductance drops to zero as there are no modes in the minigap. For the smooth potential the conductance drops to zero in almost the same way as the drop from G 0 to G 0 /2 at E F − U 0 > EM where one of the minivalleys is blocked. This is because the conductance drop can also be viewed as the barrier becoming opaque for both minivalleys.
Since ideally the barrier only transmits electrons from one minivalley/layer for a given valley, one might try to design a valley valve [33,34] with two barriers in series. This setup is illustrated by the inset of the bottom panel of Fig. S4. When the barriers are gated similarly, one minivalley can go through. For opposite polarities, one would hope that all current is blocked. Indeed, electrons that pass the first barrier, and therefore belong to the minivalley for which the second barrier is opaque, are reflected by the semiclassical mechanism in the second barrier which flips the minivalley. The reflected electrons now belong to the correct minivalley to pass the first barrier from the other side and hence all electrons reflected back. This process is illustrated in Fig. S5.
In reality, not all current is blocked due to the presence of n-p junctions on the left and right side of one of the barriers. Scattering at these junctions, which are smooth on the moiré scale, conserves the minivalley as the required momentum transfer k θ ∝ L −1 m cannot be provided. Therefore, after the electrons are reflected semiclassicaly in the second barrier, they can reflect back at the n-p junction without changing the minivalley and pass the barrier. Such processes spoil the valley valve but are suppressed when the barrier is still sharp on the scale of λ F as we discussed in the previous section. Hence, there is a trade-off between smooth junctions for which the minivalley filtering works better and less smooth junctions for which there is less reflection at the n-p junctions.
Theoretically, there is no small parameter that ensures that the reflection at a sharp n-p junction (relative to λ F ) is strongly suppressed. Nevertheless, we find that the valley valve works reasonably well as shown in the bottom panel of Fig. S4. This is again because of the fact that Dirac electrons at normal incidence and effectively many other transverse modes are passing through without reflection (Klein tunneling). In the figure, we see that when the barriers have the same polarities, the conductance is approximately given by G 0 /2 for the slow-varying potential (d = 7) and reflection at NP junction In the main text, we have considered a scattering region of L = 1 µm long and W = 3 µm wide attached to six leads in the presence of a smooth junction (see Fig. 4 of the main text), and showed the normalized transmission functions. The scattering region, along with all the six leads, are described by the effective lattice model for a θ = 2 • twisted BLG. See the insets of Fig. S6, where the original transmission functionsT p←1 from lead 1 to lead p = 1, 2, · · · , 6 are shown. These six transmission functions sum to the number of modes M 1 of the incoming lead 1, as required by the sum rule for an N -lead system where M q is the number of incoming modes in lead q. The horizontal gray line in Fig. S6 shows that the above sum rule is fulfilled, and the number of modes M 1 = 7. The normalized transmission functions shown in Fig. 4 (c) of the main text as well as in Fig. S7 are defined by  Figure S7 shows normalized transmission functions with various widths of the injection lead 1, also considering n L = 10 11 cm −2 , same as Fig. S6 as well as Fig. 4 (c) of the main text. The expectedT 2←1 ≈T 4←1 ≈ 0.5 are found to be better seen with a wider injection lead. Overall, all curves shown in Fig. S7 exhibit similar behaviors and are not sensitive to the width of the injection lead.

More on imaging of local current density
In the main text, we have shown a few local current density profiles to illustrate the valley-selective transverse current deflection. Here, we consider the same twisted BLG (L = 1 µm, W = 3 µm, θ = 2 • ) to show more images subject to different density configurations (n L = 10 11 cm −2 remains fixed) and injection lead widths. Note that we only consider the K valley here, while the K valley is obtained from K by time reversal.
In Fig. S8, the first row of panels (a)-(c), second row of panels (d)-(f), and third row of panels (g)-(i) consider W i = 500, 400, and 300 nm, respectively. Despite the decreasing incoming number of modes (see Fig. S7), the current flow behaves similarly upon different W i . The first column of Fig. S8 [panels (a), (d), and (g)] considers (n L , n R ) = (1.0, 2.3)×10 11 cm −2 , so that the twisted BLG is in the low-energy unipolar state, showing typical single-layer graphene nn behavior. The second column of Fig. S8 [panels (b), (e), and (h)] considers (n L , n R ) = (1.0, −2.3) × 10 11 cm −2 , so that the twisted BLG is in the low-energy bipolar state. In this case, the beam structure arises from the Klein collimation due to the smooth np junction, also a typical behavior of single-layer graphene [49]. Finally, the third column of Fig. S8 [panels (c), (f), and (i)] considers (n L , n R ) = (1.0, 6.7) × 10 11 cm −2 , revealing the valley-selective transverse current deflection, insensitive to the width of the injection lead. Note that Fig. S8(b), (c), and (i) are shown in the main text as Fig. 4 (d), (e), and (f), respectively.

C. Classical trajectories
To calculate the classical trajectories shown in Fig. 4 of the main text, we numerically solve the Hamilton equations, where we replaced the classical Hamiltonian H classical → E(k) + U (y), E(k) is the energy of the first conduction band, and the potential profile U (y) is determined by the relation between the electron density n and E shown in Fig. S5 (b) to obtain a linear density profile n(y) =      n L y < 0 n L + (n R − n L ) y/y 0 0 < y < y 0 n R y > y 0 .