A Tight-binding approach to electron-phonon couplings

# Rippling transition from electron-induced condensation of curvature field in graphene

## Abstract

A quantum field theory approach is applied to investigate the dynamics of flexural phonons in a metallic membrane like graphene, looking for the effects deriving from the strong interaction between the electronic excitations and elastic deformations. Relying on a self-consistent screening approximation to the phonon self-energy, we show that the theory has a critical point characterized by the vanishing of the effective bending rigidity of the membrane at low momentum. We also check that the instability in the sector of flexural phonons takes place without the development of an in-plane static distortion, which is avoided due to the significant reduction of the electron-phonon couplings for in-plane phonons at large momenta. Furthermore, we analyze the scaling properties of the many-body theory to identify the order parameter that opens up at the point of the transition. We find that the vanishing of the effective bending rigidity and the onset of a nonzero expectation value of the mean curvature are concurrent manifestations of the critical behavior. The results presented here imply that, even in the absence of tension, the theory has a critical point at which the flat geometry becomes an unstable configuration of the metallic membrane, with a condensation of the mean curvature field that may well reproduce the smooth distribution of ripples observed in free-standing graphene.

## I Introduction

During the last decade many efforts have been devoted to the investigation of graphene(1), the material consisting of a one-atom-thick layer of carbon atoms with outstanding electronic, elastic and optical properties(2). These derive to a great extent from the reduced dimensionality of the material. From the point of view of the elastic behavior, the carbon sheet has shown to be a very unusual kind of membrane, placed in a regime where anharmonic effects are very important. Regarding the electronic properties, graphene displays also unconventional features arising from the relativistic-like dispersion of the valence and conduction bands, that touch only at isolated points in the corners of the Brillouin zone.

The interplay between the electronic excitations and elastic degrees of freedom in graphene has been investigated in recent years, counting on the fact that the electron-phonon couplings should be rather large in the carbon layer(3). Most part of the studies have concentrated at first on the effect of the elastic modes on the electronic behavior. Focusing on the out-of-plane deformations (so-called flexural phonon modes), it has been found that they can have a significant impact on the electronic transport properties of the material(4); (5); (6). Considering the interaction with the substrate, it has been also shown that the corrugation induced on the carbon sheet may lead to local changes in the Fermi energy, making graphene to behave as a metallic membrane susceptible to develop charge puddles(7). There have been also other studies stressing the potential of the irregularities in the shape of the membrane to modify the band structure of the material(8).

Conversely, there have been recent works devoted to investigate the effect of the electronic degrees of freedom on the dynamics of the flexural phonons(9); (10); (11). The main aim of these investigations is the explanation of the formation of ripples in graphene. In this respect, the carbon layer is known to have corrugation, which is mainly induced in some cases by the substrate(12); (13), but ripples have been also observed in suspended samples of the material(14). In Ref. (10), it has been shown that the coupling to the electron-hole polarization may give rise to a significant renormalization of the elastic constants of the metallic membrane. For sufficiently strong deformation potential (electron-phonon coupling), it has been found that the bending rigidity undergoes a very drastic reduction, leading to a critical point and spontaneous symmetry breaking in the system of flexural phonons. A different approach has been proposed recently in Ref. (11), where the hybridization of electron-hole and flexural phonon excitations is investigated in the strong coupling regime, finding the simultaneous development of charge puddles and spontaneous breakdown of symmetry of the out-of-plane deformations.

In these studies of the dynamics of flexural phonons, the critical point is found for values of the deformation potential in the range of eV. This point marks the instability of the flat geometry of the membrane, but the relatively large values of the electron-phonon coupling oblige to check that a Kohn anomaly with the consequent in-plane static distortion does not arise before the symmetry breaking in the sector of flexural phonons. This analysis requires in turn some knowledge about the behavior of the electron-phonon couplings when going from the long-wavelength limit to the regime of large momenta.

Moreover, the symmetry-breaking pattern giving rise to ripple formation may depend on the balance between different order parameters. In Ref. (10), the condensation was found to be related to the gaussian curvature of the membrane, with a mechanism very similar to that opening a nonvanishing expectation value for a Higgs field. In that approach, the spontaneous symmetry breaking is driven by any slight amount of negative tension, being however much unclear the state of the membrane in the limit situation where that parameter vanishes. This raises the question about the dominant order parameter of the instability when the metallic membrane is not under the influence of any external perturbation.

In this paper we address these last two questions, applying self-consistent methods in the renormalization of the many-body theory of flexural phonons. The general framework goes back to the pioneering work of Ref. (15), which considered the problem of the stability of a crystalline membrane in the classical regime. For this kind of statistical system, the renormalization of elastic constants was further refined in Ref. (16) with the use of a self-consistent screening approximation (see also (17)). This computational scheme has proven to be very appropriate for the study of the many-body properties of flexural phonons(18), whose interactions account for important anharmonic effects in graphene(19); (20); (21).

