Gutzwiller Charge Phase Diagram of Cuprates, including Electron-Phonon Coupling Effects

Besides significant electronic correlations, high-temperature superconductors also show a strong coupling of electrons to a number of lattice modes. Combined with the experimental detection of electronic inhomogeneities and ordering phenomena in many high-T_c compounds, these features raise the question as to what extent phonons are involved in the associated instabilities. Here we address this problem based on the Hubbard model including a coupling to phonons in order to capture several salient features of the phase diagram of hole-doped cuprates. Charge degrees of freedom, which are suppressed by the large Hubbard U near half-filling, are found to become active at a fairly low doping level. We find that possible charge order is mainly driven by Fermi surface nesting, with competition between a near-(pi,pi) order at low doping and antinodal nesting at higher doping, very similar to the momentum structure of magnetic fluctuations. The resulting nesting vectors are generally consistent with photoemission and tunneling observations, evidence for charge density wave (CDW) order in YBa_2Cu_3O_{7-delta} including Kohn anomalies, and suggestions of competition between one- and two-q-vector nesting.


Introduction
The existence of charge density wave (CDW) order is now well established for a large class of high-temperature superconductor materials.Starting from the pioneering studies in lanthanum cuprates [1], recent resonant and hard x-ray scattering data have revealed CDW modulations also in YBCO [2,3,4,5,6] and Bi2201 [7] compounds.However, whereas in lanthanum cuprates a concomitant spin scattering with twice the period of the CDW suggests the formation of charge-spin stripes [8], there seems to be no apparent relationship between the two periodicities in the YBCO and Bi2201 materials.Instead, the analysis of resonant x-ray scattering and angle-resolved photoemission spectra indicates that the nesting properties of the underlying Fermi surface (FS) or 'Fermi arcs' are at play.This has led to proposals [9,10] in which CDW formation related to FS features is driven by magnetic interactions.In this paper we examine the simpler possibility that electron-phonon interactions could play a role in the formation of the observed CDW modulations.Earlier proposals of phonon induced CDWs [11] did not rely on nesting features.
We adopt an intermediate coupling approach in our analysis in this study.[12] This is justified by recent studies [13,14,15,16,17,18], which have indicated that correlations in the cuprates are not as strong as initially believed, and that cuprates fall, instead, in an intermediate coupling regime, with 6 ≤ U/t ≤ 9, where U is the effective Hubbard interaction and t is the nearest-neighbor hopping parameter.We have shown [12] that intermediate coupling corresponds approximately to 4 ≤ U/t ≤ 13.6 = U BR /t, where U BR is the mean-field Brinkman-Rice energy where double occupancy vanishes.[19,14] In this regime competing phase transitions often evolve from Stoner instabilities, which can be described by Hartree-Fock (HF) or, better, Gutzwiller approximation (GA) based calculations.[20,21] Furthermore, at large doping cuprates behave as Fermi liquids, so that one can hope to obtain information on ordered phases by studying the instabilities that disrupt the Fermi liquid behavior, provided the correlated physics is included in the analysis.We have shown that in the weak and intermediate coupling range, peaks in the bare susceptibility of 2D materials, which determine magnetic instabilities, form a map of the FS, and the dominant instabilities are generally related to the double nesting features [14,22], where two branches of the map cross.In particular, the T = 0 magnetic phase diagram of the cuprates was derived using the time-dependent GA (TDGA).[22] In the electron-doped cuprates, the magnetic phase diagram is dominated at all dopings by a commensurate (π, π) antiferromagnetic (AFM) order.[22] In contrast, for the hole doped cuprates, we find a wide doping range over which the magnetic order is incommensurate.
Given the large onsite Coulomb repulsion (Hubbard U), magnetic order should be favored near half-filling, but for larger hole doping the experimental evidence is more consistent with incommensurate charge order as noted above.In this paper, we apply the TDGA technique to examine the charge response, including effects of finite electron-phonon coupling of Su-Schrieffer-Heeger form [23].When phonons are included, we find charge density wave (CDW) phases with nesting vectors similar to those for the magnetic instabilities, which are controlled by a generalized Stoner criterion and the double nesting features in the susceptibility.We present the full evolution with doping of the leading charge-phonon instabilities for several families of cuprates, including La 2−x A x CuO 4+δ , A = Sr (LSCO) or Ba (LBCO) and Bi 2 Sr 2 CuO 6 (Bi2201).
The magnetic [charge] instabilities are usually determined by zeroes of the Stoner denominator, [22] 1 − [+]U eff (q)χ 0 (q, ω = 0). ( In a HF plus random-phase approximation (RPA) calculation U eff (q) would simply be the Hubbard U and χ 0 (q, ω = 0) the susceptibility for local magnetic [charge] fluctuations.But in the TDGA the situation is more complex since local and transitive fluctuations are coupled so that U eff (q) depends on the corresponding susceptibilities and the associated coupling constants.Here, by local we mean that the phonon modulates the on-site energies, as in the Holstein model, while transitive means that the hopping parameters are modified, as in the Su-Schrieffer-Heeger model.Thus, the leading HF+RPA instability is simply associated with the maximum of the bare susceptibility χ 0M = max q χ 0 (q, 0), while the leading Gutzwiller instability can be shifted by the q-dependence of U GA (q).It is clear that such instabilites cannot arise for a purely local electron-phonon coupling since local charge fluctuations are significantly suppressed in the presence of correlations.However, the situation is different for the coupling of phonons to transitive fluctuations, [24] which can induce a CDW in a system with sizeable electronic correlations for moderate values of the electron-phonon coupling.
This paper is organized as follows.Section 2 describes our model system, and focuses on the electron-phonon coupling related aspects.Section 3 presents results for the renormalized phonon dispersions and the resulting charge phase diagrams for various types of high-T c materials.In Section 4 we compare our results with experiments on the cuprates, in particular we examine evidence for a crossover between charge and magnetic instabilities in the underdoped compounds, the doping dependence of nesting vectors, and their relationship to Kohn anomalies and soft phonons.We conclude our discussion in Section 5.In Appendix A we describe our TDGA formalism for CDWs, and in Appendix B we discuss one-vs two-q nesting.Further applications of the model are briefly considered in the Supplementary Online Materials (SOM), and include extensions to photoemission and tunneling studies [Section SOM1], a search for purely electronic CDWs [Section SOM2], and a discussion on stacks of Kohn anomalies [Section SOM3].

Model and Formalism
Our investigations are based on the following Hamiltonian i,σ destroys (creates) an electron on lattice site R i and n i,σ = c † i,σ c i,σ .We incorporate band structure effects as in our earlier magnetic phase diagram calculations [22], by using for the hopping parameters t ij a one-band tight-binding fit to the local density approximation (LDA) dispersion [25] for the single-layer cuprates LSCO, Nd 2−x Ce x CuO 4 (NCCO), and Bi2201.For convenience, the band parameters are listed in Table I, and the dispersion is given by where and α is an integer.Here k z dispersion is neglected, approximating the cuprates as 2D.[26,27] Interaction effects are incorporated via the TDGA, which is used for deriving the charge susceptibility.Details of the formalism are discussed in Refs.[24,28,29], and are also summarized in Appendix A. It is obvious that phonons which couple to the local charge density have only a negligible effect on the associated electronic fluctuations which are strongly suppressed by the onsite correlation U.The situation is different for phonons which couple to transitive fluctuations as has been shown in Ref. [24].For this reason, the electronphonon coupling H el−ph in Eq. 2 is described by a generic Su-Schrieffer-Heeger [23] phonon model consisting of only longitudinal and [in-plane] transverse acoustic branches, ionic mass M, and electron-phonon coupling, which arises through a modulation of the hopping integral δt.The key ingredient is that the phonons modulate the hopping parameters with 'longitudinal' and transverse modulations.Here, longitudinal refers to modulation δa along the phonon propagation direction, and transverse to those at right angles to the propagation direction.The corresponding operator is: where u µ i denotes the displacement of the atom at site R i in direction µ, and with δr ij = (u µ j − u µ i ), and the distance between atoms i and j is r ij = R i − R j .Finally, the phonon part is given by which can be diagonalized to yield the bare phonon frequencies [Ω 0 qµ ] 2 = 2(K/M) µ (2 − cos(q x a) − cos(q y a)).Here K µ and M µ denote the effective spring constant and ionic mass for the longitudinal and in-plane transverse (µ = L [T]) acoustic mode, respectively.Since the hopping parameters are labeled t for the nearest-neighbor, t ′ for the second-nearest-neighbor, etc., we label the corresponding γ's as γ, γ ′ , and so on.The nearest-neighbor electron-phonon coupling constant thus becomes λ ep = 4 N(0)g 2 /K, [30,31,32,33,34] with average density-of-states N (0) ∼ 2/8t, and electron-phonon coupling g = γt/a, or Note that λ ep is independent of the ionic mass.
In reality, a strong modulation of δt can be produced by several phonons, including those involving motion of oxygen atoms perpendicular to the CuO 2 planes [35,36].Hence, we assume the same bare acoustic frequencies for all cuprates, adjusted to approximate oxygen modes in undoped La 2 CuO 4 [as a generic single-layer cuprate], [37] which gives (K/M) LA = 2(K/M) T A ≡ [Ω 0 L ] 2 , taking Ω 0 L = 12.4 meV, and M is the oxygen mass (where LA[TA] = longitudinal [transverse] acoustic mode).For these parameters the bare dispersion is shown in Figs. 1 and 2 (dashed lines) for a selected cut through the Brillouin zone.
Thus, the model is completely specified in terms of the coupling constant λ ep , Eq. 8, which in turn is known up to the hopping coefficient γ (where t ∼ 1/r γ ).For Bi2201, with t = 435 meV, λ ep0 = λ ep /γ 2 = 0.047.Although the γs have proven difficult to calculate, [38,39] in the large distance limit, direct wave function overlap on different atoms falls off exponentially with r, and in the cuprates, all hopping parameters except the closest Cu − O hopping t CuO are dominated by an indirect chain of hoppings, [40] yielding for the nearest-neighbor [Cu-Cu] hopping t ≃ t 2 CuO /∆, where ∆ is the on-site energy difference between Cu and O. Now t CuO ∼ r −γ CuO , with γ CuO ≃ 3 − 3.5.[38,39] The problem is, what is ∆, and how does it vary with r?We note the following: (1) ∆ remains finite for infinite separation, suggesting a weak r dependence; (2) In a strongly correlated system, the Cu on-site energy has a contribution from U. The Cu−d electrons in the anti-bonding CuO 2 band are mostly second electrons on each site, so that ∆ would be dominated by the U-term; (3) We expect U to decrease with decreasing r, due to enhanced screening effects.Thus, ∆ should probably decrease as r decreases, suggesting γ > 6 − 7. Thus, we estimate that γ lies in the range between γ CuO ∼ 3.25 and 10.In the present paper, we assume γ ′ = γ ′′ = ... = 0, unless noted otherwise.
From Eq. 8 we estimate that λ ep lies between 0.57 (using γ = 3.25) and 5.4 (for γ = 10).Note that the lower estimate is comparable to values in the literature, assuming a modest anisotropy.Recent linear response calculations have found Brillouin zone-averaged coupling strengths λ ave of ∼ 0.4 for optimally doped LSCO [41] and Ca 0.27 Sr 0.63 CuO 2 , [42] and 0.27 for YBa 2 Cu 3 O 7 .[43,44] However, λ has a strong momentum dependence, for instance, the nodal value is considerably smaller: in LSCO, λ nodal = 0.14 -0.22 at optimal doping, and 0.14 -0.20 in the overdoped regime.[41] Since nodal electrons dominate transport [15], this accounts for the smaller λs estimated from transport measurements.Further, several calculations suggest that correlation effects can enhance the anisotropy of electron-phonon coupling, generally leading to a larger value for AN nesting [28,45,46,47].Finally, using the correct density-of-states rather than the average N (0) will further increase λ ep in the doping regime near the Van Hove singularity.Thus our lower estimate for λ AN N is likely to be on the conservative side.The electron-phonon interaction renormalizes the frequencies of bare phonons as well as the charge response of the electrons.In the presence of intermediate electronic correlations, both the phonon propagator as well as the charge susceptibility can be conveniently evaluated via the time-dependent Gutzwiller approximation (TDGA) as outlined in Appendix A.

Renormalized phonon dispersions
The TDGA charge correlations induce a screening of the phonons and thus renormalize the bare dispersion Ω 0qµ according to A detailed derivation of the correction to the elastic spring constant δK µµ within the TDGA is given in Appendix A (cf. also Ref. [36]).
There is a close connection between the present results and the earlier magnetic results in that the instability is controlled by an effective Stoner criterion.To see this, we write Then, if for some q, U ef f,q χ 0q = 1, the corresponding phonon frequency will vanish, leading to an instability.Figure 1 compares the bare and renormalized LA and [in-plane] TA phonon frequencies in Bi2201.Although the modes along the (π, 0) → (π, π)-branch are mixed, these are labeled as being predominantly longitudinal or transverse.The sharp dips in the dressed frequencies are caused by peaks in the bare susceptibility associated with FS nesting.Each peak in χ 0q leads to a prominent Kohn anomaly in the phonon spectrum, which can lead to an instability if the renormalized Ω 2 qµ becomes negative.By comparing the present results with earlier calculations for magnetic stripes, [22] we find that the 0 Ω (meV) instabilities fall at nearly the same q-values for both kinds of stripes as a function of doping, being controlled by the same FS nesting.However, the relative strengths of the instabilities can be modulated by the electron-phonon coupling.Additionally, the detailed analysis of the instabilities in the charge sector reveals that they tend to favor lower-q values due to the momentum structure of Ω 0qµ .Thus, in the magnetic phase diagram an instability of the Fermi liquid towards vertical incommensurate order near (π, π) [(π − δ, π)] was found at low hole doping, x ≤ 0.16, while an instability towards a diagonal incommensurate phase near Γ [(δ, δ)] was found from x = 0.16 to x = x V HS = 0.42, where x V HS is the doping where the Fermi level crosses the Van Hove singularity (VHS).The latter instability is associated with nesting of the flat FS sections in the antinodal parts of the FS, and hence is referred to as antinodal nesting (ANN).In contrast, we will refer to the former phase as near nodal nesting (NNN).Both of these instabilities are controlled by double nesting (simultaneous nesting of two sections of FS at the same q).We find that the same instabilities dominate the charge phase diagram, but that additional single nesting charge instabilities start to become competitive because they correspond to lower Ω 0qµ .Near (π, π) there are both double nesting vertical [(π − δ, π)] (feature 2 in Fig. 1(a)) and single nesting diagonal [(π − δ, π − δ)] instabilities (feature 3).Likewise, the ANN instabilities can be either double nesting diagonal [(δ, δ)] (feature 4) or single nesting vertical [(δ, 0)] (feature 1).Proximity to the Γ-point tends to favor ANN nesting over NNN, particularly in Bi2201, where t ′ is larger and the (π, π)-plateau is weaker.Thus in Bi2201 the near-(π, π) Kohn anomalies are always subdominant, except for small U at extremely low doping.On the other hand, the diagonal charge ANN anomaly is unstable over the  full doping range x ≤ 0.4.While for the present choice of electron-phonon couplings the diagonal instability is dominant, for other choices the vertical ANN instability wins out at moderate doping (x < 0.2).Note in Fig. 1 that the q vector of the leading instability shifts toward Γ as doping increases toward x V HS , as found previously for magnetic instabilities.While there are Kohn anomalies in both phonon branches, for all dopings, the ANN instabilities are in the LA phonon branch.
Fig. 2 presents phonon dispersion maps for LSCO.Unlike the magnetic phase diagram, LSCO behaves similarly to the other cuprates, despite the smaller value of t ′ .With doping, there is a competition between a predominantly TA phonon soft mode near (π, π) and an LA soft mode near Γ.For low (x < 0.05) or very high (x ≥ 0.24) doping the vertical TA instability at (π, π − δ) is dominant, while at intermediate x the diagonal LA instability closer to Γ dominates.

Charge phase diagrams
Figs. 1 and 2 show that there is competition between Mott physics and strong electronphonon coupling.At each doping and fixed electron-phonon coupling λ ep there is a critical U c (x) such that the charge order (CO) phase exists only for values of U < U c .[For the cuprates, U ∼ 0.6U BR .]As U → U c , generally the q at the instability is sharply defined, and mainly determined by peaks in the bare susceptibility χ0q .Therefore, the corresponding momenta match the nesting curves introduced in Ref. [22] for magnetic instabilities, for which analytic formulas are available.
The phase diagram can also be described in terms of a critical electron-phonon coupling λ epc vs x for fixed U, as in Fig. 3.As doping changes, the threshold q-vector varies, and its symmetry can also change, as denoted by the different symbols in the figure.These phase diagrams show a close resemblance to the magnetic phase diagrams, including the dominant q-vectors.However, there are some notable differences.Whereas U drives magnetic instabilities (in the Stoner denominator), U acts to oppose charge instabilities, which are instead driven by the electron-phonon coupling λ ep .The dotted line in Fig. 3  We briefly comment on the overall doping dependence of the phase diagram.The minimum value of λ epc occurs when the VHS lies at the Fermi energy, x = x V HS .Near half filling, a large U leads to strong screening of charge excitations [large U GA ], which rapidly quenches the e-ph instability.At U = U BR and x → 0, λ epc → ∞.As doping increases, U GA decreases and the CDW nesting effect takes over, so that near x V HS , the e-ph instability is nearly U-independent.For x > x V HS , there is a residual nesting effect, now largest near (π, π), but this falls off rapidly for higher x.

Hubbard index
The current results, especially those of Fig. 3, provide insight into the role of U in suppressing transitive charge fluctuations.The regime U ≥ U BR is particularly important as Hartree-Fock calculations have difficulty capturing the underlying physics.As hole doping increases, electrons can hop around without causing double occupancy, so that at some point U becomes unimportant in suppressing the CO.This is not captured in HF, where the Stoner denominator remains 1 + Uχ at all dopings.In contrast, the GA solution shows a clear crossover, Fig. 3, to a regime where the critical electronphonon coupling becomes independent of U.Here we quantify this crossover in terms of a 'Hubbard index' defined by which measures how sensitive CO is to U at the Brinkman-Rice energy, and it may be approximated as for the two highest Us in Fig. 3.While this should be accurate at large x, at x = 0, λ c diverges, so that δλ c /λ c ∼ 1, while δU/U → 0, and H → ∞ -that is, the undoped cuprates become incompressible as U → U BR .In contrast, our approximate ∆U/U = 0.25, so that H ≤ 4. Bearing that in mind, we show the approximate H in Fig. 3(c) as filled circles (Bi2201) or open squares (LSCO).Remarkably, for both materials, H ≃ 0.1/x (solid line).This is quite different from the Hartree-Fock expectation, H ∼ constant, dotted line.If one looks closely at the data, both curves show hints of a plateau for 0.1 ≤ x ≤ 0.2, but for larger x, U has an anomalously small effect on CO.While U = U BR has a well defined meaning only at x = 0 as a crossover, or as a phase transition in mean-field or infinite dimensions, it still represents a regime of strong suppression of double occupancy D. We find that when U = U BR , D is approximately while the renormalization function z 0 is which is somewhat larger than its D = 0 value 2x/(1 + x).

Charge-Magnetic Crossover in Underdoped Cuprates
As noted in Sections 3.1 and 3.2, for most cuprates there should be two kinds of competing density wave orders, with a transition as a function of doping between an incommensurate phase with q near (π, π) (NNN) and an antinodal nesting (ANN) phase.This is true for both magnetic instabilities [22] and CDWs.The NNN-ANN crossover occurs near x ∼ 0.12 − 0.16 for the magnetic phases and at a somewhat lower doping for the CDWs.The question arises, which density wave has lower energy?At half-filling, the answer is clear, and in good agreement with experiment.Due to the large Hubbard U, the CDW is only marginally stable, whereas (π, π) AFM order can open a full gap over the FS, leading to a much greater energy lowering.Thus, near half filling we expect NNN magnetic order for all cuprates.As doping increases, the pseudogap gets much smaller, and most experiments see a weakening of magnetic fluctuations.At the same time, with increasing doping, the suppression of CDW order by U becomes relatively less important.Finally, CDW order tends to be commensurate with the lattice, hence belonging to the Ising universality class, which is more robust against fluctuations than magnetism which always breaks an SO(3) continuous symmetry.Given the uncertainty in γ, it will be hard to determine the exact crossover, but it is likely that the CDW will win out at higher doping.
Recently, there has been considerable experimental evidence for such a transition, between a low-doping SDW phase and a higher-doping ANN CDW phase, both in the Bi-cuprates and in YBCO.Here, we will summarize the experimental evidence for such a crossover, while in the remainder of this paper we will explore the consequences of this assumption.In deeply underdoped Bi2212, hints of this transition have been observed in recent scanning tunneling microscopy (STM) studies: Ref. [48] found that the phase we identify here as the ANN CDW seems to weaken below 1/8th doping, while in a similar doping range in CCOC Ref. [49] found that islands of a phase with very weak CO and without C 4 symmetry breaking become more prevalent with reduced dopingsuggestive of domains of predominantly magnetic order, invisible to STM.
Results for YBCO are even more interesting.A CDW phase has been found [2,3,4,5,6,50] in the doping range near x = 1/8 where quantum oscillations have been observed [51,52,53,54,55,56].In the absence of superconductivity, the CDW correlation length ξ remains finite, growing as T decreases [2], as expected for a 2D system from the Mermin-Wagner theorem [57,58].The results are reminiscent of the growth of (π, π) AFM order in electron-doped cuprates, except that in the latter case a transition to long range magnetic order occurs in underdoped samples, where superconductivity is suppressed [59].In YBCO, when an external magnetic field is used to suppress the superconducting (SC) order, the CDW correlation length is found to increase [50,3,5].
In this connection, we further note the following points: (a) The nesting vectors are very similar as a function of doping for all cuprates studied, Fig. 4, down to details of vertical vs diagonal nesting and 1-q vs 2-q order (see Subsection 4.2); (b) The nesting seems to follow the bonding FS of YBCO [4,60], just as we have found in Bi2212 (Fig. 4); (c) Ref. [6] reports that in the Ortho-II phase in YBCO 6.54 CDW order along the b-axis and SDW order along the a-axis are simultaneously present at zero field.This is similar to the two types of patches observed in extremely underdoped Bi2212 [61], which we ascribed to competing CDW and SDW orders; and, (d) In Bi2212, the CDW order is strongest near 1/8th doping, but persists well into the overdoped regime [48].In ) Nesting vectors for ANN (q, q) CO for Bi2201 (green solid line), Bi2212 with (violet dot-dashed) or without (violet dotted) bilayer splitting, compared to experimental (0, q) CO vectors for Bi2201 (green circles [66] or triangles [68]), Bi2212 (violet circles [67]), YBCO (blue squares [2,3,4,6]), and Hg1201 (silver diamond) [69].Band parameters for the theoretical calculations are from Refs.[65] and [25].Results for Hg1201 and YBCO are after Ref. [69].
YBCO, the CDW effects also peak near x = 1/8.The range has been estimated best from measurements of the Hall coefficient, which in YBCO becomes negative at low temperatures in the doping range ∼ 0.08 ≤ x ≤∼ 0.16, and this is considered to be due to FS reconstruction associated with CO [62].At the same time magnetic fluctuations are found in the region below x = 0.08, associated with the competing near-(π, π) order.[63] We thus adduce that there is considerable evidence for a crossover as a function of doping between a predominantly magnetic phase with incommensurate, near-(π, π) order to an ANN CO phase as doping is increased, both in YBCO and in Bi2201/Bi2212, consistent with our nesting model.One should keep in mind, however, that the current experimental situation remains fluid.

