{'fraction_non_alphanumeric': 0.0344151510095111, 'fraction_numerical': 0.011722009010512265, 'mean_word_length': 4.6340775558166865, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
All-strain based valley filter in graphene nanoribbons using snake states

A pseudo-magnetic field kink can be realized along a graphene nanoribbon using strain engineering. Electron transport along this kink is governed by snake states that are characterized by a single propagation direction. Those pseudo-magnetic fields point towards opposite directions in the K and K' valleys, leading to valley polarized snake states. In a graphene nanoribbon with armchair edges this effect results in a valley filter that is based only on strain engineering. We discuss how to maximize this valley filtering by adjusting the parameters that define the stress distribution along the graphene ribbon.


I. INTRODUCTION
The advent of graphene 1,2 not only represented the beginning of a new era of atomically thin materials, with potential technological applications in future electronic and photonic devices, but also brought the possibility of observing several novel phenomena due to its unique band structure, consisting of Dirac cones in points labeled as K and K ′ in its first Brillouin zone. In fact, the existence of two inequivalent cones is of special importance, since it enables a new degree of freedom to be explored in novel valley-tronic devices.
Several suggestions have been made to harvest valley polarization in graphene: Rycerz et al. 3 demonstrated that specific combinations of armchair and zigzag edges in a monolayer graphene ribbon lead to efficient valley filtering. Non-uniform substrate induced masses can also be used to obtain valley polarization, as shown in Refs. [4,5]. As for bilayer graphene, valley filtering can be obtained by specific configurations of external potentials, 6 or boundaries with monolayer graphene regions. 7,8 On the other hand, recent studies have demonstrated that pseudo-magnetic fields can be induced in graphene by specific strain configurations and, since these fields point towards opposite directions in different Dirac cones 9 , several suggestions of strain-based valley filters have been proposed in the literature. Most of these proposals involve combinations of the strain induced fields with applied magnetic and electric fields. [10][11][12][13][14] Indeed, strainbased valley filters are specially interesting, because of graphene's ability to withstand large mechanical stress. 15 Very large pseudo-magnetic fields have been experimentally observed in e.g. naturally formed bubbles in a monolayer graphene system on a Pt substrate. 16 In this paper, we propose a very different valley filter device made of a single graphene layer 17 that does not depend on specific edge configurations, 3 substrate induced masses, 4 or external magnetic fields 10,18 . Yet, it is all based on a particular kind of strain induced in a monolayer graphene nanoribbon, that provides a pseudo-magnetic kink barrier along the ribbon width.
Such strain is expected to be attainable by using flexible substrates in combination with e.g. appropriate piezos or contacts that deform when cooled down. 10,15,19,20 The valley polarization originates from a combination of (i) the uni-directional motion of snake states along the pseudo-magnetic kink, and (ii) the opposite direction of the pseudo-magnetic field felt in the K and K ′ valleys. This combination leads to electrons moving in single orbits propagating in opposite directions in the two different valleys. In order to verify the efficiency of such valley filtering device, we calculate the transmission probability of wavepackets through this structure within the tightbinding model. Our results demonstrate that a valley polarization efficiency up to 90% can be reached, provided specific conditions are met by the system, as we will discuss in what follows.

II. THEORETICAL MODEL FOR THE STRAIN
Our system consists of a monolayer graphene nanoribbon with width W y ≈ 6387Å and length W x ≈ 2214Å corresponding to 1801 × 3000 carbon atoms, that is mechanically strained in a specific configuration, as sketched in Fig. 1: along a certain region of length β, the ribbon is distorted into two circular arcs of radius R, in opposite directions. Such in-plane circular bending is obtained by defining the displacement of the atomic sites as 15 and where u x,y is the in-plane lattice distortion due to strain and the radius has its sign reversed at the y = 0 axis, i.e. R = |R| (2θ(y) − 1), with θ(y) being the step function. A sharp transition between strained and unstrained regions of the ribbon would be clearly impossible, since it would lead to unrealistically large atomic distances in the vicinity of the transition region, specially for small R. Therefore, we consider a smooth (Gaussian) variation of the curvature K R = 1/R = (1/R 0 )e (−j 2 /β 2 ) , where j is the index of the column to which a given atomic site belongs in the lattice, 21 the length is described by the dimensionless parameter β, while R 0 provides the maximum radius of the curve (namely, at the central column of atoms, where x i,j = 0 in the absence of strain).
In general, strain effects on the electronic properties of graphene can be mapped into the analogous problem of an electron under a pseudo-magnetic field distribution, 22 whose magnitude and orientation may oscillate over the space, thus making the production of a local non-zero pseudo-flux challenging. 23 However, it has been recently demonstrated 9,15 that an in-plane circular distortion, as the one proposed here, deforms the Brillouin zone, shifting the Dirac cones with respect to each other, just like when an uniform magnetic field is applied perpendicular to the graphene plane, leading to K→K+2πA/Φ 0 , where Φ 0 = e/h is the flux quantum. Such lattice distortion changes the hopping energies and thus induce an effective vector potential 24,25 where δa ij is the vector distance between the adjacent atoms i and j in the strained lattice, v F is the Fermi velocity, δτ ij is the difference between the strained and unstrained hopping energies, and the pseudo-magnetic field is given by B ps = ∇ × A. Moreover, the distortion in different directions for y > 0 and y < 0 provides a pseudo-magnetic kink barrier, with pseudo-magnetic field regions that change sign at y = 0. A schematic example of the pseudo-magnetic field distribution induced by such strain configuration is presented in Fig. 1(b) for a representative set of parameters (β, R 0 ). We point out that the sample considered in Fig. 1(b) is much smaller than the one investigated throughout this paper, since calculating and plotting a vector potential distribution along the 1801 × 3000 atomic sites of our actual sample requires high computational costs. Therefore, the pseudo-magnetic field in Fig. 1(b) is discussed here only in a qualitative way. The pseudo-magnetic field is found to be zero at input and output leads, where the lattice displacements vanish, and assume its minimum and maximum values along the ribbon width around x = 0, where the strain is maximum. Two additional kinks are also consistently observed on the left and right sides of this main central kink. They are however much smaller than the central one and, thus, do not play an important role in the valley filtering process, as we will demonstrate further on. In fact, we observe that as we increase R 0 , these additional kinks become even lower as compared to the main kink, so that their importance for the transport properties of the actual sample studied throughout the paper (with larger R 0 ) is negligible. This can be verified by comparing the bottom panels in Fig. 1(b) for two different maximum radii R 2 0 > R 1 0 assuming a fixed width for the strained region β.

III. SNAKE STATES ALONG A MAGNETIC FIELD KINK
Keeping with the analogy between this strain configuration and a magnetic field kink, let us first calculate the energy dispersion along the ribbon in the presence of such a magnetic barrier. We assume an inhomogeneous magnetic field B = Bẑ that depends only on the transversal coordinate, given by B(y) = B(θ(y)−θ(−y)). Low energy electrons in graphene exhibit a linear energy dispersion, so that they behave as massless Dirac-Weyl fermions, thus, obeying the Dirac equation: with energy E = v FĒ around the K valley (a similar analysis 2 can be made for electrons in K ′ ).
Defining the vector potential in the Landau gauge, A =Ā(y)ê x , withĀ(y) = −B(y)y, the general solution for the wavefunction with translational invariance in the x-direction is Ψ(x, y) = ψ(y)e ikx . Therefore, we obtain from Eq. (4) where A = e Ā . This leads us to decoupled equations for each component: for instance, for the upper spinor component, If we use the magnetic length l b = /e|B| as the unit of distance, we obtain The energy dispersion along the y-direction is obtained quasi-analytically by solving this equation in terms of parabolic cylinder functions 26 , D p (q). Notice the magnetic field B(y) is piecewise constant, hence, one can separate solutions for each region as where q = √ 2 [sgn(B)kl b ], and p = ǫ 2 2 − 1. The continuity of the wavefunction and its derivatives at these regions provides boundary conditions that lead to quantization of the energy of the system. Eqs. (8) and (9) represent solutions for the first region of the system, but they also can express solutions for the second region by replacing the coefficients a ± → c ± .
In order that the wavefunctions are normalizable we demand a + = c − = 0. Then, the boundary condition at y = 0 gives the equation that generates the energy quantization condition. Therefore, we obtain FIG. 2: (Color online) Energy states for a Dirac particle in the presence of a (pseudo-)magnetic kink. Results for K and K ′ valleys are the same for an external magnetic field, whereas for a pseudo-magnetic field, the K ′ spectrum (red dashed curves) differs from that from the K valley (black solid curves).
where the functions are given by Notice that these results for the wavefunction are very closely related to the one for a magnetic kink profile in a normal 2D semiconductor. 27 Numerical results for this system are illustrated in Fig. 2, where we observe an asymmetry in the energy bands along the kink with respect to the k x = 0 axis. Namely, the energy states (predominantly) monotonically decrease with k x , implying a negative velocity v = (1/ )dE/dk, that eventually converges to zero as k = k x → +∞. It is straightforward to verify that by inverting the sign of the magnetic field kink, this figure is reflected with respect to the k x = 0 axis, and the propagation velocities are then predominantly positive. The physical interpretation of this result has its basis on the existence of snake states 28,29 that propagate along the kink, which can also be inferred from a simple classical analysis of this problem, involving Lorentz force, cyclotron orbits and the right-hand rule. Nevertheless, this result is of special importance in the context of pseudomagnetic kinks discussed here: since the strain-induced pseudo-magnetic field points towards opposite directions in the different Dirac cones K and K ′ , electrons in each cone will propagate in opposite directions in the pseudomagnetic kink proposed here, thus yielding an efficient valley filtering process, as we will verify with our numerical results afterwards.

IV. WAVEPACKET PROPAGATION METHOD
In order to investigate the transport properties in our strained graphene, we will use a wavepacket propagation method. A comparison between this kind of method and those based on Green's function formalism can be found in Ref. [30]. The advantage of using this approach is in the possibility of observing the trajectories of the wave packet describing the electron propagating across the scattering (strained) region, which reveals important information about the physics behind any unusual behavior of the current through the system, being due to e.g. inter-valley scattering, skipping orbits, snake states, etc., as we will demonstrate in the following Section.
We use a Hamiltonian within the tight-binding model where the operator c † i (c i ) creates (annihilates) an electron on site i, and τ ij is the hopping energy between adjacent atoms i and j (nearest-neighbors), that depends on the distance δa ij between them according to 9 We consider an initial Gaussian wavepacket where N is a normalization factor and calculate its time evolution using the split-operator method, 18,21 in which the time-evolution operator for the Hamiltonian H = H i + H j is split as where H i(j) is the term of the tight-binding Hamiltonian H that corresponds to a horizontal (vertical) hopping between atomic sites Notice that for the horizontal term H i , one has to differentiate between hoppings to the right and left neighbouring sites, since in the honeycomb lattice, each site has only horizontal hops to one side. Hence, τ ′ ij = τ ij =⇒ τ ′′ ij = 0 and τ ′′ ij = τ ij =⇒ τ ′ ij = 0. The advantage of such splitting lies in the fact that these operators can be represented by tridiagonal matrices, that are easily handled by standard computational routines. The wavefunction after a single time step t+∆t is then obtained in three steps Each of these equations is re-written using the Cayley form for the exponentials, e.g.
and the remaining tridiagonal matrix equation is then numerically solved by standard computational routines. 31 For our study, we used a wavepacket width of σ = 300 A and its wave vector k has a modulus of k = 0.06Å −1 , unless otherwise explicitly stated in the text. Using the linear spectrum approximation for low-energy electrons in graphene, in which E = v f k, the wavepacket energy is estimated to be E = 343 meV. Besides, as we intend to demonstrate the valley polarization of the wavepacket, we place it in different valleys in reciprocal space by shifting the wave vector towards the two inequivalent Dirac points: where the positive (negative) sign refers to a displacement towards the K (K ′ ) point of the Brillouin zone, and a ≈ 1.42Å is the inter-atomic distance.
As the Gaussian wavepacket propagates, we calculate the probability of finding the electron before (P 1 ), within (P 2 ), and after (P 3 ) the strained region, as the integral of the square modulus of the wavepacket, taken within the intervals −3, 000Å ≤ x ≤ −400Å, −400Å ≤ x ≤ 400 A, and 400Å ≤ x ≤ 3, 000Å, respectively. Transmission probabilities are assumed to be the converged value of P 3 as t → ∞. Besides, we keep track of the wavepacket trajectories by calculating the average value of the position, ( x , y ), at each time step.
The armchair edges of the ribbon do not support edge states, therefore, modelling the electron propagating through the system as a wavepacket, whose tails do not reach the ribbon edges, is justified. Moreover, any improvement to come from other calculation methods, involving e.g. plane waves, scattering matrices and the Landauer-Buttiker formalism, would lead to rather quantitative corrections to our results, while the qualitative behavior of the system and the proof-of-concept of valley filtering with a pseudo-magnetic kink proposed here, which are the main goals of this work, would still hold, since they are based on more fundamental physical properties of the proposed structure, as we will discuss in what follows.

V. RESULTS AND DISCUSSION
The existence of snake states in such a strained graphene lattice, as due to the induced pseudo-magnetic fields, is confirmed by the trajectories drawn in Fig. 3 of the center-of-mass of a k = 0.06Å −1 wavepacket propagated in time through the system described by Fig. 1(a), assuming β → ∞ and R 0 = 10 4Å , as a test case. If curves are just mirror-reflected with respect to the y = 0 axis of Fig. 3, and the situation remains the same. Conversely, if the wavepacket propagates from right to left, it is the K ′ packet that draws a snake trajectory, whereas the K packet is reflected. One could think that wavepackets deflected towards the edges of the system (i.e., further away from its center) would be reflected by the ribbon edges, perform skipping orbits, and eventually pass through the strained region. However, since the ribbon has armchair edges, reflected wavepackets are scattered to the other Dirac cone, where the pseudo-magnetic field is opposite, thus the skipping orbit follows the opposite direction and the wavepacket comes back anyway. 17 This non-propagating edge state is emphasized in Fig. 3 for a wavepacket that started at the bottom-half (tophalf) of the ribbon and around K ′ (K) Dirac valley. In this way, one completely avoids the problem of having valley mixing of the snake and edge states at the end of the ribbon, which would otherwise occur in the case of propagating edge states e.g. in the presence of an external applied magnetic field. Valley mixing by scattering at the contacts can also be further suppressed by using graphene electrodes. 32 Notice, however, that the present proposal will not work very efficiently for zigzag graphene nanoribbons, where such inter-valley edge scattering does not occur and skipping orbits are allowed to propagate at the zigzag edges. Although this represents a limitation of the proposed system, fabrication techniques have been advancing fast, and armchair graphene nanoribbons with very high edge quality have already been experimentally demonstrated. 33 Results in Fig. 3, thus, allow us to conclude that electrons in K (K ′ ) cones in such a strained armchair graphene ribbon can only propagate towards the right (left). Analogously, if the strain configuration is inverted, trajectories drawn by K and K ′ packets are switched.
Such a picture of snakes states strongly suggest a valley filtering effect. In fact, if one now considers a system with finite strain region β = 900 (in units of the inter-atomic distance a 0 = 1.42Å) and R 0 = 10, 000Å, a wavepacket with k = 0.06Å −1 around the K cone passes through this region with a high probability P 3 ≈ 0.9, whereas the same packet in K ′ would have a much lower transmission probability P ′ 3 ≈ 0.3. Let us now search for an optimization of the valley polarization effect. The polarization, as defined by P = 1 − P 3 /P ′ 3 , where P ′ 3 is the transmission probability for a wavepacket in the K ′ cone, is shown in Fig. 5 for k = 0.06Å −1 , assuming different parameters β and R 0 . For a fixed strain radius R 0 , increasing the length of the strain region β increases the polarization, as shown in Fig. 5(a). Besides, results in this panel also suggest that decreasing R 0 would always improve the polarization; this would be reasonable, since smaller radii yield stronger distortions in the lattice and, consequently, larger pseudo-magnetic fields. This is however not always the case: Fig. 5(b) shows that even for β as large as 900Å, decreasing the strain radius will always lead to a maximum polarization at an intermediate value R 0 ≈ 5, 000Å, so that the polarization is reduced as the radius is further decreased. This is due to the fact that the smooth connection between the unstrained ribbon leads and the strained region might end up creating a complicated pseudo-magnetic field distribution, with regions with fields pointing to opposite directions, which would harm the polarization effect investigated here. This also explains the negative polarization observed for small β in Fig. 5(a).
So far, all results were obtained for k = 0.06Å −1 , which corresponds to a wavepacket energy E = 343 meV. It is however important to check how the polarization depends on the wavepacket energy. This is shown in Fig.  6, where we verify that the valley filtering process proposed here has an optimal range of energies. Indeed, if the energy is too low, the pseudo-magnetic Lorentz orbits would have a very small radius, so that only portions of the wavepacket that are very close to the y = 0 line would pass through the system as snake states, whereas the rest of the wavepacket readily turns back. On the other hand, if the energy is too high, orbit radii may end up being larger than the strained region length, so that the snake-like propagation that leads to valley polarization no longer occurs. Moreover, the Lorentz orbit radius is inversely proportional to the pseudo-magnetic field intensity, therefore, increasing the strain by reducing R 0 would also lead to Lorentz orbits with smaller radius. This classical picture is consistent with our numerical findings: in Fig. 6(a), the largest strain radius R 0 = 12, 000Å provides the fastest decay of polarization as the energy increases, since it yields a lower strain-induced pseudo-magnetic field and, therefore, energies slightly higher than the optimal E ≈ 125 meV already provide orbit radii larger than the strained region length. Conversely, for energies lower than E ≈ 125 meV, R 0 = 12, 000Å provides the best polarization, as its weaker pseudo-magnetic field compensates for the low energy and prevents the orbits radius of becoming too small. Also, Fig. 6(b) shows that, for higher energies, where orbit radii are larger, polarization is more efficient for larger length β, and the energy for optimal polarization increases with this parameter.
The valley filtering effect by pseudo-magnetic kinks E (meV) demonstrated here is also expected to be robust against impurity and defects scattering: as already discussed, electrons in each valley have only one possible direction of propagation (see Fig. 3), due to the monotonic behavior of all the energy states as a function of momentum (see Fig. 2), which provides a single direction for the group velocity in each valley. After scattering by impurities or defects, the electron must end up in one (or a combination) of the states in Fig. 2. If the electron is already in the valley that allows its propagation through the system (as a snake state), with positive velocity (i.e. monotonically increasing energies as a function of momentum), any component of the scattered electron wave function that ends up in the other valley must be deflected backwards, since there is simply no energy state in that valley with positive velocity as well. Thus, provided the pseudo-magnetic field kink distribution is preserved, only electrons in one of the valleys are allowed to reach the other end of the ribbon, even after scattering events.

VI. CONCLUSIONS
We have investigated the wavepacket propagation through a graphene nanoribbon with armchair edges for a specific strain distribution. The latter provides a pseudomagnetic barrier kink along the ribbon. By following the trajectory of the center-of-mass of the wavepacket, calculated by solving the time-dependent Schrödinger equation for the tight-binding Hamiltonian, one observes snake states, which have a fixed propagation direction, consistent with the pseudo-magnetic kink picture. How-ever, one can analytically verify that, by reversing the magnetic kink, the propagation direction of snake states must be reversed.
Since the pseudo-magnetic field points towards opposite direction in the different Dirac cones, wavepackets in the different cones can only have fixed opposite directions of propagation. This effect results in an efficient valley filtering process, which does not require either lattice defects, edge engineering, or any externally applied fields or potentials. Our numerical results show significant valley polarization through this system, which can be optimized by the parameters (β, R 0 ) that depend on the electron energy (i.e. the Fermi energy).
Notice that the in-plane circular deformation of a graphene nanoribon proposed here is just one particular way of inducing a kink pseudo-magnetic field barrier: any other strain distribution that produces such a pseudo-magnetic kink would lead to similar valley filtering effect, which requires only a graphene ribbon with armchair boundaries (to avoid edge propagation) and a pseudo-magnetic field that flips its direction across a line parallel to them.