Going one step beyond, we will consider the quantum statistical formulation of the many-body theory of flexural phonons, which is pertinent for a stiff membrane like graphene where quantum fluctuations may be relevant even for not very low temperatures. In this framework, we will compute the flexural phonon self-energy by means of a self-consistent screening approximation, in a theory where the effective interaction between flexural phonons encodes the coupling to the electron-hole excitations as well as to in-plane phonons. We will incorporate in particular the momentum dependence of the electron-phonon couplings for a sensible determination of the critical point of the metallic membrane, characterized by the vanishing of the effective bending rigidity at a low value of the momentum.

We will also identify the order parameter that arises naturally at the critical point. When the effective bending rigidity is drastically reduced, the flexural modes become very soft at low momenta, making the system susceptible to undergo condensation with a definite content of phonon fields. We will see that the mean curvature (the laplacian of the flexural phonon field) is the field that gets a nonvanishing expectation value at the point of the transition, playing therefore the role of order parameter. This complements the analysis of Ref. (10) showing that, even in the absence of tension, there is a natural condensation mechanism leading to the destabilization of the flat geometry of the metallic membrane. In this respect, the mean curvature is a quantity clearly constrained on the average in a planar geometry, which does not prevent however that the transition may take place with the formation of spatial domains with alternating sign of the order parameter, producing a landscape that may well reproduce the smooth distribution of ripples observed experimentally in graphene samples.

## Ii Effective action of flexural phonons

In the case of graphene considered as a metallic membrane, there is a clear hierarchy regarding the different degrees of freedom of the system. The energy scale of the electronic excitations is much higher than that of the in-plane phonons, that are found in turn at energies above those of the flexural phonons. Thus, in order to obtain an effective theory for the latter, it is pertinent to integrate out first the fields of the electron quasiparticles, carrying out the same task subsequently for the in-plane phonons.

We begin then with the system of electron quasiparticles, that we represent in terms of a set of two-component Dirac fields accounting for valley and spin degrees of freedom, coupled to the elastic deformations of the membrane. These are described by the displacement within the reference plane and the vertical shift , which build together the strain tensor

 uij=12(∂iuj+∂jui+∂ih∂jh) (1)

Assuming that the electronic correlations are weak in graphene, we may characterize the dynamics of the electron quasiparticles by means of a quadratic Dirac hamiltonian. The action including the electron-phonon coupling can be written in the form

 Se=∫dtd2rψ†n(r)(i∂t−ivF\boldmathσ⋅\boldmath∇)ψn(r)−g∫dtd2rψ†n(r)ψn(r)Tr(uij(r)) (2)

where is the Fermi velocity and stands for the so-called deformation potential (measured in energy units).

Upon integration of the Dirac fields , we get a contribution to the effective action that depends on and . While this generates a sequence of terms with increasing powers of those fields, particular attention has to be paid to the contribution proportional to the square of the strain tensor, since it gives rise to a direct renormalization of the couplings for the elastic deformations. The original action for the elastic degrees of freedom is

 Sel=12∫dtd2r(ρ(∂tu)2+ρ(∂th)2−κ0(∇2h)2−2μTr(u2ij)−λ(Tr(uij))2) (3)

where is the mass density, is the bending rigity, and stand respectively for the in-plane shear and bulk moduli. Upon integration of the fermion fields, we obtain a correction to

 ΔSel=−12g2∫dtd2r∫dt′d2r′(∂iui+12∂ih∂ih)∣∣∣r,tχ(r−r′,t−t′)(∂juj+12∂jh∂jh)∣∣∣r′,t′ (4)

where the new interaction is mediated by the electron-hole polarization .

For graphene with Dirac fermion flavors, we have the Fourier transform . As long as the typical energies of phonons (in-plane as well as out-of-plane) are much lower that those of the electronic excitations, it is a good approximation to set in for the evaluation of the effective action . This can be written then in frequency and momentum space as

 Seff = 12∫d2q(2π)2dω2π(ρω2−κ0q4)h(q,ω)h(−q,−ω) (5) +12∫d2q(2π)2dω2π(ρω2ui(q,ω)ui(−q,−ω)−2μuij(q,ω)uij(−q,−ω)−λ′(q)uii(q,ω)ujj(−q,−ω))

with , and being the Fourier transform of the strain tensor (1).

We can integrate out the in-plane phonons at this point, using a known trick that consists in separating the transverse component of the strain tensor(15). Introducing the transverse projector , we can recast the strain tensor in the form

 uij=12(∂i¯uj+∂j¯ui+P(T)ij(P(T)kl∂kh∂lh)) (6)