Nesting Vectors
For electron doped cuprates, the dominant nonsuperconducting order remains commensurate at (π, π), associated with the magnetic instability, and there has been little evidence for any secondary CO.The undoped case is special, since the full FS can be gapped and conventional nesting plays no role.Nevertheless, the optimal energy lowering is associated with q = (π, π), [14] which again is the magnetic ordering vector.Similarly, for hole-doped LSCO the primary order is magnetic, with CO forming on antiphase boundaries at a near-(π, π) incommensurate q-vector [1].However, LSCO is a special case, which will be considered elsewhere [64], and here we will limit ourselves to the discussion of CDW modulations in other hole-doped cuprates.
In contrast, in Bi2201 and Bi2212, STM studies are sensitive to charge modulations, and we find that the associated q-vector is consistent with an important role of ANN nesting in CO. Figure 4 plots the calculated doping dependence of the ANN charge nesting vectors for the Bi cuprates, [65,25] displaying a strong material dependence.Shown also are the experimental superlattice periodicities for CO found in the tunneling spectra in Bi2201 [66] (green circles) and Bi2212 [67] (violet circles), and compared to other measures of the CDW q-vector for YBCO (blue squares) [2,3,4,6], Bi2201 (green triangles) [68], and Hg1201 (silver diamond) [69].Clearly the experimental superlattices in the Bi-cuprates are close to the predicted ANN periodicities, although in Bi2201 there are hints of shifts to nearby commensurate values.In Bi2212 there are two nesting vectors, associated with bonding and antibonding combinations of the bilayersplit bands, and the experimental data fall close to the bonding band nesting vector.For both materials, the observed q-vector follows the doping dependence of ANN nesting, and is incompatible with an interpretation as a secondary CO as in LSCO, having the wrong doping dependence.
It should be noted that the experimental q-vectors represent vertical (q, 0) ANN nesting, whereas the theory predicts nesting at (q ′ , q ′ ), where q ′ is close to q in magnitude.We believe that this is a question of the evolution of the CDW ground state as U ef f is increased above threshold, in which case vertical ANN nesting is energetically favored (Appendix B).A similar competition of 1-q vs 3-q nesting has been found in NbSe 2 [70].
A closely related issue for both Bi-cuprates and YBCO is whether the CDWs are modulated along a single q vector, either along the x or the y direction (denoted 1-q nesting), [6,50] or whether there is modulation simultaneously along both the x-and y− axes (2-q nesting).[2,3,4,5] While some experiments cannot distinguish true 2 − q nesting from patches of 1 − q (x or y) CDWs, others can [5].This point is discussed further in Appendix B.

