Diamagnetic Response of Disordered Graphene to Nonuniform Magnetic Field

The diamagnetic susceptibility is calculated in disordered monolayer graphene in spatially varying magnetic field within a self-consistent Born approximation. The universality that the susceptibility is a function of q/2kF in ideal graphene is nearly valid in the weak-disorder case except that the susceptibility becomes nonzero even for slowly-varying field. In the case of strong disorder the universal behavior is completely destroyed. The mirror-image-like response to a magnet is still valid at the Dirac point as long as the disorder is sufficiently weak.


§1. Introduction
The graphene is an atomically thin monolayer of graphite, recently fabricated, and has been attracting attentions theoretically and experimentally.−12) The graphene is known to exhibit singular diamagnetic susceptibility given by a delta function at zero energy often called the Dirac point. 7,13)In the presence of disorder, the delta function is broadened into a function quite different from a simple Lorentzian. 14)n nonuniform magnetic fields, ideal graphene without disorder does not respond to magnetic fields slowly varying within the scale of the Fermi wavelength. 15)The purpose of this paper is to study effects of disorder on this peculiar response to spatially varying magnetic field.
−25) In this paper we use the self-consistent Born approximation assuming scatterers with potential range smaller than the typical electron wavelength (but larger than the lattice constant) in contrast to recent extensions to more general scatterers. 26,27)It should be mentioned that the orbital magnetism in uniform field was also studied for related materials, such as graphite intercalation compounds, 28−30) carbon nanotube, 11,31−33) graphene ribbons, 34) few-layer graphenes, 3,19,35) bulk graphite, 36) bismuth, 37−39) and organic compounds having Dirac-like spectrum. 40)he paper is organized as follows: In §2, following a brief review on the electronic states, the diamagnetic susceptibility in ideal graphene, and a self-consistent Born approximation, we shall obtain an explicit expression for the susceptibility in the presence of disorder.In §3 some examples of numerical results are presented.In §4 response to a small magnet is discussed using obtained susceptibility.A short summary is given in §5.§2.Formulation

Preliminary
In a graphene sheet the conduction and valence bands consisting of π orbitals cross at K and K' points of the Brillouin zone, where the Fermi level is located. 41,42)lectronic states of the π-bands near a K point are described by the k•p equation: ) where F is a two-component wave function, ε is energy, σ = (σ x , σ y ) is the Pauli spin matrix, k = ( kx , ky ) = −i∇ is a wave-vector operator, and γ is a band parameter, given by γ = ( √ 3/2)aγ 0 with a the lattice constant (2.46 Å) and γ 0 the nearest-neighbor hopping integral (∼ 3 eV).
The energy eigen-value, density of states, and electron concentration are given by ) respectively, where s = ±1, g v = 2 is the valley degeneracy associated with K and K' points, g s = 2 is the spin degeneracy, ε F is the Fermi energy, and sgn(ε) denotes the sign of ε.
In the presence of external magnetic field B(r) = B(q)e iq•r +c.c., induced current density j(r)= j(q)e iq•r + c.c. is given by with Due to the gauge invariance and the current conservation, the response function should be written as (2.9) According to Maxwell's equation, we have the relations between magnetic moment per unit area, m(r), perpendicular to the system and the induced current j(r), j x = c(∂m/∂y) and j y = −c(∂m/∂x), giving orbital susceptibility (2.10) In ideal graphene without disorder at zero temperature, the susceptibility has been calculated as 15) with k F being the Fermi wave vector and The response vanishes in range q < 2k F .At the Dirac point ε F = 0, particularly, we have The susceptibility of the carbon nanotube in a field perpendicular to the axis is obtained by replacing q with 2π/L with tube circumference L. 31,33) As a function of ε F at fixed q, it is nonzero only in a finite region satisfying |ε F | < γq/2 and its integral over ε F becomes constant, and in the limit of q → 0 it goes to Let us consider the system described by Hamiltonian where we have assumed scatterers with strength u and potential range smaller than the electron wavelength but larger than the lattice constant (in refs.16 and 19 such scatterers were called long-range because they do not cause intervalley scattering between the K and K' points).In a self-consistent Born approximation, shown in Fig. 1(a), the (2,2) matrix Green's function is written as 16) Ĝ with diagonal self-energy where ε k = γk (k = |k|) and we have introduced and cutoff function g(ε) which satisfies g(ε) ≈ 1 for ε < ε c and decays rapidly for ε > ε c .Further, we have introduced dimensionless parameter describing scattering strength where n i is the concentration of scatterers per unit area.
As will be clear in the following, quantities related to energy are all scaled by ε c and those related to length are scaled by k −1 with ε c = γk c .Without cutoff, the imaginary part of the self-energy converges but the real part diverges logarithmically because of the linear density of states.Due to the weak logarithmic divergence, actual value of the self-energy does not depend on details of the density of states in the high energy region, but is essentially determined by the total number of states available.Therefore, we should choose ε c ∼ 3γ 0 , corresponding to the half of the π-band width, or k c ∼ π/a.Actual value of ε c may have to be slightly modified for detailed comparison with experiments.
Because of the weak logarithmic divergence, we have This self-consistency equation has been solved as 19) where W L (z) is Lambert's W function, defined by z = W L exp(W L ), and (2.23) In the following we shall confine ourselves to such systems.