Given that the energy scale of the in-plane phonons is higher than that of the flexural phonons, we can take again a static approximation to study the dynamics of the latter in the background of in-plane phonons. Then, only the longitudinal part of the field couples to , and we get an expression for the effective action of the field

 Seff = 12∫d2q(2π)2dω2π(ρω2−κ0q4)h(q,ω)h(−q,−ω) (7) −18∫d2q(2π)2dω2πK(q)P(T)ij(q)hij(q,ω)P(T)kl(q)hkl(−q,−ω)

with

 K(q)=2μ+λ−g2|q|/4vF−(λ−g2|q|/4vF)22μ+λ−g2|q|/4vF (8)

and standing for the Fourier transform of .

The action (7) constitutes our starting point to study the effects of the interaction of the flexural phonons in the metallic membrane. It can be actually seen that it defines a scale invariant field theory, since its expression remains unmodified (except for the scaling of the electron-phonon coupling) under the change of variables

 ω′=sω,q′=√sq,h′(q′,ω′)=1s2h(q,ω) (9)

It makes sense therefore to investigate the scaling of the theory at the quantum level. This leads to characterize the effective energy dependence of the parameters in the action, whose low-energy behavior allows in turn to get information about the main elastic properties of the membrane in the long-wavelength limit.

## Iii Self-consistent screening approximation

For sufficiently large values of the deformation potential , the coupling function may reach negative values, driving the interaction into a regime of attraction for a certain range of the momentum . As we are going to see, this may lead to an instability in the metallic membrane, either from the decay of the flexural phonons or from a drastic renormalization of the bending rigidity. In order to claim the instability in the sector of flexural phonons, we have to make sure however that the pole in the denominator of is not hit for any allowed momentum, since that would imply a different instability pertaining to the in-plane phonons.

The last term at the right-hand-side of (8) represents actually the contribution to the coupling function from the exchange of in-plane phonons, and reaching the pole in their propagator would give rise to a large Kohn anomaly implying the formation of an in-plane static distortion. Before that happens as the value of is increased, it is possible however to find the instability in the sector of flexural phonons at a lower value of the deformation potential, signaling the spontaneous breakdown of symmetry in the metallic membrane but with regular behavior of the in-plane phonons.

We study these effects using a self-consistent screening approximation, which can be defined by the diagrammatic equations shown in Fig. 1. The propagator of the flexural phonon interaction is dressed in the form

 D−1(q,ωq)=K−1(q)+Π(q,ωq) (10)

with the polarization given in terms of the flexural phonon propagator as

 iΠ(q,ωq)=12∫d2k(2π)2dωk2π(k2−(q⋅k)2q2)2G(k,ωk)G(q−k,ωq−ωk) (11)

In turn, we use the dressed interaction to obtain the flexural phonon self-energy

 iΣ(q,ωq)=∫d2k(2π)2dωk2π(q2k2−(q⋅k)2)2|q−k|4D(q−k,ωq−ωk)G(k,ωk) (12)

which must give back the phonon propagator according to the equation

 G−1(q,ωq)=ρω2q−κ0q4+Σ(q,ωq) (13)

The self-energy in (12) becomes proportional to , which means that its main effect is the renormalization of the bending rigidity . Therefore, we may carry out the self-consistent resolution modifying the flexural phonon propagator to account for a dressed bending rigidity

 G(q,ωq)=1ρω2q−κ(q,ωq)q4+iη (14)

Regarding the interaction, it is crucial to keep the deformation potential below the upper limit at which the pole in the last term of in (8) is reached (which happens first for large momenta, of the order of the inverse of the lattice spacing). For a more precise evaluation of that regime, it is pertinent to discern the own momentum dependence of the electron-phonon couplings which, in the case of the in-plane phonons, turn out to be at large momenta much smaller than expected from the nominal value of the deformation potential approached at long wavelengths.

In this analysis, one is led to distinguish between the electron-phonon coupling for in-plane phonons appearing in the denominator at the right-hand-side of (8) and the electron-phonon coupling for flexural phonons making the direct quadratic contribution to . For the self-consistent resolution, we have refined this coupling function using the tight-binding expressions for both types of electron-phonon couplings discussed in the Appendix. This approach is simple enough to produce manageable expressions of the couplings, while capturing the nonlinear behavior at large momenta. Taking for instance for the elastic constants the values and (22), we can estimate in this way the limit set by the in-plane instability at a deformation potential eV.

In practice, we have resorted to a numerical computation to obtain self-consistently and . We have evaluated the expressions (11) and (12) by performing a Wick rotation and summing over a set of discrete Matsubara frequencies with , which corresponds to placing the theory at a finite temperature . In this framework, we have for the polarization

 Π(q,i¯¯¯ωq)=12∫d2k(2π)2(k2−(q⋅k)2q2)2T∑¯¯¯ωn1ρ¯¯¯ω2n+κ(k,i¯¯¯ωn)k41ρ(¯¯¯ωq−¯¯¯ωn)2+κ(q−k,i¯¯¯ωq−i¯¯¯ωn)(q−k)4 (15)

