Diffusion of Dirac fermions across a topological merging transition in two dimensions
DDiffusion of Dirac fermions across a topological merging transition in two dimensions
P. Adroguer, D. Carpentier, G. Montambaux, and E. Orignac Laboratoire de Physique, ´Ecole Normale Sup´erieure de Lyon, 46 all´ee d’Italie, 69007 Lyon, France Laboratoire de Physique des Solides, CNRS, Universit´e Paris-Sud,Universit´e Paris-Saclay, 91405 Orsay Cedex, France (Dated: November 5, 2018)A continuous deformation of a Hamiltonian possessing at low energy two Dirac points of oppositechiralities can lead to a gap opening by merging of the two Dirac points. In two dimensions, thecritical Hamiltonian possesses a semi-Dirac spectrum: linear in one direction but quadratic in theother. We study the transport properties across such a transition, from a Dirac semi-metal througha semi-Dirac phase towards a gapped phase. Using both a Boltzmann approach and a diagrammaticKubo approach, we describe the conductivity tensor within the diffusive regime. In particular, weshow that both the anisotropy of the Fermi surface and the Dirac nature of the eigenstates combineto give rise to anisotropic transport times, manifesting themselves through an unusual matrix self-energy.
I. INTRODUCTION
The discovery of graphene has triggered a lot of workon the exotic transport properties of Dirac-like particlesin solids . Indeed, the graphene electronic spectrum ismade of two sub-bands which touch at two inequiva-lent points in reciprocal space. Near the touching points,named Dirac points, the spectrum has a linear shape andthe electron dynamics is well described by a 2D Diracequation for massless particles. Due to the structure ofthe honeycomb lattice, the wave functions have two com-ponents corresponding to the two inequivalent sites of thelattice, and the Hamiltonian is a 2 × H = ± c (cid:126)p . (cid:126)σ (1)where the velocity c (cid:39) m.s − . This linearization ispossible because the energy of the saddle point separatingthe two Dirac cones (valleys) is very large ( (cid:39) I under pressure , and has been observed in artificiallyassembled nanostructures and ultracold atoms .Besides these two dimensional realizations, the existenceand properties of semi-metallic phases in three dimen-sions have recently been studied .To go beyond and in order to account for a structurewhich consists in two Dirac points separated by a saddlepoint, one needs an appropriate low energy 2 × . Thismerging is therefore a topological transition. It has beenshown that, at the transition, the electronic dispersion isquite unusual since it is quadratic in one direction andlinear in the other direction (the direction of merging). This ”semi-Dirac” spectrum has new properties inter-mediate between a Schr¨odinger and a Dirac spectrum.The vicinity of the topological transition can be describedby the following Hamiltonian in two dimensions : H = (cid:18) ∆ + p x m (cid:19) σ x + c y p y σ y . (2)It has been coined ”Universal Hamiltonian” sincethe merging scenario of two Dirac points related bytime reversal symmetry is uniquely described by thisHamiltonian . The parameter ∆ drives the transition(∆ = 0) between a semi-metallic phase (∆ <
0) with twoDirac points and a gapped phase (∆ > . (a) (b) S M (c) ² ~k k x k y ² F ² F ² F ² F ² F Figure 1. This work addresses the transport properties for anelectronic spectrum undergoing a topological merging transi-tion as depicted in this figure, and commented in more detailsin Fig. 2.
In this paper we address the evolution of the conduc-tivity tensor across the merging transition (Figs. 1,2). Afirst objective of this work is to characterize the trans-port properties as a possible signature of the evolution of a r X i v : . [ c ond - m a t . m e s - h a ll ] O c t the underlying band structure. On a more fundamen-tal perspective, an additional interest of this problemstands from two important ingredients in the descriptionof diffusive transport. First, at low energy the electronicwave functions have a spinorial structure which leads toeffective anisotropic scattering matrix elements (similarto the case of a scalar problem with anisotropic scatter-ing due to a disorder potential with finite range). Thisleads to a transport scattering time τ tr different from theelastic scattering time τ e , as in graphene for point-likeimpurities where τ tr = 2 τ e . Second, the anisotropy ofthe dispersion relation leads to an additional complex-ity: the scattering times become themselves anisotropicand depend on the direction of the applied electric field.We show that within the Green’s function formalism thisanisotropy manifests itself into a rather unusual matrixstructure of the self-energy. A comparison between aBoltzmann approach and a perturbative Green’s func-tion formalism allows for a detailed understanding of thisphysics.The outline of the paper is the following. In the nextsection, we recall the model, i.e. the Universal Hamilto-nian with coupling to impurities described by a point-likewhite noise potential. We define a directional density ofstates and derive the angular dependence of the elasticscattering time. In section III, we use the Boltzmannequation to calculate the conductivity tensor. As a re-sult of the two important ingredients mentioned above,the conductivity along a direction α is not simply propor-tional to the angular averaged squared velocity (cid:104) v α ( θ ) (cid:105) because : (i) the elastic scattering time has also anangular dependence due to the angular anisotropy ofthe spectrum, so that one should consider the average (cid:104) v α ( θ ) τ e ( θ ) (cid:105) ; (ii) since the matrix elements of the interac-tion get an angular dependence, it will lead to transporttimes different of the elastic time. These transport timesdepend on the direction α and, to obtain the conduc-tivity, we will have to consider the average (cid:104) v α ( θ ) τ tr α ( θ ) (cid:105) .These results obtained from Boltzmann equation are con-firmed by a diagrammatic calculation presented in sec-tion IV. We discuss our results in the last section. II. THE MODELA. Hamiltonian and Fermi surface parametrization
We consider the model described by the Hamiltonian H = H + V, (3)where the disorder potential V is defined and discussedin section II C and the Hamiltonian for the pure systemis defined as H = (cid:20) ∆ + p x m (cid:21) σ x + c y p y σ y . (4)In the present and the following sections (II A and II B)we start by discussing a few properties of the Hamilto- MS ( a ) ( b ) ( c ) ✏ ✏ ~k < > Figure 2. Typical energy spectrum of the model (4) for various∆ but a fixed energy (cid:15) >
0. Dirac phase with (a) ∆ < − (cid:15) , (S)∆ = − (cid:15) (saddle-point) , (b) − (cid:15) < ∆ <
0. Critical semi-Diracmetal (M) ∆ = 0. Gapped phase (c) ∆ > nian H without disorder. For ∆ > <
0, it describes twoDirac cones with opposite chiralities, hereafter named aDirac phase. Note that these Dirac cones are in generalanisotropic with respectives velocities in the x and y di-rections c x = (cid:112) | ∆ | /m and c y . The energy spectrum isgiven by (cid:15) = (cid:18) p x m + ∆ (cid:19) + ( c y p y ) . (5)We will consider only the case of positive energies (cid:15) > (cid:15) < p x parity.For each half plane p x ≶ p x m + ∆ = (cid:15) cos θ ; c y p y = (cid:15) sin θ ; η p = sign( p x ) = ± (6)where θ ∈ [ − θ , θ ] is a coordinate along the constantenergy contour. Its range depends on the topology ofthe constant energy contour, and thus on the energy (cid:15) ,see Fig. 3. Specifying the discussion to the Fermi surfaceassociated with the Fermi energy (cid:15) F , we can distinguishtwo cases :(i) Low energy metal with two disconnected Fermi sur-faces when ∆ < (cid:15) F < − ∆. In this case θ = π . This corresponds to the energy spectrum( a ) of Fig. 2.(ii) High energy metal with a single connected Fermisurface for (cid:15) F > | ∆ | . In this case cos θ = ∆ /(cid:15) F .For ∆ < θ varies from π for (cid:15) F = − ∆ to π/ (cid:15) F (cid:29) − ∆. For ∆ > θ varies from π/ (cid:15) F (cid:29) ∆ to 0 for (cid:15) F → ∆ . This corresponds to theenergy spectra ( b ) , M, ( c ) of Fig. 2.The eigenstates of positive energy corresponds to wave ✓ ✓ ✓ = ⇡ ( a )( b ) S k x k y < ✓ ✓ k x k y M ( c ) Figure 3. Left : Constant energy contours (cid:15) ( k x , k y ) for different energies and ∆ < θ parametrizing in a unambiguous way eachhalf k x > k x < (cid:126)k and energy (cid:15) . Right : Same quantities at the merging point ∆ = 0 (M) andfor (cid:15) ( k x , k y ) > ∆ > c ). functions conveniently expressed with the parametriza-tion of the constant energy contour ψ (cid:126)k ( (cid:126)r ) = 1 √ L (cid:18) e iθ (cid:126)k (cid:19) e i(cid:126)k.(cid:126)r , (7)where θ (cid:126)k is defined by inversion of Eq. (6), and (cid:126)p = (cid:126) (cid:126)k .From now on, we will set L = 1. The group velocityvaries along the constant energy contour according to : v x ( η p , θ ) = η p (cid:114) (cid:15)m cos θ √ cos θ − δ , (8a) v y ( η p , θ ) = c y sin θ , (8b)where throughout this paper we use the reduced parame-ter δ = ∆ /(cid:15) . The evolution of the velocity along constantenergy contours is shown on Fig. 4. Figure 4. Velocity (cid:126)v ( (cid:126)k ) along constant energy contours (here∆ < B. Density of States
We define a directional density of states along the con-stant energy contour parametrized by θ from the equality (cid:90) dp x dp y (2 π (cid:126) ) = (cid:90) dk x dk y (2 π ) = (cid:90) ρ ( (cid:15), θ ) d(cid:15)dθ , (9)with ρ ( (cid:15), θ ) = √ m(cid:15) (2 π (cid:126) ) c y √ cos θ − δ . (10)The density of states in then obtained by the integral ρ ( (cid:15) ) = 2 (cid:90) θ − θ dθ ρ ( (cid:15), θ ) (11)where the extra factor 2 accounts for the sign of p x . Theintegral gives ρ ( (cid:15) ) = √ m(cid:15) (2 π (cid:126) ) c y I ( δ ) , (12)with the function I ( δ ) = (cid:90) θ − θ dθ √ cos θ − δ , (13)where θ = Arccos( δ ) when | δ | < θ = π otherwise.From Eqs. (10,12), we can rewrite ρ ( (cid:15), θ ) as ρ ( (cid:15), θ ) = ρ ( (cid:15) )2 I ( δ ) √ cos θ − δ . (14) C. Disorder Potential and Elastic Scattering Time
The disorder part V of the Hamiltonian accounts forthe inhomogeneities in the system. This random poten-tial V ( (cid:126)r ) is assumed to describe a gaussian point-like un-correlated disorder, characterized by two cumulants V ( (cid:126)r ) = 0 , V ( (cid:126)r ) V ( (cid:126)r (cid:48) ) = γ δ ( (cid:126)r − (cid:126)r (cid:48) ) . (15)where the overline denotes a statistical average over re-alizations of the random potential. The presence of thisrandom potential induces a finite lifetime for the eigen-states of momentum (cid:126)k of the pure model (4), called elas-tic scattering time, and obtained from the Fermi goldenrule : (cid:126) τ e ( (cid:126)k ) = 2 π (cid:90) d (cid:126)k (cid:48) (2 π ) δ ( (cid:15) (cid:126)k − (cid:15) (cid:126)k (cid:48) ) |A ( (cid:126)k, (cid:126)k (cid:48) ) | , (16)where the scattering amplitude is defined by A ( (cid:126)k, (cid:126)k (cid:48) ) = (cid:104) ψ (cid:126)k | V | ψ (cid:126)k (cid:48) (cid:105) . (17)For uncorrelated point-like disorder, the angular depen-dence of this scattering amplitude originates from theeigenstates overlap and one has |A ( (cid:126)k, (cid:126)k (cid:48) ) | = γ (cid:0) θ (cid:126)k − θ (cid:126)k (cid:48) ) (cid:1) . (18)Defining τ e ( (cid:15), θ ) = τ e ( (cid:126)k ) where (cid:15), θ and (cid:126)k are relatedthrough Eq. (6), we can express the elastic scatteringtime as an integral (cid:126) τ e ( (cid:15), θ ) = 2 πγ (cid:90) θ − θ dθ (cid:48) ρ ( (cid:15), θ (cid:48) ) [1 + cos( θ − θ (cid:48) )] . (19)Introducing the bare scattering time τ e ( (cid:15) ) = (cid:126) πγρ ( (cid:15) ) , (20)we can rewrite (19) in the form τ e ( (cid:15), θ ) = τ e ( (cid:15) )1 + r ( δ ) cos θ , (21)where the density of states ρ ( (cid:15) ) is given by (12). Thedenominator of this expression exactly accounts for theanisotropy of the scattering time. As a convenientparametrization of this property, we have introduced theanisotropy function r ( δ ) which will be used throughoutthis paper : r ( δ ) = J ( δ ) /I ( δ ) , (22)with the function I ( δ ) defined in (13) and J ( δ ) = (cid:90) θ − θ dθ cos θ √ cos θ − δ , (23) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) ∆ r (cid:72) ∆ (cid:76) Figure 5. Function r ( δ ) parametrizing the angular depen-dence of the elastic scattering time τ e plotted as a function of δ = ∆ /(cid:15) . It has the limits r ( δ → −∞ ) (cid:39) / (4 δ ), r ( −
1) = − r (0) = 2Γ(3 / /π (cid:39) . r (1) = 1. In this figure,as in following figures, we systematically reserve the colors :blue for the Dirac phase ( δ < >
0) and red for the semi-Dirac point. where θ = Arccos( δ ) when | δ | < θ = π otherwise.The function r ( δ ) is plotted in Fig. 5. Deep in the Diracphase (∆ (cid:28) (cid:15) (cid:28) | ∆ | ), one has r ( δ ) → − | δ | = − (cid:15) | ∆ | (cid:28) , (24)so that the anisotropy can be neglected in Eq. (21) andwe recover a scattering time independent of the directionof propagation as standard for Dirac fermions. III. DIFFUSIVE REGIME FROM THEBOLTZMANN EQUATION
We now consider the transport properties of the model(3) at a fixed energy (cid:15) large enough so that the condi-tion kl e (cid:29) l e being a typical elastic meanfree path. Therefore we will not consider the close vicin-ity of a Dirac point and the associated physics of mini-mal conductivity . For a system of typical size muchlarger than this mean free path l e , this corresponds to theregime of classical diffusion. We describe this regime firstwith a standard Boltzmann equation, before turning toa complementary but equivalent diagrammatic approachbased on Kubo formula for the conductivity. The use ofthese two approaches will reveal the physics hidden be-tween the technical specificities of the diffusive transportfor the model we consider. A. Boltzmann equation
We start from the Boltzmann equation expressingthe evolution of the distribution function f ( (cid:126)k, (cid:126)r ) : dfdt + d(cid:126)rdt ∇ (cid:126)r f + d(cid:126)kdt ∇ (cid:126)k f = I [ f ] , (25)where I [ f ] is the collision integral defined below. The po-sition (cid:126)r and momentum (cid:126)k parametrizing the distributionfunction f ( (cid:126)k, (cid:126)r ) are classical variables, whose time evolu-tions entering Eq.(25) are described by the semi-classicalequations d(cid:126)rdt = (cid:126)v ( (cid:126)k ) + d(cid:126)kdt × (cid:126)F (cid:126)k (26) (cid:126) d(cid:126)kdt = − e (cid:126)E − e d(cid:126)rdt × (cid:126)B (27)with the group velocity (cid:126)v ( (cid:126)k ) = (cid:126) − ∂(cid:15) ( (cid:126)k ) /∂(cid:126)k , (cid:126)B is a lo-cal magnetic field, and (cid:126)F (cid:126)k = i ∇ (cid:126)k × (cid:104) ψ (cid:126)k |∇ (cid:126)k ψ (cid:126)k (cid:105) is theBerry curvature. In the present case, we consider the re-sponse of the distribution function f ( (cid:126)k, (cid:126)r ) due to a uni-form weak electric field (cid:126)E : we can neglect the gradient ∇ (cid:126)r f in Eq. (25) and drop the spatial dependence of f .Due to the absence of magnetic field, we deduce fromEqs. (25,27) that a stationary out-of-equilibrium distri-bution f ( (cid:126)k ) satisfies the simpler equation − e (cid:126) (cid:126)E. ∇ (cid:126)k f = I [ f ] . (28) where f is now as function of (cid:126)k and the collision integralis expressed as I [ f ] = 2 π (cid:90) d (cid:126)k (cid:48) (2 π ) δ ( (cid:15) (cid:126)k − (cid:15) (cid:126)k (cid:48) ) |A ( (cid:126)k, (cid:126)k (cid:48) ) | (cid:16) f ( (cid:126)k (cid:48) ) − f ( (cid:126)k ) (cid:17) . (29)By assuming the perturbation to be weak, we canexpand the stationary out-of-equilibrium distribution f ( (cid:126)k ) around the equilibrium Fermi distribution f ( (cid:126)k ) = n F ( (cid:15) (cid:126)k ) following the ansatz f ( (cid:126)k ) = f ( (cid:126)k ) + e ∂n F ∂(cid:15) (cid:126) Λ( (cid:126)k ) . (cid:126)E , (30)where the vector (cid:126) Λ has the dimension of a length, andits components correspond to transport lengths in thedifferent spatial directions. They are related to trans-port times through the definition Λ α ( (cid:126)k ) = v α ( (cid:126)k ) τ tr α ( (cid:126)k ).Eq. (30) can be rewritten as a shift of energies by thefield : f ( (cid:126)k ) = n F (cid:16) (cid:15) (cid:126)k + e(cid:126) Λ( (cid:126)k ) . (cid:126)E (cid:17) . In the case of anisotropic Fermi surface, we do not expect this shift todepend on the direction of application of the field (cid:126)E :in that case a unique transport time τ tr is necessaryto describe the stationary distribution . Here, for ananisotropic Fermi surface such as (5), we generically ex-pect the response of the distribution function to dependon the direction of the electric field (cid:126)E . For an elec-tric field applied in the x or y direction, this leads to thedefinition of different anisotropic transport times τ tr x , τ tr y .From Eqs. (28,29,30), one obtains (cid:126)v ( (cid:126)k ) = 1 (cid:126) ∂(cid:15)∂(cid:126)k (cid:12)(cid:12)(cid:12)(cid:12) (cid:15) = (cid:15) (cid:126)k = 2 π (cid:90) d (cid:126)k (cid:48) (2 π ) δ ( (cid:15) (cid:126)k − (cid:15) (cid:126)k (cid:48) ) |A ( (cid:126)k, (cid:126)k (cid:48) ) | (cid:16) (cid:126) Λ( (cid:126)k ) − (cid:126) Λ( (cid:126)k (cid:48) ) (cid:17) , (31)By using the parametrization (6) on the contour of constant energy (cid:15) , each component α of the velocity obeys theequation (to lighten notation, we omit the energy (cid:15) in the argument of the quantities in the next expressions) : v α ( η p , θ ) = Λ α ( η p , θ ) τ e ( θ ) − πγ (cid:126) (cid:88) η (cid:48) p = ± (cid:90) θ − θ dθ (cid:48) ρ ( θ (cid:48) ) [1 + cos( θ − θ (cid:48) )] Λ α ( η (cid:48) p , θ (cid:48) ) . (32)The transport times τ tr α ( (cid:15), θ ) are defined asΛ α ( (cid:15), η p , θ ) = v α ( (cid:15), η p , θ ) τ tr α ( (cid:15), θ ) . (33)We now assume the following ansatz, namely that the transport times and the elastic scattering time have the sameangular dependence: τ tr α ( (cid:15), θ ) = λ α ( (cid:15) ) τ e ( (cid:15), θ ) , (34)so that the parameters λ α ( (cid:15) ) are obtained from the self-consistent equation (at fixed energy (cid:15) ) v α ( η p , θ ) = λ α v α ( η p , θ ) − πγ (cid:126) λ α (cid:88) η (cid:48) p = ± (cid:90) θ − θ dθ (cid:48) ρ ( θ (cid:48) ) [1 + cos( θ − θ (cid:48) )] v α ( η (cid:48) p , θ (cid:48) ) τ e ( θ (cid:48) ) (35)where v α ( η p , θ ) is defined in Eq. (8). Then from Eq. (14) and (21), we finally get v α ( η p , θ ) = λ α v α ( η p , θ ) − λ α I ( δ ) (cid:88) η (cid:48) p = ± (cid:90) θ − θ dθ (cid:48) θ − θ (cid:48) )1 + r ( δ ) cos θ (cid:48) v α ( η (cid:48) p , θ (cid:48) ) √ cos θ (cid:48) − δ . (36)We now consider the two directions α = x, y separately. Along the x direction , since the velocity is an odd func-tion of k x , the sum over η p in Eq. (36) vanishes: we ob-tain λ x ( (cid:15) ) = 1, i.e. the transport time is equal to thescattering time τ tr x ( (cid:15), θ ) = τ e ( (cid:15), θ ). Along the y direction , where v y ( θ ) = c y sin θ indepen-dent of η p , Eq. (36) possesses a self-consistent solution,and we obtain λ y ( δ ) = 11 − I ( δ ) /I ( δ ) (37)where I ( δ ) = (cid:90) θ − θ dθ sin θ √ cos θ − δ (1 + r ( δ ) cos θ ) . (38)The function I ( δ ) is defined in (13). Note that the ex-pression (37) of the renormalization factor of the trans-port time λ y ( δ ) reflects the iterative structure of the ver-tex correction to the bare conductivity that will be ob-tained within a diagrammatic treatment in section IV B(see Eqs. 66,69). The dependence λ y ( δ ) is plotted inFig. 6. (cid:45) (cid:45) (cid:45) ∆ Λ y (cid:72) ∆ (cid:76) Figure 6. Dependence on δ = ∆ /(cid:15) of the renormalizationfactor of the transport time τ tr y with respect to the elasticscattering time : λ y ( δ ) = τ tr y ( (cid:126)k ) /τ e ( (cid:126)k ). Having obtained the transport times along the x and y directions, we now turn to the calculation of the con-ductivities. B. Conductivity
We can express the current density (cid:126)j occurring in re-sponse to the application of the electric field (cid:126)E as (cid:126)j = (cid:90) d (cid:126)k (2 π ) (cid:104) f ( (cid:126)k ) − n F ( (cid:15) (cid:126)k ) (cid:105) ( − e(cid:126)v ( (cid:126)k )) . (39)By using ∂n F /∂(cid:15) (cid:39) − δ ( (cid:15) − (cid:15) F ) and the ansatz (30) forthe distribution function f ( (cid:126)k ) we obtain (cid:126)j = e (cid:88) η p = ± (cid:90) θ − θ dθ ρ ( (cid:15) F , θ ) (cid:126)v ( (cid:15) F , η p , θ ) (cid:104) (cid:126) Λ( (cid:15) F , η p , θ ) . (cid:126)E (cid:105) . (40)The symmetries of this equation imply that off-diagonalterms of the conductivity tensor vanish ( σ α,β (cid:54) = α = 0)while the diagonal terms can be written as σ αα = 2 e (cid:90) θ − θ dθ ρ ( (cid:15) F , θ ) v α ( (cid:15) F , θ )Λ α ( (cid:15) F , θ ) . (41)where the factor 2 originates from the two possible signsof η p = ± . We end up with the Einstein relation σ αα = e ρ ( (cid:15) F ) D α , (42)with the diffusion coefficients D α = 2 λ α ( (cid:15) F ) (cid:90) θ − θ dθ ρ ( (cid:15) F , θ ) ρ ( (cid:15) F ) v α ( θ ) τ e ( θ ) (43a)= (cid:10) v α ( (cid:15) F , θ ) τ tr α ( (cid:15) F , θ ) (cid:11) θ (43b)= λ α ( (cid:15) F ) (cid:10) v α ( (cid:15) F , θ ) τ e ( (cid:15) F , θ ) (cid:11) θ . (43c)where we have defined the average along the constant en-ergy contour (cid:104)· · · (cid:105) θ = 2 (cid:82) θ − θ dθ · · · ρ ( (cid:15) F , θ ) /ρ ( (cid:15) F ). Thiscorresponds to the result announced in the introduction :the diffusion coefficients D α are obtained by an averageover the Fermi surface of v α τ tr α instead of v α τ e . Withour solution of the Boltzmann equation, this difference isaccounted for by a renormalization factor λ α ( (cid:15) F ) of thediffusion coefficients, which does not depend on the di-rection along the Fermi surface but on the direction α ofapplication of the electric field . We now specify explicitlythe conductivities along the two directions x and y . Along the x direction , there is no renormalization of thetransport time ( λ x = 1, τ tr x = τ e ) and the conductivity σ xx reads σ xx = 2 e (cid:90) θ − θ dθ ρ ( (cid:15) F , θ ) v x ( θ ) τ e ( θ )= e (cid:126) πγ (cid:15)m I (∆ /(cid:15) ) I (∆ /(cid:15) ) , (44)where we define I ( δ ) = (cid:90) θ − θ dθ cos θ √ cos θ − δ r ( δ ) cos θ , (45)and the function I ( δ ) is given in (13).For the conductivity along the y direction , the renor-malization of the transport time is given by (37) and weobtain σ yy = 2 e λ y (∆ /(cid:15) ) (cid:90) θ − θ dθ ρ ( (cid:15) F , θ ) v y ( θ ) τ e ( θ )= e (cid:126) πγ c y I (∆ /(cid:15) ) I (∆ /(cid:15) ) − I (∆ /(cid:15) ) , (46)where the functions I ( δ ) and I ( δ ) are respectively givenby Eqs. (13) and (38). Eqs. (44, 46) constitute the mainresults of this work. We discuss them in section V. In thenext section, we use a diagrammatic approach which pro-poses a complementary description of the anisotropy oftransport and allows to confirm the ansatz made to solvethe Boltzman equation and recover exactly the results ofEqs. (44, 46) . IV. DIAGRAMMATIC APPROACH
An alternative approach to describe the diffusive trans-port of electron consists in a perturbative expansion indisorder of the conductivity tensor using a diagrammatictechnique . Beyond confirming the ansatz made to solvethe Boltzmann equation described above, this method al-lows for an instructive alternative treatment of the dif-ferent transport anisotropies. In the diagrammatic ap-proach, the transport coefficients of the model are ob-tained from the Kubo formula. A perturbative expansionis then used to express the transport coefficients usingthe average single particle Green’s function. In this for-malism, the anisotropy of scattering and transport timesare cast into a unusual matrix form for the self-energyoperator Σ. Beyond the present model, such a techniqueallows to describe anisotropy of diffusion of Dirac fermionmodels due e.g. to the warping of the Fermi surface intopological insulators or anisotropic impurity scatter-ing, the study of which goes beyond the scope of thepresent paper. Nevertheless our work provides a phys-ical understanding of the technicalities naturally occur-ring in these other problems. In the next subsections,we first discuss the self-energy and the single particleGreen’s function. We then turn to the calculation of theconductivity. A. Green’s functions and self-energy
The retarded and advanced Green’s functions are de-fined by : G R/A ( (cid:126)k, (cid:126)k (cid:48) , (cid:15) F ) = (cid:104)(cid:16) ( (cid:15) F ∓ i I − H ( (cid:126)k ) (cid:17) δ ( (cid:126)k − (cid:126)k (cid:48) ) − V ( (cid:126)k, (cid:126)k (cid:48) ) I (cid:105) − (47)In the case of the model without disorder defined byEq. (4), the Green’s function is expressed as a 2 × G ( (cid:126)k, (cid:15) ) = (cid:16) (cid:15) I − H ( (cid:126)k ) (cid:17) − = (cid:15) I + (cid:16) (cid:126) k x m + ∆ (cid:17) σ x + c y (cid:126) k y σ y (cid:15) − (cid:16) (cid:126) k x m + ∆ (cid:17) − c y (cid:126) k y , (48)where I is the identity matrix. Disorder is perturbativelyincorporated in the averaged Green’s ¯ G function througha self-energy matrix Σ( (cid:126)k, (cid:15) ) such that G R/A ( (cid:126)k, (cid:15) ) = (cid:104) ( (cid:15) ∓ i I − H ( (cid:126)k ) ∓ i Im Σ( (cid:126)k, (cid:15) ) (cid:105) − . (49)The real part of the self-energy has been neglected. Theelastic scattering rates will be defined below from theimaginary part of the self-energy. To lowest order in thedisorder strength γ , this self-energy, solution of a Dysonequation, readsΣ( (cid:126)k, (cid:15) ) = (cid:90) d(cid:126)k (cid:48) (2 π ) V ( (cid:126)k (cid:48) ) V ( − (cid:126)k (cid:48) ) G ( (cid:126)k − (cid:126)k (cid:48) , (cid:15) ) . (50)Its imaginary part is then obtained as − Im Σ( (cid:15) ) = πγ (cid:90) θ − θ dθ ρ ( (cid:15), θ ) [ I − cos θ σ x ] (51)= (cid:126) τ e ( (cid:15) ) [ I + r ( δ ) σ x ] . (52)The densities of states ρ ( (cid:15), θ ) and ρ ( (cid:15) ), the bare scat-tering time τ e ( (cid:15) ) and the anisotropy factor r ( δ ) havebeen defined in Eqs. (10,12,20,22). It is worth notingthat this self-energy acquires an unusual matrix struc-ture in pseudo-spin space: this manifests within the dia-grammatic approach the anisotropy of the scattering time τ e ( (cid:15), θ ), which was described in Eq. (21) previously. In-deed, in the Green function formalism, the direction ofpropagation of eigenstates of the Hamiltonian (4) is en-coded into their spinor structure (the relative phase be-tween their components, see Eq. (7)). Hence the scatter-ing time in the corresponding direction will be obtainedas the matrix element of the above self-energy in the as-sociated spinor eigenstate. (cid:104) G R (cid:105) ab ( (cid:126)k ) a b (cid:104) G A (cid:105) cd ( (cid:126)k ) d c V ( (cid:126)k ) = γ Figure 7. Conventions for the diagrammatic representationof perturbation theory of transport. J α ( (cid:126)k ) G R ( (cid:126)k ) G A ( (cid:126)k ) j α ( (cid:126)k )Figure 8. Diagrammatic representation of the classical con-ductivity with the conventions of Fig. 7. The renormalizedcurrent operator is defined in Fig. 9. B. Conductivity
1. Kubo formula
The longitudinal conductivity can be deduced from theKubo formula ( α = x, y ) : σ αα = (cid:126) πL Tr (cid:104) j α ( (cid:126)k ) G R ( (cid:126)k, (cid:126)k (cid:48) , (cid:15) F ) j α ( (cid:126)k (cid:48) ) G A ( (cid:126)k (cid:48) , (cid:126)k, (cid:15) F ) (cid:105) , (53)where Tr corresponds to a trace over the pseudo-spinand momentum quantum numbers : Tr = tr (cid:80) (cid:126)k (cid:39) L tr (cid:82) d(cid:126)k/ (2 π ) and tr is a trace over the pseudo-spinindices only. For clarity, throughout this section on trans-port coefficients, we will omit the dependence on theFermi energy (cid:15) F of various quantities. The current den-sity operators are also operators acting on both spin andmomentum spaces. They are deduced from the Hamilto-nian (4) as: j x ( (cid:126)k ) = − em (cid:126) k x σ x ; j y ( (cid:126)k ) = − ec y σ y . (54)Note that j x is linear in momentum while j y dependsonly on spin quantum numbers. Perturbation in the dis-order amplitude of the conductivity (53) is obtained byexpanding the Green’s function in the disorder potential V before averaging over the gaussian distribution. In theclassical diffusive limit, the dominant terms which deter-mine the averaged classical conductivity are representeddiagrammatically on Fig. 8 and lead to σ αα = (cid:126) πL Tr (cid:104) J α G R j α G A (cid:105) (55)where J α is the renormalized current density operator.The discrepancy between J α and the bare current opera-tor j α accounts for the appearance of transport time τ tr α in the Boltzmann approach due to the anisotropy ofscattering. This renormalized current operator is easierto define diagrammatically, as shown on Fig. 9.
2. Conductivity along x In this direction, the current operator is linear in k x ,while the averaged Green’s functions G R ( (cid:126)k ) , G A ( (cid:126)k ) areeven functions of k x . Hence all the terms in the expres-sion of the renormalized current J x with at least a Green’sfunction vanish by k x → − k x symmetry, and J x ( (cid:126)k ) = j x ( (cid:126)k ) = − em (cid:126) k x σ x . (56)There is no renormalization of the current operator, inagreement with the result τ tr x = τ e from the Boltzmannequation approach. In the x direction, the expression(55) reduces to σ xx = (cid:18) (cid:126) em (cid:19) (cid:126) πL Tr (cid:104) k x σ x G R ( (cid:126)k ) σ x G A ( (cid:126)k ) (cid:105) . (57)Using L − (cid:80) (cid:126)k (cid:39) (cid:82) ρ ( (cid:15), θ ) d(cid:15)dθ and the parametrizationdefined in Eq. (6) of the contours of constant energy (cid:15) we perform the integration over energy to obtain σ xx = e τ e (cid:15) F m (cid:90) + θ − θ dθ ρ ( (cid:15), θ ) (cos θ − δ )1 + r ( δ ) cos θ × tr (cid:20) σ x [ I + cos θσ x + sin θσ y ] σ x [ I + cos θσ x + sin θσ y ] (cid:21) . Performing the spin trace first, we obtain σ xx = 4 e τ e (cid:15) F m (cid:90) + θ − θ dθ ρ ( (cid:15), θ ) (cos θ − δ )1 + r ( δ ) cos θ cos θ . (58)By using eq. (14) for the directional density of states werecover exactly the integral expression for the result (44)of Boltzmann approach: σ xx = e (cid:126) πγ (cid:15)m I (∆ /(cid:15) ) I (∆ /(cid:15) ) . (59)
3. Renormalized current operator along y In the y direction, the current operator is renormal-ized : the bare current operator j y is independent of themomentum (cid:126)k and the symmetry argument used for the x direction does not hold anymore. This renormalizedcurrent operator satisfies a Bethe-Salpeter equation rep-resented in Fig. 9: J y = j y + J y Π γ (60) J α ( (cid:126)k ) = + + + + · · · J α ( (cid:126)k ) ab = j α ( (cid:126)k ) ab + J α ( (cid:126)k (cid:48) ) G R ( (cid:126)k (cid:48) ) G A ( (cid:126)k (cid:48) ) γabcd Figure 9. Schematic representation of renormalized current operator [ J α ] ab ( (cid:126)k ) as the infinite sum of vertex corrections to thebare current operator (top), and corresponding recursive equation satisfied by J α (bottom). where tensor product in spin space are assumed andΠ( (cid:15), ∆) = (cid:90) d(cid:126)k (2 π ) G R ( (cid:126)k, (cid:15) ) ⊗ G A ( (cid:126)k, (cid:15) ) T . (61)Due to the spinorial structure of the wave functions,this propagator is here an operator acting as the ten-sor product of two spin spaces. The notation · · · T corresponds to a transposition of spin matrices. Usingthe parametrization defined in Eq. (6) of the contours ofconstant energy (cid:15) we perform the integration over energyto obtain for Π( (cid:15), ∆) ≡ Π( δ = ∆ /(cid:15) ):Π( δ ) = πτ e (cid:126) (cid:90) + θ − θ dθ ρ ( (cid:15), θ )1 + r ( δ ) cos θ × [ I + cos θσ x + sin θσ y ] ⊗ [ I + cos θσ x − sin θσ y ] . (62)The expression (14) for the directional density of statesallows to rewrite it asΠ( δ ) = 12 γI ( δ ) (cid:20) I ( δ ) I ⊗ I + ( I ( δ ) − I ( δ )) σ x ⊗ σ x − I ( δ ) σ y ⊗ σ y + J ( δ )( I ⊗ σ x + σ x ⊗ I ) (cid:21) , (63)where we introduced the functions: I ( δ ) = (cid:90) θ − θ dθ √ cos θ − δ (1 + r ( δ ) cos θ ) , (64) J ( δ ) = (cid:90) θ − θ dθ cos θ √ cos θ − δ (1 + r ( δ ) cos θ ) , (65)whereas I and I are defined in Eqs. (13,38).The inversion of the Bethe-Salpeter equation (60) isdone in the appendix A and we find J y = j y ( I ⊗ I − γ Π( δ )) − = (cid:18) − I ( δ ) I ( δ ) (cid:19) − j y . (66)
4. Conductivity along y Following the formula (55), the average conductivityalong y is expressed as σ yy = (cid:126) π tr [ J y . Π( δ ) .j y ] . (67)From the eq. (60), we express J y . Π( δ ) = γ − ( J y − j y ) toobtain from (67): σ yy = (cid:126) πγ tr [( J y − j y ) j y ] . (68)The expression for the renormalized current operator (66)leads to the final result σ yy = e (cid:126) πγ c y I ( δ ) I ( δ ) − I ( δ ) , (69)which is precisely the result (46) obtained within theBoltzmann equation approach.This concludes the derivation of the conductivity ten-sor within the diagrammatic approach. In doing so, wehave identified the encoding of the anisotropic scatteringrates through the matrix self-energy (52), while the cor-responding transport times are hidden into the renormal-ization of vertex operators (56,66). Comparison with theBoltzmann approach allows to unveil the physical mean-ing of these technical structures, which we believe to beapplicable to other situations of anisotropic transport ofDirac-like states. V. RESULTS AND DISCUSSION
We now turn to a discussion of our results for varioussituations corresponding to energy spectra represented inFig. 2.0 A. ∆ = 0 : Semi-Dirac spectrum Focusing first on the merging point (∆ = 0), wefind that the conductivities are expressed, from Eqs.(44,46,59,69) as : σ xx ( (cid:15) ) = e (cid:126) πγ (cid:15)m I (0) I (0) (cid:39) . e (cid:126) πγ (cid:15)m (70a) σ yy ( (cid:15) ) = e (cid:126) πγ c y I (0) I (0) − I (0) (cid:39) . e (cid:126) πγ c y . (70b)This case, which corresponds to an hybrid dispersionrelation, linear in one direction and quadratic in theother direction, has been previously studied in Ref. 34.However these authors have neglected both the spino-rial structure of the wave function and the angular de-pendence of the elastic scattering time caused by theanisotropic dispersion, whose importance is emphasizedin the present paper. Using the numerical values of theintegrals given in appendix B, we find λ y (0) (cid:39) . : σ xx = I (0) I (0) σ Bxx (cid:39) . σ Bxx (71a) σ yy = λ y (0) I (0) I (0) σ Byy (cid:39) . σ Byy , (71b)where σ Bxx and σ Byy are the values obtained in Ref. [34].The energy dependence of the conductivities (70) arisesfrom the energy dependence of the average squared ve-locities. It is therefore independent of the energy alongthe y direction, while it increases linearly with energyalong the x direction. B. ∆ > : gapped spectrum When ∆ >
0, the energy spectrum exhibits a gap andwe study here the conductivity above this gap at energies (cid:15) > ∆. Along the x direction, the renormalization factor λ x = 1 so that τ trx ( θ ) = τ e ( θ ). The energy dependencearises mainly from the energy dependence of the averagesquared velocity. Therefore we expect a roughly lineardependence in energy : σ xx ( (cid:15) ) = e (cid:126) πγ (cid:15)m I (∆ /(cid:15) ) I (∆ /(cid:15) ) = e (cid:126) πγ c x (cid:15) ∆ I (∆ /(cid:15) ) I (∆ /(cid:15) ) (72) ≈ . e (cid:126) πγ c x (cid:15) − ∆∆ (73)where c x = (cid:112) /m is the velocity along x of the mas-sive Dirac equation describing the spectrum for small mo-menta.The dependence in energy of the conductivity σ yy ismainly due to the energy dependence of the renormal-ization factor λ y between transport time and relaxation time : σ yy ( (cid:15) ) = e (cid:126) πγ c y I (∆ /(cid:15) ) I (∆ /(cid:15) ) − I (∆ /(cid:15) ) (74) ≈ e (cid:126) πγ c y I (0) I (0) λ y (∆ /(cid:15) ) . (75)The energy dependence of the conductivities σ xx and σ yy is plotted in Fig. 10. Figure 10. (Conductivities σ αα in units of e (cid:126) c α /γ for α = x, y ) as functions of (cid:15)/ ∆ in the gapped phase (∆ > σ xx arises from the energy dependenceof the velocity along x , while for σ yy it comes mainly fromthe energy dependence of the renormalization factor λ y . C. ∆ < : Dirac spectrum Figure 11. Conductivities σ αα in units of e (cid:126) c α /γ for α = x, y ) as functions of (cid:15)/ | ∆ | in the Dirac phase (∆ < (cid:15) = | ∆ | ). Thevicinity of the saddle point should be treated with a self-consistent Born approximation (see text). In this phase, we have two regimes separated by thesaddle point energy | ∆ | . At high energy above the sad-1dle point, (cid:15) (cid:29) | ∆ | , the energy dependence of the con-ductivities are still given by Eqs. (73, 75). In the lowenergy limit (cid:15) (cid:28) | ∆ | , expanding these expressions usingEq. (B17), we recover the conductivities associated to aconic dispersion of characteristic velocities c x and c y : σ xx ( (cid:15) →
0) = e (cid:126) πγ c x (76) σ yy ( (cid:15) →
0) = 2 e (cid:126) πγ c y . (77)Note however the factor 2 between the two expressions.This is due to the fact that τ tr y = 2 τ e like in graphenewhile τ tr x = τ e . It is instructive to compare with thislimit with the case of graphene, where it is known that τ tr = 2 τ e is all directions , and where Einstein relation σ αα = e ( c α τ tr / ρ ( (cid:15) F ) together with the Fermi goldenrule τ tr = 2 τ e = 2 (cid:126) / ( πρ ( (cid:15) F ) γ ) leads to σ (graphene) = e (cid:126) πγ v F . (78)Using the fact that c x = c y = v F /
2, we find the sameresult for σ yy but the conductivity is twice smaller alongthe x direction. The difference by a factor 2 between σ xx and σ yy results from intervalley scattering taking placealong the x direction (see Ref. 33 or a related discussionof diffusion within graphene with different intervalley andintravalley disorder rates).It is important to note that our calculations predictvanishing conductivities at the saddle point (cid:15) = | ∆ | .That is the result of the logarithmic divergence of thedensity of states producing a vanishing elastic scatteringtime in Eq. (21). However, in such a limit, k F l e →
0, soour approximations are no longer valid. To describe cor-rectly the behavior of the scattering time in the vicinity ofthe saddle point, it is necessary to go beyond second orderperturbation theory, using for instance the self-consistentBorn approximation , in order to obtain a finite densityof states and a non-zero elastic scattering time. Quali-tatively, we expect that the zero of the conductivity willbe replaced by a minimum for (cid:15) (cid:39) | ∆ | . D. Evolution of conductivites across the transition
We are now in position to discuss the evolution ofthe conductivity at fixed energy (cid:15) F , as a function ofthe parameter ∆ as we cross the merging transition.Such evolution, derived for eq. (44,46) is represented onFig. 12, where we have plotted σ αα in units of e (cid:126) ˜ c α /γ for α = x, y and ˜ c x = (cid:112) (cid:15) F /m and ˜ c y = c y . Be-low the saddle point for ∆ < − (cid:15) F , σ yy is nearly con-stant, while σ xx decreases almost linearly with ∆. Atthe saddle point ∆ = − (cid:15) F where the topology of theFermi surface changes, a dip in both σ xx and σ yy isvisible, down to minimal values not quantitatively cap-tured by the present approach. Past the saddle point, while σ xx remains linearly decreasing with ∆, albeit moreslowly, σ yy is first increasing, presenting a maximum for∆ /(cid:15) F (cid:39) − .
39 and then decreases to zero. No signatureof the underlying transition at ∆ = 0 is manifest in thetransport at high Fermi energy (cid:15) F . Figure 12. Conductivities (in units of e (cid:126) ˜ c α /γ for α = x, y ) asa function of the parameter ∆ for a fixed chemical potential (cid:15) . The conductivity σ xx decreases monotonically with ∆, inalmost linear fashion with a change of slope at the saddle pointS. The conductivity σ yy is almost constant below the saddlepoint. Above the saddle point, the behavior of σ yy becomesnon-monotonous with ∆. At the metal-insulator transition,both σ xx and σ yy vanish linearly. Symbols (a), (b), (c), Mand S refer to the regions presented in Fig. 2. VI. CONCLUSION
We have studied the behavior of the conductivity inboth the Dirac phase, the critical semi-Dirac and abovethe gapped phase. Using the complementary Boltzmannand diagrammatic techniques we have identified the dif-ferent nature of anisotropy of the elastic scattering timesand transport times. Indeed the transport is inherentlyanisotropic due both to the spinorial structure of theeigenstates and the anisotropy of the dispersion rela-tion. The approaches developed in this paper can begeneralized to study the diffusive transport in other semi-metallic phases, including the various three dimensionalspecies recently identified.
ACKNOWLEDGMENTS
We thank the hospitality of the Institut Henri Poincar´ewhere part of this work was completed. This workwas supported by the French Agence Nationale de la2Recherche (ANR) under grants SemiTopo (ANR-12-BS04-0007), IsoTop (ANR-10-BLAN-0419).
Appendix A: Current vertex renormalization1. Inverse of tensor product
We use the notation M abcd = A ab ⊗ B cd for the coef-ficients of a tensor product M = A ⊗ B . The inverse(for the outer product) N = A − ⊗ B − of M satisfiesthe relation M abcd N bedf = δ ae δ cf . In the obtention of thediffuson propagator, we need to invert a tensor productof the form M = a I ⊗ I + b σ x ⊗ σ x + c σ y ⊗ σ y + d ( I ⊗ σ x + σ x ⊗ I ) . (A1)Its inverse M − can be parametrized as M − = ∆ − (cid:20) A I ⊗ I + B σ x ⊗ σ x + C σ y ⊗ σ y + D ( I ⊗ σ x + σ x ⊗ I ) + E σ z ⊗ σ z (cid:21) , (A2)with A = a + 2 bd − a ( b + c + 2 d ) (A3a) B = − b ( a − b + c ) + 2( a − b ) d (A3b) C = c ( − a − b + c + 2 d ) (A3c) D = − d (cid:2) ( a − b ) − c (cid:3) (A3d) E = 2 c ( − ab + d ) (A3e)∆ = (cid:2) ( a − b ) − c (cid:3) (cid:2) ( a + b ) − c − d (cid:3) . (A3f)Let us now focus on the following contraction (cid:2) M − (cid:3) abcd ( σ y ) bd = A − B − C − E ∆ ( σ y ) ac (A4)= 1 a − b − c ( σ y ) ac . (A5)irrespective of d and hence of ∆. In particular eq. (A5)is valid even if ∆ vanishes.
2. Current renormalization
Let us now use the above parametrization (A1) for thetensor M = I ⊗ I − γ Π( δ ). From the parametrization ofeq. (63), we obtain the following identification of coeffi-cients a = 1 − I ( δ )2 I ( δ ) , b = I ( δ ) − I ( δ )2 I ( δ ) , (A6a) c = I ( δ )2 I ( δ ) , d = − J ( δ )2 I ( δ ) . (A6b) Then the equation (A5) provides the expression for therenormalized current operator: J y = j y ( I ⊗ I − γ Π( δ )) − = (cid:18) − I ( δ ) I ( δ ) (cid:19) − j y . (A7)A word of caution is necessary at this stage :(1 − γ Π( (cid:126)q )) − is the structure factor which encodes thepropagation of the diffuson modes . In the symplecticclass which we consider, there is one such mode whichis diffusive : 1 − γ Π( (cid:126)q ) possesses a vanishing eigenvalue ∝ ( Dq ). Hence in the limit q → − γ Π is no longer invertible. In principle, we shouldhave kept a finite momentum q during the calculation,and taken the limit q → M : this isa classical contribution, which depends on these massivediffuson modes, while the diffusive long distance modesenter the quantum correction not discussed in this paper. Appendix B: Special Functions
In this appendix, we discuss a few useful results ofthe various integrals entering the expressions of transportcoefficients in the paper and arising as integrals alongthe constant energy contours of the model. Let us firstconsider the integrals : I ( δ ) = (cid:90) θ − θ dθ √ cos θ − δ , J ( δ ) = (cid:90) θ − θ dθ cos θ √ cos θ − δ , with cos θ = δ if | δ | < θ = π otherwise. Defining X = (cid:112) / (1 − δ ), we find : I ( δ ) = 4 √ − δ K ( X ) δ < − I ( δ ) = 2 √ K (1 /X ) − < δ < J ( δ ) =4 √ − δ [(1 − δ ) E ( X ) + δK ( X )] for δ < − , (B3)and J ( δ ) =2 √ E (1 /X ) − K (1 /X )] for − < δ < . (B4)The anisotropy factor r ( δ ) = J ( δ ) /I ( δ ) reads : r ( δ ) = (1 − δ ) E ( X ) /K ( X ) + δ for δ < − r ( δ ) = 2 E (1 /X ) /K (1 /X ) − − < δ < . (B6)3The functions I ( δ and J ( δ ) are plotted in Fig. 13, andthe dependence on δ of r ( δ ) is shown on Fig. 5. (cid:45) (cid:45) (cid:45) (cid:45) ∆ I (cid:72) ∆ (cid:76) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) ∆ J (cid:72) ∆ (cid:76) Figure 13. Functions I ( δ ) and J ( δ ) (cid:45) (cid:45) (cid:45) (cid:45) ∆ (cid:73) (cid:72) ∆ (cid:76) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) ∆ J (cid:72) ∆ (cid:76) (cid:45) (cid:45) (cid:45) (cid:45) ∆ (cid:73) (cid:72) ∆ (cid:76) (cid:45) (cid:45) (cid:45) (cid:45) ∆ (cid:73) (cid:72) ∆ (cid:76) Figure 14. Functions I ( δ ), J ( δ ), I ( δ ) and I ( δ ) Let us now consider the four integrals I ( r, δ ) = (cid:90) θ − θ dθ (1 + r cos θ ) √ cos θ − δ (B7) J ( r, δ ) = (cid:90) θ − θ dθ cos θ (1 + r cos θ ) √ cos θ − δ (B8) I ( r, δ ) = (cid:90) θ − θ dθ sin θ (1 + r cos θ ) √ cos θ − δ (B9) I ( r, δ ) = (cid:90) θ − θ dθ cos θ √ cos θ − δ (1 + r cos θ ) (B10)We first focus on first the integral I ( r, δ ). For δ < − I ( r, δ ) = 4(1 + r ) √ − δ (cid:90) π dθ (cid:16) − r r sin θ (cid:17) (cid:113) − − δ sin θ (B11)= 4(1 + r ) √ − δ Π (cid:32) π , r r , (cid:114) − δ (cid:33) , (B12) where Π is an elliptic integral of the third kind. For | δ | <
1, we use the change of variable sin θ = sin θ sin ϕ to rewrite this function as: I ( r, δ ) = 2 √
21 + r (cid:90) π dϕ (cid:16) − r r sin θ sin ϕ (cid:17) (cid:113) − sin θ sin ϕ = 2 √
21 + r Π (cid:18) π , r r sin θ , sin θ (cid:19) . (B13)Moreover the integrals J ( r, δ ), I ( r, δ ) and I ( r, δ ) canbe expressed in terms of I ( δ ), J ( δ ) and I ( r, δ ): J ( r, δ ) = 1 r [ I ( δ ) − I ( r, δ )] , (B14) I ( r, δ ) = 1 r I ( δ ) − r J ( δ ) + (cid:18) − r (cid:19) I ( r, δ ) , (B15) I ( r, δ ) = 1 r (cid:18)
13 + δr + 1 r (cid:19) I ( δ ) − r (cid:18) δ r (cid:19) J ( δ ) − r (cid:18) δ + 1 r (cid:19) I ( r, δ ) . (B16)Finally, the integrals used in the text are I ( δ ) = I [ r ( δ ) , δ ] J ( δ ) = J [ r ( δ ) , δ ] I ( δ ) = I [ r ( δ ) , δ ] I ( δ ) = I [ r ( δ ) , δ ] . The following special values are of particular interest forthe expressions in the text : I (0) = 2 √ K (cid:18) √ (cid:19) (cid:39) . J (0) = π √ /K (cid:18) √ (cid:19) (cid:39) . r (0) = J (0) /I (0) (cid:39) . I (0) = 43 √ K (1 / √ (cid:39) . I (0) (cid:39) . I (0) (cid:39) . I (0) (cid:39) . . as well as the limits when δ → −∞ : I ( δ ) I ( δ ) − I ( δ ) → I ( δ ) I ( δ ) → − δ . (B17)4 A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S.Novoselov, and A. K. Geim, “The electronic propertiesof graphene,” Rev. Mod. Phys. , 109 (2009). Akito Kobayashi, Shinya Katayama, Yoshikazu Suzumura,and Hidetoshi Fukuyama, “Massless fermions in organicconductor,” Journal of the Physical Society of Japan ,034711 (2007). M. O. Goerbig, J.-N. Fuchs, G. Montambaux, andF. Pi´echon, “Tilted anisotropic dirac cones in quinoid-type graphene and α -(BEDT-TTF) I ,” Phys. Rev. B ,045415 (2008). Yoshikazu Suzumura and Akito Kobayashi, “Berry curva-ture of the dirac particle in α -(BEDT-TTF) I ,” Journalof the Physical Society of Japan , 104701 (2011). Takeshi Choji, Akito Kobayashi, and Yoshikazu Suzu-mura, “Zero-gap state in organic conductor α -(BEDT-TTF) NH Hg(SCN) ,” Journal of the Physical Society ofJapan , 074712 (2011). M. Polini, F. Guinea, M. Lewenstein, H. C. Manoha-ran, and V. Pellegrini, “Artificial honeycomb lattices forelectrons, atoms and photons,” Nature Nanotech. , 625(2013). S. Bittner, B. Dietz, M. Miski-Oglu, P. Oria Iriarte, andF. Richter, A. Sch¨afer, “Observation of a dirac point inmicrowave experiments with a photonic crystal modelinggraphene,” Phys. Rev. B , 014301 (2010). M. Bellec, U. Kuhl, G. Montambaux, and F. Mortessagne,“Topological transition of dirac points in a microwave ex-periment,” Phys. Rev. Lett. , 033902 (2010). M.C. Rechtsman, J.M. Zeuner, A. T¨unnermann, S. Nolte,M. Segev, and A. Szameit, “Strain-induced pseudomag-netic feld and photonic landau levels in dielectric struc-tures,” Nature Photon. , 153 (2013). T Jacqmin, I. Carusotto, I. Sagnes, M. Abbarchi, D. D. Sol-nyshkov, G. Malpuech, E. Galopin, A. Lemaˆıtre, J. Bloch,and A. Amo, “Manipulation of edge states in microwaveartificial graphene,” Phys. Rev. Lett. , 116402 (2014). Leticia Tarruell, Daniel Greif, Thomas Uehlinger, GregorJotzu, and Tilman Esslinger, “Creating, moving and merg-ing dirac points with a fermi gas in a tunable honeycomblattice,” Nature , 302 (2012). L.-K. Lim, J.-N. Fuchs, and G. Montambaux, “Bloch-zener oscillations across a merging transition of diracpoints,” Phys. Rev. Lett. (2012). Oskar Vafek and Ashvin Vishwanath, “Dirac fermions insolids - from high tc cuprates and graphene to topologicalinsulators and weyl semimetals,” Annual Review of Con-densed Matter Physics , 83–112 (2014). Pavan Hosur and Xiaoliang Qi, “Recent developments intransport phenomena in weyl semimetals,” Comptes Ren-dus Physique , 857 (2013). Y. Hasegawa, R. Konno, H. Nakano, and M. Kohmoto,“Zero modes of tight-binding electrons on the honeycomblattice,” Phys. Rev. B , 033413 (2006). B. Wunsch, F. Guinea, and F. Sols, “Dirac-point engineer-ing and topological phase transitions in honeycomb opticallattices,” New. J. Phys. , 103027 (2008). V. M. Pereira, A. H. Castro Neto, and N. M. R. Peres,“Tight-binding approach to uniaxial strain in graphene,”Phys. Rev. B , 045401 (2009). G. Montambaux, F. Pi´echon, J.-N. Fuchs, and M.O. Goer-big, “Merging of dirac points in a two-dimensional crystal,”Phys. Rev. B , 153412 (2009). G. Montambaux, F. Pi´echon, J.-N. Fuchs, and M.O. Go-erbig, “A universal hamiltonian for the motion and themerging of dirac cones in a two-dimensional crystal,” Eur.Phys. J. B , 509 (2009). S. Banerjee, R. R. P. Singh, V. Pardo, and W. E. Pickett,“Tight-binding modeling and low-energy behavior of thesemi-dirac point,” Phys. Rev. Lett. , 016402 (2009). M. Katsnelson, “Zitterbewegung, chirality, and minimalconductivity in graphene,” Eur. Phys. J. B , 157–160(2006). J. Tworzydlo, B. Trauzettel, M. Titov, A. Rycerz,and C.W.J Beenakker, “Sub-poissonian shot noise ingraphene,” Phys. Rev. Lett. , 246802 (2006). A.A. Abrikosov,
Fundamentals of the Theory of Metals (North Holland, Amsterdam, Netherlands, 1988). J. M. Ziman,
Principles of the Theory of Solids (Cam-bridge University Press, Cambridge, 1972). Di Xiao, Ming-Che Chang, and Qian Niu, “Berry phaseeffects on electronic properties,” Rev. Mod. Phys. , 1959(2010). D. T. Son and B. Z. Spivak, “Chiral anomaly and classicalnegative magnetoresistance of weyl metals,” Phys. Rev. B , 104412 (2013). E.H. Sondheimer, “The boltzmann equation for anisotropicmetals,” Proc. R. Soc. London, Ser. A , 100 (1962). R.S. Sorbello, “On the anisotropic relaxation time,” J.Phys. F , 503 (1974). R.S. Sorbello, “Effects of anisotropic scattering on elec-tronic transport properties,” Phys. Cond. Matter , 303(1975). This is not true in graphene where λ x = 2. Note thatthis peculiar result ( λ x = 1) is due the fact that the ma-trix elements of the disorder potential are supposed hereto have no momentum dependence. Assuming an oppositelimit where the disorder would not couple valleys, then inthe Dirac limit 0 < (cid:15) (cid:28) − ∆, one would recover λ x = 1and τ tr x = 2 τ e . E. Akkermans and G. Montambaux,
Mesoscopic Physicsof Electrons and Photons (Cambridge Univ. Press, 2007). P. Adroguer, D. Carpentier, J. Cayssol, and E. Orignac,“Diffusion at the surface of topological insulators,” New J.Phys. , 103027 (2012). E. McCann, K. Kechedzhi, V. I. Fal’ko, H. Suzuura,T. Ando, and B. L. Altshuler, “Weak-localization mag-netoresistance and valley symmetry in graphene,” Phys.Rev. Lett. , 146805 (2006). S. Banerjee and W. E. Pickett, “Phenomenology of asemi-dirac semi-weyl semimetal,” Phys. Rev. B , 075124(2012). Note the correspondance between our notations I k ( x ) , ( k =1 , ,
3) and those I k of Ref. 34 : I k (0) = 4 I k . More precisely, for (cid:15) (cid:28) ∆ we have : σ xx ( (cid:15) ) ≈ . e (cid:126) / ( πγ ) c x ( (cid:15)/ ∆ − . D. Carpentier, A. A. Fedorenko, and E. Orignac, “Effectof disorder on 2D topological merging transition from adirac semi-metal to a normal insulator,” Europhys. Lett. , 67010 (2013). We use the definition of the elliptic integrals from Grad-shteyn and Ryzhik . They differ from those used inMathematica : K Grad. ( x ) = K Math. ( x ), E Grad. ( x ) = E Math. ( x ) and Π Grad. ( φ, n, x ) = Π Math. ( n, φ, x ). A. Gradshteyn and R. Ryzhik,