Scattering Theory of Graphene Grain Boundaries

The implementation of graphene-based electronics requires fabrication processes that are able to cover large device areas, since the exfoliation method is not compatible with industrial applications. The chemical vapor deposition of large-area graphene represents a suitable solution; however, it has an important drawback of producing polycrystalline graphene with the formation of grain boundaries, which are responsible for the limitation of the device’s performance. With these motivations, we formulate a theoretical model of a single-layer graphene grain boundary by generalizing the graphene Dirac Hamiltonian model. The model only includes the long-wavelength regime of the charge carrier transport, which provides the main contribution to the device conductance. Using symmetry-based arguments deduced from the current conservation law, we derive unconventional boundary conditions characterizing the grain boundary physics and analyze their implications on the transport properties of the system. Angle resolved quantities, such as the transmission probability, are studied within the scattering matrix approach. The conditions for the existence of preferential transmission directions are studied in relation with the grain boundary properties. The proposed theory provides a phenomenological model to study grain boundary physics within the scattering approach, and represents per se an important enrichment of the scattering theory of polycrystalline graphene. Moreover, the outcomes of the theory can contribute to understanding and limiting the detrimental effects of graphene grain boundaries, while also providing a benchmark for more elaborate techniques.


Introduction
Graphene (G) is a two-dimensional honeycomb lattice that is constituted by carbon atoms. Its reciprocal lattice determines a hexagonal Brillouin zone that has six corners (K/K' points) where the low-energy part of the band's structure is well described by a linear energy-momentum dispersion relation, defining the so-called Dirac cone. The existence of Dirac cones in the graphene band's structure can be understood by using a tight-binding model. Consequently, the electron motion in the graphene lattice follows the Dirac equation [1], with the latter being the manifestation of an emergent ultra-relativistic behavior in a many-body system that was initially described by the Schrödinger equation. Due to its unique band structure, in the past few years, graphene has attracted much attention and subsequently its intriguing transport properties, such as Klein tunneling [2], Zitterbewegung effect [3], antilocalization [4], anomalous quantum Hall effect [5], and Veselago focusing effect [6], have been suggested and, in some cases, experimentally proven. Apart from its theoretical interest, graphene is a two-dimensional chemical homogeneous system that is characterized by very high electrical mobility [7] and extraordinary mechanical properties [8], making it appealing in nanoelectronics and for flexible electronics implementations [9]. Due to its single-atom-thick nature, graphene appears to be the ideal candidate to study field effect devices [10]. However, common problems with graphene-based field effect transistors (FET) are: device performances are usually limited by the electrode-graphene contact resistance [11][12][13][14]; the position of the Dirac point is strongly affected by a random doping induced in the fabrication process; a low on-off currents ratio in G-FET compared to the commercial silicon-based FET, due to the lack of a bandgap; graphene channel contamination [15,16] occurs quite easily due to chemical impurities, ambient conditions (e.g., humidity), etc.; Dirac point modulation along the graphene channel due to the chemical doping induced by the electrodes and controlled by the difference of the work functions at the metal-graphene interface.
Moreover, the implementation of graphene-based electronics requires the chemical vapor deposition technique, which is able to cover large device areas and is compatible with industrial fabrication processes. However, an important drawback of this technique is represented by the polycrystalline nature of the samples. Polycrystalline graphene presents grain boundaries [17][18][19], which are responsible for the limitation of the device's performance. For this reason, the morphology [20,21] and the electrical properties [22,23] of grain boundaries in polycrystalline graphene have received great attention in the last decade. In particular, it has been shown that one-dimensional defects in graphene accumulate electrostatic charge via a self-doping mechanism [24]. Self-doping takes place because one-dimensional defects in the sp 2 carbon structure accommodate localized states, acting as trapping centers for the charge carriers. The effects of charged line defects in graphene have been studied in the literature and are recognized as an important factor for mobility degradation in chemically vapor-deposited graphene [25,26].
From a theoretical point of view, linear defects in graphene have been studied by J.N.B. Rodrigues and coworkers in Refs. [27][28][29][30], where specific boundary conditions have been derived from the tight-binding model of the extended defect. The mentioned approach provides a continuous model to study transport properties in the presence of selected linear defects, whose atomic arrangement determines the boundary conditions of the scattering problem.
In this work, we provide an alternative method and develop a continuous model of polycrystalline graphene in which general boundary conditions are derived without specifying the microscopic structure of the linear defect.
With these motivations, we formulate a theoretical model of the graphene grain boundary by generalizing the graphene Dirac Hamiltonian model. A continuous approach is used that is justified within the long wavelength limit assumed throughout this work. Using symmetry-based arguments, we derive unconventional boundary conditions characterizing the grain boundary physics, and analyze their implications on the transport properties of the system. The scattering matrix approach [31] is used to derive the conduction properties of the grain boundary interface, and angle-resolved quantities, such as the transmission probability, are presented. The conditions for the existence of preferential transmission directions are studied in connection with the grain boundary properties.
The work is organized as follows. In Section 2, we derive the Dirac Hamiltonian within a rotated reference frame, current density conservation, and boundary conditions at a grain boundary. The matching matrix method is introduced and adapted to the scattering problem at the grain boundary. In Section 3, we formulate a grain boundary Hamiltonian model with the position-dependent rotation angle θ(x). The model is studied by using space-dependent unitary transformation, which allows the use of conventional boundary conditions in the scattering problem. In Section 4, the results of the scattering theory are reported. Conclusions are given in Section 5. Appendix A is included to discuss K point displacement in momentum space induced by off-diagonal potentials in sublattice representation.

Dirac Hamiltonian within a Rotated Reference Frame, Current Density Conservation, and Boundary Conditions at a Grain Boundary
In this section, we provide a first approach to the grain boundary problem in monolayer graphene, while studying the grain boundaries in multilayer graphene is beyond the purposes of this work. We deal with a theory of spinless particles described by the continuous Dirac Hamiltonian, and consider a single valley degree of freedom. Neglecting the electron spin is appropriate in the absence of magnetic fields or magnetized regions, which is the case considered in this work. On the other hand, describing graphene electrons adopting a single-valley perspective requires an appropriate justification. Indeed, in principle, the grain boundary region can be the source of valley-flipping scattering events. The existence of such kind of phenomena has inspired the emerging field of valleytronics [32], aiming at manipulating the valley degree of freedom as it is done with the spin in spintronics. However, valleytronics manipulations require engineered scattering centers, which are quite challenging to implement using current technologies. In view of this circumstance, we can argue that valley-mixing scattering events are quite rare in spontaneously formed grain boundaries. Moreover, to a good approximation, electrons that have different valley degrees of freedom can be described as two independent quantum fluids. The latter argument justifies the single-valley treatment that is adopted. Generalization of the ideas exposed hereafter to spinful particles having two-valley degrees of freedom is in principle immediate.

Dirac Hamiltonian within a Rotated Reference Frame
In this subsection, we describe the continuous Dirac Hamiltonian of a single graphene sheet within a rotated reference frame. The continuous approach considered in this work is appropriate for describing the system conduction properties (e.g., the differential conductance) when the electron wavelength, which is associated with the particle momentum, is much longer than the interatomic distance of the honeycomb lattice. While the latter assumption is clearly invalidated in graphene nanoribbons devices, it is appropriate for a large variety of graphene-based devices where micrometric graphene sheets are employed as conduction channels.
Let us describe the Dirac Hamiltonian of the graphene sheet GS2 in Figure 1, in terms of the reference frame RF1 of the graphene sheet GS1. The Hamiltonian of GS2 written in terms of reference frame RF2 takes the usual form: where ξ = ±1 represents the valley quantum number and v F is the Fermi velocity, whilep x = −i ∂ x andp y = −i ∂ y are the quantum mechanical operators associated with the linear momentum components of quasiparticles. The matrix structure of the Hamiltonian originates from the presence of two atoms inside the unit cell of the graphene Bravais lattice. Accordingly, the wave function Ψ(x , y ) = (Ψ A , Ψ B ) T describing the quantum state of the charge carriers is a two-component spinor whose components are related to the probability of finding the particle on atom A or B of the unit cell. For the reasons explained above, we neglect the valley quantum number and set ξ = 1.
The correspondence between the wave functions expressed in RF1 or RF2 is established, observing that the α-component of the wave function written in RF1, i.e., Φ α (x, y), is related to the homologous component in RF2 by the equation Φ α (x, y) = Ψ α (x (x, y), y (x, y)), with α ∈ {A, B}. The relation between the distinct reference frames is simply given by a two-dimensional rotation: with θ being an appropriate rotation angle. From the above observations, one obtains useful relations linking the partial derivatives of the wave functions in the form: Using Equation (3) in evaluating the quantity H + Ψ(x , y ), one easily gets the equality: where: represents the Hamiltonian of GS2 written in terms of RF1. Solving the eigenvalues problem H θ ψ = E ψ, with the usual planewave ansatz ψ = (α, β) T e i(k x x + k y y) , one finds eigenstates: with associated eigenvalues E ν = ν v F |k|. Here, ν = ±1 represents the band index (ν = 1 for conduction band and ν = −1 for valence band), → k = |k|(cos(ϕ), sin(ϕ)) is the particles wavevector, and → r = (x, y) is the coordinates vector. The above results show that the physical properties of the system, e.g., the energy spectrum, are reference-frame independent and are clearly reminiscent of the rotational invariance of a bulk (monocrystalline) graphene sheet. Up to now, we have limited our attention to the problem of a single graphene sheet and its description within a rotated reference frame. We have found that the rotated Hamiltonian depends on a phase factor e ±iθ , which is not observable. Hereafter, we study the junction depicted in Figure 1, where graphene sheets with different reticular orientations are connected by means of a grain boundary region. Setting the RF1 as a global reference frame, the Hamiltonian of the entire system H gb takes the piecewise form: The line x = 0 defines the grain boundary region in the limit W → 0 , which is an appropriate approximation when the grain boundary extension W along the x-direction is negligible. This requirement is always met by real grain boundaries, which are defected regions where W covers at least a few lattice sites. The latter observation suggests that the scattering properties of the grain boundary region are adequately captured by appropriate boundary conditions within the framework of a continuous model. Such boundary conditions are non-standard, and are the object of the subsequent analyses. Indeed, when the grain boundary problem formalized so far is treated using ordinary boundary conditions, a violation of the current conservation is found, which is a clear indication that modified boundary conditions are required. Thus, the solution of the problem first requires the derivation of the current density conservation law.

Current Density Conservation
Before treating the grain boundary problem, let us consider the current density conservation law originated by the rotated Hamiltonian (4). A continuity equation in the form ∂ t ρ + → ∇· → J = 0 is obtained by taking the time derivative of the charge density ρ = Ψ + Ψ = |Ψ A | 2 + |Ψ B | 2 (the electric charge is omitted) and using in the computation the rotated Dirac equation in the form ∂ t Ψ = (i ) −1 H θ Ψ. After straightforward computation, the current density components, namely J x = Ψ +Ĵ x Ψ and J y = Ψ +Ĵ y Ψ, are easily recognized. Here, the first quantized current density operators take the following form:Ĵ whose structure explicitly depends on the rotation angle θ. It is worth mentioning here that the expressions for the θ = 0 case coincide with the usual relationsĴ x = v Fσx andĴ y = v Fσy , whereσ x andσ y represent ordinary Pauli matrices. In the following discussion, we denote the quantities in Equations (7) and (8) asĴ θ x,y in order to stress the dependence on the rotation angle θ, while notation J 0 x,y will be adopted to indicate the same quantities when θ = 0.

The Mathematical Problem of Boundary Conditions at a Grain Boundary
We are now ready to treat the problem of a grain boundary. This situation is schematized using the Hamiltonian model given in Equation (6). Accordingly, the current density operator does not admit a global definition, while the following piecewise definition: is required. The form of Equation (9) fully explains the ultimate reasons of the failure of usual boundary conditions in describing the scattering problem defined by the Hamiltonian in Equation (6).
Hereafter, we provide a careful explanation of this point. First of all, we explicitly observe that under translational invariance along the y-axis, the y-component of the current density, namely J y = Ψ +Ĵ y Ψ, does not depend on the coordinate y, and thus the quantity ∂ y J y = 0 does not play any role in the continuity equation. Therefore, current density conservation entirely depends on the conservation of the x-component of the current density. The current density conservation problem is related to the invariance of the current density operator under appropriate transformations. This point can be easily understood by observing that, within the transfer matrix formalism, boundary conditions for the spinorial wave function at the x = 0 interface can be written as Ψ(0 + ) = MΨ(0 − ), with M being a 2 × 2 matching matrix. Here, we have introduced the notation 0 + and 0 − , meaning a spatial position close to x = 0 and belonging to the right or left side of the junction, respectively. Within this formalism, the current density on the right (left) side of the interface, namely J R (J L ), takes the form of J R = Ψ + (0 + )Ĵ x (0 + )Ψ(0 + ) (J L = Ψ + (0 − )Ĵ x (0 − )Ψ(0 − )), while current density conservation requires the condition J R = J L . Observing that J R = Ψ + (0 − )M +Ĵ x (0 + )MΨ(0 − ) and using the continuity of the current density at the interface, we get the important relation M +Ĵ x (0 + )M =Ĵ x (0 − ). Here, we explicitly notice that the current density operator in the absence of lattice mismatching at the interface (i.e., θ = 0) is globally defined so thatĴ x (0 + ) =Ĵ x (0 − ), and the usual boundary conditions are recovered. This is not the case for the problem under study whereĴ x (0 + ) =Ĵ x (0 − ), as evident by direct inspection of Equation (9). According to the above arguments, proper boundary conditions for the scattering problem at the grain boundary interface are assigned once an opportune matching matrix M has been identified. Based on the current density conservation law, a matching matrix M has to respect the following matrix equation: for arbitrary choices of the rotation angle, θ. The matching matrix M is completely determined by the properties of the grain boundary local potential U(x) to be added to the Dirac Hamiltonian (6). Such potential is a 2 × 2 Hermitian operator, which is different from zero only inside the grain boundary region (the GB region in Figure 1b), while its specific form depends on the realization of the interface between grains at the atomic level. The latter information is clearly not available within the framework of a long wavelength (continuous) model, and thus, U(x) has to be meant as a phenomenological potential that is able to reproduce the transport properties of the interface. A useful simplification that will be used in the following discussion parameterizes the interface potential as U(x) = Bδ(x) with δ(x) being the Dirac delta function and B as a 2 × 2 Hermitian operator. Once the structure of U(x) = Bδ(x) is known, the matching matrix M can be determined as described in Section 2.5. Figure 1. (a) Schematic of the region close to a grain boundary (GB region). Apart from the grain boundary region, where lattice distortions and vacancies affect the honeycomb lattice structure of graphene, the bulk of the system is described by a regular atomic arrangement in which the lattice orientation is subject to a rotation of a definite angle going from the left to the right side of the junction. (b) The long wavelength limit of the grain boundary junction depicted in panel (a). The graphene sheets 1 and 2, namely GS1 and GS2, represent the two sides of the junction, while the GB region is represented by a region of negligible extension W, whose effects on the scattering properties of the junction are captured by appropriate boundary conditions. The boundary conditions are a direct manifestation of the presence of a grain boundary effective potential U(x). The presence of the GB imposes a lattice rotation of the GS2 compared to the GS1. This effect is easily described using distinct local reference frames, RF1 and RF2, in describing the two sides of the system. Translational invariance along the y-direction of the RF1 is assumed.
An alternative and interesting approach is completely based on the algebraic properties of the matching matrix M. The algebraic method allows the identification of all of the possible matching matrices in the absence of any information on the scattering potential U(x). The generality of this approach, which will be presented in Section 2.4, is quite appealing in describing real systems where precise information on the interface are not available.

Algebraic Classification of the Matching Matrix M for the Grain Boundary Problem
In this subsection, we study the general problem of a grain boundary formed between two graphene sheets whose reticular axes are rotated by the angles γ and θ, respectively. Under this general condition, the grain boundary Dirac Hamiltonian takes the following piecewise form: Moreover, we assume that an unknown grain boundary potential U(x) is present at the interface located at x = 0. Our purpose is deriving the general structure of all of the possible matching matrices M and disregarding the unknown properties of the grain boundary potential. First of all, we do observe that current density conservation takes the general form M +Ĵ θ x M =Ĵ γ x , which is insensitive to the transformation M → e iϕ M , with ϕ being an arbitrary phase factor. The latter property will be taken into account in the following presentation by omitting the arbitrary phase factor.
In order to follow our program, we consider two relevant classes of matching matrices: (i) matching matrices belonging to the special unitary group SU(2) constituted by 2 × 2 unitary matrices with a determinant equal to 1; (ii) matching matrices belonging to the special linear group SL(2,C) constituted by 2 × 2 matrices with a determinant equal to 1 defined over the complex number field.
Let us start with the case of M ∈ SU(2). This class of matching matrices contains the identity matrix M = 1 2×2 which is involved in the matching condition Ψ(0 + ) = Ψ(0 − ) (continuity of the wave function) used, for instance, in describing the Klein tunneling phenomenon through a rectangular electrostatic potential. We solve the matrix equation coming from the current density conservation by substituting the following general parameterization of M ∈ SU(2): where α, β, λ are real parameters to be determined. The solution of the matrix equation allows the identification of the unknown parameters that appear in Equation (12). After straightforward computation, we obtain that a general SU(2) matching matrix takes the form: As expected, M 1 only depends on the relative angle θ − γ = 0 formed by the crystallographic axes of the graphene sheets forming the two sides of the grain boundary junction. Interestingly, M 1 cannot be equal to the identity matrix within the SU(2) parametrization, and thus the usual boundary condition Ψ(0 + ) = Ψ(0 − ) is not compatible with the current density conservation in the presence of a grain boundary. In particular, while the off-diagonal contribution in M 1 can be eliminated by setting λ = 0, diagonal terms are always different from one when θ − γ = 0, with the latter being the condition of interest in the grain boundary case. Ordinary boundary conditions, namely Ψ(0 + ) = Ψ(0 − ), are correctly recovered by setting θ − γ = 0. The physical meaning of the parameter λ will be clarified in Section 2.5.
A similar procedure can be applied to the case of M ∈ SL(2, C), where the unitary condition is relaxed compared to the previous case. Under this assumption, we obtain the following matching matrix parameterizations: with a, b, c, d ≥ 0, |ad − bc| = 1 and σ = sign(ad − bc). An additional matching matrices set, which cannot be obtained from Equation (14) under suitable assumptions, is given by: with a, b, c, d ≥ 0 and ad + bc = 1. Interestingly, the only diagonal matching matrix such that M ∈ SL(2, C) takes the form , with a > 0. Moreover, matching matrices belonging to SL(2, C) admit a description with three parameters, while SU(2) parameterization requires only one parameter.
The matching matrices M 1,2,3 define a rather general set of physically relevant matching matrices whose structure is fully defined by algebraic properties that are the direct manifestation of the conservation laws of the problem. The matching matrix parameterizations given in Equations (13)- (15) are one of the main results of this work. The physical consequences of the matching matrix structure and its algebraic properties will be discussed later on.

Direct Derivation of the Matching Matrix M from the Scattering Potential U(x) = Bδ(x)
In this section, we provide a link between the features of the scattering potential U(x) and the parameters a, b, c, d, λ that were introduced in Equations (13)- (15). To follow this program, we set aside the problem of grain boundary for a moment and focus our attention on the scattering problem described by the Hamiltonian H = H θ→0 + Bδ(x), with B + = B being an Hermitian operator in the sublattice space. As will be discussed in Section 3, this problem is related to the grain boundary problem by a local unitary transformation, and thus it is not collateral in our discussion. First of all, we do observe that deriving proper boundary conditions for a Dirac delta potential in a Dirac equation requires some care [33]. Indeed, the usual procedure that is followed in the case of the Schrödinger equation for massive particles does not work when the quantum operator associated with the particle kinetic energy is represented by a first-order derivative with respect to the space variable. Under this circumstance, the particle wave function Ψ(x) is not constrained to be continuous in , but this conclusion is incorrect and leads to contradictions. To avoid these difficulties, we provide an alternative derivation. Let us start with the stationary Dirac equation HΨ(x, y) = EΨ(x, y), with Ψ(x, y) = Ψ(x)e i k y y due to the translational invariance along the y-direction. Since we are interested in the wave function behavior in close vicinity of x = 0, we only retain the diverging potential in the Dirac equation. This procedure leads to the first-order differential equation , whose formal solution shares similarities with that of a time-dependent Schrodinger problem in which time-ordered operators are involved. Similarly, the formal solution of our problem involves the definition of a space-ordering operator P x [. . .] instead of the usual time-ordering operator. In particular, the formal solution can be derived by integration with the following result: An infinite iteration of Equation (16) leads to the formal solution Ψ(x) = P x exp x x 0T (y)dy Ψ(x 0 ), where the action of the space-ordering operator is defined as P represents the Heaviside step function. An important observation for our derivation is that when the commutator [A(x 1 ), B(x 2 )] is a vanishing quantity, space ordering simply gives P x [A(x 1 )B(x 2 )] = A(x 1 )B(x 2 ). When the above arguments are applied to our original problem using a limit procedure (i.e., x 0 → − and x → , with → 0 ), we immediately get: A direct inspection of Equation (17) x B and its dependence on the scattering potential. In the following, we provide relevant examples of matching matrices based on the analysis of Equation (17). We use the general parameterization for the Hermitian matrix ( v F ) −1 B reported below: with w a,b,X and ϕ being real parameters. Once the scattering potential has been characterized using the parameterization given in Equation (18), the matching matrix can be computed by direct exponentiation of the matrix: Hereafter, we provide explicit examples of the implementation of this procedure and its outcomes.
This case corresponds to consider an interface potential of the form: Following the procedure described before, we obtain the matching matrix: which is an SU(2) matrix (i.e., det(R) = 1 and R + R = 1) parameterized by the dimensionless quantity Λ. The comparison of Equation (21) with Equation (13) shows that the identification R = M is possible by taking −Λ = λ and θ = γ = 0. The above comparison makes a link between the dimensionless potential amplitude Λ in Equation (20) and the parameter λ, which is the free parameter coming from the algebraic analysis carried out in Section 2.4. In this way, the physical meaning of the parameter λ can be recognized.

Case 2: w a =Λ, w b =−Λ and w X =0 (mass term)
This case corresponds to consider a mass term potential of the form: The matching matrix corresponding to this interface potential takes the form: which is not unitary (R + R = 1) and belongs to the SL(2, C), implying that det(R) = 1. It is easy to prove that Equation (23) can be obtained by properly fixing the free parameters that appear in Equation (14).

Case 3: w a = w b =0 and w X =Λ
This case corresponds to consider an interface potential of the form: The physical meaning of this kind of interface potential is explained in Appendix A; here, we just comment that the introduction of an off-diagonal potential in the Dirac equation produces a displacement of the K, K' points in momentum space, with the latter effect being usually the result of lattice deformations, e.g., induced by mechanical strain. The matching matrix corresponding to the interface potential given in Equation (24) takes the form: which is defined by a pure phase factor e −i Λcos(ϕ) multiplied by a determinant one-diagonal matrix belonging to SL(2, C). According to our previous discussion, this kind of matching matrices preserves the current density at the interface; thus, they are admissible matching matrices despite the determinant being different from one. A direct proof of this statement is provided here for the sake of completeness. The Jacobi formula for a generic complex coefficient square matrix C implies that det(exp(C)) = exp(Tr(C)), with Tr(C) being the trace of C. When C is substituted with A in Equation (19), we get det(M) = e −2iw X cos(ϕ) , with M = e A being the matching matrix of the problem. In this way, we have proven that the determinant of a matching matrix is fixed to one only for diagonal interface potentials, while it is a pure phase factor when off-diagonal terms affect the structure of the interface potential (i.e., when w X = 0).

Grain Boundary Hamiltonian Model with Position-Dependent Rotation Angle θ(x)
Before characterizing the scattering properties of a grain boundary, we would like to understand if a phase gradient alone can generate a non-vanishing scattering potential at the grain boundary junction.
To answer this question, we need a convenient model in which the phase profile at the interface is a smooth function of the position. Thus, a regularization of the Hamiltonian model presented in Equation (11) is needed. A regularized Hamiltonian model is obtained by taking the rotation angle θ as a position-dependent function θ(x) defining a smeared step function; it provides a gentle matching between the lattice rotation angle θ L on the left of the junction and θ R , characterizing the right side. In this way, a finite phase gradient is present inside the grain boundary region (i.e., the region close to x = 0). The resulting Hamiltonian model cannot be obtained just by making the substitution θ → θ(x) in Equation (4). This is because, according to quantum mechanics prescriptions, terms such as e iθ(x)p x have to be symmetrized to avoid non-Hermitian contributions to the Hamiltonian. After the symmetrization procedure has been implemented, the following regularized grain boundary Hamiltonian is obtained: where we have introduced the anticommutator Â ,B =ÂB +BÂ between generic quantum operatorsÂ andB. It is easy to demonstrate that the dependence on θ(x) in Equation (26) can be eliminated using the local unitary transformation: In particular, as anticipated in Section 2.5, we obtain U , with the latter result being independent of the specific functional form of θ(x). Moreover, terms in the Hamiltonian Equation (26) depending on the phase gradient ∂ x θ(x) are exactly cancelled by the unitary transformation. According to these arguments, the mismatch of the crystallographic axes at the grain boundary of a Dirac material does not provide disturbance for particles transmission, at least in the framework of the present model. This conclusion crucially depends on the first-order nature of the differential operator representing the kinetic energy of the Dirac particles.
When the Hamiltonian presented in Equation (26) includes a generic interface potential of the form: with w a , w b , and w X being dimensionless real parameters. The transformed potential W explicitly contains a dependence on the phase value θ(0) being assumed at the grain boundary position x = 0. This behavior is specific to the Dirac Delta potential, which is appropriate to describe the grain boundary interface in the long wavelength limit. Here, we notice that the relation between θ(0) and θ R,L depends on the details of the θ(x) dependence. To be specific, a simple model of θ(x) can be expressed, e.g., by the following piecewise function: according to which θ(0) = (θ L + θ R )/2. In Equation (30), where the limit d → 0 has to be meant, x ∈ − d 2 , d 2 represents the grain boundary region. As a final comment, we argue that the scattering problem related to the transformed Hamiltonian H θ(x)=0 + W can be treated using ordinary boundary conditions, as described in Section 2.5.

Transmission Properties of a n/n' Grain Boundary Junction with Unconventional Boundary Conditions
The purpose of this section is to provide simple examples of the use of unconventional matching matrices, and create a link between the transmittance properties of the problem and the specific matching matrix that is adopted. A complete analysis of all of the conduction regimes of the junction (e.g., the p/n case) is beyond the scopes of this section. Thus, we study the scattering problem of a n/n' grain boundary junction described by the Hamiltonian Equation (6). The peculiar properties of the transmittance are derived using matching matrices with a simplified structure compared to those described in Section 2.4.
The scattering states in the left and right part of the junction are expanded in terms of local eigenstates, taking into account the translational invariance of the scattering region along the y-direction. According to this procedure, the scattering states, which are normalized by the absolute value of the group velocity v L,R x along the propagation direction, can be written in the following general form: Here, Ψ L (x, y) describes the reflection of a particle belonging to the conduction band (E > 0) impinging from the left on the grain boundary interface with a group velocity of v L x = v F cosφ and the linear momentum k x = E v F cosφ, φ ∈ [−π/2, π/2], in which the incidence angle is formed by the particle trajectory with the x-axis. The expression of v L x is related to the dispersion relation on the left side of the junction, namely E = v F k 2 x + k 2 y , by the usual relation v L x = −1 ∂ k x E. In writing the reflected wave contribution (weighted by the reflection coefficient R), we have explicitly taken into account that the reflection angle φ r is related to the incidence angle φ by the relation φ r = π − φ. In the right part of the junction (x > 0), the wave function Ψ R (x, y) describes a transmitted electron with group velocity v R x = v F cosφ t and linear momentum q x = E v F cosφ t . Since we are assuming the alignment of the Dirac points of the left and right part of the junction, the dispersion relation of the transmitted electron is given by E = v F q 2 x + q 2 y . Here, it is worth mentioning that the spinor of the transmitted particle explicitly depends on the rotation angle θ characterizing the rotation of the cristallographic axes of the right part of the junction. In the elastic scattering of a translational invariant system along the y-direction, both the energy E and the linear momentum along the y-direction are conserved quantities. This observation implies that k y = q y . Observing that k y = E v F sinφ and q y = E v F sinφ t , we conclude that φ = φ t . The scattering states in Equation (31) are used to derive the scattering problem equation Ψ R (x = 0, y) = MΨ L (x = 0, y), whose structure is fixed once the matching matrix M has been specified. Once the scattering equation has been solved, the scattering coefficients R and T are determined; moreover, current density conservation implies |R| 2 + |T | 2 = 1. The transmittance |T | 2 is an angle-resolved quantity that is directly related to the scattering properties of the interface and determines the differential conductance of the system. In the following, several relevant cases are treated, while the transmittance curves as a function of the incidence angle φ are presented in Figure 2a-c.  (32) and setting the model parameters as specified in the legend. The matching matrix that is used belongs to SU(2), and accordingly, full transmission is observed for φ = 0 incidence, with the latter result being insensitive to the parameter choice and reminiscent of the Klein tunneling phenomenon. (b) Polar plot of the transmittance (transmission probability) |T (φ)| 2 versus the incidence angle φ computed according to Equation (36) and setting the model parameters as specified in the legend. Diagonal SL(2, C) matching matrix (see Equation (35)) has been used. A reduction of the transmittance is observed for arbitrary values of the incidence angle as the parameter λ is increased. (c) Polar plot of the transmittance (transmission probability) |T (φ)| 2 versus the incidence angle φ computed according to Equation (38) and setting the model parameters as specified in the legend. A non-diagonal SL(2, C) matching matrix (see Equation (37)) has been used. An anisotropic reduction of the transmittance is observed for the arbitrary values of the incidence angle as the parameter µ is increased. Preferential transport directions are obtained.

Case 1: SU(2) matching matrix (Equation (13))
We have solved the scattering problem using the SU(2) matching matrix given in Equation (13) by fixing γ = 0. The transmittance of the problem has been derived and takes the following (energy-independent) form: The main property of Equation (32) is that |T (φ = 0)| 2 = 1 for any choice of the parameter λ (see Figure 2a). This behavior, which is reminiscent of the Klein tunneling, is clearly inappropriate to describe the grain boundary physics for which a reduction of the transmittance is expected when the strength of the scattering potential (controlled by λ) is increased. No dependence on θ is present in the transmittance, since the transmission coefficient T (φ) presents an irrelevant prefactor e iθ/2 .

Case 2: Matching matrix belonging to SL(2,C); example 1
We have solved the scattering problem using the SL(2, C) matching matrix: which is a special case of the matching matrix given in Equation (14), and it is appropriate to describe the interface potentials proportional toσ z (mass term). The solution of the scattering problem provides the scattering coefficients: implying an angle-independent transmittance given by |T (φ)| 2 = 1/ 1 + λ 2 , which is a decreasing function of the scattering strength parameter λ. This kind of matching matrix seems to be appropriate to describe the conductance reduction induced by a grain boundary region, even though no dependence on θ is detected. The above findings confirm the confining properties of a potential proportional toσ z .

Case 3: Matching matrix belonging to SL(2,C); example 2
We have solved the scattering problem using the diagonal SL(2, C) matching matrix: The transmittance of the problem has been derived, and takes the following (energy-independent) form: A reduction of the transmittance (see Figure 2b) is observed when the scattering strength parameter λ is increased, with the latter behavior being appropriate to describe the grain boundary physics.

Case 4: Matching matrix belonging to SL(2,C); example 3
We have solved the scattering problem using the following SL(2, C) matching matrix: which depends on the dimensionless parameters λ and µ. The transmittance of the problem has been derived, and takes the following (energy-independent) form: A reduction of the transmittance is observed when the scattering strength parameter λ is increased, with the latter behavior being appropriate to describe the grain boundary physics. Moreover, depending on the choice of parameters, preferential transport directions are obtained (see Figure 2c).
The analysis of cases 1-4 shows that the grain boundary transmittance can be independent of θ, despite the matching matrix and the spinorial wave function (see Equation (31)) presenting a dependence on this parameter. Furthermore, direct computation (not reported here) shows that θ-independent transmittance is not a peculiarity of the single barrier case, but rather survives in the double barrier case even when the different regions of the conduction channel are described by a Dirac Hamiltonian rotated by different angles. The above behavior suggests that the algebraic deduction of the matching matrix structure (see Section 2.4) is not able to capture the hidden dependence on the rotation angle θ, which in principle can be important.
To clarify this point, in the subsequent section, we analyze a simple model of the n/n' grain boundary junction based on the Hamiltonian treatment proposed in Section 3. The physical conditions under which the system presents a θ-dependent differential conductance are studied.

Grain Boundary Junction with θ-Dependent Differential Conductance
We now formulate a theory of the n/n' grain boundary junction based on the Hamiltonian model presented in Section 3. Let us consider the grain boundary junction described by the Hamiltonian mimicking a potential profile induced by charge transfer at the grain boundary. Such potential is assumed to be diagonal in the sublattice representation, I 2×2 = diag(1, 1) representing the identity in the sublattice space. Using the unitary transformation in Equation (27), the Hamiltonian problem can be written in the equivalent form is given in Equation (29), while V s (x) is invariant under unitary transformation. The junction described by H can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set  (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set ( ) = ( ) = . The latter relation allow us to deduce the moduli of the linear momenta on the different junction sides, namely | | = /ℏ and | | = ( − 0 )/ℏ with − 0 > 0 , since we are considering an n/n' junction. The  (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set  (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set  (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set  (29), while ( ) is invariant under unitary transformation. The junction described by ℋ can be treated using ordinary boundary conditions, which are implemented using the matching matrix formulation given in Section 2.5. Since we are interested in verifying the existence of preferential transport directions, we focus our analysis on the interface potential: which is obtained from Equation (29) by setting = 0. The matching matrix associated with Equation (40) is given by: which defines a (2, ) boundary condition with parameters ℊ = (0) , ℊ * = − (0) , and = . Equation (41) presents an explicit dependence on the rotation angle (0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely ( / ) , are different, and in particular we have that ( ) = ℏ √ 2 + 2 = ℏ | | for < 0 and ( ) = ℏ √ 2 + 2 + 0 = ℏ | | + 0 for > 0. The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy of a particle is conserved during the scattering event, and thus, we can set ( ) = ( ) = . The latter relation allow us to deduce the moduli of the linear momenta on the different junction sides, namely | | = /ℏ and | | = ( − 0 )/ℏ with − 0 > 0 , since we are considering an n/n' junction. The y-component of the linear momentum on different sides of the junction is thus given by = * = w X e −iθ(0) , and w a = αw X . Equation (41) presents an explicit dependence on the rotation angle θ(0) at the grain boundary junction, and thus a dependence of the transport properties of the system on this parameter is expected.
In order to present the results of the model, we focus our attention on the n/n' junction described in Figure 3a,b. The scattering problem related to the conduction properties of the junction is solved using the scattering states given in Equation (31) with θ = 0. However, differently from the case treated in Section 4.1, the dispersion relation in the distinct sides of the junction, namely E (L/R) k , are different, and in particular we have that E The group velocities are the same as considered before, and remain unaffected by the presence of a potential step at the interface. Energy E of a particle is conserved during the scattering event, and thus, we can set E The latter relation allow us to deduce the moduli of the linear momenta on the different junction sides, namely |k| = E/ v F and |q| = (E − V 0 )/ v F with E − V 0 > 0, since we are considering an n/n' junction. The y-component of the linear momentum on different sides of the junction is thus given by k y = |k|sinφ and q y = |q|sinφ t , while translational invariance implies k y = q y . Conservation of the y-component of the linear momentum implies the following relation between the incidence and the transmission angle: Figure 3. (a) Bands alignment of the n/n' junction model considered in the main text. Grain boundary interface is located in x = 0. The system Fermi level is located above the Dirac points, and thus, a n-type conduction regime is established. (b) In the translational invariant system along the y-direction, the corresponding linear momentum is a conserved quantum number in the scattering events. When the incidence angle φ is greater than the critical angle φ c , the y-component of the linear momentum cannot be conserved, and no current can flow through the interface.
Since the quantity E E−V 0 > 1, there exists a critical angle of incidence for which Equation (42) cannot be satisfied, and no current can flow through the interface. Such a critical value takes the following form: Once the scattering problem has been solved using the matching matrix given in Equation (41), the angle-resolved transmittance |T (φ)| 2 is derived. The zero-temperature differential conductance of the junction G nn is related to the transmittance evaluated at the Fermi level, E = E F by the following relation: with g s g v = 4 being a degeneracy factor coming from the spin (g s = 2) and the valley (g v = 2) degeneracy, L being the transverse dimension of the junction assumed to be a macroscopic quantity, and k F = E F /( v F ) being the modulus of the Fermi wave vector. Here is worth mentioning that the Fermi level can be tuned by using a back gate with the twofold effect of modulating the number of transverse channels N (E F ) = g s g v k F L 2π involved in the conduction and changing the value of the critical angle φ c , which is an energy-sensitive quantity.
Before treating the general case of V 0 = 0, we study the solution of the scattering problem under the simplifying assumption V 0 = 0, which allow us to obtain the following simple expression for the angle-resolved transmittance: A direct inspection of Equation (45) evidences a dependence on θ(0), which is the relic of the cristallographic axes mismatching between the two sides of the grain boundary junction. In the following analysis, we use the model for θ(x) given in Equation (30) with θ L = 0. Accordingly, we set θ(0) = θ R /2. Furthermore, the analysis of the formation energy of a graphene grain boundary as a function of the misorientation angle θ R shows an M-shaped behavior with a local minimum at θ R ∼ = 30 • and absolute minima in θ R ∼ = 0 • and θ R ∼ = 60 • [19]. In Figure 4a-c, we have studied the behavior of Equation (45) as a function of the relevant junction parameters. A reduction of the transmission probability is observed as the misorientation angle at the grain boundary θ R is increased. Preferential transport directions are clearly related to a non-vanishing value of the parameter α. The general expression of the transmittance pertaining to the V 0 = 0 case is quite lengthy, and thus, we only present a numerical evaluation of it. This analysis is performed in Figure 5a-c, where the effect on the transmittance of a finite potential step at the interface is presented, setting the values of the remaining junction parameters as done in Figure 4b. Here, we observe that the critical angle φ c of the scattering problem is an energy-sensitive quantity that is reduced when the potential step intensity is increased. Figure 5 clearly shows the effect of the critical angle reduction, and suppression of the junction transmission at high incidence angles is clearly visible.
In Figure 6a-d, we report the differential conductance given in Equation (44) normalized to the quantity g s g v e 2 h k F L 2π ≡ G 0 , with the resulting dimensionless quantity being g nn = G nn /G 0 . In particular, in Figure 6a,b, we show g nn as a function of the misorientation angle θ R for different values of the potential step and of the matching matrix parameters. A suppression of the junction conductance as a function of θ R is observed, in which the latter behavior is more pronounced for small values of the potential step intensity. In Figure 6c,d we show g nn as a function of the potential step intensity V 0 . Different curves are related with distinct choices of the misoriantation angle of the grain boundary junction. An almost linear lowering of the g nn versus V 0 curves is obtained, which is related to the critical angle reduction causing the progressive closure of the conduction channels. The general conclusion associated with the analysis performed in Figure 6 is that grain boundary junctions with large misorientation angles manifest the tendency to be more opaque compared to those with a better lattice matching. Moreover, charge transfer at the grain boundary interface induces a potential step whose effect is relevant in determining the transparency of the interface.

Conclusions
We have proposed a continuous model to study the grain boundary effects in graphene. The model provides a description of the grain boundary based on a Dirac Hamiltonian written in a rotated side-dependent reference frame describing crystallographic axes mismatching at a grain boundary junction. We have shown that the scattering problem related to the transmission properties of a grain boundary junction requires modified boundary conditions, which can be implemented by using the matching matrix method in the scattering problem. We have characterized the algebraic properties of all of the possible matching matrices showing that the SU(2) and SL(2, C) matching matrices are admissible. In particular, we have proven that the SU(2) matching matrices support Klein tunneling, and are not adequate in describing the conductance lowering that is associated with a grain boundary physics, which is instead captured by SL(2,C) matching matrices. We have studied specific matching matrix examples and the associated transmittance properties. It has been demonstrated that under opportune assumptions, preferential transport directions are supported by the interface microscopic properties, which are encoded in our approach by the matching matrix structure.
Moreover, using a space-dependent unitary transformation, we have formulated a grain boundary model where usual boundary conditions can be used. Conditions to observe scattering properties related to the crystallographic axes mismatching at the grain boundary interface are studied. Reduction of the junction conductance has been demonstrated as the effect of the misorientation angle between the two junction sides.
Compared to alternative theoretical methods [27][28][29][30], the present approach does not require a tight-binding model of the linear defect (grain boundary) and boundary conditions, which are expressed via the matching matrix method, are derived using the algebraic properties of the scattering problem.
The proposed theory provides a phenomenological model to study grain boundary physics within the scattering approach, and represents per se an important enrichment of the scattering theory of graphene, having the potential to stimulate new experiments and further theoretical investigations.
The spectrum of the problem is obtained by solving the equation det[H(k) − 1E] = 0 with respect to E. Straightforward computation shows that the dispersion relations of the conduction and valence bands take the following form: which clearly shows a dependence on the rotation angle θ which is absent when the off-diagonal potential is switched off (V = 0). Moreover, the Dirac point, which is defined as the point in the momentum space where E + = E − , is displaced from the original K point in momentum space by the vector . In real systems, Dirac cone displacement from the unperturbed K point is related to strain-induced change in the graphene electronic structure [34]. Since grain boundary regions present defective lattice accompanied by mechanical deformations, the off-diagonal potential introduced in (A1) can well emulate strain-induced modifications of the electronic properties. Clearly, in modeling these effects non-vanishing off-diagonal potential is only admissible in close vicinity of the grain boundary region, while it is expected to be irrelevant when an ordered (without defects and vacancies) graphene lattice is restored. The dependence of the spectrum in (A2) on the rotation angle Accordingly, the spectrum of H θ (k) does not depend on θ. The same conclusion is valid in the presence of a diagonal potential. To demonstrate this point, let us apply the same argument to the Hamiltonian H = H θ (k) + W, with We observe that the transformed Hamiltonian takes the form U + HU = H θ=0 (k) + W since W is invariant under unitary transformation, i.e., U + WU = W. The above observation implies that the energy spectrum remains insensitive to the rotation angle even in the presence of a diagonal single-particle potential W.
On the other hand, when the above procedure is applied to the Hamiltonian H(k) = H θ (k) +V considered in (A1), we obtain U + H(k)U = H θ=0 (k) + U +V U with the latter being the transformed potential. Since U +V U retains a dependence on the rotation angle, a θ-dependent energy spectrum is obtained for the Hamiltonian in (A1). The above results suggest that the relative angle formed by the crystallographic axes of the two sides of a grain boundary junction can affect the scattering properties of the system in the presence of an off-diagonal potential in the sublattice indices. This conclusion is clearly supported by the arguments presented in Section 4.2.