On the other hand, the equation for the phonon self-energy becomes a condition for the self-consistent renormalization of the bending rigidity

 κ(q,i¯¯¯ωq)=κ0+∫d2k(2π)2(k2−(^q⋅k)2)2|q−k|4T∑¯¯¯ωnD(q−k,i¯¯¯ωq−i¯¯¯ωn)ρ¯¯¯ω2n+κ(k,i¯¯¯ωn)k4 (16)

with .

The remaining integrals in Eqs. (15) and (16) can be approximated by summing over the points of a grid covering the momentum space, up to the boundary of the Brillouin zone. Usually, taking 200 divisions in both the angular variable and the modulus of the momentum leads to results with reasonable precision, assuring also the convergence of the recursive resolution of the integral equations within iterations. The behavior of the bending rigidity obtained with such a procedure is represented in Fig. 2(a), for a temperature such that meV. It becomes clear the downward renormalization of the bending rigidity for increasing values of the deformation potential, leading to significant deviations from the resolution without electron-phonon interaction for eV.

Following the decreasing trend of , there is however a value of the deformation potential at which the iterative resolution of Eqs. (11) and (16) becomes unstable for a purely real function . The convergence of the procedure can be restored by assuming that the effective bending rigidity has an imaginary part, which starts having a nonvanishing value for a deformation potential eV (for meV). This is shown in Fig. 2(b), which represents the real and imaginary parts of for different values of the deformation potential at meV.

These results reflect that, before a regime with unphysical negative values of is reached, the flexural phonons have a more conventional instability characterized by a nonzero imaginary part of the self-energy. This is the signature of a finite lifetime of the phonon modes, which can decay into collective excitations of the strongly interacting system. Indeed, it can be shown already at a perturbative level that single-phonon states are unstable towards the emission of flexural phonon pairs or more exotic hybrid excitations. We stress anyhow that the decay process described here has nonperturbative character, since it arises even after performing the Wick rotation devised to circumvent the pole in the phonon propagator. This explains that the phonon decay rate may have an unusual behavior in the limit , as observed from the plot of Fig. 2(b). It can be shown that the ratio converges however to a finite value in that limit, implying that the decay rate of the low-energy phonon modes is proportional to their energy.

The growth of the imaginary part of at low momenta is actually a consequence of the crossover to the regime dominated by thermal effects. At zero temperature, the effective momentum dependence of the bending rigidity can be obtained as a renormalization effect characteristic of the quantum theory, which leads to the low-energy behavior (neglecting the electron-phonon interaction)(23)

 κ∼|log(q)|4/7 (17)

As , however, there is always some regime, corresponding to phonon energies much smaller than the thermal energy, where the thermal effects take over, changing the behavior of the bending rigidity to a power-law

 κ∼1|q|η (18)

The divergence of at vanishing momentum can be taken as a consequence of the asymptotic behavior (18), as long as the ratio between the real and the imaginary part of the bending rigidity has a finite limit as .

Finally, for a sufficiently large value of the deformation potential, there is a point where the real part of the bending rigidity vanishes, marking the threshold beyond which it becomes negative for an increasingly wide range of frequencies and momenta. The vanishing of starts taking place at a frequency eV, and it points at the transition to a different phase of the system as one cannot make physical sense of negative values of the bending rigidity. The approach to the critical point is shown in Fig. 3 for two different values of the temperature. The critical value of the deformation potential depends slightly on , and it is placed always close but below the upper limit marking the instability of the in-plane elastic deformations.

Negative values of the real part of point at the transition to a phase that cannot be cast in terms of the original phonon fields. The vanishing of the bending rigidity at the critical point implies the appearance of very soft modes which, in turn, make the system more prone to develop some sort of spontaneous breakdown of symmetry with condensation of flexural phonons. The investigation of this effect requires however a different approach to discern the possible divergence of any of the susceptibilities of the field theory, a question that we address specifically in the next section.

## Iv Condensation of flexural phonons

In the model considered above, the effective bending rigidity has a bounce as a function of the momentum for large electron-phonon coupling. This is a consequence of the nontrivial behavior of the coupling function , which has different attractive and repulsive regimes. The attraction between flexural phonons is anyhow the actual reason for the drastic decrease observed in the bending rigidity. For the sake of clarifying the nature of the corresponding instability, it is then convenient to carry out the investigation of a simpler model where the interaction does not change its character, having an attractive potential

 K(q)=−G (19)