Kohn Anomalies and Soft Phonons
To compare the Kohn anomaly to experiment, it must be kept in mind that our model involves effective acoustic modes of a copper plane, which must be embedded into the full phonon band dispersions.That means the Kohn anomaly will show anticrossing behavior with the different bands.Here we provide one example of this.In Fig. 5, we plot an expanded view of the Kohn anomaly in our model.In particular, the anomaly along the Γ → (π, 0) longitudinal branch in LSCO, dark blue curve in Fig. 2(a) or Fig. 2(d), can be compared with the experimental results of Ref. [71].While the shape and position are qualitatively correct, some differences can be expected due to our oversimplified phonon model.In particular, the experimental Kohn anomaly is in an LO branch, whereas our model only has acoustic branches.Since the bare LO branch has a maximum at Γ and decreases towards (π, 0), this reverses the left-right asymmetry of the Kohn anomaly.To see this more clearly, we have replotted the data of Fig. 2 in Fig. 5(b) with a reversed horizontal axis.A second difference is that the experimental Kohn anomaly softens, but does not go unstable, unlike the theory.This is again an anticrssing effect, and is discussed further in SOM Section III.

'Non-BCS' CDWs
CDWs in cuprates seem anomalous when compared to a BCS-like mean-field picture, which predicts a second order transition with a diverging correlation length, and the ratio of the T = 0 CDW gap to T CDW of η CDW ≡ 2∆ CDW (0)/k B T CDW = 3.53.In practice, however, such a description hardly ever works even for conventional CDWs in the correlation length does not diverge and the gap ratio η CDW is typically >> 3.53.This was originally discussed in 2H-TaSe 2 by McMillan, [72] who suggested that the short correlation length means that the electron couples to phonons with a wide range of q-values, and that the transition is therefore controlled by the phononic and not electronic entropy.He then modeled the transition as a double transition, first to short-range order at a high T SRO consistent with the BCS ratio, and then to long-range order at a much lower temperature.The phonon entropy effect is now understood as a breakdown of the RPA due to mode coupling effects.[31,32] For a 2D system, modecoupling effects account for the Mermin-Wagner physics, suppressing the transition to T = 0. [57,32,58] 2D systems are also sensitive to impurity effects, which can further limit correlation lengths.[73] If the cuprates display a 'conventional' CDW instability, we should expect: (a) the transition is characterized by phonon softening; (b) the phonons will soften over a range of q-values, [74] where the range is related to the electron correlation length-perhaps associated with an order-disorder transition and a central peak [75]; (c) the electronic correlation length will not diverge at the transition; and (d) the gap ratio should be anomalously large, with the magnitude of the anomaly also related to the correlation length.In both Bi2212 and YBCO the CDW appears to have an anomalously small correlation length.[76,4]