Diamagnetic susceptibility
In terms of Greens' function, the response function is written as where • • • means average over configurations of scatterers and Tr stands for the summation over all states.This is rewritten as This becomes where k ± = k ± (q/2), Ĵν (q, ε) is the current vertex part, and Tr now stands for the summation over two pseudo-spin states.Figure 1(b) shows diagramatic rep-resentation of I μν (q, ε).
In the self-consistent Born approximation, we have (2.28) as shown in Fig. 1(c).For convenience, we shall define Ĵz (q, ε) by the same expression for spin matrix σ z .In order to solve this, we consider (2.29) Explicitly, we have where θ q is the direction of q and , (2.33) with θ being the direction of k.Note that A, B, and C are functions of q and ε and independent of θ q .Therefore, we can set Then, the Bethe-Salpeter type equation (2.28) is solved as giving Because the integral giving C logarithmically diverges while those giving A and B do not, |Ξ zz | becomes much larger than other elements of Ξ.Let us rewrite matrix Ξ as Ξ(q, ε) = Ξ 0 (q, ε) + Ξ 1 (q, ε), (2.38)   with Then, we make expansion with respect to Ξ 1 as follows: (2.41) Terms containing Ξ 1 can be neglected because of divergent Ξ zz and the final expression becomes with μ, ν = x, y.This shows that only the x and y elements of Ξ(q, ε) are required and other elements are all irrelevant.
Let us define (2.44) Then, we have Therefore, we immediately have At first glance, it does not seem that the above result satisfies the gauge invariance and the current conservation condition because Λ = 0.However, we should note (2.47)Because the integral converges, the result does not change if we introduce cutoff function g(ε k ) in the integrand in the left hand side.Then, we have for energies and wave vectors satisfying |ε| ε c and q k c , where we have used the fact that the integral is dominated by wave vector close to cutoff k c .Then, we have (2.49) Because Λ(q, ε) is real, the first term proportional to δ μν of eq.(2.46) does not contribute to the response function and the resulting K μν (q) satisfies eq.(2.9).With the use of W 1, the final expression for I μν (q, ε) becomes (2.50) It should be noted that the first term in the expression of I μν (q, ε) is present even in ideal graphene without disorder.
The susceptibility is given by (2.51) In the long-wavelength limit, q → 0, we have (2.52) Therefore, we have (2.53) in agreement with the expression obtained previously. 19)Numerical Results Figure 2 shows some examples of the density of states.In the presence of disorder, it is considerably enhanced in comparison with the density of states in ideal graphene and effects become larger with increasing W .The same is applicable to the electron concentration for a fixed value of ε F .In the following, therefore, we shall define effective Fermi wave vector k F and Fermi energy ε F for given electron concentration n s by Figure 2 also shows some examples of ε F for several values of k F .For γ/h ≈ c/300, we have roughly n s = 6.2 × 10 11 cm −2 for k F /k c = 0.01, 1.4 × 10 12 cm −2 for 0.015, and 5.6×10 12 cm −2 for 0.03. Figure 3 shows calculated diamagnetic susceptibility versus wave vector q for k F /k c = 0.01, 0.015, and 0.03 and for W = 0, 0.02, 0.05, and 0.1.The susceptibility is measured in units of −χ c with χ c = g s g v e 2 γ/(h 2 c 2 k c ).In ideal graphene W = 0, the susceptibility vanishes up to q = 2k F , suddenly increases, takes a maximum at q slightly larger than 2k F , and then decreases roughly in proportion to q −1 .In the presence of disorder, the susceptibility becomes nonzero even for q < 2k F and the maximum value decreases.In the case of strong disorder W = 0.1, the susceptibility in the region q <2k F becomes almost the same as the maximum, showing that the special feature that graphene does not respond to slowly-varying magnetic field is completely destroyed.
Figure 4 shows a re-plot of χ(q) shown in the previous figure as a function of q/k F .The unit of χ(q) is χ F = g s g v e 2 γ/(h 2 c 2 k F ).In this figure, we can clearly see that the universal behavior of χ(q), i.e., a function of q/2k F , is gradually destroyed with the increase of disorder.In the case of weak disorder W = 0.02 and 0.05, the deviation is apparent only for q > 2k F .In the case of strong disorder W = 0.1, on the other hand, the universal feature is completely destroyed.The universal behavior and the nonvanishing χ(0) for weak disorder manifest themselves in the magnetic field induced by a small magnet, as will be shown in the next section.
Figure 5 shows susceptibility as a function of k F /q for q/k c = 0.02 and 0.05 and for W = 0, 0.02, 0.05, and 0.1.The susceptibility is measured in units of −χ q with χ q = g s g v e 2 γ/(h 2 c 2 q).In the ideal graphene, χ(q) is nonzero only for k F < q/2 and totally vanishes for k F > q/2 as has been discussed above.With the increase of disorder, χ(q) is reduced for small k F and becomes nonzero exhibiting long tails for large k F .This disorder effect is larger for small values of q as is expected.§4.Discussion With the use of Maxwell's equation, we can easily calculate magnetic field associated with magnetization m(q)= χ(q)B ext (q) of the graphene due to external field B ext (q).We choose the graphene as the xy plane and introduce coordinate z in the perpendicular direction.We then have (4.1)This shows that the induced field vanishes for uniform magnetic field q → 0. It is consistent with the fact that no magnetic field comes out of a thin flat-plane magnet.
Let us consider, for example, a situation where a magnetic charge (mono-pole) ρ m is located at distance d (d > 0) and r = 0 above the graphene.This is nearly realized when we put a very long and thin magnet vertically on top of the graphene.Then, the field component perpendicular to the graphene becomes ( The induced field in perpendicular direction at the magnetic charge is calculated as In ideal graphene with the Fermi level at the Dirac point for which χ(q) is given by eq.(2.13), in particular, we have with This corresponds to the field of a mirror-image magnet discussed in ref. 15.With the typical value γ/h ≈ c/300 with c being the light velocity, α is estimated as ≈ 4 × 10 −5 , showing that the counter field is much smaller than the original.This mirroring holds whenever χ(q) is proportional to 1/q, but this field is expected to be reduced due to modified χ(q) in the presence of disorder.Figure 6 shows induced magnetic field at the point of the magnetic charge as a function of distance d for different values of k F and for W = 0 and 0.05.At the Dirac point, the field is reduced by disorder, but the reduction is not so considerable for sufficiently small W .When the Fermi level is away from the Dirac point, the field is reduced in comparison with that in ideal graphene, but remains larger for larger distance in contrast to the ideal case.This result is to be expected because the susceptibility in the long-wavelength region remains nonzero in the presence of disorder and may be used for determination of W , if being observed experimentally.Both in ideal and disordered graphenes (as long as disorder is sufficiently small), d 2 × B ind is a universal function of dk F corresponding to the universal behavior χ(q, 2k F ) = χ(q/2k F ), as discussed in the previous sections.
The diamagnetic susceptibility of graphene is much larger than the Pauli paramagnetism determined by the electron mass in vacuum.−46) §5.Summary In summary, we have calculated the diamagnetic susceptibility of disordered monolayer graphene in spatially varying magnetic field in the self-consistent Born approximation.The universality that the susceptibility is a function of q/2k F in ideal graphene is nearly valid in the weak disorder case except that the susceptibility becomes nonzero for slowly varying field.In the case of strong disorder the universal behavior is completely destroyed.The mirror-image-like response is still valid at the Dirac point as long as the disorder is sufficiently weak.
Fig. 5 (Color online) Calculated diamagnetic susceptibility versus wave vector k F /q for q/k c = 0.02 and 0.05 and for W = 0, 0.02, 0.05, and 0.1.The susceptibility is measured in units of −χ q with χ q = g s g v e 2 γ/(h 2 c 2 q).In ideal graphene, the result is universal in this plot.
Fig. 1 Diagramatic representation of (a) self-energy, (b) response function K μν (q), and (c) current vertex Ĵν (q, ε) in the self-consistent Born approximation.Fermi Wave Vector (units of q) Susceptibility (units of −χ q ) g s g v e 2 γ h 2 c 2 q _ χ q = Fig. 5 (Color online) Calculated diamagnetic susceptibility versus wave vector k F /q for q/k c = 0.02 and 0.05 and for W = 0, 0.02, 0.05, and 0.1.The susceptibility is measured in units of −χ q with χ q = g s g v e 2 γ/(h 2 c 2 q).In ideal graphene, the result is universal in this plot.