with constant . We note that one can make physical sense of this model as an approximate description of a graphene layer under suitably large doping. In that situation, the electron-hole susceptibility may have a leading term proportional to the chemical potential (measured from the Dirac point), leading to . It is feasible then that, at sufficiently large doping levels (with eV), the coupling obtained from (8) may be dominated by a negative term in the relevant range of strong electron-phonon coupling.

If one applies now the self-consistent screening approximation, it can be seen that the model with the coupling (19) has a critical point which is again characterized by the vanishing of the effective bending rigidity, taking place now at momentum . The self-consistent resolution (keeping the phonon interaction unrenormalized at this point) leads to a critical value of the coupling that vanishes for large values of the cutoff in the momentum integrals as

 Gc(Λ)∼√ρκ3/20log(Λ) (20)

This behavior leads to a vanishing in the continuum limit, which is consistent with the renormalization of the model that we discuss in what follows.

The simplicity of the interaction (19) makes possible to use powerful scaling methods for the characterization of the critical point. As observed at the end of Sec. II, the field theory of flexural phonons has a scale invariance that allows to encode the most important quantum corrections into the renormalization of a finite number of parameters (the bending rigidity and the coupling in our case). In the scaling approach, a key distinction is made between the original bare theory, written with the set of unrenormalized parameters and integrals running down to a microscopic cutoff, and the renormalized theory, in which one has to be able to make predictions independent of the value of the cutoff(24). The passage is done through the redefinition of the parameters

 κ0 = ZκκR (21) K = ZKKR (22)

If the theory is well-behaved (renormalizable), it must be possible to absorb all the dependences on the microscopic cutoff into the factors and , allowing to express the observables of the theory just in terms of cutoff-independent parameters and .

To see how this program works in our model to large orders in perturbation theory, it is convenient to implement a dimensional regularization of the integrals instead of using a more conventional high-energy cutoff in momentum space. We start then by writing the action with a slight deviation of the spatial dimension from to , in such a way that (7) becomes

 Seff = 12∫dDq(2π)Ddω2π(ρω2−ZκκRq4)h(q,ω)h(−q,−ω) (23) −18ZKKRμϵ∫dDq(2π)Ddω2πP(T)ij(q)hij(q,ω)P(T)kl(q)hkl(−q,−ω)

being an auxiliary momentum scale. In this approach, the factors arising in the theory regularized with a momentum cutoff are replaced by poles, which are supposed to be absorbed into an appropriate definition of the factors and .

The critical point of the theory has in particular a reflection in the renormalization of the bending rigidity . This can be analyzed from the phonon self-energy corrections or, equivalently, studying the renormalization of the composite operator

 Φ(r,t)=∇2h(r,t)∇2h(r,t) (24)

The vertex corresponding to the operator (24) is given by

 Γ(q,ωq;k,ωk)=⟨Φ(q,ωq)h(k−q,ωk−ωq)h(−k,−ωk)⟩1PI (25)

where 1PI denotes that we take the one-particle irreducible part of the correlator. It can be easily seen, from inspection of the respective diagrammatic contributions, that the vertex (25) is related to the phonon self-energy through the equation

 ∂Σ(k,ωk)∂κ0=Γ(0,0;k,ωk) (26)

In order to study the scaling of the bending rigidity, we can therefore focus on , which is in general easier to analyze in comparison to .

The bare vertex (25) is not a finite quantity in the physical limit , as it can be seen from (26) that it inherits the cutoff dependence of . We can rewrite that equation in the form

 ∂Σ(k,ωk)∂κR=∂κ0∂κRΓ(0,0;k,ωk) (27)

The left-hand-side of Eq. (27) must be a finite quantity in the limit , since it can be written as a function only of the cutoff-independent bending rigidity . This means that it must be possible to absorb all the poles of into a multiplicative renormalization of the vertex

 Γren=ZΓΓ (28)

with

 ZΓ=∂κ0∂κR (29)

leading to a finite vertex in the physical limit .

The identification of a possible critical point in the renormalized vertex requires the adoption of a nonperturbative approach to the many-body theory. For this purpose, we may resort to a most affordable ladder approximation for the computation of , which is encoded in the self-consistent diagrammatic equation represented in Fig. 4. This can be expressed as

 iΓ(0,0;k,ωk)=ik4+Gμϵ∫dDp(2π)2dωp2π(k2p2−(k⋅p)2)2|k−p|4Γ(0,0;p,ωp)1(ρω2p−κRp4+iη)2 (30)

At this stage, the ladder approximation amounts to a sum of ladder diagrams that do not contain themselves any phonon self-energy loop. Therefore, it is safe to take the renormalized bending rigidity in the phonon propagators in (30). The solution of the self-consistent equation can be expressed as a power series in the dimensionless coupling

 Γ(0,0;k,ωk)=k4+k4∞∑n=1˜Gnμnϵsn|k|nϵ (31)