Purely Electronic CDWs
We comment on recent papers on purely electronic CDW models [77,78,79,80,81,9,10], showing how our work relates to these papers.One paper [82] noted, "because the Q = 0 modulations exhibit wave vectors generated by scattering regions ('hot spots') moving along the k-space lines (±π, 0) → (0, ±π), FS nesting provides an inadequate explanation for the cuprate density waves."Our analysis, however, indicates otherwise in that diagonal hot spot nesting was predicted from a nesting model and its origin carefully described in Ref. [22].The issues raised in Ref. [22] concerning vertical vs diagonal CDWs and one-or two-q order seem to arise in all the current CDW models.In order to stress that the 'hot-spot' nesting is a FS effect unrelated to the (π, π) magnetic order, a different argument is presented in Appendix B.
All the preceding models are based on an assumed (π, π)-dominated spin susceptibility.Within a Stoner-type framework, this would suggest the presence of susceptibility peaks at (π, π) at all dopings, a quantum critical point x c when Uχ (π,π) (T = 0) ∼ 1, and strong commensurate fluctuations for x > x c .However, this is not found to be the case in most cuprates [22]: as doping increases, there is a crossover to a regime where the 'hotspot' susceptibility is the largest, and the near-(π, π) fluctuations are cut off.Ironically, only in LSCO is the hotspot susceptibility weak and the near-(π, π) fluctuations dominate for all dopings.However, in LSCO the CDW seems to be absent, and the high-doping regime is consistent with a spin-density wave with q-vector given by the incommensurate (π, π − δ) nesting vector.[22] While the above models all describe purely electronic CDWs, the observed CDW couples strongly to phonons.Indeed, the x-ray diffraction intensities of Ref. [4] could only be explained by assuming that the CDW was accompanied by a conventional Peierls distortion, which increased the x-ray intensity by a factor of ∼600.Hence, an improved model should incorporate effects of both the electronic CDW and the accompanying lattice distortion.

