Dynamical polarization and plasmons in a two-dimensional system with merging Dirac points

1.9 MB·

by  P. K. Pyatkovskiy

Dynamical polarization and plasmons in a two-dimensional system with merging
Dirac points
P. K. Pyatkovskiy and Tapash Chakraborty
Department of Physics and Astronomy, University of Manitoba, Winnipeg, Canada R3T 2N2
(Dated: March 1, 2016)
We have studied the dynamical polarization and collective excitations in an anisotropic two-
dimensionalsystemundergoingaquantumphasetransitionwithmergingoftwoDiracpoints. Ana-
lyticalresultsfortheone-looppolarizationfunctionareobtainedatthefinitemomentum,frequency,
andchemicalpotential. Theevolutionoftheplasmondispersionacrossthephasetransitionisthen
6 analyzed within the random phase approximation. We derive analytically the long-wavelength dis-
1 persionoftheundampedanisotropiccollectivemodeandfindthatitevolvessmoothlyatthecritical
0 merging point. The effects of the van Hove singularity on the plasmon excitations are explored in
2 detail.

b
e
F I. INTRODUCTION sition is governed by a single parameter ∆ that changes
itssignacrossthecriticalpoint(Fig.1). Thismodelpro-
9
vides a universal description [17] for a two-dimensional
2 For more than a decade, we have been witnessing the
system in the vicinity of the phase transition with two

] riseofaplethoraofevernewtwo-dimensional(2D)mate- merging Dirac points related by time-reversal symmetry.
l rials displaying their unique electronic properties, which Previous results related to our present study include the
l
a has initiated major activities in those systems. Leading long-wavelength plasmon dispersion at the critical point
h the pack was, of course, monolayer and bilayer graphene
(∆=0) obtained in Ref. [35], spectrum of collective ex-
- displayingthebehaviorof“Diracfermions”ofthecharge
s citations in a single- and few-layer BP [36, 37] (only the
e carriers with “Dirac points” where the two energy bands
conduction band or the valence band was taken into ac-
m meet [1, 2] with linear dispersions in the vicinity that
count due to the large value of the gap), and the spec-
. forms the characteristic “Dirac cones” [3]. Their many trum of plasmons across the phase transition obtained
t
a exoticphysicalproperties,inparticular,inastrongmag- numerically within a tight-binding model for bilayer BP
m netic field, have been well documented, and range from
[38].
- the magnetic field effects in the extreme quantum limit Wecalculatetheone-loopdynamicalpolarizationfunc-
d [4] to Hofstadter butterflies [5, 6]. These were then fol-
tionatzerotemperatureforarbitraryvaluesoftheFermi
n
lowed by other graphenelike systems, such as silicene
o energy and the gap. In general, we are able to perform
andgermanene[7–9],the2Dversionofblackphosphorus
c one momentum integration and derive an expression in
[ (BP)[10,11],and,mostrecently,theplanarelectronsys- terms of a single integral valid for arbitrary complex
tems in ZnO heterojunctions [12, 13]. Interestingly, an frequencies. This expression is used to study numeri-
2
anisotropic two-dimensional system can undergo a tran-
v cally the evolution of the plasmon dispersion across the
sitionbetweenaninsulatingstatewithgappedspectrum
3 phasetransitionwithintherandomphaseapproximation
and a semimetal state with two Dirac cones separated in
1 (RPA). The imaginary part of the vacuum polarization
0 the momentum space. The possibility of such a quan- function and the long-wavelength spectrum of collective
1 tum phase transition has been considered theoretically
excitations are evaluated analytically.
0 in honeycomb lattice models [14–18] and few-layer black
1. phosphorus where the band inversion can be induced by
0 anexternalperpendicularelectricfield[19–21]orbydop- II. POLARIZATION FUNCTION
6 ing [22]. The gapless spectrum at the phase transition
1 pointmayariseintheTiO /VO nanostructures[23,24].
2 2 We use a universal low-energy two-band Hamiltonian
:
v Experimentally, merging or creation of Dirac points has [17] describing the merging transition,
i been observed in systems of ultracold atoms [25], pho-
X
tonic crystals [26], microwave analog of graphene [27], H =(∆+ak2)σ +vk σ , (1)
x x y y
r and, more recently, in a potassium-doped few-layer BP
a where the Pauli matrices σ , σ act on the two-
[28]. x y
component wave functions. The spin-orbit coupling is
Various properties of a system undergoing this phase neglectedandthepresenceoftwospinstatesisaccounted
transition have been reported in the literature, which for by the degeneracy factor g = 2. The energy eigen-
s
include the Landau levels and the Hofstadter spectrum values are given by
[29], the Hall conductivity [30], effects of disorder [31],
(cid:113)
thequantumcriticalbehavior[32,33], andthetransport Eλ =λ (ak2+∆)2+v2k2, λ=±. (2)
k x y
characteristics [34]. In this paper, we consider the dy-
namical polarization and collective excitations utilizing The Hamiltonian (1) can also be used to describe the
themodelintroducedinRef.[17]inwhichthephasetran- single-layer BP (∆≈0.8 eV) [39] when the difference in
2
ture T is given by
E
HaL HbL HcL ky Π(iω ,q)=g T (cid:88)∞ (cid:90) d2k Tr(cid:2)G(iΩ ,k)
m s (2π)2 n
n=−∞
kx (cid:3)
×G(iΩ +iω ,k+q) , (5)
n m
D>0 D=0 D<0
where ω =2πmT, Ω =(2n+1)πT aretheMatsubara
m n
FIG. 1. Evolution of the electron energy spectrum at the frequencies, and the Green’s function is
phase transition.
iΩ +µ+(ak2+∆)σ +vk σ
G(iΩ ,k)= n x x y y. (6)
n (iΩ +µ)2−(ak2+∆)2−v2k2
n x y
the effective masses of the positive- and negative-energy
Evaluating the trace and the sum over n yields
bandsisneglected. Inthecasewhenthechemicalpoten-
tlaiarlgeµ(0lie<s µwi−th∆in(cid:28)the∆c)o,nwdeuccatannnceegbleacntdthaendcotnhteribguatpioins Π(iω ,q)=g (cid:90) d2k (cid:88)Fλ,λ(cid:48) nF(Ekλ)−nF(Ekλ+(cid:48) q),
m s (2π)2 k,k+q Eλ−Eλ(cid:48) +iω
from the negative-energy band for energies close to the λ,λ(cid:48)=± k k+q m
Fermi level and, at small momenta, approximate Eq. (2) (7)
by where n (x)=[e(x−µ)/T +1]−1 and
F
E+−µ≈ kx2 + ky2 −µ (3) Fkλ,,kλ(cid:48)(cid:48) = 21(cid:20)1+ (∆+akx2)(∆E+kλEakλk(cid:48)(cid:48)x(cid:48)2)+v2kyky(cid:48) (cid:21) (8)
k 2m 2m 0
x y
is the wave-function overlap factor. In the following, we
consider only the case of T =0, when n (x)→θ(µ−x)
F
with m = 1/(2a) and m = ∆/v2 being the effective and the polarization function can be written as the sum
x y
masses in the x and y directions, respectively, and the of two terms,
chemicalpotentialµ =µ−∆measuredfromthebottom
0
of the positive-energy band. At the critical point, ∆=0 Π(iω,q)=Π0(iω,q)+Π1(iω,q), (9)
[Fig.1(b)],thespectrumislinearintheydirection,while
where Π (iω,q) is the “vacuum” polarization at µ = 0
quadraticinthexdirection(withthesameeffectivemass 0
andΠ (iω,q)givesanadditionalcontributionwhenµ>
m ). Such a system is often referred to in the literature 1
x ∆ (we choose µ (cid:62) 0, and the case µ < 0 is equivalent
as the “semi-Dirac” system.
becauseoftheelectron-holesymmetry). Thesetwoterms
In the case of ∆ < 0, the spectrum has two Dirac are given by
cones [Fig. 1(c)] located at k = (±K ,0) with K =
x x
(cid:112)−∆/a. In the vicinity of these points, the linearized Π (iω,q)=−χ−(iω,q),
0 ∞
Hamiltonian (1) reads Π (iω,q)=χ+(iω,q)+χ−(iω,q), (10)
1 µ µ
where
H (cid:39)±v (k ∓K )σ +vk σ , (4)
x x x x y y
χ±(ω,q)=g (cid:90) d2k (cid:88) θ(µ−Ek+)Fk+,k,±+q. (11)
√ µ s (2π)2 E+−E± +σω
where the velocities v =2 −a∆ in the x direction and σ=± k k+q
x
vintheydirectionaredifferentingeneral. Thespectrum
The k integration in the above equation can be per-
hassaddlepointsE± =±|∆|withdivergentdensityof y
k=0 formed analytically for an arbitrary complex frequency
states (the van Hove singularity).
ω away from the real axis (see Appendix A). The result-
The one-loop polarization function at finite tempera- ing expressions are
∞ (cid:34)
Π (ω,q)=− gs (cid:90) dk v2q2(v2q2−ω2)+ηξ(v2q2+ω2)Re(cid:2)arctanh(ξ/η)(cid:3)
0 2π2v(v2q2−ω2)2 x y y y
y
−∞ (12)
(cid:18) √α−ω2 (cid:112)β−ω2(cid:19) (cid:112)min(α,β)−ω2(cid:35)
− ξ2ω2 +v2q2η2√ arctanh ,
(cid:112)β−ω2 y α−ω2 (cid:112)max(α,β)−ω2
3
√
(µ−∆)/a
g θ(µ−∆) (cid:90) (cid:20) µ
Π (ω,q)= s dk θ(µ+∆+ak2) ηξ(v2q2+ω2)arctanh (cid:101)
1 4π2v(v2q2−ω2)2 x x y µ
y √
− (µ−∆)/a
√
(cid:18) α−ω2 (cid:112)β−ω2(cid:19) µ(ηξ+v2q2−ω2)−2vq (ak2+∆)2 (13)
− ξ2ω2 +v2q2η2√ arctanh (cid:101) √y y x
(cid:112)β−ω2 y α−ω2 µ α−ω2(cid:112)β−ω2
ηξ+v2q2−ω2+2vq µ (cid:21)
−2vq ηξωarctanh y y(cid:101) +(q →−q ) ,
y 2µω y y
where thevanHovesingularityappearsingraphene[40,41]. In
contrast to the case of graphene, in our model this sin-
ξ =aq (q +2k ), gularity occurs only at a single point in the momentum
x x x
η =ξ+2(ak2+∆), space between the two Dirac cones. Because of this, the
x
divergent term is proportional to q and vanishes for the
α=η2+v2q2, (14) y
y momentum directed along the x axis.
β =ξ2+v2q2, For ∆ = 0, some limiting cases of vacuum polariza-
y
(cid:112) tion (12) can be evaluated analytically,
µ= µ2−(ak2+∆)2,
(cid:101) x
Γ(5/4) g vq2
and the retarded polarization on the real ω axis is ob- Π(ω,q =0,q )=− √ s y , (18)
x y Γ(3/4)6 πa(v2q2−ω2)3/4
tained using the prescription ω → ω + i0. We use y
Eqs. (12) and (13) in our numerical calculations of the
collectiveexcitationspectrumandalsotoanalyticallyob-
g |q |
tain some important limits. Π(ω =0,q ,q =0)=− s x , (19)
x y 16v
The imaginary part of the vacuum term can be cal-
culated analytically (see Appendix B). If ∆ < 0 and inagreementwiththepreviouslyreportedresults[32,33].
q < 2K , it has a logarithmic singularity at ω = ±ω,
x x (cid:101)
where
(cid:113) III. PLASMONS
ω = v2q2+a2(2K2−q2/2)2. (15)
(cid:101) y x x
Plasmondispersionω (q)intheRPAisobtainedfrom
p
In the vicinity of this singularity, ImΠ0(ω,q) is given by zeros of the dielectric function
(including terms finite at ω =±ω)
(cid:101)
(cid:15)(ω,q)=1−V(q)Π(ω,q), (20)
(cid:34)
∓g v2q2 512a2q6
ImΠ (ω,q ,q )(cid:39) √s y ln (cid:101)
0 x y 64π av a3/2q3 ω|q2−q2||ω∓ω| where Π(ω,q) is the one-loop polarization function and
(cid:101) (cid:101) x (cid:101) (cid:101) V(q) = 2πe2/(κq) is the Coulomb potential screened
− ω(cid:101)2(qx2−3q(cid:101)2)+v2qy2(qx2+q(cid:101)2) (16) onlybythesubstratewiththecorrespondingbackground
a3/2q2q3 dielectric constant κ.
x(cid:101)
(cid:35) We numerically found the real solutions ω = ω (q)
ω2−8a2(q4+2q2q2−q4) |q −q| p
+ (cid:101) x x(cid:101) (cid:101) ln x (cid:101) , of Eq. (20) in the regions where ImΠ(ω,q) = 0,
a3/2qx3 qx+q(cid:101) i.e., the Landau damping is absent. For the SPE
regions where the imaginary part of the polarization
(cid:112)
where q = K2/2−q2/8. For q =0, this simplifies to function is nonzero, we calculate the energy-loss func-
(cid:101) x x x
tion −Im[1/(cid:15)(ω,q)], the peaks of which represent the
∓g vq2 (cid:18) 128∆2 8(cid:19) damped plasmons. Our approach assumes the strictly
ImΠ (ω,0,q )(cid:39) √ s y ln −
0 y π a(−8∆)3/2 ω |ω∓ω | 3 two-dimensional system and does not take into account
(cid:101)0 (cid:101)0 the charge distribution in the perpendicular direction.
(17)
(cid:113) Nevertheless,ourresultswillbevalidifthecharacteristic
with ω = 4∆2+v2q2. This logarithmic divergence
(cid:101)0 y lengthofthisdistributionl (e.g.,theinterlayerdistance
z
for ∆ < 0 is due to the van Hove singularity which re- in the case of bilayer BP) is much smaller then 1/|q| for
sults in the saddle point in the interband single-particle the considered wave vectors q.
excitation (SPE) energy E+−E− [the real frequency In the case of µ = 0, there are no real solutions
k k+q
corresponding to the pole of the integrand in Eq. (7)] as of Eq. (20) and the energy-loss function is shown in
a function of k for a given external wave vector q. A Fig. 2. For ∆ < 0, the logarithmic divergence due
similar divergence of ImΠ (ω,q) due to the presence of to the van Hove singularity manifests itself as a dip in
0
4
-Im@1(cid:144)ΕHΩ,qLD -Im@1(cid:144)ΕHΩ,qLD -Im@1(cid:144)ΕHΩ,qLD
0.0 0.4 0.8 0 1 2 0.0 0.5 1.0 1.5 2+
2.5
D=-L D=-L HaL D(cid:144)Μ=-3 HbL D(cid:144)Μ=-3
4 2.0
2
3 1.5
L Μ
Ω(cid:144) 2 (cid:144)Ω 1.0
1
1 0.5
HaL HbL
0.0
0 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4
D=0 D=0 2.5
4 HcL D(cid:144)Μ=-1 HdL D(cid:144)Μ=-1
2 2.0
3
L 1.5
Ω(cid:144) 1 2 (cid:144)ΩΜ 1.0
1
0.5
HcL HdL
0 0 0.0
8 D=L D=L 0.0 0.2 0.4 0.6 0.8 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4
4 2.5
HeL D(cid:144)Μ=-0.5 HfL D(cid:144)Μ=-0.5
6
3 2.0
L
Ω(cid:144) 4 2 Μ 1.5
2 1 (cid:144)Ω 1.0
HeL HfL 0.5
0 0
0 1 2 3 0 1 2 3 4 0.0
qx a(cid:144)L qyΥ(cid:144)L 2.50.0 0.2 0.4 0.6 0.8 0.0 0.2 0.4 0.6 0.8 1.0
HgL D(cid:144)Μ=0 HhL D(cid:144)Μ=0
2.0
FIG. 2. Energy-loss function −Im[1/(cid:15)(ω,q)] at µ = 0 for Μ 1.5
positive, zero, and negative ∆. Left: qy = 0, κv =√10−3c, (cid:144)Ω 1.0
and Λ is an arbitrary energy scale. Right: q =0, 2κ aµ=
x 0.5
10−3c, and κ2Λ=4×10−3e2/a. The boundaries of the SPE
0.0
regions are marked by dotted lines. The insets in the right 0.0 0.2 0.4 0.6 0.8 0.0 0.2 0.4 0.6 0.8 1.0
panels schematically show the band structure and filling of HiL D(cid:144)Μ=0.5 HjL D(cid:144)Μ=0.5
the bands. 3.0
Μ 2.0
(cid:144)Ω
−Im[1/(cid:15)(ω,q)] for q (cid:54)=0 and ω =ω followed by a peak 1.0
y (cid:101)
at a larger energy [Fig. 2(b)]. Analogous behavior has
0.0
also been reported in graphene [41]. 0.0 0.2 0.4 0.6 0.8 1.0 1.2 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5
2.5
The evolution of the plasmon spectrum across the HkL D(cid:144)Μ=0.9 HlL D(cid:144)Μ=0.9
phasetransitionatnonzerochemicalpotentialisshownin 2.0
Fig. 3. The momentum is chosen to be aligned with one 1.5
Μ
of the principal axes. In the case of q = 0, the dielec- (cid:144)Ω 1.0
y
tric function (20) expressed in terms of the dimension-
0.5
(cid:112)
less momentum q a/µ and energy ω/µ depends only
x 0.0
on a single adjustable parameter κv. Similarly, (cid:15)(ω,q) 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.5 1.0 1.5 2.0
at q = 0 can be represented as a function of the di- qx a(cid:144)Μ qyΥ(cid:144)Μ
x
mensionless momentum q v/µ and energy ω/µ, which
y FIG. 3. Energy-loss function −Im[1/(cid:15)(ω,q)] at q = 0,
depends on a single parameter κ2aµ. In our numeri- √ y
κv = 10−3c (left) and q = 0, 2κ aµ = 10−3c (right) for
cal calculations, we choose the values of the parameters x
different values of ∆/µ. The undamped plasmon mode and
κv =10−3c≈3×105 m/s, aµ=v2/4. The latter choice its long-wavelength asymptote (23) are shown by solid and
corresponds,e.g.,tom =1/(2a)=m andµ≈0.26eV, dashed lines, respectively. The boundaries of the SPE re-
x e
where m is the bare electron mass. gions are marked by dotted lines. The insets in the right
e √
In the regime aq ,vq (cid:28) ω (cid:28) µ, i.e., relevant for panels schematically show the band structure and filling of
x y
the long-wavelength plasmons, the asymptotic behavior the bands.
of the polarization function is
√
g µ
Π(ω,q)= s√ (cid:2)µaq2f (∆/µ)+v2q2f (∆/µ)(cid:3),
4π2v aω2 x x y y
(21)
5
20 When the Fermi level crosses the van Hove singularity
fxHD(cid:144)ΜL at |µ| = −∆ and the Fermi surfaces of the two Dirac
15 cones merge into a single one [the plasmon spectrum at
L
(cid:144)DΜ 10 this transition is shown in Figs. 3(c) and 3(d)], the plas-
H mon frequency in the long wavelengths changes continu-
fi
ously as a function of ∆/µ but has a logarithmic singu-
5 f HD(cid:144)ΜL
y larity of its derivative, as seen in Fig. 4. The functions
0 (22) in the vicinity of this crossing are
-3 -2 -1 0 1
√
D(cid:144)Μ f (−1+ε)=2 2(cid:0)56 +εln|ε|(cid:1)+O(ε),
x √ 15 (26)
f (−1+ε)= 2(cid:0)8 −εln|ε|(cid:1)+O(ε).
FIG. 4. Factors fx(∆/µ) and fy(∆/µ) determining the long- y 3
wavelength plasmon dispersion (23).
For ∆ < 0, there also exists an additional damped plas-
mon mode in the x direction for momenta 0<q <2K
x x
with its maximum at q ∼K , which lies entirely in the
where the functions f (δ) are defined as x x
x,y
intraband SPE region [Figs. 3(a), 3(c), and 3(e)].
√
(cid:90) 1 t2 t−δ At the crossing of the critical point ∆ = 0 [Figs. 3(g)
f (δ)=8 dt√ , and 3(h)] the spectrum of the undamped plasmon
x 1−t2
t changes smoothly and we have
0 √ (22)
(cid:90) 1 1−t2 √
fy(δ)=2 t dt √t−δ , fx(0)=3√πΓ(3/4)/Γ(9/4)≈5.751, (27)
0
f (0)=2 πΓ(5/4)/Γ(7/4)≈3.496,
y
witht =max(δ,−1)andshowninFig.4. Theanalytical
0
expressions for f (δ) in terms of the complete elliptic in agreement [43] with Ref. [35].
x,y
integrals are given in Appendix C. The effect of a small positive ∆/µ is similar to that
UsingEq.(21),weobtainthelong-wavelengthplasmon in the case of gapped graphene [44]: the plasmon mode
dispersion, becomes extended to larger values of momenta due to
the opening of the gap between the interband and intra-
(cid:114)g e2µq(cid:20)√µa band SPE continua [Figs. 3(i) and 3(j)]. In the regime
ω (q)(cid:39) s f (∆/µ)cos2θ
p 2πκ v x ∆→µ[Figs.3(k)and3(l)],electronshaveapproximately
parabolic anisotropic dispersion (3). Using the asymp-
v (cid:21)1/2
+ √ f (∆/µ)sin2θ , (23) totic behavior of (22) for δ →1,
µa y
√
f (1−ε)=2 2πε+O(ε2),
x
wherecosθ =q /q,sinθ =q /q. Ithastheusualsquare- √ (28)
root dependencxe on the momyentum just as expected for fy(1−ε)= 2πε+O(ε2), 0(cid:54)ε(cid:28)1,
a 2D system, with the anisotropy of the spectrum be-
√ we recover from Eq. (23) in this limit,
ing fully determined by the dimensionless factor µa/v
and the two functions f (∆/µ). This undamped plas-
mon mode lies in the gapx,ybetween the interband and the ω (q)(cid:39)(cid:114)gse2µ0q(cid:20)(cid:114)my cos2θ+(cid:114)mx sin2θ(cid:21)1/2,
intraband SPE regions (Fig. 3). p κ m m
x y
For −∆ (cid:29) µ, i.e., for well separated Dirac cones (29)
[Figs. 3(a) and 3(b)], using the asymptotics as previously reported in Ref. [37] for the monolayer BP.
f (δ)=4π|δ|1/2+O(|δ|−1/2),
x
(24)
f (δ)=π|δ|−1/2+O(|δ|−3/2), −δ (cid:29)1, IV. CONCLUSIONS
y
one can obtain the plasmon dispersion from Eq. (23) in We have evaluated the polarization function and the
this limit spectrumofcollectiveexcitationsinthetwo-dimensional
system undergoing a topological phase transition with
(cid:114)g e2µq(cid:20)v v (cid:21)1/2 two merging Dirac points. A single integral represen-
ω (q)(cid:39) s x cos2θ+ sin2θ . (25)
p κ v v tation for Π(ω,q) has been derived which is suitable
x
for calculations on both real and imaginary frequency
ThisresultcorrespondstothelinearizedHamiltonian(4) axes. An analytic expression was obtained for the imag-
and provides a generalization of the long-wavelength inary part of the vacuum polarization and its asymp-
plasmon spectrum in a single-layer graphene [42] to the toticbehaviornearthelogarithmicdivergenceduetothe
case of the different Fermi velocities in the x and y di- van Hove singularity. We analytically found the long-
rection. wavelength plasmon dispersion and numerically studied
6
the spectrum of collective excitations for arbitrary mo-
menta for both zero and nonzero values of the chem- 4 HaL ImP0=0 1 HbL ImP0=0
iwctvtcihiaoeatenanlhteHpaiaoxvoconiertsvsoetefensoanstusticinatnehlzdg.eeourbflomBoattryethhiretegymeuivgnnpaiagndelpuartptaamrhtteapuiednnresgeeddilteatiaacomhnntnerddp.oeednsdTnateuhsmprpdeglpeaiyepces-dtdrmlreouspotsemhsnlnaesfcimulmreenoaoeocddvnftesoitoelhutaxnoet--, Ψ ----864202 IDmP0=0 AC B (cid:144)-21222LHDΩ-Υqy --210 B A C ImP0=D0
zero chemical potential in the semimetal phase. At fi-
0 1 2 3 4 5 6 7 0.0 0.2 0.4 0.6 0.8 1.0
niteµ,thereisoneundampedanisotropiccollectivemode
AA jAA aq2HΩ2-Υ2q2L-1(cid:144)2
withthesquare-rootdispersion,whichliesinthegapbe- AA AA x y
tween the interband and intraband SPE regions. In the
gapped phase (∆ > 0), this undamped mode is gener- FIG. 5. Regions (B4) plotted for parameters (a) φ and ψ
ically extended to larger values of momenta due to the or (b) q, ω, and ∆. At the boundaries, ImΠ0(ω,q) has a
enhanced separation between the two SPE continua. At logarithmic singularity (solid line), logarithmically divergent
the critical point (∆ = 0), the undamped plasmon dis- derivative(dashedline)orjumpdiscontinuityofthederivative
(dotted line).
persion changes smoothly, while an additional damped
and strongly anisotropic mode emerges at ∆ < 0 in the
interbandSPEcontinuum. ThecrossingofvanHovesin-
The function F (ω,q,k) does not have any singulari-
gularity by the Fermi level manifests itself in a divergent ±
ties for Imω (cid:54)=0 and the multivalued functions taken on
derivative of the long-wavelength plasmon frequency.
their principal branches. Therefore, the definite k inte-
y
gral of (A1) is obtained straightforwardly by evaluating
F (ω,q,k) in the integration limits, which, after some
±
ACKNOWLEDGMENTS algebra, yields Eqs. (12) and (13).
TheworkhasbeensupportedbytheCanadaResearch
Chairs Program of the Government of Canada. Appendix B: Analytic expression for ImΠ (ω,q)
0
The nonzero imaginary part of the expressions (12)
Appendix A: Calculation of the polarization function and(13)forthepolarizationfunctionoriginatesfromthe
regionswheretheargumentofinversehyperbolictangent
is real and larger than unity,
The integrand in Eq. (11) can be written as
(cid:88) Fk+,k,±+q = ∂ F (ω,q,k), (A1) Im(cid:2)arctanh(x±i0)(cid:3)=±π2θ(x2−1). (B1)
σ=±Ek+−Ek±+q+σω ∂ky ± Thestepfunctionabovedeterminestheintegrationlimits
which, for the nonvacuum term (13), involves the roots
where of the higher-order polynomials that cannot be written
in a closed form in the general case (for nonzero q , q ,
(cid:20) x y
1
F (ω,q,k)= −v2q γ(v2q2−ω2) and ∆). For the vacuum term, on the other hand, the
± 2v(v2q2−ω2)2 y y imaginarypartcanbeevaluatedintermsofthecomplete
y
2k −q elliptic integrals [45]:
+ηξ(v2q2+ω2)arctanh y y
y γ (cid:112)
g |ω|θ(1−τ)θ(2−ψ)θ(ν)
−2vq ηξωarctanhηξ−ω2+2v2qyky ImΠ0(ω,qx,qy)=− s 128π√av(1−τ)3/4
y vγω (cid:26)
(cid:18) √α−ω2 (cid:112)β−ω2(cid:19) ×sgn(ω) 4(2ψ+3ϕ−τϕ+2τψ)hiE(ti)
− ξ2ω2 +v2q2η2√ (B2)
(cid:112)β−ω2 y α−ω2 +(1−b)h−1(cid:104)r K(t )+s (cid:0)(16−4ψϕ)(1−τ)
i i i i
(2k −q )(ηξ+v2q2−ω2)−q (η−ξ)2
×arctanh y y γ√α−ω2y(cid:112)β−ω2 y −ϕ2(3−τ)(cid:1)Π(π/2,ρ ,t )(cid:105)(cid:27),
i i
(cid:21)
(cid:0) (cid:1)
∓ qy →−qy, ξ →−ξ , (A2) where
d=4∆(ω2−v2q2)−1/2, τ =v2q2/ω2,
variables η, ξ, α, β are defined in Eq. (14) and y y
ϕ=(ω2−v2q2)1/2/(aq2), ψ =d+1/ϕ, (B3)
y x
(cid:113)
γ = (η−ξ)2/v2+(2k −q )2. (A3) ν =(ψ+ϕ+2)/4, b=(2−ψ)/ϕ,
y y
7
andthesubscripti=A,B,C,Ddeterminestheregionin For q = 0, only regions A and B survive and Eq. (B2)
x
the (ϕ,ψ) space (see Fig. 5): simplifies to
(cid:112)
A: −2<ψ <2, ϕ>2−ψ, g τ |ω|θ(1−τ)θ(2−d)sgn(ω)
ImΠ (ω,0,q )=− s √
B: ψ <−2, ϕ>2−ψ, 0 y 12π av(1−τ)3/4
(B4)
C: ψ <−2, −2−ψ <ϕ<2−ψ, (1−d)K(u)+2dE(u), |d|<2,
 (cid:101) (cid:101)
D: −2<ψ <2, 0<ϕ<2−ψ, × 2+d2 d(2−d)
 K(1/u)+ E(1/u), d<−2,
2u (cid:101) 2u (cid:101)
and (cid:101) (cid:101)
(B6)
ρ =ρ−1 =b, ρ =ρ−1 =ν, √
A D C √B where u= 2−d/2.
(cid:101)
t =t =t−1 =t−1 = bν,
A C B D
s =s =1, s =s =−1,
C D A B
√
h =h =1, h =h = bν,
A C B D
r =−ϕ(cid:2)12+4ψ+3ϕ−τ(4+4ψ+ϕ)(cid:3), (B5)
A
(cid:2) (cid:3)
r =4ν 2ψ−3ϕ+τ(2ψ+ϕ) ,
B
(cid:2) (cid:3)
r =−4 4+3ϕ−τ(4+ϕ) ,
C
r =(ψ−2)(8+2ψ+3ϕ)
D
+τ[ψ(4+2ψ−ϕ)+2(8+ϕ)].
Appendix C: Expressions for f (∆/µ)
x,y
The functions (22) determining the long-wavelength plasmon frequency can be written in terms of the complete
elliptic integrals as
8√2(cid:40)2(9−2δ2)E(u)−(1+δ)(9−2δ)K(u), |δ|<1,
f (δ)= (C1)
x 15 2u(cid:2)(9−2δ2)E(1/u)+2δ(1+δ)K(1/u)(cid:3), δ <−1,
√ (cid:40)
4 2 (1+δ)K(u)−2δE(u), |δ|<1,
f (δ)= (C2)
y 3 2u(cid:2)(1+δ)K(1/u)−δE(1/u)(cid:3), δ <−1,
(cid:112)
where u= (1−δ)/2. According to the notations used in Ref. [34],
f (δ)=4I (0,δ), f (δ)=I (0,δ). (C3)
x 3 y 2
[1] Physics of Graphene,editedbyH.AokiandM.S.Dres- [5] D. Langbein, Phys. Rev. 180, 633 (1969); D. R. Hof-
selhaus (Springer, New York, 2014). stadter, Phys. Rev. B 14, 2239 (1976); L. A. Pono-
[2] D.S.L.Abergel,V.Apalkov,J.Berashevich,K.Ziegler, marenko, R. V. Gorbachev, G. L. Yu, D. C. Elias, R.
and T. Chakraborty, Adv. Phys. 59, 261 (2010). Jalil, A. A. Patel, A. Mishchenko, A. S. Mayorov, C. R.
[3] J. Wang, S. Deng, Z. Liu and Z. Liu, Nat. Sci. Rev. 2, Woods, J. R. Wallbank, M. Mucha-Kruczynski, B. A.
22 (2015) [arXiv:1410.5895]. Piot, M. Potemski, I. V. Grigorieva, K. S. Novoselov,
[4] T. Chakraborty and V. Apalkov, in [1], Chap. 8; Solid F. Guinea, V. I. Fal’ko, and A. K. Geim, Nature 497,
StateCommun.175–176,123(2013)[arXiv:1302.2940]; 594 (2013) [arXiv:1212.5012]; B. Hunt, J. D. Sanchez-
T. Chakraborty and P. Pietila¨inen, Phys. Rev. Lett. 76, Yamagishi, A. F. Young, M. Yankowitz, B. J. LeRoy,
4018 (1996); Phys. Rev. B 39, 7971 (1989); D.S.L. K. Watanabe, T. Taniguchi, P. Moon, M. Koshino, P.
Abergel and T. Chakraborty, Phys. Rev. Lett. 102, Jarillo-Herrero, and R. C. Ashoori, Science 340, 1427
056807 (2009) [arXiv: 0811.0574]; V. M. Apalkov and (2013) [arXiv:1303.6942]; C. R. Dean, L. Wang, P. Ma-
T. Chakraborty, ibid. 97, 126801 (2006) [arXiv: cond- her, C. Forsythe, F. Ghahari, Y. Gao, J. Katoch, M.
mat/0606037]; 105, 036801 (2010) [arXiv: 1007.4188]; Ishigami,P.Moon,M.Koshino,T.Taniguchi,K.Watan-
107, 186803 (2011) [arXiv:1104.2020]. abe, K. L. Shepard, J. Hone, and P. Kim, Nature 497,
8
598 (2013) [arXiv:1212.4783]. 103, 016402 (2009) [arXiv:0906.1564].
[6] V. M. Apalkov and T. Chakraborty, Phys. Rev. Lett. [24] H. Huang, Z. Liu, H. Zhang, W. Duan, and D. Vander-
112, 176401 (2014) [arXiv: 1311.4935]; A. Ghazaryan bilt,Phys.Rev.B92,161115(2015)[arXiv:1507.01742].
and T. Chakraborty, Phys. Rev. B 91, 125131 (2015) [25] L. Tarruell, D. Greif, T. Uehlinger, G. Jotzu, and T.
[arXiv:1412.2792];A.Ghazaryan,T.Chakraborty,andP. Esslinger,Nature483,302(2012)[arXiv:1111.5020];L.-
Pietila¨inen,J.Phys.: Condens.Matter27,185301(2015) K. Lim, J.-N. Fuchs, and G. Montambaux, Phys. Rev.
[arXiv: 1408.3424]; A. Ghazaryan and T. Chakraborty, Lett. 108, 175303 (2012) [arXiv:1201.1479].
Phys. Rev. B 92, 235404 (2015) [arXiv:1502.01373]. [26] M. C. Rechtsman, Y. Plotnik, J.M. Zeuner, D. Song, Z.
[7] P. De Padova, C. Ottaviani, C. Quaresima, B. Olivieri, Chen, A. Szameit, and M. Segev, Phys. Rev. Lett. 111,
P. Imperatori, E. Salomon, T. Angot, L. Quagliano, C. 103901 (2013) [arXiv:1211.5683].
Romano, A. Vona, M. Muniz-Miranda, A. Generosi, B. [27] M. Bellec, U. Kuhl, G. Montambaux, and F. Mortes-
Paci,andG.LeLay,2DMater.1,021003(2014);J.Sone, sagne, Phys. Rev. Lett. 110, 033902 (2013) [arXiv:
T. Yamagami, Y. Aoki, K. Nakatsuji, and H. Hirayama, 1210.4642].
New J. Phys. 16, 095004 (2014); M. E. Da´vila, L. Xian, [28] J. Kim, Seung S. Baik, S. H. Ryu, Y. Sohn, S. Park, B.-
S. Cahangirov, A. Rubio and G. Le Lay, New J. Phys. G.Park,J.Denlinger,Y.Yi,H.J.Choi,andK.S.Kim,
16, 095002 (2014) [arXiv:1406.2488]. Science 349, 723 (2015) [arXiv:1508.05523].
[8] V. M. Apalkov and T. Chakraborty, Phys. Rev. B 90, [29] P. Dietl, F. Pi´echon, and G. Montambaux, Phys. Rev.
245108 (2014) [arXiv: 1407.4441]; 91, 235447 (2015) Lett. 100, 236405 (2008) [arXiv:0707.0219]; P. Delplace
[arXiv:1504.01749]. and G. Montambaux, Phys. Rev. B 82, 035438 (2010)
[9] W. Luo and T. Chakraborty, Phys. Rev. B 92, 155123 [arXiv:1001.3952].
(2015) [arXiv:1506.07896]. [30] S. Yuan, M. I. Katsnelson, and R. Rolda´n, arXiv:
[10] Xi Ling, H. Wang, S. Huang, F. Xia, and M. S. Dres- 1512.06345.
selhaus, PNAS 112, 4523 (2015) [arXiv:1503.08367]; H. [31] D. Carpentier, A. A. Fedorenko, and E. Orignac, Euro-
Liu, Y. Du, Y. Deng and P. D. Ye, Chem. Soc. Rev. 44, phys. Lett. 102, 67010 (2013) [arXiv:1303.4320].
2732 (2015) [arXiv:1411.0056]. [32] H. Isobe, B.-J. Yang, A. Chubukov, J. Schmalian, and
[11] A. Ghazaryan and T. Chakraborty, Phys. Rev. B 92, N.Nagaosa,Phys.Rev.Lett.116,076803(2016)[arXiv:
165409 (2015) [arXiv:1507.04709]. 1508.03781].
[12] A. Tsukazaki, A. Ohtomo, T. Kita, Y. Ohno, H. [33] G. Y. Cho and E.-G. Moon, Sci. Rep. 6, 19198 (2016)
Ohno, and M. Kawasaki, Science 315, 1388 (2007); A. [arXiv:1508.03777].
Tsukazaki,S.Akasaka,K.Nakahara,Y.Ohno,H.Ohno, [34] P. Adroguer, D. Carpentier, G. Montambaux, and E.
D.Maryenko,A.Ohtomo,andM.Kawasaki,Nat.Mater. Orignac, arXiv:1511.00036.
9, 889 (2010); J. Falson, D. Maryenko, B. Friess, D. [35] S.BanerjeeandW.E.Pickett,Phys.Rev.B86,075124
Zhang, Y. Kozuka, A. Tsukazaki, J.H. Smet and M. (2012) [arXiv:1207.5760].
Kawasaki, Nat. Phys. 11, 347 (2015). [36] T. Low, R. Rolda´n, H. Wang, F. Xia, P. Avouris, L. M.
[13] W. Luo and T. Chakraborty, arXiv:1512.05815. Moreno, and F. Guinea, Phys. Rev. Lett. 113, 106802
[14] Y. Hasegawa, R. Konno, H. Nakano, and M. Kohmoto, (2014) [arXiv:1404.4035].
Phys. Rev. B 74, 033413 (2006) [arXiv: cond- [37] A. S. Rodin and A. H. Castro Neto, Phys. Rev. B 91,
mat/0604433]. 075422 (2015) [arXiv:1412.0811].
[15] S.-L. Zhu, B. Wang, and L.-M. Duan, Phys. Rev. Lett. [38] F. Jin, R. Rolda´n, M. I. Katsnelson, and S. Yuan, Phys.
98, 260402 (2007) [arXiv:cond-mat/0703454]. Rev. B 92, 115440 (2015) [arXiv:1507.07013].
[16] B. Wunsch, F. Guinea, and F. Sols, New J. Phys. 10, [39] A. N. Rudenko and M. I. Katsnelson, Phys. Rev. B 89,
103027 (2008) [arXiv:0807.4245]. 201408 (2014) [arXiv:1404.0618].
[17] G.Montambaux,F.Pi´echon,J.-N.Fuchs,andM.O.Go- [40] T. Stauber, J. Schliemann, and N. M. R. Peres, Phys.
erbig,Phys.Rev.B80,153412(2009)[arXiv:0904.2117]; Rev. B 81, 085409 (2010) [arXiv:0910.1242].
Eur. Phys. J. B 72, 509 (2009) [arXiv:0907.0500]. [41] S.Yuan,R.Rolda´n,andM.I.Katsnelson,Phys.Rev.B
[18] V. M. Pereira, A. H. Castro Neto, and N. M. R. Peres, 84, 035439 (2011) [arXiv:1103.5350].
Phys. Rev. B 80, 045401 (2009) [arXiv:0811.4396]. [42] Kenneth W.-K. Shung, Phys. Rev. B 34, 979 (1986);
[19] Q. Liu, X. Zhang, L. B. Abdalla, A. Fazzio, and A. E.H. Hwang and S. Das Sarma, ibid. 75, 205418 (2007)
Zunger, Nano Lett. 15, 1222 (2015) [arXiv:1411.3932]. [arXiv:cond-mat/0610561]; B. Wunsch, T. Stauber, F.
[20] K. Dolui and S. Y. Quek, Sci. Rep. 5, 11699 (2015) Sols, and F. Guinea, New J. Phys. 8, 318 (2006) [arXiv:
[arXiv:1503.03647]. cond-mat/0610630].
[21] J. M. Pereira and M. I. Katsnelson, Phys. Rev. B 92, [43] According to the notations of Ref. [35], f (0) = 16I ,
x 3
075437 (2015) [arXiv: 1504.02452]; A. N. Rudenko, S. f (0)=4I .
y 2
Yuan, and M. I. Katsnelson, Phys. Rev. B 92, 085419 [44] X.-F. Wang and T. Chakraborty, Phys. Rev. B 75,
(2015) [arXiv:1506.01954]. 033408 (2007) [arXiv: cond-mat/0605498]; 75, 041404
[22] S.S.Baik,K.S.Kim,Y.Yi,andH.J.Choi,NanoLett. (2007) [arXiv: cond-mat/0611635]; 81, 081402 (2010)
15, 7788 (2015) [arXiv:1508.04932]. [arXiv:1001.2762];P.K.Pyatkovskiy,J.Phys.: Condens.
[23] V. Pardo and W. E. Pickett, Phys. Rev. Lett. 102, Matter 21, 025506 (2009) [arXiv:0808.0931].
166803 (2009) [arXiv: 0903.4820]; Phys. Rev. B 81, [45] I. Gradshteyn and I. Ryzhik, Table of Integrals, Series,
035111 (2010) [arXiv: 0910.4411]; S. Banerjee, R.R.P. and Products, 7th ed. (Academic, New York, 2007).
Singh, V. Pardo, and W. E. Pickett, Phys. Rev. Lett.