The coefficients diverge in the limit with powers of , which can be absorbed into the multiplicative redefinition (28) with a renormalization factor that has the general structure

 ZΓ=1+∞∑i=1ci(˜G)ϵi (32)

The factor contains important information about the scaling of the renormalized vertex, which develops an anomalous dependence on the momentum scale leading to the long-wavelength behavior

 Γren∼(μ|k|)γ (33)

The original vertex , as any quantity in the bare theory, does not know about the auxiliary scale , so that the anomalous exponent can be obtained as(24)

 γ=μZΓ∂ZΓ∂μ (34)

On the other hand, the own independence of (31) on requires the implicit dependence of on the auxiliary scale given by

 μddμ(˜Gμϵ)=0 (35)

Thus we get

 γ=−ϵ˜GZΓdZΓd˜G (36)

In principle, multiple poles from may contribute to the right-hand-side of (36), but they can cancel out provided that the hierarchy of equations

 dd˜Gci+1(˜G)−ci(˜G)dd˜Gc1(˜G)=0 (37)

is satisfied. In that case, the anomalous dimension is given by(25)

 γ=−˜Gdd˜Gc1(˜G) (38)

making it only dependent on the coupling .

It can be shown that the field theory of flexural phonons fulfills indeed the conditions (37) in the ladder approximation. Successive terms of the power series in (31) can be obtained recursively by using the formula

 Gμϵ∫dDp(2π)2d¯¯¯ωp2π(k2p2−(k⋅p)2)2|k−p|4p4μnϵ|p|nϵ1(ρ¯¯¯ω2p+κRp4)2=k4˜Gμ(n+1)ϵAn(ϵ)|k|(n+1)ϵ (39)

with

 An(ϵ)=364π(4π)ϵ/2Γ(n+12ϵ)Γ(2−n+12ϵ)Γ(1−12ϵ)Γ(1+n2ϵ)Γ(3−n+22ϵ) (40)

We get in this way the relation

 sn+1=An(ϵ)sn (41)

The different poles arising from the coefficients can be then canceled out in the renormalized vertex by adjusting appropriately the residues in (32). We can compute analytically for instance the first perturbative orders

 c1(˜G) = −λ˜G−18λ2˜G2−124λ3˜G3−5256λ4˜G4−7640λ5˜G5+… (42) c2(˜G) = 12λ2˜G2+18λ3˜G3+19384λ4˜G4+19768λ5˜G5+… (43) c3(˜G) = −16λ3˜G3−116λ4˜G4−11384λ5˜G5+… (44) c4(˜G) = 124λ4˜G4+148λ5˜G5+… (45) c5(˜G) = −1120λ5˜G5+… (46)

with . It can be seen that the analytic expressions in (42)-(46) satisfy the conditions (37). A most important property is also that the residues do not contain any nonlocal dependence (in fact any dependence) on the momentum of the vertex, which guarantees that the renormalization only involves the redefinition of a finite number of local operators in the field theory.

Furthermore, the computation of the residues can be carried out numerically to much larger orders in perturbation theory, which allows to inspect the behavior of in the nonperturbative regime. Thus, we have been able to evaluate the residues to order , finding that the terms in each power series in the coupling approach a geometric progression. We have paid particular attention to the residue , which has an expansion

 c1(˜G)=∞∑n=1c(n)1λn˜Gn (47)

with coefficients that approach a constant value at large . The result of their numerical computation is represented in Fig. 5. A very precise fit of the dependence of the coefficients on the order can be indeed achieved by assuming the scaling

 c(n)1c(n−1)1=c+c′n+c′′n2+c′′′n3+… (48)

Using the results to order , we find that , and values of that are at most of order . We may conclude therefore that the series (47) has a finite radius of convergence given by

 ˜Gc=1cλ=64π3 (49)

The singularity of the residue at has important consequences from the physical point of view, since it implies the divergence of the anomalous dimension of the vertex according to (38). In fact, any expectation value involving the operator (24) needs to be made cutoff-independent multiplying it by the factor , which means that the anomalous exponent has to contribute to the scaling dimension of the correlator. The divergence of the susceptibility of the field implies in particular that the composite operator must have a nonvanishing expectation value at , that is

 ⟨|∇2h(r,t)|2⟩≠0 (50)

It becomes clear that represents a critical point of the many-body theory, with the pattern of symmetry breaking dictated by the condensation of the field .

It is now easy to see that the critical coupling must correspond to the singular behavior characterized by the vanishing of the effective bending rigidity in the self-consistent screening approximation. From Eq. (29) we have

 Zκ+∂Zκ∂κRκR=ZΓ (51)

We infer that must also diverge at the critical coupling . If we think of the relation between and as giving the correspondence between the values of the bending rigidity at short and long wavelengths, we conclude that any finite value of at the microscopic scale must lead to a vanishing effective bending rigidity when the coupling is set at . This argument highlights again that the origin of the critical behavior is the divergence of at the critical interaction strength.

