Table Of ContentFunctional renormalization group and variational Monte Carlo studies of the
electronic instabilities in graphene near 1/4 doping
Wan-Sheng Wang,1 Yuan-Yuan Xiang,1 Qiang-Hua Wang,1 Fa Wang,2 Fan Yang,3 and Dung-Hai Lee4,5
1National Lab of Solid State Microstructures, Nanjing University, Nanjing, 210093, China
2Department of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA
3Department of Physics, Beijing Institute of Technology, Beijing 100081, China
4Department of Physics, University of California at Berkeley, Berkeley, CA 94720, USA
5Materials Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA
2
1
Westudytheelectronicinstabilitiesofnear1/4electrondopedgrapheneusingthesingular-mode
0
functional renormalization group, with a self-adaptive k-mesh to improve the treatment of the van
2
Hovesingularities, and variational Monte-Carlo method. At 1/4 doping the system is a chiral spin
n density wavestate exhibiting theanomalous quantized Hall effect. When thedoping deviates from
a 1/4,thedx2−y2+idxy Cooperpairingbecomestheleadinginstability. Ourresultssuggestthatnear
J 1/4electron-orhole-doping(awayfromtheneutralpoint)thegrapheneiseitheraCherninsulator
3 or a topoligical superconductor.
1
PACSnumbers: 74.20.-z,74.20.Rp,74.70.Wz, 81.05.ue,71.27.+a
]
n
o
INTRODUCTION willbe trulyexciting,wefeelitismeaningfulto examine
c
- this problemusing the more realisticband structureand
r interaction parameters.
p Graphene, a single atomic layer of graphite, has been
u afocusofinterestsincethepioneeringworkofNovoselov InviewoftheheavydopingweusetheHubbardinter-
s action to model the screened Coulomb interaction. We
and Geim[1]. At the fundamental level the past research
.
t activities on graphene mostly focused on exploring the perform singular-mode functional renormalizationgroup
a
(SM-FRG)[14] and variational Monte Carlo (VMC) cal-
m consequences of the unique Dirac-like bandstructure[2].
culations to address the possible electronic instabilities.
Onthe experimentalside,few exceptionsinclude the ob-
-
d servation of the fractional quantum Hall effect[3, 4], the Since the interaction strength is estimated to be a frac-
n tionofthe bandwidth, we believe SM-FRGshouldyield
detection of the Fermi velocity renormalization[5], and
o qualitativelycorrectanswer. TheVMCisusedtofurther
thepossibleobservationof“plasmaron”inangle-resovled
c
confirm such belief. The main results are summarized
[ photosemission[6]. In general the effects of electron-
as follows. At 1/4 electron doping and for interaction
electron interaction on the properties of graphene re-
2 strength appropriate for graphene we found the chiral
main a frontier of this field. Previously based on the
v
4 resonating-valence-bond[7] concept Pathak et al.[8] and spin density wave (SDW) state is the dominating insta-
bility. When the doping level slightly deviates from 1/4
8 Black-Schaffer and Doniach [9] proposed that doped
8 graphene should be a high temperature superconduc- we find the d+id′ pairing instability surpasses that of
3 tor with d+id′ pairing symmetry. (Henceforth d and the chiral SDW. We propose a schematic phase diagram
9. d′ are used to denote interchangeably dx2−y2 and dxy inFig.6(b). As inpnictides andoverdopedcuprates[15],
0 symmetries, respectively.) In particular, the possibility thepairingmechanismisduetoastrongscatteringchan-
1 nel shared by the SDW and pairing.
of unusual superconductivity and other orders in doped
1
graphene with van Hove singularities at (or near) the
:
v Fermilevelbecomesahotissue.[10,11]Mostrecentlyby
i MODEL
X aperturbativerenormalizationgroupcalculationNandk-
ishoreet al. concludedthatthevanHovesingularitieson
r
a the Fermi surface drive chiral d+id′ superconductivity The real-space hamiltonian we used is given by
in the limit of vanishing interaction strength[12].
H = (c† t c +h.c.) µN +U n n
Ona differentfrontTaoLirecently proposedthatdue − iσ ij jσ − e i↑ i↓
X X
to the existence of Fermi surface nesting 1/4 electron (ij)σ i
doped Hubbard model on honeycomb lattice favors the 1
+ V n n , (1)
i i+δ
formationof a magnetic insulating state which possesses 2 X
iδ
nonzero spin chirality and exhibit the anomalous quan-
tized Hall effect, hence is a Chern insulator[13]. Thus where (ij) denotes bonds connecting sites i and j, σ is
near quarter doping graphene suddenly becomes a play- the spin polarity, µ is the chemical potential, N is the
e
ing ground where either a Chern insulator or a topolog- totalelectronnumberoperator,theU-termistheon-site
ical superconductor can potentially be realized. Because Hubbard interaction and V is the Coulomb interaction
therealizationofeitherphasesinheavilydopedgraphene onnearest-neighborbondsδ. Thehoneycomblatticehas
2
two sublattices which we denote as A and B henceforth.
(a) (b)
AssuggestedinRef.[2]wetaket1 =2.8eV,t2 =0.1eV, σ σ σ σ
and t = 0.07 eV for hoppings between the first, second 1 4 1 4
3 k+q p+q
the third neighbors, respectively, and set U = 3.6t1. As (m) (n)
for V, we expect V < U in doped graphene, and take
k p
V =t1 asatypicalupperbound. Theoretically,enriched 2 τ τ 3 2 −τ −τ 3
phases mayappear forevenlargervaluesof V/U.[11,16]
(c) (d)
σ σ
1 4 1 3
k+q p+q σ k+q p+q τ
(m) (n) (m) (n)
METHOD k p σ k p τ
3 2 4 2
τ τ
The SM-FRG method[14] we used is a modification
FIG. 1: A generic 4-point vertex (a) is rearranged into the
of the FRG method[17] applied to the cuprates[18] and
pairing (b),crossing (c) and direct (d)channels. Herek,q,p
pnictides[19]. Fig.1 (a) shows a generic 4-point ver-
aremomenta,σ andτ denotespinswhich areconserveddur-
tex function Γ which appears in the interaction
1234 ingfermionpropagation,andm,ndenotetheformfactor(see
c†1c†2(−Γ1234)c3c4. Here1,2,3,4representmomentum(or thetextfordetails). Theopenarrowsindicatecollectiveprop-
realspaceposition)andsublattice. Thespinσ andτ are agation.
conserved during fermion propagation, and will be sup-
pressed henceforth. Figs.1(b)-(d) are rearrangements of
tion can be written as, in matrix form,
(a) into the pairing (P), the crossing (C) and the direct
(D)channelsinsuchawaythatacollectivemomentumq ∂P/∂Λ=Pχ′ P,
pp
can be associated and the other momentum dependence
∂C/∂Λ=Cχ′ C,
can be decomposed as, ph
∂D/∂Λ=(C D)χ′ D+Dχ′ (C D), (3)
− ph ph −
Γ f∗(k)P (q)f (p), where the collective momentum q is left implicit for
k+q,−k,−p,p+q→X m mn n brevity, χ′ are loop integrations projected by the
mn pp/ph
form factors (See Appendix for details), and Λ is the
Γ f∗(k)C (q)f (p),
k+q,p,k,p+q → m mn n running cutoff energy. Integrating Eq. (3) with respect
X
mn
to Λ yields the ladder approximation. However, since
Γ f∗(k)D (q)f (p). (2)
k+q,p,p+q,k → m mn n ∂P, ∂C and ∂D add up to the full change dΓ, the full
X
mn flow equations for P,C and D should be given by
dP/dΛ=∂P/∂Λ+Pˆ(∂C/∂Λ+∂D/∂Λ),
Here f is a set of orthonormal lattice form fac-
m
tors. {For}honeycomb lattice the form factor label m dC/dΛ=∂C/∂Λ+Cˆ(∂P/∂Λ+∂D/∂Λ),
also includes a sublattice label, m = (m,a) with a = dD/dΛ=∂D/∂Λ+Dˆ(∂P/∂Λ+∂C/∂Λ), (4)
A/B, within our choice of C point group with re-
3v where the Pˆ,Cˆ and Dˆ are the projection operators in
spect to the atomic site. (See Appendix.) The de-
the senseofEq.(2), andweusedthe factthatKˆ(∂K)=
composition into each channel would be exact if the
∂K for K = P,C,D. In Eq. (4) the terms preceded
set is complete. In practice, however, a set of a few
by the projection operatorsrepresentthe overlapsof the
form factors is often sufficient to capture the symme-
three different channels. It is those terms which allow
try of the order parameters associated with leading
pairing to be induced by virtual particle-hole scattering
instabilities.[14] The momentum space form factors are
processes[15].
the Fourier transform of the real-space ones: 1) on-site,
It can be shown that the effective interaction in the
f (r) = 1; 2) 1st neighbor bonds, f (r) = 1/3,
s0 s1
f (r) = 2/3cos(lθ ) and f (r) = 2/3sinp(lθ ), superconducting (SC), spin density wave (SDW) and
d1 r d1′ r chargedensity wave(CDW) channels are givenby V =
where l =p2 and θ is the azimuthal angple of r;[20] 3) sc
r
P, V = C, and V = C 2D, respectively. By
2nd neighbor bonds, fs2(r) = 1/6, fp2 = 1/3cosθr, −singulasrdwvalue decompocsdiwtion, we−determine the leading
f (r) = 1/3sinθ , f (r) =p 1/3cos2θp, f (r) =
p2′ r d2 r d2′ instability in each channel,
1/3sin2pθ , f (r) = 1/6cos3pθ . These form factors
r f2 r
(pcombinedwithsublattpicelabels)areusedinallchannels Vmn(q )= Sαφα(m)ψα(n), (5)
X X X X X
on equal footing. We have tested that further neighbor Xα
form factors do not change the results qualitatively. where X = sc,sdw,cdw, Sα is the singular value of the
X
The one-loopcorrectiontothe flowof the vertexfunc- α-th singular mode, φα and ψα are the right and left
X X
3
the outer hexagon). The flow of the most negative sin-
gular values (denoted as S) in the SC, SDW and CDW
channels are shownin Fig.2(b). Clearly the SDW (green
dot-dashed line) is the leading instability. This is be-
cause the SDW scattering is already attractive at high
energies and is further enhanced by the Fermi surface
nesting shown in Fig.2(a) down to the lowest energies.
This SDW singularmode φ (m) has a dominantvalue
sdw
for m = (s ,A/B), showing that the magnetic order-
0
ing moment is site-local. It is also identifiable in the
renormalizedinteraction Vmm(q) for m=(s ,A) shown
sdw 0
in Fig.2(c), which has strong peaks at three indepen-
dent momenta Q = (0,2π/√3), Q = ( π,π/√3),
1 2
−
Q = (π,π/√3) and their symmetric images. They
3
define the possible ordering vectors q for the SDW
sdw
order. In contrast, the scattering associated with pair-
ing (blue solid line) is initially repulsive in all channels,
and only becomes attractive in the d-wave channel af-
ter the SDW scatteringgrowsstrong. RegardingCooper
pairing we find two degenerate leading form factors: the
FIG. 2: (Color online) Results for δ = 1/4 and V = 0. (a) dx2−y2 and dxy doublets. One of these is used to gener-
The Fermi surface. The brightness on the surface represents atethemomentumspacegapfunctionshowninFig.2(a).
themomentumspacegapfunctionassociated withoneofthe The singularpairingmode φ (m)is nonzeroforthe 2nd
sc
degenerate pairing modes. (b) FRG flow of (the inverse of)
neighbor bonds, but the amplitude is about one order of
the most negative singular values S in the SC (blue solid
magnitude smaller than that for the 1st neighbor bonds,
line), SDW (green dot-dashed line) and CDW (red dashed
hence justifying the cutoff in the real-space range of the
line) channels. (c) and (d) are the renormalized interaction
Vsmdwm form=(s0,A),andVsmcm form=(d1,A),asfunctions form factors. The renormalizedVsmcm(q) for m=(d1,A)
of the collective momentum q. The hexagons in (a), (c) and in Fig.2(d) has negative(but weak) peaks at q=0, con-
(d) indicate theBrillouine zone boundary. firming the Cooper instability at this momentum. The
CDW channel (red dashed line) also becomes weakly at-
tractive from Fig.2(b). We checked that the singular
eigenvectorsofVX,respectively. We fixthephaseofthe mode φcdw(m) has significant values for m = (s1,A/B)
eigen vectors by requiring Re[ mφαX(m)ψXα(m)] > 0 so and m = (d1,1′,A/B), showing that it is an extended
that SXα < 0 corresponds to aPn attractive mode in the CDW.Themixtureofs1 andd1,1′ is duetothe factthat
X-channel. In the pairing channel qsc = 0 addresses theCDWwavevectorq=(2π/3,π)(oritssymmetryim-
the Cooper instability. The pairing function in the sub- ages) is not invariant under the point group operations.
lattice basis can be constructed from φα, and a fur- However, the CDW channel remains weak in our case
sc
ther unitary transformis needed to get the pairing func- hencewillnotbediscussedintherestofthepaper. (The
tion in the band basis. (See the Appendix for details.) merging of SC and CDW lines in Fig.2(b) is induced by
The ordering wave vector in the SDW/CDW channel the diverging SDW channel via the overlapping among
q = qsdw/cdw is chosen at which Vsdw/cdw(q) is max- these channels.)
| |
imal. We note that such a vector has symmetry-related
The above results indicate three independent and de-
images,andmaychangeduringtheFRGflowbeforeset-
generate SDW momenta (apart from the global spin
tlingdowntofixedvalues. TheRGflowisstoppedifany
SU(2) symmetry). A calculation that keeps the
of P , C , or D becomes roughly 10 times
max max max symmetry-breaking self-energy flow is needed to de-
| | | | | |
of the bandwidth.[21]
cide which combination of them is realized in the or-
More technical details can be found in the Appendix.
dered state, but this is beyond the scope of the present
work. Alternatively, one may resort to a Landau the-
ory or mean field theory. Indeed, according to the mean
RESULTS AND DISCUSSION field theory of Ref. [13], a particular linear combination
[13] S = M eiQ3·R + M eiQ1·R + M eiQ2·R and
R,A 3 1 2
h i
We define the doping level by δ = n 1 where n S = M eiQ3·R M eiQ1·R M eiQ2·R, gives the
e e R,B 3 1 2
− h i − −
is the number of electrons per site. We first discuss the most energetically favorable spin structure. Here R la-
results for δ =1/4 and V =0. Fig.2(a) shows the Fermi bels unit cell and M are three mutually orthogonal
1,2,3
surface in this case, which is well nested and close to and equal length vectors. The handedness of the M
1,2,3
the van Hove singularities (the mid points of edges of triad breaks time-reversal and spatial reflection symme-
4
4
(a) (b) d+id′, 12 × 12
2
d+id′, 18 × 18
0 chiral SDW, 18 x18
)1
4t
−
10 −2
E (
∆
−4
−6
FIG. 3: (Color online) Chiral SDW order on the honey-
comb lattice. (a) The spins on the black,red,green and blue −8
0 1 2 3 4 5 6
sublattices (of different gray scales) are M1 + M2 + M3, ∆ (10−2t )
−M1−M2+M3,M1−M2−M3,−M1+M2−M3 respec- 1
tively. (b) A 3D perspectiveview of thechiral SDW order.
FIG. 4: (Color online) Variational Monte Carlo results for
the energy gain per site, ∆E, versus the variational order
′
try. The resulting four-sublattice chiral SDW order is parameters ∆ for the d+id-SC (circles and triangles) and
the chiral SDW states (squares) at the doping level δ =1/4.
shown in Fig. 3.
Thelattice sizesare givenin thelegends. Linesare guidesto
Since our result is at odd with that of Ref. [12] which
theeye.
applies in the limit of vanishing interaction strength,
we further check the above conclusion by a variational
Monte-Carlo (VMC) calculation using exactly the same
parameters as for Figs.2. We adopted the partially-
projectedmean-fieldwave-functions[22]asourtrialwave-
functions,withthetheformfactorsguidedbythepresent
SM-FRG result. Fig. 4 shows the energy gain per site
due to d+id′ pairing on 12 12 (circles) and 18 18
× ×
(triangles) lattices, showing negligible size dependence.
We then compare to the energy gainassociatedwith the
chiral SDW (squares). It is clear that the SDW state is
far more energetically favorable than the chiral d+id′-
SC state at this doping level. The reason that the SDW
state is realized in our lattice model lies in the fact that
the perfect Fermi surface nesting is as important as the
inter-saddle scattering addressed in Ref.[12].
Below1/4doping,thetheFermilevelmovesawayfrom
the van Hove points and the Fermi surface nesting wors-
ens. This is shown in Fig. 5(a) for δ = 0.211 as an ex-
ample. The bare interactions are still set as U = 3.6t
1
and V = 0. The SDW scattering is still attractive at
FIG. 5: (Color online) The same plots as in Fig. 2 but for
high energies. As seen in Eq.(4), this relatively strong δ =0.211. Note the splitting of the SDW peaks in panel (c)
SDW channel causes attraction in the SC channel via signifies theincommensurate SDW instability.
overlap between these channels (terms with the projec-
tion operator in Eq. (4)). At even lower energy scales
the pairing channel attraction (with q=0) continues to verging superconducting pairing channel again find the
grow due to the Cooper instability, while the enhance- degenerate dx2−y2 and dxy doublets, with dominant am-
ment of the SDW scattering is saturated due to the lack plitudes for m = (d1,1′,A/B). The momentum space
of Fermi surface nesting. As the result the pairing insta- gap function of one of the pairing modes is shown in
bility surpassesthe SDWinstability atthe lowestenergy Fig. 5(a). Fig. 5(c) shows the renormalized Vmm(q) for
sdw
scale. This is shown by the flow of the singular values m = (s0,A), which shows weak peaks at six indepen-
in Fig. 5(b). It is worth to mention that precisely the dentandincommensuratemomenta. Fig.5(d)showsthe
same phenomenon was observed in the FRG studies of renormalized Vsmcm(q) for m = (d1,A), which shows a
the cuprates and pnictides[15, 18]. A close inspection strong negative peak at q=0.
of the eigenvectors φ (m) associated with the most di- The above results imply degenerate d-wave pairing
sc
5
The phase diagram for V = t is also drawn in Fig.6(a)
1
(a) VV==0t1 (b) (dashedline). Thecriticalscaleisslightlyhigherthanthe
case of V = 0. In the SC region this is due to the slight
10−2
Λ / tc1 SC SDW SC T CSDhiWral Uennhlaiknecetmhaetntcloafimd-ewdavine Rpaeifr.[i2n3g],onin2anldl cnaesigeshbsoturdbieodndisn.
this paperthe f-wavepairingis nota leadinginstability.
10−4
SC SC We end by presenting Fig. 6(b) as a schematic phase
0 .15 0.2 0.25 0.3 1/4 diagram in the temperature-doping plane. In reality
δ δ
when the doping level slightly deviates from1/4,the ex-
tra charges will be localized by the presence of disorder,
FIG. 6: (Color online) (a) The FRG diverging energy scale which enables the system to stay in the Chern insulator
Λc plottedasafunction ofdopinglevelnearδ =1/4. crosses state fora finite doping interval. In the transitionregion
and open circles represent Λc associated with the SC and marked by gray, where the doped charges delocalize, in-
SDW channel, respectively. V = 0,t1 for solid and dashed commensurate SDW states with wave vectors shown in
lines,respectively. (b)Aschematictemperature-dopingphase
Fig.5(c) will emerge.
diagramnearδ =1/4inlinearscales. Thegreyregiondenotes
thetransition between SC and SDW.
SUMMARY
instabilities. As for the degenerate SDW instabilities,
additional analysis, such as the mean field theory or
Insummary,we haveperformedSM-FRGcalculations
Ginzburg-Landautheory,isneededtofixthestructureof
for parameters suitable for graphene. Our results indi-
the pairing function in the ordered state. To a large ex-
cate that near 1/4 electron- or hole- doping, graphene is
tent,thiskindofanalysishasbeenperformedinRef.[12].
either a Chern insulator or a chirald-wavesuperconduc-
We have also performed simple mean field calculations
tors. Both phases are topological in nature, and deserve
using the renormalized pairing interaction. The result is
experimental efforts in searching them.
that a time-reversalbreaking dx2−y2±idxy-wavepairing We thank Hong Yao for helpful discussions, and
is always more favorable. This could have been antici-
are grateful to Tao Xiang for computing resources.
patedsincebothdx2−y2 anddxy formfactorshavenodes QHW acknowledges the support by NSFC (under
on the Fermi surface, a natural way to gain energy is to
grant No.10974086 and 10734120) and the Ministry
formtheabovechirald-wavepairing,whichgapsoutthe
of Science and Technology of China (under grant
entire Fermi surface.
No.2011CBA00108 and 2011CB922101). FY acknowl-
We have also checkednearby doping and analyzedthe
edges the support by NSFC (under grant No.10704008).
competition between the incommensurate SDW and the
DHL acknowledges the support by the DOE grant num-
SC state. In Fig. 6(a) we plot the higher diverging scale
ber DE-AC02-05CH11231.
among these two competing orders as a function of dop-
ing (solid line). We found a similar phase diagram near
-1/4 (hole) doping (not shown), mirroring that of elec-
APPENDIX
tron doping. (Notice that the particle-hole symmetry is
not exact in the presence of hopping integrals t ).
2,3
We now discuss briefly the effect of the nearest neigh- Toillustratetheideaofthemethod,wefirstignorethe
borinteractionV. Asatypicalexample,wesetU =3.6t sublattice index, and return to it at a later stage.
1
and V = t , and perform the FRG calculations. We The one-loop contributions to the flow of the irre-
1
findthatthe resultsarequalitativelysimilarto the cases ducible 4-point vertex function is shown in Fig.7, where
with V = 0, except that in the leading pairing singu- (a) and (b) lead to the flow of P and C, respectively,
lar mode, φsc(d2,2′,A/B) becomes slightly stronger, but and (c)-(e) lead to the flow of D. The internal Greens
still smaller than φsc(d1,1′,A/B) by a factor of 4 6. functions are convoluted with the form factors, hence
∼
6
∂ dω d2p
(χ′ ) = n f (p)G(p+q,iω )G( p, iω )f∗(p)θ(ω Λ)
pp mn ∂ΛZ 2π Z S m n − − n n | n|−
BZ
1 d2p
= f (p)G(p+q,iΛ)G( p, iΛ)f∗(p) +(Λ Λ),
−2π Z S m − − n →−
BZ
∂ dω d2p
(χ′ ) = n f (p)G(p+q,iω )G(p,iω )f∗(p)θ(ω Λ)
ph mn ∂ΛZ 2π Z S m n n n | n|−
BZ
1 d2p
= f (p)G(p+q,iΛ)G(p,iΛ)f∗(p) +(Λ Λ), (6)
−2π Z S m n →−
BZ
(a) (b)
(a) (b)
1 1
P P C C 0.5 0.5
m m′ n′ n m m′ n′ n πk / y 0 πq / y 0
−0.5 −0.5
−1 −1
(c)
−1 0 1 −1 0 1
k / π q / π
D D x x
m m n n
′ ′
FIG. 8: (Color online) Examples of self-adaptive k-mesh (a)
and q-mesh (b) used for loop integrations and collective mo-
(d) (e) menta for interactions, respectively. Notice that in (a) the
mesh is too densenear the Fermi surface to be differentiated
D C C D bynaked eyes.
m m n n m m n n
′ ′ ′ ′
1,2,3,4 in Figs.1, so in principle the label m in the
FIG. 7: One loop diagrams contributing to the flow of the
form factor f (r) also includes sublattice indices. How-
m
the 4-point vertex function in the pairing channel (a), cross-
′ ′ ever, once r is fixed they are not independent. We ab-
ing channel (b), and direct channel (c)-(e). Here m,mn,n
sorb an independent index into the form factor label,
denote form factors, while the momentum and spin indices
areleftimplicit. Theopenarrowsindicatetheflowofthecol- m (m,a), where a labels fermion field 1 or 4 in
→
lectivemomentum. Theslashed linesaresingle-scale fermion Fig.1(b), 1 or 4 in (c), and 1 or 3 in (d). This is rea-
propagators. The slash can be placed on either internal lines sonable because the point group operations do not mix
associated with theloop. the sublattice when the center is chosen to be atomic
sites. 2)In the presence ofsublattice indices the Green’s
functions are matrices. 3) In orderto ensure that in mo-
where G is the bare fermion propagator,S is the area
BZ mentumspaceP, C andD transformexactlyasproduct
of the Brillouine zone. Here Λ > 0 is the infrared cut-
offormfactors,caremustbetakeninchoosingthephase
off of the Matsubara frequency ω .[24] As in usual FRG
n of the Bloch states for complex unit cells.
implementation,theselfenergycorrectionandfrequency
The pairing function is determined as follows. A
dependence of the vertex function are ignored.
singular mode φα corresponds to a pairing operator
Ingeneral,theformfactorf (k)= f (r)exp( ik sc
m r m − · c† (k)φα(m)f (k)∗c† ( k) in the momen-
r), where fm(r) transforms accordingPto an irreducible m=(m,a) a↑ sc m am↓ −
representation of the point group, and r is the relative Ptum space, where am is determined by m = (m,a) (for
position vector between the two fermion fields on each allassociatedvectorsr). Theparityofthepairingmatrix
sideofthediagramsinFig.1(b)-(c). Fortwotypesofdi- function under space inversiondetermines automatically
agramstooverlap,allofthefourfermionfieldssitwithin whether itis a spinsingletorspin triplet. Afurther uni-
the range set by the form factors. Hence the projections tary transform can be used to get the momentum space
in Eq.(4) are all preformed in real space. gap function.
We now return to the honeycomb lattice with two In the current implementation of SM-FRG the sam-
sublattices. The necessary modifications are as follows: plingofmomentumspace(kandq)isperformedonself-
1) The sublattice label can be absorbed into the labels adaptive meshes. As an example, Fig.8(a) shows the k-
7
mesh, which is progressively denser in approaching the Ostler,T.Seyller,L.Moreschini,Y.J.Chang,M.Polini,
FermisurfaceorvanHovepoints (iftheyareclosetothe R.Asgari,A.H.MacDonald,K.Horn,andE.Rotenberg,
Fermi level). This is important for our problem because Phys. Rev.B 84, 085410 (2011).
[7] P. W. Anderson,Science 235, 1196 (1987).
of the rapid variation of the Fermi velocity near the van
[8] S.Pathak,V.B.Shenoy,andG.Baskaran,Phys.Rev.B
Hove singularities. The k-mesh is obtained as follows.
81, 085431 (2010); arXiv:0809.0244.
First define an energy scale Ω (of the order of the band-
[9] A. M. Black-Schaffer, and S. Doniach, Phys. Rev. B 75,
width), and begin with 6 equal-area triangles spanning 134512 (2007).
the BZ. Breaka specific triangle into four smaller equal- [10] J. Gonz´alez, Phys. Rev. B 78, 205431 (2008); C. A.
area ones if any eigen energy ǫ Ω in the originaltri- Lamas, D. C. Cabra, and N. Grandi, Phys. Rev. B 80,
k
angle. Then lower the energy|sc|a≤le as Ω Ω/b (b > 1) 075108(2009);D.Makogon,R.vanGelderen,R.Rold´an,
→ and C. M. Smith, Phys.Rev. B 84, 125404 (2011).
and repeat the above process recursively. The centers of
[11] E. V. Castro, A. G. Grushin, B. Valenzuela, M. A. H.
thetrianglesformthemeshpoints. Combiningthetrian-
Vozmediano,A.Cortijo,andF.deJuan,Phys.Rev.Lett.
gleareas,theyareusedintheloopintegrationsinEq.(6).
107, 106402 (2011).
In our implementation, we perform the above processes
[12] R. Nandkishore, L. Levitov, and A. Chubukov,
8-10 times, so that the last generation of triangles has a arXiv:1107.1903.
linear size oforder2−8π 2−10π, andthe centerof such [13] Tao Li, arXiv:1103.2420.
∼
triangles are sufficiently close to (or accidentally on) the [14] C. Husemann, and M. Salmhofer, Phys. Rev. B 79,
Fermi surface. The q-mesh as in Fig.8(b) used for the 195125 (2009).
[15] H. Zhai, F. Wang, and D.-H. Lee, Phys. Rev. B 80,
interactions includes all important scattering momenta:
064517 (2009).
theoriginandthehighsymmetrynestingvectors. Wede-
[16] C. Honerkamp,Phys.Rev. Lett.100, 146404 (2008).
vise a function η such that it is zero atthose important
q [17] M. Salmhofer and C. Honerkamp, Prog. Theor. Phys.
scattering momenta, and generate the mesh in a similar 105, 1 (2001).
fashionasfork,exceptthatΩbecomesanartificialscale [18] C. Honerkamp, M. Salmhofer, N. Furukawa, and T. M.
and η is used in place of ǫ . Similar q-mesh already Rice, Phys.Rev.B 63, 035109 (2001).
q k
appears in Ref.[14]. [19] F. Wang, H. Zhai, Y. Ran, A. Vishwanath, and D.-H.
Lee, Phys.Rev.Lett. 102, 047005 (2009).
[20] Wenotebysettingl=1abovedoesnotproduceindepen-
dent spatial form factors because the above three form
factors are complete for the three 1st neighbor bonds.
Hence by using f we are able to investigate, e.g.,
[1] K.S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, d1,d1′
bothp-andd-wavepairing(ofcoursetheirtotalformfac-
Y. Zhang, S. V. Dubonos, I. V. Gregorieva, and A. A.
tors have completely different spin and sublattice struc-
Firsov, Science 306, 666 (2004).
tures).
[2] ForareviewseeA.H.CastroNeto,F.Guinea,N.M.R.
[21] Translating P, C and D back to Γ1234, the divergence
Peres,K.S.NovoselovandA.K.Geim,Rev.Mod.Phys.
criterion is smaller since two form factors must be ac-
81, 109 (2009).
counted for.
[3] X.Du,I.Skachko,F.Duerr,A.Luican,andE.Y.Andrei,
[22] FanYang,HuiZhai,FaWang,Dung-HaiLee,Phys.Rev.
Nature462, 192 (2009).
B 83, 134502 (2011)
[4] K.I.Bolotin,F.Ghahari,M.D.Shulman,H.L.Stormer,
[23] M. Kiesel, C. Platt, W. Hanke, D. A. Abanin, R.
P.Kim, Nature 462, 196 (2009).
Thomale,arXiv:1109.2953.Inthispapertheelectronfill-
[5] D. A. Siegel, Cheol-Hwan Park, C. Hwang, J. Deslippe,
ing f is defined as the number of electrons per site per
A.V.Fedorov,S.G.Louie,andA.Lanzara,PNAS108, spin, which is related to our doping level as δ=2f −1.
11365 (2011); D. C. Elias, R. V. Gorbachev, A. s. May-
[24] C. Karrasch, R. Hedden, R. Peters, Th. Pruschke, K.
orov,S.V.Morozov,A.A.Zhukov,P.Blake,L.A.Pono-
Sch¨onhammer, and V.Meden, J. Phys.Cond. Matt. 20,
marenko, I. V. Grigorieva, K. S. Novoselov, F. Guinea,
345205 (2008).
and A. K.Geim, NaturePhys. 7, 701 (2011).
[6] A. L. Walter, A. Bostwick, K.-J. Jeon, F. Speck, M.
