Gap generation and phase diagram in strained graphene in a magnetic field
aa r X i v : . [ c ond - m a t . s t r- e l ] J a n Gap generation and phase diagram in strained graphene in a magnetic field
D.O. Rybalka, E.V. Gorbar,
1, 2 and V.P. Gusynin Department of Physics, Taras Shevchenko National Kiev University, 03022, Kiev, Ukraine Bogolyubov Institute for Theoretical Physics, 03680, Kiev, Ukraine (Dated: November 5, 2018)The gap equation for Dirac quasiparticles in monolayer graphene in constant magnetic and pseu-domagnetic fields, where the latter is due to strain, is studied in a low-energy effective model withcontact interactions. Analyzing solutions of the gap equation, the phase diagram of the systemin the plane of pseudomagnetic and parallel magnetic fields is obtained in the approximation ofthe lowest Landau level. The three quantum Hall states, ferromagnetic, antiferromagnetic, andcanted antiferromagnetic, are realized in different regions of the phase diagram. It is found that thestructure of the phase diagram is sensitive to signs and values of certain four-fermion interactioncouplings which break the approximate spin-value SU (4) symmetry of the model. I. INTRODUCTION
Among many remarkable properties of graphene its response to strain can be singled out as one of efficient meansto change and control the characteristics of the electronic states in graphene. Since graphene is only one atom thick, itis easily subjected to mechanical deformations. Various proposals to engineer strain in graphene were discussed in theliterature . It is known that strain induces effective gauge fields with the corresponding effective ”magnetic” fieldsof opposite sign in valleys K and K ′ that means that elastic deformations, unlike real magnetic field, preserve timereversal symmetry . Since time reversal symmetry is unbroken, strain induced fields are known as pseudomagneticfields in the literature (for a review of gauge fields in strained graphene, see Refs.[8,9]). It was proposed in Refs.[10,11]that a designed strain may induce uniform pseudomagnetic field, which can easily reach values exceeding 10 T.The observation of anomalous integer quantum Hall (QH) effect with the filling factors ν = ± | n | + 1 /
2) ( n is theLandau level index) in graphene in a magnetic field , in accordance with theoretical studies in Refs.[13,14], was amilestone in graphene research as it became a direct experimental proof of the existence of gapless Dirac quasiparticlesin graphene. The four-fold degeneracy of the Landau levels in graphene is due to the SU (4) symmetry connectedwith valley and spin. Later the plateaus ν = 0 , ± , ± in a strongmagnetic field B ≥ T . These plateaus are connected with the magnetic field induced splitting of the n = 0 and n = 1 Landau levels and the degeneracy of the lowest Landau level (LLL) is thus completely resolved.The Landau levels related to pseudomagnetic fields were observed in spectroscopic measurements. It waspointed out that pseudomagnetic fields due to strain can interfere in many ways with real magnetic fields. Forexample, the interplay of pseudomagnetic and magnetic fields in the quantum Hall regime causes backscattering inthe chiral edge channels that can destroy the quantized conductance plateaus and gives rise to unconventional QHeffect in strained graphene with oscillating Hall conductivity.The gap generation in graphene in the presence of a pseudomagnetic field was studied in Ref.[20]. Interestingly,it was found that unlike magnetic field which catalyses the generation of the time reversal invariant Dirac mass,pseudomagnetic field catalyses the generation of time-reversal symmetry breaking Haldane mass. Various competingground states in monolayer graphene in pseudomagnetic fields were recently studied in Ref.[21]. Finally, we wouldlike to add that very strong 50 − T pseudomagnetic fields may be realized in molecular graphene .The interplay between different possible ground states in strained graphene in a magnetic field represents animportant unsolved problem at the moment. In the present paper, we study a gap generation for quasiparticles inmonolayer graphene in the presence of both constant magnetic and pseudomagnetic fields and using the model withlocal four-fermion interactions considered in Refs.[23,24]. Local four-fermion terms in the Hamiltonian are remnantsof the interactions on the atomic scale, and in spite being much smaller than Coulomb interaction, they play animportant role in deciding how the SU (4) symmetry is broken in monolayer graphene as well as bilayer graphene.This especially concerns the nature of the QH state with half-filled zero-energy Landau level. We obtain the phasediagram for competing quantum Hall states in the LLL approximation when the chemical potential is tuned to thecharge neutrality point, i.e., the state with the zero filling factor.The paper is organized as follows. We begin by presenting in Sec.II the model describing low-energy quasiparticlesexcitations in strained monolayer graphene in an external magnetic field and in the presence of local four-fermioninteractions. The derivation of the gap equation is given in Sec.III and its solutions are presented in Sec.IV. The phasediagram of the system is derived and discussed in Sec.V. The main results are summarized in Sec.VI. Appendices atthe end of the paper contain technical details and derivations used to supplement the presentation in the main text. II. HAMILTONIAN OF THE MODEL
The low-energy quasiparticles excitations in graphene can be described in terms of a four-component Dirac spinorΨ Tα = ( ψ KAα , ψ
KBα , ψ K ′ Bα , ψ K ′ Aα ) which combines the Bloch states with spin index α = 1 , A, B ) and with momenta near the two nonequivalent valley points (
K, K ′ ) of the Brillouin zone. The freequasiparticle Hamiltonian has a relativistic-like form with the Fermi velocity v F = 10 m/s playing the role of thespeed of light H = Z d r (cid:2) v F ¯Ψ( γ π x + γ π y )Ψ + ǫ Z Ψ † σ z Ψ (cid:3) , (1)where ¯Ψ = Ψ † γ is the Dirac conjugated spinor and r = ( x, y ). The matrices γ ν with ν = 0 , , × { γ µ , γ ν } = 2 g µν , where g µν = diag(1 , − , −
1) and µ, ν = 0 , ,
2. These matrices belong to a reducible representation of the Dirac algebra γ ν = ˜ τ z ⊗ ( τ z , iτ y , − iτ x ), wherethe Pauli matrices ˜ τ i and τ i with i = x, y, z act in the subspaces of the valley ( K, K ′ ) and sublattices ( A, B ) indices,respectively.The canonical momentum π = − i ~ ∇ + e A /c + γ γ e A /c includes the vector potential in the Landau gauge A = (0 , B ⊥ x ) corresponding to the component B ⊥ of an external magnetic field B orthogonal to the plane ofgraphene, and A is the vector potential describing the strain induced gauge fields . In the representation ofthe Dirac matrices that we use, γ and γ = iγ γ γ γ matrices equal γ = i ˜ τ y ⊗ I and γ = ˜ τ x ⊗ I , where I is the2 × γ γ = ˜ τ z ⊗ I is a matrix diagonal in thesubspace of valleys that ensures that the term in the canonical momentum with the vector potential A takes oppositesigns in K and K ′ valleys. In what follows, we will consider only the case of constant magnetic B and pseudomagnetic B fields. The pseudomagnetic field B points always in the direction perpendicular to the plane of graphene and,therefore, is described by the vector potential A = (0 , B x ), where B = | B | .The last term in the free Hamiltonian (1) is the Zeeman interaction ǫ Z = µ B B with µ B = e ~ / (2 mc ) being the Bohrmagneton and σ z is Pauli spin matrix whose eigenstates describe spin states directed along or against the magneticfield B . Here B = q B ⊥ + B || is the strength of the magnetic field and B || is its component parallel to the plane ofgraphene. We note that the standard Zeeman interaction µ B B σ can be reduced to this form using a rotation in spinspace.The Coulomb interaction between electrons is described by the following Hamiltonian: H C = 12 Z d rd r ′ ¯Ψ( r ) γ Ψ( r ) U C ( r − r ′ ) ¯Ψ( r ′ ) γ Ψ( r ′ ) , where U C ( r ) is the Coulomb potential. In order to simplify the analysis, we follow the approach of Ref.[27] and replacethe Coulomb interaction U C ( r ) by the contact interaction G int δ ( r ). The Hamiltonian H + H C in the absence of theZeeman term possesses a global SU (4) symmetry connected with valley and spin degrees of freedom.Although the Coulomb interaction is the strongest interaction between electrons in graphene, local four-fermioninteractions play a crucial role too. Although these interactions are much smaller than the Coulomb one, theybreak, in general, the SU (4) symmetry and crucially affect the selection of the ground state of the system. A setof local valley and sublattice asymmetric four-fermion interactions was introduced in Ref.[23]. The ν = 0 quantumHall state was studied and it was shown that the phase diagram, obtained in the presence of generic valley andsublattice anisotropy and the Zeeman interaction, consists of four phases: ferromagnetic, canted antiferromagnetic(CAF), charge density wave, and Kekule distortion. The Hamiltonian of generic local four-fermion interactions reads H contact = 12 Z d r X j,k g jk [ ¯Ψ( r ) γ T jk Ψ( r )] , T jk = ˜ τ j ⊗ τ k , (2)where j, k = x, y, z . We do not include in H contact the term with g as it corresponds to the local Coulomb interaction,which has already been taken into account by G int . In addition, we dot not include in our model the terms with g k and g j , which vanish in the first order in the Coulomb interactions and arise only in the second order due to virtualtransitions to other bands . The coupling constants g jk are not all independent. As shown in Ref.[23], symmetryand other considerations lead to the following equalities for nonzero constants: g ⊥⊥ = g xx = g xy = g yx = g yy , g ⊥ z = g xz = g yz , g z ⊥ = g zx = g zy . (3)Thus, totally we have four interaction coupling constants, G int , g ⊥⊥ , g ⊥ z , g z ⊥ , in the considered model. Finally, letus present T jk in terms of the γ -matrices T xx = − iγ γ , T xy = iγ γ , T xz = − γ γ , T yx = − γ γ , T yy = γ γ , T yz = iγ γ , T zx = iγ , T zy = − iγ , T zz = γ . (4)All these matrices are normalized as T ij = 1 . This presentation is useful for the derivation of the gap equation in thenext section.
III. GAP EQUATION
We will solve the gap equation in the Hartree-Fock (mean-field) approximation which is conventional andappropriate in this case. In the subsection III A, we will derive the gap equation in the case where only real magneticfield is present. In the next subsection, we will generalize the gap equation to the case where both magnetic andpseudomagnetic fields are present.
A. Magnetic field
At zero temperature and in the clean limit (no impurities), the Schwinger-Dyson equation for the quasiparticlepropagator G ( u, u ′ ) = ~ − h | T Ψ( u ) ¯Ψ( u ′ ) | i in graphene in the mean-field approximation takes the form iG − ( u, u ′ ) = iS − ( u, u ′ ) − ~ G int γ G ( u, u ) γ δ ( u − u ′ ) + ~ G int γ tr (cid:2) γ G ( u, u ) (cid:3) δ ( u − u ′ ) − ~ X j,k g jk { γ T jk G ( u, u ) γ T jk − γ T jk tr (cid:2) γ T jk G ( u, u ) (cid:3) } δ ( u − u ′ ) , (5)where u = ( t, r ). In this subsection, we will derive the gap equation in graphene in a magnetic field. The generalizationto the case of both magnetic and pseudomagnetic fields is rather straightforward and will be considered in the nextsubsection.The inverse free propagator in the case under consideration is given by iS − ( u, u ′ ) = [( i ~ ∂ t − ǫ Z σ z ) γ − v F ( π · γ )] δ ( u − u ′ ) . (6)For the full quasiparticle propagator, we will use an ansatz which is a generalization of the ansatz used in the previouswork by two of us iG − ( u, u ′ ) = [ i ~ ∂ t γ + µγ + ˜ µγ γ γ − v F ( π · γ ) − ˜∆ + ∆ γ γ ] δ ( u − u ′ ) , (7)where matrices µ, ˜ µ, ∆ , ˜∆ are defined as µ = µ ν σ ν , ˜ µ = ˜ µ ν σ ν , ∆ = ∆ ν σ ν , ˜∆ = ˜∆ ν σ ν , and index ν runs the values ν = 0 , x, z with σ x and σ z being Pauli spin matrices and σ the unit 2 × σ y matrix is consistent with subsequent analysis of a gap equation]. In what follows, we consider twelve dynamicallygenerated parameters µ ν , ˜ µ ν , ∆ ν , and ˜∆ ν as constant that is consistent with our mean-field analysis of the presentmodel with contact interactions.The parameters µ j and ˜ µ j with j = x, z are generalized chemical potentials connected with the QHferromagnetism . On the other hand, ∆ j and ˜∆ j are related to the magnetic catalysis scenario andare Haldane and Dirac masses, respectively, and correspond to excitonic condensates (for a brief review of the QHferromagnetism and magnetic catalysis scenario, see Refs.[39–41]). Actually, it was shown in Ref.[27] that the QHferromagnetism and magnetic catalysis scenario order parameters necessarily coexist. The physics underlying theircoexistence is specific for the systems with relativistic-like quasiparticle spectrum that makes the quantum Hall dynam-ics of the SU (4) breakdown in graphene to be quite different from that in conventional systems with non-relativisticquasiparticle spectrum.According to Eq.(7), the full propagator G ( u, u ′ ) can be written in the form G ( u, u ′ ) = i h u | h ( i ~ ∂ t + µ ) γ − v F ( π · γ ) + i ˜ µγ γ + i ∆ γ γ γ − ˜∆ i − | u ′ i , (8)where the states | u i are eigenstates of the time-position operator ˆ u : ˆ u | u i = u | u i , h u | u ′ i = δ ( u − u ′ ). In Appendix A,we derive an explicit expression for the propagator G ( u, u ′ ) in the form of a sum over Landau levels.The symmetry-breaking generalized chemical potentials and gaps µ ν , ˜ µ ν , ∆ ν , ˜∆ ν are related to the correspondingorder parameters through the following relationship: h ¯Ψ O ν Ψ i = − ~ tr[ O ν G ( u, u )] , (9)where 8 × O ν = γ σ ν , γ γ γ σ ν , γ γ σ ν , σ ν , respectively. Compared to our previous analysis, weincluded the spin matrix σ x in order to be able to describe the canted antiferromagnetic state.Since the right-hand side of the gap equation (5) contains the full propagator at the coincidence limit u ′ = u , weshould calculate G ( u, u ′ ) | u = u ′ = ¯ G ( u, u ), where ¯ G is the translation invariant part of the full propagator defined inthe mixed frequency-momentum representation in Eq.(A15). By making use of Eqs. (A19) and (A20), we find that G ( u, u ) = ∞ Z −∞ dωd k (2 π ) ~ ¯ G ( ω, k ) = i πl ∞ X n =0 ∞ Z −∞ dω π ~ W [ P − + P + θ ( n − M − nǫ B , (10)where l = p ~ c/ | eB ⊥ | is the magnetic length, ǫ B = p ~ v F | eB ⊥ | /c ≃ p | B ⊥ | [T]K is the Landau energy scale, P ± are projectors given by Eq.(A11), and W , M are matrices expressed through µ , ˜ µ , ∆, ˜∆ and defined in Eqs.(A2) and(A3). We note that the filling factor ν = 2 πl ρ is related to the carrier imbalance ρ = n e − n h , where n e and n h arethe densities of electrons and holes, respectively, and ρ is determined through the Green’s function as ρ = h | Ψ † ( u )Ψ( u ) | i = − ~ tr[ γ G ( u, u )] . (11)Since W and M contain only γ and γ γ Dirac matrices, it is convenient to work with eigenvectors of thesematrices. The equality ( γ γ ) = − | s i of the matrix γ γ are purely imaginary γ γ | s i = is | s i , s = ± . (12)Similarly, since ( γ ) = 1, the eigenvectors of the matrix γ are real and given by γ | s i = s | s i , s = ± . (13)Furthermore, since γ and γ γ commute, we can consider states | s s i which are simultaneously eigenvectors of γ and γ γ with eigenvalues s and is , respectively. The vectors | s s i form a complete basis, and since W and M contain only γ and γ γ matrices, the propagator G ( u, u ) is diagonal in the basis of | s s i vectors and is given by G ( u, u ; s , s ) = i πl Z ∞−∞ dw π ~ ∞ X n =0 [ s ( ω + m ) + d ][ c + c j σ j ] − × (cid:16) s s ⊥ + [1 − s s ⊥ ] θ ( n − / (cid:17) , j = x, y, z. (14)Here s ⊥ = sgn( eB ⊥ ), the matrices m and d are defined in Appendix A, and the coefficients c , c j are c = ω + 2 ωm − nǫ B + m ν − d ν , c y = 2 is ( d x m z − d z m x ) , (15) c x = 2( ωm x + m m x − d d x ) , c z = 2( ωm z + m m z − d d z ) , (16)where m ν = µ ν − s s ˜ µ ν , d ν = ˜∆ ν + s s ∆ ν ( ν = 0 , x, z ) and summation over dummy index ν is meant. For strongmagnetic fields, we write G ( u, u ) = G LLL ( u, u ) + G hLL ( u, u ) , (17)where we separated the contributions of the zero Landau level, G LLL ( u, u ) with n = 0, and higher Landau levels, G hLL ( u, u ) with n ≥
1, in Eq.(14).Let us calculate first the lowest Landau level propagator G LLL ( u, u ). In order to integrate over ω in Eq.(14), werewrite the integrand by using the relation c + c j σ j = [ s ( ω + m ) − d ][ s ( ω + m ) + d ] (18)valid for n = 0 and assume as usual that ω is replaced by ω + iǫ sgn ω and ǫ → + . Hence we obtain that Eq.(14)implies the following propagator at the limit of coinciding points in the LLL approximation: G LLL ( u, u ; s , s ) = i πl s s ⊥ Z ∞−∞ dω π ~ ( s ω + s m ν σ ν − d ν σ ν ) − = − s π ~ l s s ⊥ (cid:20) ( m i − s d i ) σ i E θ ( E − | µ L | ) + sgn( µ L ) θ ( | µ L | − E ) (cid:21) , i = x, z (19)where the factor (1 + s s ⊥ ) / s ⊥ = sgn( eB ⊥ ) reflects the presence of the spin projector P − = (1 − is ⊥ γ γ ) / E = p ( m x − s d x ) + ( m z − s d z ) , (20)and we introduced the notation µ L = m − s d = µ − s s ˜ µ − s ˜∆ − s ∆ (21)for an ”effective chemical potential” in the lowest Landau level.By integrating over ω , it is not difficult to check that the higher Landau level contribution G hLL ( u, u ) diverges as P ∞ n =1 n − / . Indeed, making the change of the variable ω → √ nǫ B ω and taking into account that all the dynamicallygenerated parameters are much less than the scale ǫ B , we find that the leading contribution at large nǫ B is given by G hLL ( u, u ; s , s ) ≃ π ~ l ∞ X n =1 √ n d ν σ ν ǫ B . (22)This means that the right-hand side of the gap equation (5) diverges too. This result is the well-known artefact ofusing a model with local four-fermion interactions. For a long-range interaction like, for example, the Coulomb oneconsidered in Ref.[42], such a divergence is absent because the gap equation contains the quasiparticle propagatorat different points u and u ′ . To proceed further, we regularize the divergence in the model under considerationintroducing a cutoff n max in the sum over Landau levels (a slightly different approach was used in Ref.[27]), which isconnected with the ultraviolet (UV) cut-off in energy Λ (band width) according to the relation n max = Λ /ǫ B . Byusing the regularization described above and retaining only the leading contribution, we find that the higher Landaulevels contribution to the propagator at the limit of coinciding points is given by G hLL ( u, u ; s , s ) = Λ4 π ~ v F d ν σ ν . (23)By combining Eqs.(23) and (19), the gap equation (5) takes the following final form: s m ν σ ν − d ν σ ν + s s ⊥ ǫ Z σ z = − ( G int + g zz ) G ( s , s ) − g ⊥⊥ G ( s , − s ) + 2 g z ⊥ G ( − s , − s )+ 2 g ⊥ z G ( − s , s ) + g zz X s ′ ,s ′ tr G ( s ′ , s ′ ) , (24)where G ( s , s ) ≡ G LLL ( u, u ; s , s ) + G hLL ( u, u ; s , s )= − s π ~ l (cid:20) ( m i − s d i ) σ i E θ ( E − | µ L | ) + sgn ( µ L ) θ ( | µ L | − E ) (cid:21) s s ⊥ π ~ v F d ν σ ν , (25)and trace in the last term in Eq.(24) is taken over the Pauli spin matrices [note that the quantities E and µ L dependon s , s according to Eqs.(20),(21) and m ν , d ν depend on them too]. In deriving Eq.(24), we omitted the third termon the right-hand side of Eq.(5) which defines the Hartree contribution due to the charge density of carriers. Thepoint is that there are other contributions due to the charges of ions in graphene and the charges in the substrate andgates. In view of the overall neutrality of the system, all these contributions should cancel exactly (the Gauss law).Finally, we would like to note that it is advantageous in deriving Eq.(24) to use T jk matrices given by Eq.(4) andutilize the Dirac algebra in order to calculate the contribution to the gap equation due to the last term in Eq.(5).[The matrices T jk have simple commutation relations with the Green’s function G ( u, u ) which contains only γ and γ γ matrices (see Eq.(10)).] B. Magnetic and pseudomagnetic fields
In this subsection, we will derive the gap equation for quasiparticles in graphene in the case where both constantmagnetic and pseudomagnetic fields are present. Then the canonical momentum has the form π = − i ~ ∇ + e A /c + γ γ e A /c , where A = (0 , B ⊥ x ) and A = (0 , B x ). Repeating the same computations as in the previous subsection,one can show that the only difference between the former and present cases is that the magnetic field B ⊥ is nowreplaced by the effective field B ⊥ + ( iγ γ γ ) B or B ⊥ − s s B in the eigenstate basis of the matrices γ and γ γ .Consequently, the pseudomagnetic field has opposite signs in the K and K ′ valleys. Therefore, in order to take intoaccount pseudomagnetic field, we should simply make the replacement B ⊥ → B ⊥ − s s B in the correspondingequations of Subsec.III A except the Zeeman energy, where ǫ Z = µ B q B ⊥ + B k , which includes the component ofthe magnetic field parallel to the plane of graphene. We find it convenient in the analysis below to use the notation b k = B k /B ⊥ and b = B /B ⊥ .Taking into account all Landau levels contributions, the gap equation for quasiparticles in graphene in the presenceof constant magnetic and pseudomagnetic fields is given by s m ν σ ν − d ν σ ν + s s ⊥ ǫ Z σ z = − ( G int + g zz ) G (5) ( s , s ) − g ⊥⊥ G (5) ( s , − s )+ 2 g z ⊥ G (5) ( − s , − s ) + 2 g ⊥ z G (5) ( − s , s ) + g zz X s ′ ,s ′ tr G (5) ( s ′ , s ′ ) , (26)where G (5) ( s , s ) = − s | eB ⊥ − s s eB | π ~ c (cid:20) ( m i − s d i ) σ i E θ ( E − | µ L | ) + sgn ( µ L ) θ ( | µ L | − E ) (cid:21) × s sgn( eB ⊥ − s s eB )2 + Λ4 π ~ v F d ν σ ν . (27)The first term on the right-hand side of Eq.(27) is the LLL contribution and the last one describes the contributiondue to higher Landau levels.We have found solutions of the gap equation taking into account the contributions due to all Landau levels. Itturned out that the higher Landau levels contribution in the weak coupling regime does not qualitatively change theresults obtained in strong magnetic fields when the LLL approximation is valid. On the other hand, the higher Landaulevels contribution essentially enlarges formulas and makes them very complicated. Therefore, in what follows, we willsolve the gap equation and present our analysis in the LLL approximation omitting the contribution due to higherLandau levels. IV. SOLUTIONS OF GAP EQUATION IN THE LLL APPROXIMATION
In this section, we consider solutions of the gap equation (26) retaining only the LLL contribution. Propagator(27) in the LLL approximation contains different projectors depending on which field B ⊥ or B is stronger. We willconsider both possibilities separately. A. | B ⊥ | > | B | Let us find solutions of the gap equation (26) in the LLL approximation in the case where magnetic field is strongerthan pseudomagnetic field. We will use the notation of Refs.[23,24] u = G int πl , u z = g zz πl , u ⊥ = g ⊥ z πl . (28)The gap equation in the LLL approximation is obtained from Eq.(26) by replacing the full fermion Green’s function G (5) with G (5) LLL and multiplying both sides by the LLL projector (1+ s s ⊥ ) / s ⊥ = +). It is easy to see that in the LLL approximation twelve parameters µ ν , ˜ µ ν and ∆ ν , ˜∆ ν enter in combinations s ( µ ν − ∆ ν ) − (˜ µ ν − ˜∆ ν ). Therefore, we have only six independent variables µ ν − ∆ ν and ˜ µ ν − ˜∆ ν . Without loss ofgenerality we can put µ ν = ˜ µ ν = 0 so that m ν = 0 and we are left only with parameters d ν . Then the gap equationtakes the following form: − d ν ( s ) σ ν + s ǫ Z σ z = − (1 − s b ) 12 ( u + u z ) (cid:20) d i ( s ) σ i E ( s ) θ ( E ( s ) − | d ( s ) | ) + sgn( d ( s )) θ ( | d ( s ) | − E ( s )) (cid:21) + (1 + s b ) u ⊥ (cid:20) d i ( − s ) σ i E ( − s ) θ ( E ( − s ) − | d ( − s ) | ) + sgn( d ( − s )) θ ( | d ( − s ) | − E ( − s )) (cid:21) + u z X s ′ (1 − s ′ b ) sgn( d ( s ′ )) θ ( | d ( s ′ ) | − E ( s ′ )) , (29)where magnetic length l = p ~ c/ | eB ⊥ | and b = B /B ⊥ . We find the following solutions of the gap equation (29)(solutions in the case of purely magnetic field are considered in Appendix B 1):(i) Ferromagnetic (F) solution:˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = ∓ b u + u z − u ⊥ ) , ∆ z = ǫ Z ±
12 ( u + u z + 2 u ⊥ ) . (30)This solution exists for | ˜∆ z | < | ∆ z | . Using Eqs.(9) and (27), it is easy to check that, for b →
0, this solution ischaracterized by the unique order parameter h Ψ † σ z Ψ i which defines uniform magnetization. According to Eq.(7),the corresponding term ∆ z σ z γ γ in the inverse full propagator describes the Haldane-type mass antisymmetricin spin.(ii) Antiferromagnetic (AF) solution:˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = ±
12 ( u + u z − u ⊥ ) , ∆ z = ǫ Z ∓ b u + u z + 2 u ⊥ ) . (31)This solution exists for | ˜∆ z | > | ∆ z | . If we neglect ∆ z compared to ˜∆ z , the unique order parameter characterizingthis solution is h ¯Ψ σ z Ψ i which defines staggered magnetization. Obviously, the term ˜∆ z σ z in the inverse fullpropagator (7) corresponds to the Dirac-type mass antisymmetric in spin.(iii) Canted antiferromagnetic (CAF) solution:˜∆ = ∆ = ∆ x = 0 , ˜∆ x = ± s(cid:18) − b cos θ (cid:19) ( u + u z − u ⊥ ) sin θ, ∆ z = cos θ u + u z − u ⊥ ) , ˜∆ z = − b θ ( u + u z − u ⊥ ) , cos θ = − ǫ Z u ⊥ . (32)This solution exists when ˜∆ x is real. For b →
0, this solution is characterized by two nonzero order parameters.They are the staggered magnetization in the x direction h ¯Ψ σ x Ψ i and the uniform magnetization h Ψ † σ z Ψ i in the z direction. It is important that the order parameter h ¯Ψ σ x Ψ i can not be transformed into h ¯Ψ σ z Ψ i by means ofa rotation in spin space without inducing the uniform magnetization in the y direction.(iv) Charge density wave (CDW) solution:˜∆ x = ∆ x = ˜∆ z = 0 , ∆ z = ǫ Z , ˜∆ = ±
12 ( u − u z − u ⊥ ) , ∆ = ∓ b u + u z + 2 u ⊥ ) . (33)For b →
0, the order parameter which characterizes this solution is h ¯ΨΨ i . The corresponding term ˜∆ in theinverse full propagator (7) describes the Dirac mass.It is easy to check that for b = 0 pseudomagnetic field induces additional staggered magnetization h ¯Ψ σ z Ψ i ∼ b inthe F solution, and vice versa, additional uniform magnetization h Ψ † σ z Ψ i ∼ b in the AF solution. As to the CAFsolution, pseudomagnetic field produces additional staggered magnetization in the z direction. Pseudomagnetic fieldleads also to the generation of the Haldane mass in the CDW solution in addition to the Dirac mass.The carrier density for strained graphene in the LLL approximation for | B ⊥ | > | B | is given by ρ = 12 πl X s (1 − s b ) sgn( µ L ( s )) θ ( | µ L ( s ) | − E ( s )) , (34)where µ L ( s ) = − ( s ˜∆ + ∆ ) and E ( s ) = q ( s ˜∆ x + ∆ x ) + ( s ˜∆ z + ∆ z ) . One can check that the carrier density ρ equals zero for the F, AF, and CAF solutions, whereas for the CDW state ρ is nonzero, thus prohibiting it. It is easyto show that for B → ρ CDW ∼ B . Consequently, for B = 0, the CDW solution becomes admissibleand coincides with the corresponding solution in the case of purely magnetic field obtained in Appendix B 1. Clearly,the F, AF, CAF solutions for B = 0 also reduce to those found in Appendix B 1. B. | B ⊥ | < | B | Similarly to the previous subsection it is more convenient in the case | b | > G int , g zz , and g ⊥⊥ : u = G int πl , u z = g zz πl , u = g ⊥⊥ πl , (35)where l = p ~ c/ | eB | and for specificity we take eB >
0. The gap equation in the LLL approximation is obtainedfrom Eq.(26) by replacing the full fermion Green’s function G (5) with G (5) LLL and multiplying both sides by the LLLprojector (1 − s ) /
2. This time twelve parameters µ ν , ˜ µ ν and ∆ ν , ˜∆ ν enter in combinations ( µ ν + ˜∆ ν ) + s (˜ µ ν − ∆ ν ).Therefore, we have only six independent variables µ ν + ˜∆ ν and ˜ µ ν − ∆ ν . Without loss of generality we can put µ ν = ˜ µ ν = 0 so that m ν = 0, and we are left only with parameters d ν .By making use of m ν = 0 and constricting the propagator on the s = − subspace, we obtain the following equationfor parameters d ν : − d ν ( s ) σ ν − ǫ Z σ z = − (1 + s b − ) 12 ( u + u z ) (cid:20) d i ( s ) σ i E ( s ) θ ( E ( s ) − | d ( s ) | ) + sgn( d ( s )) θ ( | d ( s ) | − E ( s )) (cid:21) − (1 − s b − )2 u (cid:20) d i ( − s ) σ i E ( − s ) θ ( E ( − s ) − | d ( − s ) | ) + sgn( d ( − s )) θ ( | d ( − s ) | − E ( − s )) (cid:21) + u z X s ′ (1 + s ′ b − ) sgn( d ( s ′ )) θ ( | d ( s ′ ) | − E ( s ′ )) , (36)where d ν = ˜∆ ν − s ∆ ν and E ( s ) = p d x + d z .We find the following solutions of the above gap equation (solutions in the case of purely pseudomagnetic field areconsidered in Appendix B 2):(i) F solution:˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = − ǫ Z ∓ b ( u + u z + 4 u ) , ∆ z = ±
12 ( u + u z − u ) . (37)This solution exists for | ˜∆ z | < | ∆ z | .(ii) AF solution:˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = − ǫ Z ±
12 ( u + u z + 4 u ) , ∆ z = ∓ b ( u + u z − u ) . (38)This solution exists for | ˜∆ z | > | ∆ z | .(iii) CAF solution: ˜∆ = ∆ = ˜∆ x = 0 , ∆ x = ± s (1 − z ) (cid:18) − b z (cid:19)
12 ( u + u z − u ) , ˜∆ z = z u + u z − u ) , ∆ z = − b z ( u + u z − u ) , z = ǫ Z u . (39)This solution exists only when ∆ x is real. According to the analysis performed in Appendix B 2, the CAFsolution is absent in the case where only pseudomagnetic field is present. Therefore, the question arises whathappens with the CAF solution found here as B || = 0 and B ⊥ →
0. Obviously, ˜∆ z vanishes in this limit because z is proportional to ǫ Z and tends to zero. On the other hand, b z does not depend on B ⊥ and, therefore, isnot sensitive to the limit B ⊥ →
0. Furthermore, the Hamiltonian of the system in the absence of magneticfield is invariant with respect to SU (2) spin rotations because the Zeeman term vanishes in this case. Thensolutions with different components ∆ i ( i = x, y, z ) in ∆ = ∆ x σ x + ∆ y σ y + ∆ z σ z in ansatz (7) but with thesame q ∆ x + ∆ y + ∆ z are physically equivalent because they can be easily connected by means of appropriate SU (2) spin rotations. It is not difficult to check that p ∆ x + ∆ z for the CAF solution (39) equals | ∆ z | of theF solution (37) as well as | ∆ z | of the F solution (B11) in the case where only pseudomagnetic field is present.Thus, the F and CAF solutions found in this subsection are physically equivalent if magnetic field is absent.(iv) QAH solution:˜∆ x = ∆ x = ∆ z = 0 , ˜∆ z = − ǫ Z , ˜∆ = ∓ b ( u − u z + 4 u ) , ∆ = ±
12 ( u + u z − u ) . (40)Using Eq.(27), we find that the carrier density for strained graphene in the LLL approximation equals ρ = 12 πl X s (1 + s b − ) sgn( µ L ( s )) θ ( | µ L ( s ) | − E ( s )) , (41)where µ L ( s ) = ˜∆ − s ∆ and E ( s ) = q ( ˜∆ x − s ∆ x ) + ( ˜∆ z − s ∆ z ) . For the F, AF, and CAF solutions,the carrier density equals zero while it is nonzero for the QAH solution. Thus, the QAH state is not realized for B ⊥ = 0. For B ⊥ → ρ QAH ≃ − sgn∆ π ~ c eB ⊥ . (42)Therefore, for B ⊥ = 0, the QAH solution becomes admissible and coincides with the corresponding solution inAppendix B 2. Other solutions also reduce to those in Appendix B 2 if B || = 0 and B ⊥ = 0. V. PHASE DIAGRAM OF ν = 0 QH STATES IN MAGNETIC AND PSEUDOMAGNETIC FIELDS
We found above several solutions of the gap equation for quasiparticles in graphene in magnetic and pseudomagneticfields. In order to determine which of these solutions is the ground state, we should calculate their energy densitiesand then find out the phase diagram of the system. We used the Baym-Kadanoff formalism in order to calculatethe energy density Ω of the system:Ω = − π ~ c X s s | eB ⊥ − s s eB | (cid:20)(cid:18) E + ǫ Z s d z − m z E (cid:19) θ ( E − | µ L | ) + | d − s m | θ ( | µ L | − E ) (cid:21) × s sgn( eB ⊥ − s s eB )2 , (43)where m ,z , d ,z , µ L , and E are functions of discrete variables s , s (see Eqs.(20),(21) and definitions of m ,z , d ,z after Eq.(16)). The derivation of the energy density Ω is given in Appendix C and we retained in Eq.(43) only thecontribution due to the lowest Landau level. Therefore, our analysis is valid when the magnitude of dynamicallygenerated gaps is much less than the maximum of Landau gaps p ~ v F | eB ⊥ | /c or p ~ v F | eB | /c .By making use of the energy density (43) and the solutions found in Sec.IV, we easily calculate the following energydensities for the solutions in the case | b | < F = | u ⊥ | ∓ ǫ Z πl − C < , Ω AF = | u ⊥ | b ± b ǫ Z πl − C < , Ω CAF = − ǫ Z | u ⊥ | πl − C < , C < = (1 + b )( u + u z − u ⊥ )4 πl , (44)where signs in the expressions for Ω F and Ω AF correlate with the corresponding ones in solutions (30) and (31). Theseenergy densities in the nonstrained limit equal up to a constant to the corresponding energies found in Refs.[23,24].For | b | >
1, we haveΩ F = − u b − ± ǫ Z b − πl − C > , Ω AF = − u ∓ ǫ Z πl − C > , Ω CAF = ǫ Z u πl − C > , C > = (1 + b − )( u + u z − u )4 πl , (45)0 CAFF AF - - - - - - b b þ FIG. 1: (Color online) The phase diagram of the ν = 0 QH states in graphene in the plane ( b || , b ) for g ⊥⊥ >
0. The bluecolor represents the CAF state, green - AF state, and yellow - F state. This diagram is obtained for B ⊥ = 20 T , u ⊥ = − K ,2 | e | g ⊥⊥ / ( µ B π ~ c ) = 2 .
8, and ǫ Z = 13 . q b || K . FCAF CAFAF - - - - - - b b þ AFCAF CAFF - - - - - - b b þ FIG. 2: (Color online) The phase diagram of the ν = 0 QH states in graphene in the plane ( b || , b ) for − µ B π ~ c/ (2 | e | ) < g ⊥⊥ < g ⊥⊥ < − µ B π ~ c/ (2 | e | ) (right panel). The blue color represents the CAF state, green - AF state, and yellow -F state. These diagrams similarly to Fig.1 are obtained for B ⊥ = 20 T , u ⊥ = − K , and ǫ Z = 13 . q b || K . For the leftpanel, 2 | e | g ⊥⊥ / ( µ B π ~ c ) = − .
8, whereas 2 | e | g ⊥⊥ / ( µ B π ~ c ) = − . where signs in the expressions for Ω F and Ω AF correlate with the corresponding ones in solutions (37) and (38). Theterms C < and C > describe a common shift in energy density for three solutions, therefore, they are not important fordetermining the ground state. Hence the energy density depends only on two essential four-fermion coupling constants g ⊥ z and g ⊥⊥ (related to u ⊥ and u , respectively).Using the energy densities Eqs.(44) and (45), we find the phase diagram of the ν = 0 QH states in graphene inconstant magnetic and pseudomagnetic fields. More precisely, we fix the value of perpendicular magnetic field B ⊥ and study the phase diagram in the plane ( b || , b ), where b || = B || /B ⊥ and b = B /B ⊥ . In our analysis, followingKharitonov , we assume that u ⊥ <
0. Then for 2 | u ⊥ | > ǫ Z (recall that ǫ Z = µ B | B ⊥ | q b || ), we find that thephase diagram depends on coupling constant g ⊥⊥ .The simplest phase diagram is realized for positive g ⊥⊥ , which we plot in Fig.1. The CAF solution is found to bethe ground state of the system in the region defined by | b | < ǫ Z ( b || ) / | u ⊥ | <
1, where this solution exists. The sideborders of the region with the CAF solution are determined from the condition | b | = cos θ = ǫ Z ( b || ) / (2 | u ⊥ | ). Thismakes the central blue regions, where the CAF state is realized, look like biconcave lens. The top and bottom bordersof the region with the CAF solution are determined from the condition cos θ = ǫ Z ( b || ) / (2 | u ⊥ | ) = 1.1The ferromagnetic solution with the minus sign in Eq.(44) is the ground state for | b || | > b cr || and | b | <
1, where b cr || is determined by ǫ Z ( b cr || ) = 2 | u ⊥ | , and the AF solution is the most preferable solution for the regions | b | > | b || | < b cr || , 1 > | b | > ǫ Z ( b || ) / | u ⊥ | . As to the AF solution, plus and minus sign for it in Eq.(44) is preferable for b < b >
0, respectively. The regions of the CAF, F, and AF states are shown in Fig.1 by blue, yellow, andgreen colors, respectively.The phase diagram of the system is somewhat more complicated for negative g ⊥⊥ . There are two qualitativelydifferent cases, namely, − µ B π ~ c/ (2 | e | ) < g ⊥⊥ < g ⊥⊥ < − µ B π ~ c/ (2 | e | ). The corresponding phase diagramsare plotted in the left and right panels of Fig.2, respectively, where like in Fig.1 the blue color represents the CAFstate, yellow - F state, and green - AF state. Note the symmetry of the phase diagrams with respect to the change b → − b and b || → − b || . The main difference of the phase diagram in the left panel of Fig.2 compared to that inFig.1 is the appearance of two additional regions, where the CAF state is realized as the ground state of the system.The reason for this is that the terms with u in the energy densities of the AF and CAF solutions (45) for | b | > g ⊥⊥ >
0, the energy density of the AF state is smaller than that of the CAFstate. The situation changes for g ⊥⊥ <
0, where the energy density of the CAF state is smaller than that of the AFstate for sufficiently small b || . As b || increases, the CAF state ceases to exist for ǫ Z > | u | and transforms into theAF state. The side borders of two regions with the CAF solution in the left panel in Fig.2 are determined by thecondition ( µ B π ~ c/ (2 | e || g ⊥⊥ | )) q b || = | b | .For g ⊥⊥ < − µ B π ~ c/ (2 | e | ), the phase diagram becomes even more complicated. Though in this case for | b | > µ B π ~ c/ (2 | e || g ⊥⊥ | )) q b || = | b | . The CAF-Fborders are constant in b || and determined from the condition (2 | e || g ⊥⊥ | / ( µ B π ~ c )) = q b || . The F solution ismore preferable than the AF solution for small b || and | b | > b = 0 in our phase diagrams. According to Figs.1 and 2,there is the continuous phase transition from the CAF state to the F state as parallel magnetic field B || increases thatagrees with Kharitonov’s findings. VI. SUMMARY AND DISCUSSIONS
In the present work, we studied the ν = 0 quantum Hall states in strained monolayer graphene under tilted externalmagnetic field in the presence of local Coulomb and other local four-fermion interactions, which, in general, breakthe approximate spin-valley SU(4) symmetry of the low-energy electron Hamiltonian of graphene. Solving the gapequation in the LLL approximation, we found a rich phase diagram where different QH states, the ferromagnetic,antiferromagnetic, and canted antiferromagnetic, compete when we change parallel magnetic field B || and pseudo-magnetic field B (we keep perpendicular magnetic field B ⊥ fixed and rather large). For zero strain, these states arecharacterized by the nontrivial order parameters Ψ † σ z Ψ, ¯Ψ σ z Ψ, and ¯Ψ σ x Ψ and Ψ † σ z Ψ, respectively. In the presenceof pseudomagnetic field these order parameters gain an admixture of uniform staggered magnetization or uniformmagnetization in the z direction.Assuming that the strength of the Coulomb interaction is larger than all other local interactions, we found anessential dependence of the phase diagram only on two four-fermion couplings g ⊥ z and g ⊥⊥ , which appear in thelow-energy effective Hamiltonian (2). Our main results are accumulated in Figs.1 and 2. For negative coupling g ⊥ z we found that the canted antiferromagnetic state is always preferable in the center of all figures, i.e., for not too largevalues of B || and B . When B || and B increase, the ferromagnetic, antiferromagnetic, and canted antiferromagneticstates are realized as the ground state of the system depending on the values of B || , B .The particular case of purely magnetic field corresponds to the line B = 0 in our Figs.1 and 2 where there is a phasetransition from the CAF state to the F state as parallel magnetic field B || increases that agrees with Kharitonov‘sfindings. . In the case of purely pseudomagnetic field considered in Appendix B 2, we find that the F and QAHsolutions for g ⊥⊥ < g ⊥⊥ >
0, the AF state with Dirac-type mass is the most preferable.These results and those obtained in Sec.V show that the structure of the phase diagrams in Figs. 1 and 2 is sensitiveto the sign and the strength of four-fermion coupling g ⊥⊥ as well as coupling g ⊥ z .In our analysis, we ignored boundaries of graphene. However, it is known that they may play an essential role instrained graphene. A novel magnetic ground state, where the Neel and ferromagnetic orders coexist, was reported forthe Hubbard Hamiltonian in strained graphene in Ref.[44]. Whereas the Neel order takes the same sign through the2entire system, the magnetization at the boundary takes the opposite sign from the bulk. Since the total magnetizationvanishes, the magnetic ground state is edge-compensated antiferromagnet.In conclusion, strained graphene in a magnetic field provides a unique opportunity to observe various symmetrybreaking phases. Further progress in achieving strain induced pseudomagnetic fields, especially uniform pseudomag-netic fields, might allow one to probe experimentally the obtained phase diagram. As for the future, it would beinteresting to extend the results of the present paper beyond the neutral point with the filling factor ν = 0 anddescribe the quantum Hall states with other filling factors. VII. ACKNOWLEDGMENTS
We are grateful to S.G. Sharapov for useful suggestions and discussions. The work of E.V.G. and V.P.G. wassupported partially by the European IRSES Grant No. SIMTECH No. 246937 and by the Program of FundamentalResearch of the Physics and Astronomy Division of the NAS of Ukraine.
Appendix A: Quasiparticle propagator: Expansion over LLs
For the Fourier transform in time of the full propagator G ( u, u ′ ) in Eq. (8) we can write G ( ω, r , r ′ ) = i h r | (cid:2) ( ω + m ) γ − v F ( π · γ ) − d (cid:3) − | r ′ i = i h r | (cid:2) ( ω + m ) γ − v F ( π · γ ) + d (cid:3) × (cid:2)(cid:0) ( ω + m ) γ − v F ( π · γ ) − d (cid:1) (cid:0) ( ω + m ) γ − v F ( π · γ ) + d (cid:1)(cid:3) − | r ′ i = i [ W − v F ( π r · γγγ )] h r | (cid:0) M − v F π − ieB ⊥ ( ~ v F /c ) γ γ (cid:1) − | r ′ i , (A1)where m = µ + iγ γ γ ˜ µ , d = ˜∆ − iγ γ γ ∆, the matrices µ, ˜ µ, ∆ , ˜∆ are defined after Eq.(7), and matrices W and M are W = ( ω + m ) γ + d, m = m ν σ ν , d = d ν σ ν , ν = 0 , x, z, (A2) M = c + c i σ i , i = x, y, z, (A3)with c = ω + 2 ωm + m ν − d ν , c x = 2( ωm x + m m x − d d x ) ,c y = 2 iγ ( d x m z − d z m x ) , c z = 2( ωm z + m m z − d d z ) . (A4)Our aim is to find an expression for the propagator (A1) as an expansion over LLs (we follow below the considerationin Appendix A in Ref.[27]). The operator π has well known eigenvalues (2 n + 1) ~ | eB ⊥ | /c with n = 0 , , , . . . andits normalized wave functions in the Landau gauge A = (0 , B ⊥ x ) are ψ np ( r ) = 1 √ πl p n n ! √ π H n (cid:16) xl + pl (cid:17) e − l ( x + pl ) e ipy , (A5)where H n ( x ) are the Hermite polynomials and l = p ~ c/ | eB ⊥ | is the magnetic length. These wave functions satisfythe conditions of normalizability Z d rψ ∗ np ( r ) ψ n ′ p ′ ( r ) = δ nn ′ δ ( p − p ′ ) , (A6)and completeness ∞ X n =0 ∞ Z −∞ dpψ ∗ np ( r ) ψ np ( r ′ ) = δ ( r − r ′ ) . (A7)Using the spectral expansion of the unit operator (A7), we obtain h r | (cid:0) M − v F π − ieB ⊥ ( ~ v F /c ) γ γ (cid:1) − | r ′ i = 12 πl exp (cid:18) − ( r − r ′ ) l − i ( x + x ′ )( y − y ′ )2 l (cid:19) × ∞ X n =0 M − (2 n + 1)( ~ v F /c ) | eB ⊥ | − i ( ~ v F /c ) eB ⊥ γ γ L n (cid:18) ( r − r ′ ) l (cid:19) , (A8)3where we integrated over p by using ∞ Z −∞ e − x H m ( x + y ) H n ( x + z ) dx = 2 n π / m ! z n − m L n − mm ( − yz ) , (A9)assuming m ≤ n . Here L αn are the generalized Laguerre polynomials, and L n ≡ L n . The matrix i ( ~ v F /c ) eB ⊥ γ γ has eigenvalues ± ~ v F | eB ⊥ | /c = ± ǫ B /
2. Therefore, we have L n ( ξ ) M − (2 n + 1)( ~ v F /c ) | eB ⊥ | − i ( ~ v F /c ) eB ⊥ γ γ = P − L n ( ξ ) M − nǫ B + P + L n ( ξ ) M − ( n + 1) ǫ B , (A10)where ξ = ( r − r ′ ) / (2 l ) and the projectors P ± are P ± = 12 (cid:2) ± iγ γ sign( eB ⊥ ) (cid:3) . (A11)By redefining n → n − h r | [ M − v F π − ieB ⊥ ( ~ v F /c ) γ γ ] − | r ′ i = 12 πl e i Φ( r , r ′ ) e − ξ/ ∞ X n =0 P − L n ( ξ ) + P + L n − ( ξ ) M − nǫ B , (A12)where L − ≡ r , r ′ ) = − ( x + x ′ )( y − y ′ )2 l = − e ~ c r Z r ′ dz i A i ( z ) (A13)is the Schwinger phase which appears due to the noncommutative character of magnetic translations . Since π x e i Φ = e i Φ ~ (cid:18) − i∂ x − y − y ′ l (cid:19) , π y e i Φ = e i Φ ~ (cid:18) − i∂ y + x − x ′ l (cid:19) , (A14)propagator (A1) can be presented as a product of the phase factor and a translation invariant part ¯ G ( ω ; r − r ′ ), G ( ω ; r , r ′ ) = e i Φ( r , r ′ ) ¯ G ( ω ; r − r ′ ) , (A15)where¯ G ( ω ; r − r ′ ) = i (cid:20) W − ~ v F γ (cid:18) − i∂ x − y − y ′ l (cid:19) − ~ v F γ (cid:18) − i∂ y + x − x ′ l (cid:19)(cid:21) e − ξ/ πl ∞ X n =0 P − L n ( ξ ) + P + L n − ( ξ ) M − nǫ B . (A16)The Fourier transform of the translation invariant part of propagator (A16) can be evaluated by performing theintegration over the angle, π Z dθe ikr cos θ = 2 πJ ( kr ) , (A17)where J ( x ) is the Bessel function, and then using the following formula: Z ∞ xe − αx L n (cid:18) βx (cid:19) J ( xy ) dx = ( α − β ) n α n +1 e − α y L n (cid:18) βy α ( β − α ) (cid:19) , (A18)valid for y > α >
0. We obtain¯ G ( ω, k ) = ie − k l ∞ X n =0 ( − n D n ( ω, k ) M − nǫ B , (A19)with D n ( ω, k ) = 2 W (cid:2) P − L n (cid:0) k l (cid:1) − P + L n − (cid:0) k l (cid:1)(cid:3) + 4 ~ v F ( k · γ ) L n − (cid:0) k l (cid:1) , L α − ≡ , (A20)describing the n th Landau level contribution.4 Appendix B: Solutions of gap equation in purely magnetic or pseudomagnetic field1. Purely magnetic field
Let us find solutions of the gap equation (24) in the LLL approximation in the case where only magnetic field ispresent. In this case, the gap equation takes the form: − d ν ( s ) σ ν + s ǫ Z σ z = −
12 ( u + u z ) (cid:20) d i ( s ) σ i E ( s ) θ ( E ( s ) − | d ( s ) | ) + sgn( d ( s )) θ ( | d ( s ) | − E ( s )) (cid:21) + u ⊥ (cid:20) d i ( − s ) σ i E ( − s ) θ ( E ( − s ) − | d ( − s ) | ) + sgn( d ( − s )) θ ( | d ( − s ) | − E ( − s )) (cid:21) + u z X s ′ sgn( d ( s ′ )) θ ( | d ( s ′ ) | − E ( s ′ )) , (B1)where d ν = ˜∆ ν + s ∆ ν and E = p d x + d z . Projecting on the Pauli matrices, we get d ( s ) = 12 ( u − u z ) sgn( d ( s )) θ ( | d ( s ) | − E ( s )) − ( u ⊥ + u z ) sgn( d ( − s )) θ ( | d ( − s ) | − E ( − s )) , (B2) d x ( s ) = 12 ( u + u z ) d x ( s ) E ( s ) θ ( E ( s ) − | d ( s ) | ) − u ⊥ d x ( − s ) E ( − s ) θ ( E ( − s ) − | d ( − s ) | ) , (B3) d z ( s ) − s ǫ Z = 12 ( u + u z ) d z ( s ) E ( s ) θ ( E ( s ) − | d ( s ) | ) − u ⊥ d z ( − s ) E ( − s ) θ ( E ( − s ) − | d ( − s ) | ) . (B4)These equations define a system of non-linear equations because E is a non-linear function of d x and d z . Since G int approximates the Coulomb interaction, which is the strongest electron-electron interaction in graphene, we assume inwhat follows that u ≫ ǫ Z , u z , u ⊥ . The system of equations (B2)-(B4) has the following solutions which are consistentwith the charge neutrality condition for both ± signs:(i) F solution: ˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = 0 , ∆ z = ǫ Z ±
12 ( u + u z + 2 u ⊥ ) . (B5)(ii) AF solution: ˜∆ = ∆ = ˜∆ x = ∆ x = 0 , ˜∆ z = ±
12 ( u + u z − u ⊥ ) , ∆ z = ǫ Z , (B6)which exists for | ˜∆ z | > ǫ Z that is satisfied automatically since we assumed u ≫ ǫ Z , u z , u ⊥ .(iii) CAF solution: ˜∆ = ∆ = 0 , ∆ x = ˜∆ z = 0 , cos θ = − ǫ Z u ⊥ , ˜∆ x = ±
12 ( u + u z − u ⊥ ) sin θ, ∆ z = 12 ( u + u z − u ⊥ ) cos θ. (B7)The CAF solution exists only for ǫ Z < | u ⊥ | .(iv) CDW solution: ∆ = ∆ x = ˜∆ x = ˜∆ z = 0 , ∆ z = ǫ Z , ˜∆ = ±
12 ( u − u z − u ⊥ ) . (B8)5The ferromagnetic and CAF solutions reproduce solutions obtained in Refs.[23,24]. Note that here like in Sec.V wefollow Kharitonov and assume that u ⊥ <
0. In this case, the CAF solution is the ground state of the system forsufficiently small parallel magnetic field.Since in the present paper we consider the state with zero filling factor, ν = 0, we need to control the carrier density ρ which according to Eqs.(11),(25) in the LLL approximation is given by the expression ρ = 12 πl X s sgn ( µ L ( s )) θ ( | µ L ( s ) | − E ( s )) , (B9)where µ L ( s ) = − ( s ˜∆ + ∆ ) and E ( s ) = q ( ˜∆ x + s ∆ x ) + ( ˜∆ z + s ∆ z ) . It is easy to check that for all obtainedsolutions the carrier density ρ vanishes so that the filling ν = 0 is realized. The ground state is determined by thesolution with the lowest free energy density. The case of purely magnetic field corresponds to the line B = 0 in Figs.1and 2 in Sec.V. One can see that there is a phase transition from the CAF state to the F state as parallel magneticfield B || increases that agrees with the results obtained by Kharitonov in Refs.[23,24].
2. Purely pseudomagnetic field
In this subsection, we consider solutions of the gap equation (26) in the case where only pseudomagnetic field ispresent ( B ⊥ = B || = 0) retaining only the LLL contribution. Without loss of generality we can put µ ν = ˜ µ ν = 0 sothat m ν = 0, and we are left only with parameters d ν . Moreover, due to the absence of the Zeeman term the gapequation possesses the SU(2) spin symmetry, so that we can take d x = 0. Therefore, we get the following equationfor the parameters d and d z : − d z ( s ) σ z − d ( s ) = −
12 ( u + u z ) (cid:20) sgn( d z ( s )) θ ( E ( s ) − | d ( s ) | ) σ z + sgn( d ( s )) θ ( | d ( s ) | − E ( s )) (cid:21) − u (cid:20) sgn( d z ( − s )) θ ( E ( − s ) − | d ( − s ) | ) σ z + sgn( d ( − s )) θ ( | d ( − s ) | − E ( − s )) (cid:21) + u z X s ′ sgn( d ( s ′ )) θ ( | d ( s ′ ) | − E ( s ′ )) , (B10)where d ,z = ˜∆ ,z − s ∆ ,z and E = | d z | .Assuming that u >> u z , u , we find the following solutions:(i) F solution: ˜∆ = ∆ = 0 , ˜∆ z = 0 , ∆ z = ±
12 ( u + u z − u ); (B11)(ii) AF solution: ˜∆ = ∆ = 0 , ˜∆ z = ±
12 ( u + u z + 4 u ) , ∆ z = 0; (B12)(iii) QAH solution: ˜∆ z = ∆ z = 0 , ˜∆ = 0 , ∆ = ±
12 ( u + u z − u ) . (B13)Note that the CAF solution is no longer present due to the absence of the Zeeman term. For the carrier density weobtain the following expression: ρ = 12 πl X s sgn (cid:16) ˜∆ − s ∆ (cid:17) θ (cid:16) | ˜∆ − s ∆ | − | ˜∆ z − s ∆ z | (cid:17) . (B14)One can check that the condition ρ = 0 is satisfied for all obtained solutions. The ground state of the system isdetermined as the solution with the lowest free energy density (see Sec.V). By making use of the energy density(43) for B ⊥ = ǫ Z = 0, we find that the F and QAH solutions for g ⊥⊥ < g ⊥⊥ > g ⊥⊥ in the selectionof the ground state of the system in a strong pseudomagnetic field.6 Appendix C: Free-energy density
The energy density is defined by the following expression Ω = − Γ /T V , where T V is a space-time volume. Byutilizing the Baym-Kadanoff-Jackiw-Tomboulis) formalism, we find the effective action Γ at its extrema in the mean-field approximation (for details see Ref.[27]):Γ = − i Tr (cid:20) ln G − + 12 ( S − G − (cid:21) , (C1)where the trace, logarithm, and product S − G are taken in the functional sense, G = diag( G + , G − ), 1 is the unitoperator in both matrix and coordinate sense and the expressions for the free and full propagator are given by Eqs.(6)and (7). Performing the Fourier transform in time, integrating by parts the logarithm term, and omitting the irrelevantsurface term, we arrive at the expression (for simplicity in this section we put constants ~ = c = 1):Γ = − iT Z ∞−∞ dw π Tr (cid:20) − w ∂G − ( w ) ∂w G ( w ) + 12 ( S − ( w ) G ( w ) − (cid:21) (C2)with ∂G − ( w ) ∂w = − iγ δ ( r − r ′ ) . The multiplier T came from the functional trace in time, and now Tr contains only spatial integration. By substitutingthe expression for the Green function (A15) into Ω = − Γ /T V , one can see that Schwinger phase Φ goes away andafter the Fourier transformation we get for the energy densityΩ = i Z ∞−∞ dw π Z d k (2 π ) tr (cid:20) iγ w ¯ G ( w, k ) + 12 [ − i { ( w − ǫ Z σ Z ) γ − v F ( k · γ ) } ¯ G ( w ; k ) − (cid:21) = − Z ∞−∞ dw π Z d k (2 π ) tr { [( w + ǫ Z σ Z ) γ + v F ( k · γ )] ¯ G ( w ; k ) + i } . (C3)By making use of the explicit form of the propagator, we calculate the following integrals which contribute to theenergy density: Z d k (2 π ) γ ¯ G s ( w ; k ) = i | eB | π ∞ X n =0 γ W P − + P + θ ( n − / M − nǫ B , (C4) Z d k (2 π ) v F ( k · γ ) ¯ G s ( w ; k ) = i | eB | π ∞ X n =0 nǫ B M − nǫ B , (C5)where in strained case we use the effective magnetic field B = B ⊥ − s s B in the Landau energy ǫ B and projectors P ± .By dropping an infinite divergent term independent of the physical parameters and normalizing Ω by subtractingits value at m µ = d µ = ǫ Z = 0, we obtain the following expression:Ω = − i | eB | π Z ∞−∞ dw π tr ∞ X n =0 (cid:20) ( w + ǫ Z σ Z ) γ W ( P − + P + θ ( n − / nǫ B M − nǫ B − w ( P − + P + θ ( n − / nǫ B w − nǫ B (cid:21) . (C6)After integrating this expression over frequency and taking trace, we finally arrive at Eq.(43). X. Li, X. Wang, L. Zhang, S. Lee, and H. Dai, Science , 1229 (2008). T. Low and F. Guinea, Nano Letters , 3551 (2010). J. Lu, A.C. Neto, and K.P. Loh, Nature Communications , 823 (2012). H. Suzuura and T. Ando, Phys. Rev. B , 235412 (2002). J.L. Manes, Phys. Rev. B , 045430 (2007). S.V. Morozov, K.S. Novoselov, M.I. Katsnelson, F. Schedin, D. Jiang, and A.K. Geim, Phys. Rev. Lett. , 016801 (2006). T.O. Wehling, A.V. Balatsky, A.M. Tsvelik, M.I. Katsnelson, and A.I. Lichtenstein, Europhys. Lett. , 17003 (2008). M.A. Vozmediano, M. Katsnelson, and F. Guinea, Physics Reports , 109 (2010). M.I. Katsnelson,
Graphene: Carbon in two dimensions (Cambridge University Press, Cambridge, 2012). F. Guinea, M.I. Katsnelson, and A.K. Geim, Nature Physics , 30 (2010). F. Guinea, A.K. Geim, M.I. Katsnelson, and K.S. Novoselov, Phys. Rev. B , 035408 (2010). K.S. Novoselov, A.K. Geim, S.V. Morozov, D. Jaing, M.I. Katsnelson, I.V. Grigorieva, S.V. Dubonos, and A.A. Firsov,Nature , 197 (2005); Y. Zhang, Y.-W. Tan, H.L. St¨ormer, and P. Kim, Nature , 201 (2005). V.P. Gusynin and S.G. Sharapov, Phys. Rev. Lett. , 146801 (2005). N.M.R. Peres, F. Guinea, and A.H. Castro Neto, Phys. Rev. B , 1254111 (2006). Y. Zhang, Z. Jiang, J.P. Small, M.S. Purewal, Y.-W. Tan, M. Fazlollahi, J.D. Chudow, J.A. Jaszczak, H.L. Stormer, andP. Kim, Phys. Rev. Lett. , 136806 (2006). Z. Jiang, Y. Zhang, H.L. Stormer, and P. Kim, Phys. Rev. Lett. , 106802 (2007). N. Levy, S.A. Burke, K.L. Meaker, M. Panlasigui, A.Zettl, F. Guinea, A.H. Castro Neto, and M.F. Crommie, Science ,544 (2010). E. Prada, P. San-Jose, G. Leon, M.M. Fogler, and F. Guinea, Phys. Rev. B. , 161402(R) (2010). B. Roy, Z.-X. Hu, and K. Yang, Phys. Rev. B , 121408(R) (2013). I.F. Herbut, Phys. Rev. B , 205433 (2008). B. Roy and J.D. Sau, Phys. Rev. B , 075427 (2014). K.K. Gomes, W. Mar, W. Ko, F. Guinea, and H. Manoharan, Nature , 306 (2012). M. Kharitonov, Phys. Rev. B , 155439 (2012). M. Kharitonov, Phys. Rev. B , 075450 (2012). R. Jackiw and S.-Y. Pi, Phys. Rev. Lett. , 266402 (2007). B. Roy, Phys. Rev. B , 035458 (2012). E.V. Gorbar, V.P. Gusynin, V.A. Miransky, and I.A. Shovkovy, Phys. Rev. B , 085437 (2008). J. Alicea and M.P.A. Fisher, Phys. Rev. B , 075422 (2006). I. L. Aleiner, D. E. Kharzeev, and A. M. Tsvelik, Phys. Rev. B6 , 195415 (2007). D.V. Khveshchenko, Phys. Rev. Lett. , 206401 (2001); , 246802 (2001). E.V. Gorbar, V.P. Gusynin, V.A. Miransky, and I.A. Shovkovy, Phys. Rev. B , 045108 (2002). K. Nomura and A.H. MacDonald, Phys. Rev. Lett. , 256602 (2006); K. Yang, S. Das Sarma, and A.H. MacDonald, Phys.Rev. B , 075423 (2006). M.O. Goerbig, R. Moessner, and B. Doucot, Phys. Rev. B , 161407(R) (2006). L. Sheng, D.N. Sheng, F.D.M. Haldane, and L. Balents, Phys. Rev. Lett. , 196802 (2007). V.P. Gusynin, V.A. Miransky, S.G. Sharapov, and I.A. Shovkovy, Phys. Rev. B , 195429 (2006). I.F. Herbut, Phys. Rev. B , 165411 (2007). J.-N. Fuchs and P. Lederer, Phys. Rev. Lett. , 016803 (2007). M. Ezawa, J. Phys. Soc. Japan , 094701 (2007). K. Yang, Solid State Commun. , 27 (2007). E.V. Gorbar, V.P. Gusynin, and V.A. Miransky, Low Temp. Phys. , 790 (2008). M.O. Goerbig, Rev. Mod. Phys. , 1193 (2011). E.V. Gorbar, V.P. Gusynin, V.A. Miransky, and I.A. Shovkovy, Phys. Scr.
T 146 , 014018 (2012). G. Baym and L.P. Kadanoff, Phys. Rev. , 287 (1961); J.M. Cornwall, R. Jackiw, and E. Tomboulis, Phys. Rev. D ,2428 (1974). B. Roy, F.F. Assaad, and I.F. Herbut, Phys. Rev. X , 021042 (2014). J. Schwinger, Phys. Rev. , 664 (1951). J. Zak, Phys. Rev.134