Regarding the precise location of the critical point, it has to be noted that the value given in (49) does not admit a direct comparison with the result that is obtained in the self-consistent screening approximation. This includes the iterated effect of phonon-loop corrections in the own phonon self-energy, while these corrections are absent in our renormalization of the vertex. This can be certainly refined to incorporate the effect of phonon self-energy insertions in the ladder approximation. It can be shown that, including for instance one-loop phonon corrections in our computation of the vertex, the value of in (48) increases from 1 to 1.8. In general, iterating more and more phonon-loop corrections to the vertex leads to a significant reduction of the critical coupling with respect to the value quoted in (49). The overall picture is then that the critical point is anyhow preserved after taking into account the phonon self-energy corrections, and that the behavior (20) may be consistent with a nonvanishing critical value for the dimensionless renormalized coupling, as this is defined in terms of the effective bending rigidity that must vanish at the critical point.

## V Discussion

We have studied the many-body theory of flexural phonons in a metallic membrane like graphene, where the interaction between the electronic excitations and elastic deformations may induce important effects. We have shown that, for sufficiently strong electron-phonon interaction, the theory has a critical point characterized by the vanishing of the effective bending rigidity at long wavelengths. In this investigation, we have relied on a self-consistent screening approximation to evaluate the flexural phonon self-energy, while refining the analysis to incorporate the main features of the momentum dependence of the electron-phonon couplings.

We have seen that the effective bending rigidity has in general a bounce as a function of momentum, which is a consequence of the interplay between the electron-induced attraction and the natural repulsive interaction between flexural phonons. This investigation has clarified that the instability in the sector of flexural phonons takes place without the development of an in-plane static distortion of the lattice, which is avoided due to the significant reduction of the electron-phonon couplings for in-plane phonons at large momenta.

We have applied a renormalization group analysis to identify the order parameter characterizing the new phase marked by the vanishing of the effective bending rigidity. In a model with constant attractive interaction between flexural phonons, we have shown that there is a close relation between the renormalization factors of the bending rigidity and the mean curvature field . This implies that the vanishing of and the onset of a nonzero expectation value for the mean curvature can be viewed as different manifestations of the same effect. We may therefore conclude that the critical point we have studied corresponds to a rippling transition at which the flat geometry becomes an unstable configuration of the metallic membrane.

As at the critical point, the weight of the kinetic energy of the membrane vanishes, which makes it more susceptible to an external influence and, in particular, to tension exerted by external forces. This is the route to symmetry breaking that was actually investigated in Ref. (10), finding that any slight amount of negative tension may turn the flat geometry unstable, inducing a non-vanishing expectation value of the field . In that analysis, there was however no conclusive evidence of the spontaneous out-of-plane distortion of the membrane in the limit situation with vanishing tension. The results presented here imply that, in the absence of tension, the instability at vanishing induces anyhow the condensation of the mean curvature field . The same conclusion applies to the case of very small positive or negative tension, when leading to modulations of with a period larger than the critical length scale marking the vanishing of the bending rigidity, since the presence of such a small perturbation does not invalidate in that case the results of Sec. III.

It has to be also pointed out that the vanishing of leads to a regime of strong coupling which could trigger the divergence of a different susceptibility apart from that of the mean curvature . In this respect, an alternative possibility corresponds to the other scale-invariant operator in the action (23), built from the field . The condensation of this field falls again in the route to symmetry breaking pursued in Ref. (10), and it has been also recently investigated in Ref. (11). As , the propagator of can develop a pole in a RPA framework, which has to be interpreted as the signal of an out-of-plane static distortion given by the modulation of that field.

The strong-coupling character of the regime with does not allow our renormalization group method to discern the prevalence of any of the two symmetry-breaking patterns driven respectively by and . We note anyhow that the effects deriving from the condensation of each of these fields have to be very different in the two cases. The condensation of actually implies a non-vanishing expectation value of which, being a vector, may induce vorticity in the geometry of the membrane. On the other hand, the condensation of the field will lead in general to a smooth configuration of the out-of-plane distortion. In a real material both patterns of symmetry breaking may be present. A careful analysis of experimental data of the curvature distributions may be however required to disentangle the two contributions, for the purpose of assessing the relative significance of each of them in the geometry of the real metallic membrane.

## Acknowledgments

I thank F. Guinea and P. San-José for very useful discussions. The financial support from MICINN (Spain) through grant FIS2011-23713 is gratefully acknowledged.

*

## Appendix A Tight-binding approach to electron-phonon couplings

In the tight-binding approach, the expression of the electron-phonon couplings can be obtained from the variation of the potential energy of the electron system under the displacement of the atoms in the underlying lattice. We take as starting point the expectation value of the periodic potential formed by superposition of atomic potentials centered at the different lattice sites , labelled by unit cell and the position within the unit cell,

 V(r)=∑nαv(r−Rnα) (52)