Conclusions
We have examined the charge response of the cuprates within the framework of the Hubbard model using the time-dependent Gutzwiller approximation, where effects of a finite electron-phonon coupling are included for the first time.The resulting ANN CDWs provide a good model for the higher-doping regime of the pseudogap phase in most hole doped cuprates.Specifically, the ANN phase in the Bi cuprates captures the experimentally observed doping dependence of incommensurability, and the predicted FS is found to be consistent with that seen in QO measurements.A secondary magnetic order (SOM Section I) enhances the resemblance to a conventional stripe interaction effects are contained in the second order contributions, which due to translational invariance will now be evaluated in momentum space.
For the Hubbard term one finds where various coefficients are defined in Appendix A.2 and the relevant fluctuation modes are the local density fluctuations δρ q = 1 N k,σ δρ k+q,k,σ , the intersite charge fluctuations δT q = 1 N k,σ (ǫ 0 kσ + ǫ 0 k+q,σ )δρ k+q,k,σ , and the double occupancy fluctuations δD q .We also define δρ k+q,k = σ δρ k+q,k,σ with the density matrix ρ kk ′ ,σ = c † kσ c k ′ σ .The fluctuation contribution for the electron-phonon interaction takes the form k,k+q,µ δρ k,k+q + f (2) ′ q,µ δρ −q + h q,µ δD −q ], (A.7) where k,k+q,µ Q µ q is the Fourier transform of f i,j,µ = −t i,j γ i,j (u µ j − u µ i )/r i,j , and f (0) q,µ Q µ q is the Fourier transform of f i,µ = j f i,j,µ .Explicit expressions for f (0) k,k+q,µ and f (0) q,µ are given in Appendix A.3.Finally, the double occupancy fluctuations can be eliminated by an antiadiabatic approximation [29], where one assumes that the fluctuations are faster than other degrees of freedom.Since we will be concerned with the static limit in the present paper it is always justified to take this antiadiabatic limit, ∂E (2)  ∂D −q = 0, (A.11) which allows us to express the double occupancy via the density fluctuations δD q = −(L q δρ q + z 0 z ′ D δT q + µ h qµ Q µ q )/U q .(A.12) Inserting Eq. (A.12) into Eq.(A.5) yields an energy functional δρ q δT q T A q B q B q C q δρ −q δT −q , (A. 14) the FS and q = 2k F .This is the conventional nesting, leading to a contribution to the susceptibility inversely proportional to the local FS curvature, or loosely speaking, to the length L of FS that nests in a 2D case.Let us apply this to a very simple model of AN nesting.Let the antinodal part of the FS be flat over a length L in k-space, with the two flat sections separated by q AN .If one FS is shifted vertically with respect to the other by q AN , the two FSs will nest over a length L, so that χ 0 (q AN , 0) ∼ L. If the FS is shifted diagonally, the surfaces will nest over only L − q AN , but both the x-and y-antinodal regions will be nested (double nesting), yielding χ 0 (q AN , q AN ) ∼ 2(L − q AN ).Hence, for L > 2q AN , diagonal nesting becomes unstable first.On the other hand, if there are CDWs along both x and y, we would nest substantially more FS χ 0 ∼ 2L, and hence 2-q nesting would dominate.In reality, FS sections are almost never exactly parallel, so diagonal nesting wins at threshold.However, when U ef f is larger than the threshold value U c , more of the FS will be gapped, so the above agrument holds for 'nearly-nesting' segments, and 2-q should ultimately win out.
Next, we give an additional argument that diagonal AN nesting is a form of 'hotspot' nesting, but it has its origin in a band structure effect completely unrelated to any underlying (π, π) AF order.[22] Because the FS has a mirror symmetry about the Γ → (π, π) line of the Brillouin zone, the folded FS maps q = 2k F → (2π − q x , q y ) or (q x , 2π − q y ) will overlap only when q x = q y along the zone diagonal.But when unfolded into the doubled zone, the diagonal becomes the lines (2π, 0) → (π, π) and (π, π) → (0, 2π), which are the q = 2k F image of the AF zone boundary.
Finally, we would like to clarify a point of terminology, which is often confusing in the literature.For instance, the 'packed golf ball' motif shown in the inset of Fig. 1(c) of Ref. [66] is characterized in that paper as checkerboard order.However, a checkerboard implies a density modulation (high-low-high-low) along the Cu-O bond direction.This kind of CO would induce peaks in the Fourier transformed spectra along the diagonals of the Brillouin zone, whereas the experimental data [66] clearly have the maxima lying along the reciprocal Cu-O bond direction.The observed pattern is in fact more properly termed 'crossed stripes', having three different charge densities: high in the regions where (charge) stripes cross, low on the sites which are not occupied by the (charge) stripes, and intermediate on other sites where only one stripe is occupied.For crossed 6×6 stripes (SOM Section 1) this is exactly the 2D pattern observed in Ref. [66], which also has the dominant Fourier peaks along the reciprocal Cu-O bond direction.
corresponds to our lower estimate for the experimental λ ep = 0.57; the upper limit, λ ep =5.4, is off the scale of the figure.Incorporation of long-range Coulomb interaction might shift the critical value of λ ep to larger values, but is unlikely to introduce any qualitative change in the phase diagram.In comparison to experimental CDW phases, a value of λ ep near our lower estimate seems likely.

Table 1 .
I. Band Parameter Sets e denotes the Hubbard model, H e−ph the coupling between electrons and phonons and H ph is the bare phonon part.In the Hubbard model