Fig. 2 (
Fig. 2 (Color online) Calculated density of states and Fermi energy ε F for several values of Fermi wave vector k F corresponding to fixed electron concentrations.W = 0 (ideal), 0.02, 0.05, and 0.1.

Fig. 4 (
Fig.4(Color online) Calculated diamagnetic susceptibility shown in Fig.4is re-plotted as a function of q/k F .The unit of χ(q) is χ F = g s g v e 2 γ/(h 2 c 2 k F ).

Fig. 6 (
Fig. 6 (Color online) Calculated induced magnetic field for a magnetic charge ρ m as a function of distance d from graphene for different values of k F and for W = 0 and 0.05.We use γ/h = c/300 and γk c = 9 eV to show the actual distance. n

2 )Fig. 2 (of k c ) g s g v e 2 γ h 2 c 2 Fig. 3 (
Fig. 2 (Color online) Calculated density of states and Fermi energy ε F for several values of Fermi wave vector k F corresponding to fixed electron concentrations.W = 0 (ideal), 0.02, 0.05, and 0.1.

Fig. 4 (
Fig.4(Color online) Calculated diamagnetic susceptibility shown in Fig.4is re-plotted as a function of q/k F .The unit of χ(q) is χ F = g s g v e 2 γ/(h 2 c 2 k F ).

kFig. 6 (
Fig. 6 (Color online) Calculated induced magnetic field for a magnetic charge ρ m as a function of distance d from graphene for different values of k F and for W = 0 and 0.05.We use γ/h = c/300 and γk c = 9 eV to show the actual distance.