We can then evaluate the potential energy for electronic states with well-defined momentum , represented by a Bloch superposition of wavefunctions for the atomic orbitals

 Ψk(r)=1√N∑αψα(k)∑neik⋅Rnαϕ(r−Rnα) (53)

We end up in this way with the matrix element

 Ek′,k≡⟨Ψk′|V|Ψk⟩=1N∑α,α′ψ∗α′(k′)ψα(k)∑n,n′ei(−k′⋅Rn′α′+k⋅Rnα)v(n′α′,nα) (54)

written in terms of the atomic matrix element

 v(n′α′,nα)=∫d2rϕ(r−Rn′α′)∑n′′α′′v(r−Rn′′α′′)ϕ(r−Rnα) (55)

Under a vibration of the lattice, the displacement of the atoms from their equilibrium positions gives rise to a so-called deformation potential. An approximate expression for the variation of can be obtained by truncating the sums in (54) and (55) to the terms in which the points and have some degree of proximity. In this respect, one usually makes a distinction between the on-site and the off-site deformation potential, corresponding respectively to taking on one hand, and either or on the other. In general, the component of the electron-phonon coupling arising from the on-site deformation potential is much larger than that from the off-site counterpart. For this reason, we concentrate in what follows on the on-site contribution.

Focusing then on the case with , we can write the matrix element (55) as

 v(nα,nα)=∫d2rϕ(r)∑n′′α′′v(r−Rn′′α′′+Rnα)ϕ(r) (56)

The deviation of the atoms from the equilibrium position, , leads to the on-site deformation potential

 Δvon = −∫d2rϕ(r)∑n′′α′′\boldmath∇v(r−Rn′′α′′+Rnα)⋅(d(Rn′′α′′)−d(Rnα))ϕ(r) (57) +12∫d2rϕ(r)∑n′′α′′(d(Rn′′α′′)−d(Rnα))⋅\boldmath∇\boldmath∇v(r−Rn′′α′′+Rnα)⋅(d(Rn′′α′′)−d(Rnα))ϕ(r)+…

We are interested in particular in the situation where the displacement arises from a phonon of the lattice. Then, such a field can be represented in terms of the polarization of the given phonon mode as

 d(Rnα)=eiq⋅Rnαeα(q) (58)

Inserting this expression into (57) and the atomic matrix element back into (54), we arrive at the electron-phonon couplings for the respective interactions with one and two phonon fields:

 g(1)α,α(k;q) = ψ∗α(k+q)ψα(k)∑n′α′(eiq⋅un′α′eα′(q)−eα(q))⋅n(un′α′) (59) g(2)α,α(k;q,q′) = ψ∗α(k+q+q′)ψα(k) (60) ∑n′α′(eiq⋅un′α′eα′(q)−eα(q))m(un′α′)(eiq′⋅un′α′eα′(q′)−eα(q′))

where and

 na(u) = ∫d2rϕ(r)∇av(r−u)ϕ(r) mab(u) = ∫d2rϕ(r)∇a∇bv(r−u)ϕ(r) (61)

We can now specialize the expressions (59) and (60) to the instances which are relevant for our study. First, there is the case of the in-plane acoustic longitudinal phonons, for which the polarization is simply given by

 eα(q)=q|q| (62)

One can easily see that the vector has to point in the same direction of the vector (26). Then, we have for this type of phonons

 g(1)α,α(k;q)∼ψ∗α(k+q)ψα(k)∑n′α′(eiq⋅un′α′−1)q⋅un′α′s(|un′α′|) (63)

In practice, one has to restrict the sum in (63) to a certain set of neighbors of the point when computing the electron-phonon coupling. At small , has a leading linear dependence on the momentum, which is in agreement with the form of the coupling in the continuum limit shown in (2). As the momentum increases, however, the electron-phonon coupling deviates significantly from the linear behavior, leading to a drastic reduction from what one could expect according to the linear small- dependence. This is illustrated in Fig. 6(a), where we have plotted the coupling (63) including up to third nearest-neighbor contributions. The decreasing trend of at large is consistent with the very small values predicted by detailed calculations of the coupling for in-plane acoustic phonons at the point of the Brillouin zone(27).

We shift now to the case of the out-of-plane acoustic phonons. These have a polarization

 eα(q)=z (64)

where stands for the unit vector in the direction transverse to the reference plane of the lattice. Then, it is clear that vanishes in the case of a free-standing membrane, as is orthogonal to , and we have to consider the two-phonon coupling . The tensor has certainly a nonvanishing component, so that the evaluation of (60) leads in this case to

 g(2)α,α(k;q,q′)∼ψ∗α(k+